4
$\begingroup$

I try to understand details of the Needleman–Wunsch algorithm and use an example from the book [Nello Cristianini, Matthew W. Hahn. Introduction to Computational Genomics. A Case Studies Approach, www.computational-genomics.net].

To my opinion, the algorithm is described not very neatly, so the question arises.

Consider two sequences:

I: VIVALASVEGAS

II: VIVADAVIS

The first step is clear:

enter image description here

Also, the calculation of the score for diagonal neighbor is also clear,

enter image description here

because this is an alignment without a gap:

I: V

II: V

Unclear moment is the alignment for left (or above) neighbor. Say, the left neighbor

enter image description here

corresponds to this

I: V

II: -

or this

I: -

II: V

alignment?

$\endgroup$
1
  • $\begingroup$ This is really off-topic here, but you need to read very carefully how this algorithm works, how the scoring matrix is constructed, how the traceback is done, and how that is finallly converted into an alignment. It takes me a 50 min lecture to explain it. I do not know the book you mention but have found the best explanation in Durbin et al. Biological Sequence Anlaysis. $\endgroup$
    – David
    Dec 1, 2017 at 17:31

1 Answer 1

7
$\begingroup$

This is a classic application of dynamic programming in bioinformatics, and is one of the first dynamic programs that I was exposed to in university! Let's begin..


Overview of algorithm

Step 1. Conceive of a scoring system.
Step 2. Create a table using the two character sequences.
Step 3. Iteratively fill in each cell, row by row, left to right. This is achieved by selecting the maximum (most positive) value coming from either:

(1) the cell that's located diagonally (top-left) to the current cell, plus the "cost" to align the two characters in consideration (cost value is obtained from scoring system).

(2) the cell to left of the current cell, plus the gap penalty (score issued for insertion/deletion mutation); looking at your example, the gap penalty is -1.

(3) the cell just above the current cell, plus the gap penalty.

Again, the greatest (most positive) of these three [sum] values will be placed in the current cell.

..

If you want the recurrence relation for this, then, if $T$ is the table being filling in, and $T(i, \ j) \equiv T(row, \ column)$, then

$$T(i, \ j) \ = maximum \begin{pmatrix} T (i-1, \ j-1) + [alignment \ cost],\\ T (i, \ j-1) \ + [gap \ penalty],\\ T (i-1, \ j) + [gap \ penalty] \end{pmatrix}$$


Scoring system

For simplicity sake, and to remain true to your example, I'll define and use the following:

  • Character match score: 0
  • Character mismatch score: -1
  • Gap penalty: -1

The character scores are potential alignment costs, and which to use is determined by if the two characters in consideration (row vs. column) are the same or not.


Now, to fill in table $T$..

Let sequence 1 run along the rows of $T$, and sequence 2 run along the columns, as you already have. If $T$ is initialized with a value of zero at location (0, 0), then, in order to fill in the first row [going from cell (0, 1) to (0, 9)], since there are no diagonal cells that exist yet, nor is there yet a row above the current row, for all cells in the first row we will use the second option that was mentioned above to obtain each cell value. So then, with the gap penalty being defined as -1, the first row will range from 0 to -9.

enter image description here

Make sure to keep track of which adjacent cell was used to obtain the value for each cell. This descriptive will later be used to trace $T$ once it's completely filled out, in order to obtain the optimal solution. If you're writing this on paper, people usually just draw arrows (from current to adjacent cell), however, if it's on a computer, an adjacency matrix is created & filled in.

In cases where more than one adjacent cell yields the greatest value, draw an arrow to each one. This is to say that there can be more than one optimal solution.


Next, we start with the first cell of the second row. Since a left & diagonal cell doesn't exist, that omits options 1 & 2, and leaves only option 3, which yields a value of -1 (only gap penalty). And, as you maybe have noticed by now, this behavior will always be true for the first row, and for the first column, and so, this is why you're able to so immediately fill out the table as you already have.

Continuing, to compute cell (1, 1), we now have to compare all three potential values and choose the maximum. For options 2 & 3, they both yield a value of $-1 + (-1) = -2$, however, option 1 yields $0 + [alignment \ cost]$, and, since both characters are the same (V and V), cell (1, 1) will have a value of $0 + (0) = 0$.

If we're to then fast forward a little and calculate cell (2, 1) ...

enter image description here

Considering option 1, we have $-1 + (-1) = -2$, since the characters being compared are I and V (and are a mismatch). Option 2 yields a value of $-2 + (-1) = -3$, and option 3 gives $0 + (-1) = -1$, hence, cell (2, 1) = -1.

Continue doing this for all cells until $T$ is filled.


Obtaining an optimal solution..

As mentioned, multiple optimal solutions may exist. In order to obtain an optimal solution, you must trace $T$, using the arrows you've been drawing (or the adjacency matrix), as a guide.

To obtain the path for best alignment, simply begin with the bottom-right most cell of $T$, i.e. the most recent/last computed cell, and trace your arrows to create a path back to the first cell (the top-left most).

To then construct the alignment, you will traverse this alignment path, beginning with the top-left most cell, and going to the bottom-right most cell, and, when considering the orientation of each arrow when going from cell to cell (diagonal, horizontal, vertical), you will do one of the following:

Diagonal arrow

For the first alignment sequence (top), simply add the character that's in the column of the cell you just moved to; for the second alignment sequence (bottom), use the character that's in that row. So then, when considering your example, and for both sequences, the first characters of an optimal solution when going from cell (0, 0) = 0 to cell (1, 1) would be:

V

V

Vertical arrow

When this happens, it means that there is a gap in the [initial] sequence that runs along the columns, since there are two rows for one column. Because of this, for the first alignment sequence, you will use a - to denote a gap. For alignment sequence two, use the character that's in the row.

Within the context of your sequences, for the optimal solution, a gap isn't encountered until the 7th character, at which point, the partial optimal solution will look like:

V I V A D A -

V I V A L A S

Horizontal arrow

This is the same concept as the vertical arrow rule, except it implies there is a gap in the row-oriented sequence, and so, a - will be placed in the second alignment sequence; the first alignment sequence will use the character in that column.


After all is said and done, your optimal solution is:

V I V A D A - V I - - S

V I V A L A S V E G A S

enter image description here

An interactive version of this algorithm can be found here.

$\endgroup$
2
  • $\begingroup$ Thank you very much for so detailed description. But I wished a very little detail of the algorithm. Now, due to the interactive version, all becomes clear. $\endgroup$
    – Konstantin
    Dec 4, 2017 at 9:44
  • $\begingroup$ @Konstantin I provided the entire algorithm because you mentioned the book not being clear to you in instruction, and, because the determination of alignment follows a specific set of rules, just as filling out the table does. I'm glad you were able to find your answer. :) Please mark this as the correct answer if this was helpful and accurate. $\endgroup$
    – user22020
    Dec 4, 2017 at 12:48

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.