C# ISequence.GetSubSequence方法代码示例

        /// <summary>
        /// Removes gaps that are inserted by overlap algorithm at beginning or end of sequence.
        /// </summary>
        /// <param name="inputSequence">input sequence</param>
        /// <returns>Sequence without gaps</returns>
        private static ISequence SequenceWithoutTerminalGaps(ISequence inputSequence)
        {
            //string input = inputSequence.ToString();
            long start = 0;
            while (inputSequence[start] == '-')
            {
                ++start;
            }

            long len = inputSequence.Count - start;
            while (inputSequence[len - 1] == '-')
            {
                --len;
            }

            ISequence seq = inputSequence.GetSubSequence(start, len);
            seq.ID = inputSequence.ID;
            return seq;
        }

        /// <summary>
        /// Returns the sequence for the specified start and end positions.
        /// If the accession is null or empty then the source sequence is used to construct the output sequence,
        /// otherwise appropriate sequence from the referred sequence is used to construct output sequence.
        /// </summary>
        /// <param name="start">Start position.</param>
        /// <param name="end">End position.</param>
        /// <param name="accession">Accession number.</param>
        /// <param name="sepataror">Start and End separator.</param>
        /// <param name="source">Source sequence.</param>
        /// <param name="referredSequences">Referred Sequences.</param>
        private static ISequence GetSubSequence(int start, int end, string accession, string sepataror, ISequence source, Dictionary<string, ISequence> referredSequences)
        {
            if (string.Compare(sepataror, "^", StringComparison.OrdinalIgnoreCase) == 0)
            {
                return new Sequence(source.Alphabet, string.Empty);
            }

            if (string.Compare(sepataror, "..", StringComparison.OrdinalIgnoreCase) != 0 &&
                string.Compare(sepataror, ".", StringComparison.OrdinalIgnoreCase) != 0 &&
                !string.IsNullOrEmpty(sepataror))
            {
                string str = string.Format(CultureInfo.CurrentCulture, Properties.Resource.InvalidSeparator, sepataror);
                throw new ArgumentException(str);
            }

            if (!string.IsNullOrEmpty(accession) && (referredSequences == null || !referredSequences.ContainsKey(accession)))
            {
                string str = string.Format(CultureInfo.CurrentCulture, Properties.Resource.AccessionSequenceNotFound, accession);
                throw new ArgumentException(str);
            }

            if (!string.IsNullOrEmpty(accession))
            {
                if (source.Alphabet != referredSequences[accession].Alphabet)
                {
                    string str = string.Format(CultureInfo.CurrentCulture, Properties.Resource.InvalidReferredAlphabet, accession);
                    throw new ArgumentException(str);
                }

                source = referredSequences[accession];
            }

            // as location.start is one based where as Range accepts zero based index.
            start = start - 1;
            int length = end - start;

            if (string.IsNullOrEmpty(sepataror) || string.Compare(sepataror, ".", StringComparison.OrdinalIgnoreCase) == 0)
            {
                length = 1;
            }

            ISequence newSequence = source.GetSubSequence(start, length);
            byte[] seqData = new byte[newSequence.Count];
            for (long i = 0; i < newSequence.Count; i++)
            {
                seqData[i] = newSequence[i];
            }

            return new Sequence(source.Alphabet, seqData);
        }

        /// <summary>
        /// For input sequence, constructs k-mers by sliding 
        /// a frame of size kmerLength along the input sequence.
        /// Track positions of occurance for each kmer in sequence.
        /// Constructs KmersOfSequence for sequence and associated k-mers.
        /// </summary>
        /// <param name="sequence">Input sequence.</param>
        /// <param name="kmerLength">K-mer length.</param>
        /// <returns>KmersOfSequence constructed from sequence and associated k-mers.</returns>
        public static KmerPositionDictionary BuildKmerDictionary(ISequence sequence, int kmerLength)
        {
            if (sequence == null)
            {
                throw new ArgumentNullException("sequence");
            }

            if (kmerLength > sequence.Count)
            {
                throw new ArgumentException(Properties.Resource.KmerLengthIsTooLong);
            }

            // kmers maintains the map between k-mer strings to list of positions in sequence.
            KmerPositionDictionary kmers = new KmerPositionDictionary();

            // Sequence 'kmer' stores the k-mer in each window.
            // Construct each k-mer using range from sequence.
            for (long i = 0; i <= sequence.Count - kmerLength; ++i)
            {
                ISequence kmerString = sequence.GetSubSequence(i, kmerLength);
                if (kmers.ContainsKey(kmerString))
                {
                    kmers[kmerString].Add(i);
                }
                else
                {
                    kmers[kmerString] = new List<long>() { i };
                }
            }

            return kmers;
        }

        /// <summary>
        /// Align the Gap by executing pairwise alignment.
        /// </summary>
        /// <param name="referenceSequence">Reference sequence.</param>
        /// <param name="querySequence">Query Sequence.</param>
        /// <param name="sequenceResult1">Editable sequence containing alignment first result.</param>
        /// <param name="sequenceResult2">Editable sequence containing alignment second result.</param>
        /// <param name="consensusResult">Editable sequence containing consensus sequence.</param>
        /// <param name="mum1">First MUM of Gap.</param>
        /// <param name="mum2">Second MUM of Gap.</param>
        /// <param name="insertions">Insertions made to the aligned sequences.</param>
        /// <returns>Score of alignment.</returns>
        private long AlignGap(
                ISequence referenceSequence,
                ISequence querySequence,
                List<byte> sequenceResult1,
                List<byte> sequenceResult2,
                List<byte> consensusResult,
                Match mum1,
                Match mum2,
                out List<long> insertions)
        {
            long score = 0;
            ISequence sequence1 = null;
            ISequence sequence2 = null;
            IList<IPairwiseSequenceAlignment> sequenceAlignment = null;
            byte[] mum1String;
            byte[] mum2String;

            insertions = new List<long>(2);
            insertions.Add(0);
            insertions.Add(0);

            long mum1ReferenceStartIndex = 0;
            long mum1QueryStartIndex = 0;
            long mum1Length = 0;
            long mum2ReferenceStartIndex = 0;
            long mum2QueryStartIndex = 0;
            long mum2Length = 0;

            if (mum1.Length != 0)
            {
                mum1ReferenceStartIndex = mum1.ReferenceSequenceOffset;
                mum1QueryStartIndex = mum1.QuerySequenceOffset;
                mum1Length = mum1.Length;
            }

            if (mum2.Length != 0)
            {
                mum2ReferenceStartIndex = mum2.ReferenceSequenceOffset;
                mum2QueryStartIndex = mum2.QuerySequenceOffset;
                mum2Length = mum2.Length;
            }
            else
            {
                mum2ReferenceStartIndex = referenceSequence.Count;
                mum2QueryStartIndex = querySequence.Count;
            }

            long referenceGapStartIndex = mum1ReferenceStartIndex + mum1Length;
            long queryGapStartIndex = mum1QueryStartIndex + mum1Length;

            if (mum2ReferenceStartIndex > referenceGapStartIndex
                && mum2QueryStartIndex > queryGapStartIndex)
            {
                sequence1 = referenceSequence.GetSubSequence(
                    referenceGapStartIndex,
                    mum2ReferenceStartIndex - referenceGapStartIndex);
                sequence2 = querySequence.GetSubSequence(
                    queryGapStartIndex,
                    mum2QueryStartIndex - queryGapStartIndex);

                sequenceAlignment = this.RunPairWise(sequence1, sequence2);

                if (sequenceAlignment != null)
                {
                    foreach (IPairwiseSequenceAlignment pairwiseAlignment in sequenceAlignment)
                    {
                        foreach (PairwiseAlignedSequence alignment in pairwiseAlignment.PairwiseAlignedSequences)
                        {
                            sequenceResult1.InsertRange(
                                    sequenceResult1.Count,
                                    alignment.FirstSequence);
                            sequenceResult2.InsertRange(
                                    sequenceResult2.Count,
                                    alignment.SecondSequence);
                            consensusResult.InsertRange(
                                consensusResult.Count,
                                    alignment.Consensus);

                            score += alignment.Score;

                            if (alignment.Metadata.ContainsKey("Insertions"))
                            {
                                List<int> gapinsertions = alignment.Metadata["Insertions"] as List<int>;
                                if (gapinsertions != null)
                                {
                                    if (gapinsertions.Count > 0)
                                    {
                                        insertions[0] += gapinsertions[0];
//.........这里部分代码省略.........

        /// <summary>
        /// For input sequence, constructs k-mers by sliding 
        /// a frame of size kmerLength along the input sequence.
        /// Track positions of occurrence for each kmer in sequence.
        /// Constructs KmersOfSequence for sequence and associated k-mers.
        /// </summary>
        /// <param name="sequence">Input sequence.</param>
        /// <param name="kmerLength">K-mer length.</param>
        /// <returns>KmersOfSequence constructed from sequence and associated k-mers.</returns>
        public KmersOfSequence Build(ISequence sequence, int kmerLength)
        {
            if (sequence == null)
            {
                throw new ArgumentNullException(Properties.Resource.SequenceCannotBeNull);
            }

            if (kmerLength <= 0)
            {
                throw new ArgumentException(Properties.Resource.KmerLengthShouldBePositive);
            }

            if (kmerLength > sequence.Count)
            {
                throw new ArgumentException(Properties.Resource.KmerLengthIsTooLong);
            }

            // kmers maintains the map between k-mer strings to list of positions in sequence.
            Dictionary<string, List<long>> kmers = new Dictionary<string, List<long>>();

            // Sequence 'kmer' stores the k-mer in each window.
            // Construct each k-mer using range from sequence.
            for (long i = 0; i <= sequence.Count - kmerLength; ++i)
            {
                string kmerString = new string(sequence.GetSubSequence(i, kmerLength).Select(a => (char)a).ToArray());
                if (kmers.ContainsKey(kmerString))
                {
                    kmers[kmerString].Add(i);
                }
                else
                {
                    kmers[kmerString] = new List<long>() { i };
                }
            }

            return new KmersOfSequence(
                sequence,
                kmerLength,
                new HashSet<KmersOfSequence.KmerPositions>(kmers.Values.Select(l => new KmersOfSequence.KmerPositions(l))));
        }

        /// <summary>
        /// Gets the set of kmer strings that occur in given sequences.
        /// </summary>
        /// <param name="sequence">Source Sequence.</param>
        /// <param name="kmerLength">Kmer Length.</param>
        /// <returns>Set of kmer strings.</returns>
        public static IEnumerable<ISequence> GetKmerSequences(ISequence sequence, int kmerLength)
        {
            if (sequence == null)
            {
                throw new ArgumentNullException("sequence");
            }

            if (kmerLength <= 0)
            {
                throw new ArgumentException(Properties.Resource.KmerLengthShouldBePositive);
            }

            if (kmerLength > sequence.Count)
            {
                throw new ArgumentException(Properties.Resource.KmerLengthIsTooLong);
            }

            IList<ISequence> kmers = new List<ISequence>();
            for (int i = 0; i <= sequence.Count - kmerLength; ++i)
            {
                ISequence kmer = sequence.GetSubSequence(i, kmerLength);
                kmers.Add(kmer);
            }

            return kmers;
        }

        /// <summary>
        /// Returns true if the given mismatch is an inversion, that is, whether the query subsequence matches the complement of the reference subsequence.
        /// </summary>
        /// <param name="mismatch">The Mismatch to test</param>
        /// <param name="querySequence">The sequence that is being compared to the reference</param>
        /// <returns>True if the complement of the subsequence in the reference is equal to the subsequence in the given query.</returns>
        private bool isInversion(Mismatch mismatch, ISequence querySequence)
        {
            var complement = querySequence.GetSubSequence(mismatch.QuerySequenceOffset, mismatch.QuerySequenceLength)
                                .GetComplementedSequence();

            return ConvertToString(complement, 0, complement.Count).Equals(
                ConvertToString(ReferenceSequence, mismatch.ReferenceSequenceOffset, mismatch.ReferenceSequenceLength)
            );
        }

        /// <summary>
        /// Align the Gap by executing pairwise alignment.
        /// </summary>
        /// <param name="referenceSequence">Reference sequence.</param>
        /// <param name="querySequence">Query Sequence.</param>
        /// <param name="sequenceResult1">Editable sequence containing alignment first result.</param>
        /// <param name="sequenceResult2">Editable sequence containing alignment second result.</param>
        /// <param name="consensusResult">Editable sequence containing consensus sequence.</param>
        /// <param name="mum1">First MUM of Gap.</param>
        /// <param name="mum2">Second MUM of Gap.</param>
        /// <param name="insertions">Insertions made to the aligned sequences.</param>
        /// <returns>Score of alignment.</returns>
        private long AlignGap(
                ISequence referenceSequence,
                ISequence querySequence,
                List<byte> sequenceResult1,
                List<byte> sequenceResult2,
                List<byte> consensusResult,
                Match mum1,
                Match mum2,
                out List<long> insertions)
        {
            long score = 0;
            ISequence sequence1 = null;
            ISequence sequence2 = null;
            byte[] mum1String;
            byte[] mum2String;

            insertions = new List<long>(2);
            insertions.Add(0);
            insertions.Add(0);

            long mum1ReferenceStartIndex = 0;
            long mum1QueryStartIndex = 0;
            long mum1Length = 0;
            long mum2ReferenceStartIndex = 0;
            long mum2QueryStartIndex = 0;
            long mum2Length = 0;

            if (mum1.Length != 0)
            {
                mum1ReferenceStartIndex = mum1.ReferenceSequenceOffset;
                mum1QueryStartIndex = mum1.QuerySequenceOffset;
                mum1Length = mum1.Length;
            }

            if (mum2.Length != 0)
            {
                mum2ReferenceStartIndex = mum2.ReferenceSequenceOffset;
                mum2QueryStartIndex = mum2.QuerySequenceOffset;
                mum2Length = mum2.Length;
            }
            else
            {
                mum2ReferenceStartIndex = referenceSequence.Count;
                mum2QueryStartIndex = querySequence.Count;
            }

            long referenceGapStartIndex = mum1ReferenceStartIndex + mum1Length;
            long queryGapStartIndex = mum1QueryStartIndex + mum1Length;

            /* Stich the exact matches together according to if both sequences have data
             * in the gap (in which case use a global alignment) or if only one does 
             * (in which case just insert gaps).
             */
            if (mum2ReferenceStartIndex > referenceGapStartIndex
                && mum2QueryStartIndex > queryGapStartIndex) // Both sequences have data in the gap.
            {                
                // Get the sequences in between
                sequence1 = referenceSequence.GetSubSequence(
                    referenceGapStartIndex,
                    mum2ReferenceStartIndex - referenceGapStartIndex);
                sequence2 = querySequence.GetSubSequence(
                    queryGapStartIndex,
                    mum2QueryStartIndex - queryGapStartIndex);
                
                // Do a pairwise alignment (must be needleman wunsh)
                var alignment = this.RunPairWiseReturnJustAlignment(sequence1, sequence2);
                sequenceResult1.AddRange(alignment.FirstSequence);
                sequenceResult2.AddRange(alignment.SecondSequence);
                consensusResult.AddRange(alignment.Consensus);

                score += alignment.Score;

                if (!alignment.Metadata.ContainsKey ("Insertions")) {
                    // Should never happen - can remove later.
                    throw new Exception ("NeedlemanWunsch alignment did not have an insertion entry");
                }
                List<long> gapinsertions = alignment.Metadata ["Insertions"] as List<long>;
                if (gapinsertions == null || gapinsertions.Count != 2) {
                    // Should never happen - can remove later
                    throw new Exception("Alignment Insertions were not available as a size 2 list");
                }
                insertions [0] += gapinsertions [0];
                insertions [1] += gapinsertions [1];               
            }
            else if (mum2ReferenceStartIndex > referenceGapStartIndex) // Only the reference has data, insert gaps for the query
            {
                sequence1 = referenceSequence.GetSubSequence(
                    referenceGapStartIndex,
//.........这里部分代码省略.........

        /// <summary>
        ///     Writes the sequence to the specified writer.
        /// </summary>
        /// <param name="sequence"></param>
        /// <param name="type"></param>
        /// <param name="writer"></param>
        private void WriteSeqData(ISequence sequence, string type, TextWriter writer)
        {
            if (sequence.Count > 0)
            {
                byte[] TempSeqData = null;
                type = this.GetGenericTypeString(sequence.Alphabet);

                this.WriteHeaderLine(writer, type, sequence.ID);

                for (long lineStart = 0; lineStart < sequence.Count; lineStart += MaxSequenceSymbolsPerLine)
                {
                    long length = Math.Min(MaxSequenceSymbolsPerLine, sequence.Count - lineStart);
                    ISequence subSequence = sequence.GetSubSequence(lineStart, length);
                    TempSeqData = new byte[length];
                    for (int i = 0; i < length; i++)
                    {
                        TempSeqData[i] = subSequence[i];
                    }
                    string key = Encoding.UTF8.GetString(TempSeqData, 0, TempSeqData.Length);

                    this.WriteHeaderLine(writer, key);
                }

                this.WriteHeaderLine(writer, "end-" + type);
            }
        }