4
\$\begingroup\$

I split my implementation of this sequence alignment algorithm in three methods. Where NeedlemanWunsch-method makes use of the ScoringFunction and the Traceback methods. Further I decided to go with two matices. One matrix is for the scoring, the second contains data to make the traceback easier.

public static int ScoreFunction(char a,char b,int matchScore,int mismatchScore) 
            {
                if (a == b)
                {
                    return matchScore;
                }
                else 
                {
                    return mismatchScore;
                }
            }

public static string[] NeedlemanWunsch(string sequenceA, string sequenceB, int gapPenalty, int matchScore, int mismatchScore)
                {

                #region Initialize
                int[,] matrix = new int[sequenceA.Length + 1, sequenceB.Length + 1];
                char[,] tracebackMatrix = new char[sequenceA.Length + 1, sequenceB.Length + 1];
                matrix[0, 0] = 0;

                for (int i = 1; i < sequenceA.Length + 1;i++)
                {
                    matrix[i,0] = matrix[i-1,0] + gapPenalty;
                    tracebackMatrix[i,0] = 'L';
                }

                for (int i = 1; i < sequenceB.Length + 1; i++)
                {
                    matrix[0, i] = matrix[0 , i - 1] + gapPenalty;
                    tracebackMatrix[0,i] = 'U';
                }
                #endregion

                #region Scoring
                for (int i = 1; i < sequenceA.Length + 1;i++)
                {
                    for (int j = 1; j < sequenceB.Length + 1;j++)
                    {
                        int diagonal = matrix[i - 1, j - 1] + ScoreFunction(sequenceA[i-1],sequenceB[j-1],matchScore,mismatchScore);
                        int links = matrix[i - 1, j] + gapPenalty;
                        int oben = matrix[i, j - 1] + gapPenalty;
                        matrix[i, j] = Math.Max(oben,Math.Max(links, diagonal));                           

                        if(matrix[i,j] == diagonal)
                        {
                            tracebackMatrix[i, j] = 'D';
                        }
                        else if (matrix[i, j] == links)
                        {
                            tracebackMatrix[i, j] = 'L';
                        }
                        else if (matrix[i, j] == oben)
                        {
                            tracebackMatrix[i, j] = 'U';
                        }
                    }
                }
                #endregion

                #region Traceback
                return TraceBack(tracebackMatrix,sequenceA,sequenceB);
                #endregion

            }

public static string[] TraceBack(char[,] tracebackMatrix,string sequenzA, string sequenzB) 
        {

            int i = tracebackMatrix.GetLength(0) - 1;
            int j = tracebackMatrix.GetLength(1) - 1;

            string alignedSeqA = "";
            string alignedSeqB = "";


            while(i != 0 && j != 0)
            {
                switch (tracebackMatrix[i, j])
                {
                    case 'D':
                        alignedSeqA += sequenzA[i - 1];
                        alignedSeqB += sequenzB[j - 1];                        
                        i--;
                        j--;
                        break;
                    case 'U':                        
                        alignedSeqA += "-";
                        alignedSeqB += sequenzB[j - 1];
                        j--;
                        break;
                    case 'L':
                        alignedSeqA += sequenzA[i - 1];
                        alignedSeqB += "-";
                        i--;
                        break;

                }
            }

        string[] alignments = new string[2];

        alignedSeqA = new string(alignedSeqA.Reverse().ToArray());
        alignedSeqB = new string(alignedSeqB.Reverse().ToArray());
        alignments[0] = alignedSeqA;
        alignments[1] = alignedSeqB;

        return alignments;

    }

I am especially concerned about the traceback method and whether it is a good idea to use two matrices.

\$\endgroup\$
5
\$\begingroup\$

Some general notes:

  • The TraceBack method does string concatenation in a tight loop. This is a big f***ing nono. Use a StringBuilder instead, and you will notice a considerable speedup. You have to understand that .net Strings are immutable. This means that a statement like strA += strB does not update strA, it creates a third string and updates the reference to strA. This means concatenation in a tight loop is very expansive.

  • If you do not have a compelling reason to use 2d Arrays (var bla = new int[5,5]), use jagged arrays (var bla = new int[5][5]). Most of the time, the .net JIT compiler can drop the bounds checks for jagged arrays, but not for 2D arrays.

  • It may be worthwhile to convert tracebackMatrix from chars (16 bits) to bytes (8 bits)

\$\endgroup\$
  • \$\begingroup\$ All three points make sense to me. But I don't really understand how the StringBuilder would change the method. That wouldn't make the while-loop unnecessary, would it? Or should I just replace the concatenation "" += "" by a call to a StringBuilder method? \$\endgroup\$ – チーズパン Feb 15 '14 at 15:15
  • 2
    \$\begingroup\$ I've edited the description but yes, you would replace strA += strB with strBuilderA.Append(strB) and at the end you would add a single statement like result = strBuilderA.ToString() \$\endgroup\$ – thomasch Feb 16 '14 at 11:25
3
\$\begingroup\$

First of all, your method "TraceBack" does not return anything (what should it return, both strings?). You can also simplify your first method to use this instead:

public static int ScoreFunction(char a, char b, int matchScore, int mismatchScore)
{
    return a == b ? matchScore : mismatchScore;
}

Rest does not look bad to me, however, I am not that familiar with the algorithm...

\$\endgroup\$
  • \$\begingroup\$ I somehow missed to copy/paste the return statement into the function. So edit is done. \$\endgroup\$ – チーズパン Feb 15 '14 at 12:12

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.