Needleman-Wunsch-Algorithmus

Der Needleman-Wunsch Algorithmus ist eine Methode, globale Alignments von zwei Sequenzen a und b zu berechnen. Der Algorithmus bedient sich dazu der dynamischen Programmierung. Gaps werden zugelassen. (Siehe Sequenzalignment)

Das Verfahren

Für den Needleman-Wunsch-Algorithmus benötigt man eine Matrix $D$ der Größe $(n+1)*(m+1)$, wobei n die Sequenzlänge von Sequenz a und m die Sequenzlänge von b ist. Die Matrix ist genau um eine Spalte und Zeile größer als die Sequenzen. Damit werden Gaps am Anfang der Sequenz ermöglicht.

Jeder Matrixeintrag $D(i,j)$ wird so interpretiert, dass er die beste Bewertung eines Teilalignments von $a_{1},a_{2},\dots ,a_{i}$ und $b_{1},b_{2},\dots ,b_{j}$ darstellt. Für $i=n$ und $j=m$ erhält man also die Bewertung des Alignments der gesamten Sequenzen. Diese steht an der Position $D(n,m)$.

Außerdem benötigt man eine Funktion $w(x,y)$, welche eine Bewertung der eingegebenen Buchstaben x und y berechnet oder nachschlägt. Die einfachste Bewertungsfunktion, die man sich vorstellen kann, gibt bei Gleichheit 1 zurück, ansonsten wird 0 zurückgegeben. Dabei werden also Substitutionen, Insertionen oder Deletionen geduldet und nicht bestraft.

Die Bewertungen ergeben sich aus dem Maximum folgender drei Fälle, wobei wir davon ausgehen, dass schon ein Teilalignment bis zur Position i und j vorliegt und eine Position weitergegangen wird um zu entscheiden wie die nächsten beiden Buchstaben aligniert werden. Alle Fälle werden in der Reihenfolge aufgeführt, wie sie in der Rekursionsgleichung auftreten

1. Fall: Man aligniert $a_{i}$ mit $b_{j}$. Dies entspricht in der Matrix einem Schritt diagonal nach unten rechts. Die Bewertung für $D(i,j)$ setzt sich nun aus dem Wert des vorhergehenden Eintrags $D(i-1,j-1)$ und dem Ergebnis von $w(a_{i},b_{j})$ zusammen.
2. Fall: Man aligniert $a_{i}$ mit einem Gap. Dies entspricht in der Matrix einem Schritt nach unten. Die Bewertung für $D(i,j)$ setzt sich nun aus dem Wert des Eintrages $D(i-1,j)$ und $w(a_{i},-)$ zusammen.
3. Fall: Man aligniert $b_{j}$ mit einem Gap. Dies entspricht in der Matrix einem Schritt nach rechts. Die Bewertung für $D(i,j)$ setzt sich nun aus dem Wert des Eintrages $D(i,j-1)$ und $w(-,b_{j})$ zusammen.

Der Eintrag in $D(i,j)$ wird aus dem Maximum dieser drei Fälle ermittelt. Daraus ergibt sich folgende Rekursionsgleichung.

$D(i,j)=\max {\begin{cases}D(i-1,j-1)&+\ w(a_{i},b_{j})\\D(i-1,j)&+\ w(a_{i},-)\\D(i,j-1)&+\ w(-,b_{j})\\\end{cases}}$

Der Algorithmus kann sowohl auf der Basis von Ähnlichkeit als auch von Distanz angewendet werden. Es muss die Bewertungsfunktion entsprechend gewählt werden (siehe unten) und in der Rekursion entsprechend das Maximum (Ähnlichkeit) oder das Minimum (Distanz) ermitteln werden.

Initialisierung

$D(0,0)=0$

Initialisierung der zusätzlichen (ersten) Spalte und (ersten) Zeile:

$D(i,0)=D(i-1,0)+w(a_{i},-);\ i=1,2,\dots ,n-1$

$D(0,j)=D(0,j-1)+w(-,b_{j});\ j=1,2,\dots ,m-1$

Nun kann zeilenweise ab $i=1$ und $j=1$ die Rekursion angewandt werden.

Beispiel

Anhand eines kleinen Beispiels werden hier die Schritte des Algorithmus' vorgestellt.

Als Bewertungsfunktion wird die folgende Funktion benutzt:

$w(x,y)={\begin{cases}1&,x=y\\-1&,{\textrm {sonst}}\\\end{cases}}$

a = ACGTC und b = AGTC

Zum besseren Verständnis kann man sich vorstellen, dass die Zeilen mit den Buchstaben aus Sequenz a gelabelt sind und die Spalten mit den Buchstaben aus Sequenz b. Mathematisch gesehen macht dies innerhalb der Matrix keinen Sinn, deshalb ist dies hier nur zur Anschauung.

$D={\begin{pmatrix}&-&A&G&T&C\\-&0&-1&-2&-3&-4\\A&-1&0&0&0&0\\C&-2&0&0&0&0\\G&-3&0&0&0&0\\T&-4&0&0&0&0\\C&-5&0&0&0&0\\\end{pmatrix}}$

0. Schritt: Initialisierung

Die Einträge der Matrix $D(i,j)$ für die erste Zeile und die erste Spalte wird wie oben beschrieben gefüllt. Die Bewertung für den Eintrag $D(1,0)$ wird berechnet aus der darüberliegenden Bewertung $D(i-1,j)=D(0,0)=0$ und dem Score an der Stelle $w(a_{i},b_{i})=w(a_{1},-)=w(A,-)=-1$. Also $D(1,0)=0+(-1)=-1$ die anderen Werte werden nun analog berechnet.

$D={\begin{pmatrix}0&-1&-2&-3&-4\\-1&0&0&0&0\\-2&0&0&0&0\\-3&0&0&0&0\\-4&0&0&0&0\\-5&0&0&0&0\\\end{pmatrix}}$

1. Schritt: Berechnung von $D(1,1)$: ${\begin{cases}D(0,0)+w(A,A)\Rightarrow 0+1\\D(0,1)+w(A,-)\Rightarrow -1+(-1)\\D(1,0)+w(-,A)\Rightarrow -1+(-1)\\\end{cases}}$

→ Das Maximum entsteht aus dem ersten Fall, d.h A wird mit A aligniert.

$D={\begin{pmatrix}0&-1&-2&-3&-4\\-1&1&0&0&0\\-2&0&0&0&0\\-3&0&0&0&0\\-4&0&0&0&0\\-5&0&0&0&0\\\end{pmatrix}}$

→ Erhöhung von j um 1, i bleibt gleich

2. Schritt: Berechnung von $D(1,2)$: ${\begin{cases}D(0,1)+w(A,G)\Rightarrow -1+(-1)\\D(0,2)+w(A,-)\Rightarrow -2+(-1)\\D(1,1)+w(-,G)\Rightarrow 1+(-1)\\\end{cases}}$

→ Das Maximum entsteht aus dem dritten Fall, da hier das Maximum der Berechnung, nämlich 0 entsteht, d.h ein Gap(-) würde mit G aligniert.

$D={\begin{pmatrix}0&-1&-2&-3&-4\\-1&1&0&0&0\\-2&0&0&0&0\\-3&0&0&0&0\\-4&0&0&0&0\\-5&0&0&0&0\\\end{pmatrix}}$

$\vdots$

Die gefüllte Matrix sieht nach vollständiger Ausführung der o.a. Schritte folgendermaßen aus:

$D={\begin{pmatrix}0&-1&-2&-3&-4\\-1&1&0&-1&-2\\-2&0&0&-1&0\\-3&-1&1&0&-1\\-4&-2&0&2&1\\-5&-3&-1&1&3\\\end{pmatrix}}$

Die Bewertung dieses Alignments ist 3.

Das dazugehörige Alignment sieht so aus:

${\begin{matrix}{\textrm {Sequenz}}\ a:&A&C&G&T&C\\{\textrm {Sequenz}}\ b:&A&-&G&T&C\end{matrix}}$

Berechnet wird es durch ein Traceback.

Traceback

Das Traceback wird benötigt um aus der Matrix D das Alignment zu generieren. Man startet dazu in $D(n,m)$ und berechnet rückwärts bis zum Eintrag $D(0,0)$ woher die Bewertung, die in der aktuellen Zelle steht gekommen ist. Je nach Richtung bzw. Fall baut man das Alignment von hinten entsprechend auf, d.h man aligniert zwei Buchstaben oder fügt ein gap ein.

Eine andere Möglichkeit das Traceback durchzuführen besteht darin, sich während der Berechnung der Matrix in einer zweiten Matrix T, der Dimension $n*m$, den Fall zu merken aus dem die entsprechende Bewertung entstanden ist. Für das Traceback wird jetzt nur noch der Pfad, beginnend in $D(n,m)$, in T zurückverfolgt und parallel "notiert" welche Fälle auftreten und damit wie das Alignment aussieht. Möchte man nicht nur ein optimales Alignment, sondern alle optimalen Alignments bekommen (es kann durchaus mehrere Alignments mit der gleichen Bewertung geben) muss man alle möglichen Wege in T rekursiv abarbeiten und die Alignments speichern.

Wahl der Bewertungsfunktion

Die Wahl der Bewertungsfunktion hat einen großen Einfluss auf die Ergebnisse, die man durch den Needleman-Wunsch-Algorithmus erhält. Eine einfach Bewertungsfunktion wie oben gewählt spiegelt keinesfalls den biologischen Hintergrund eines Alignments wieder und ist deshalb für praktische Zwecke eher ungeeignet. Die im Moment gebräuchlichsten Bewertungsfunktionen lesen die Bewertung aus eine sogenannte Scoring Matrix aus. Für Proteine kann man die PAM- oder Blosum-Matrizen benutzen. Diese Matrizen mit der Größe 20*20 (bzw. 24*24, wenn noch einige Sonderfälle beachtet werden) enthalten Bewertungen (sogenannte log-odds) dafür, dass eine Aminosäure durch eine andere substituiert wird. Die log-odds basieren auf Wahrscheinlichkeiten. Berechnet werden diese Scoring-Matrizen ebenfalls aus Sequenzalignments.

Die oben verwendete Bewertungsfunktion wird benutzt um die Ähnlichkeit zweier Sequenzen zu bestimmen. Um nun die Distanz bestimmen zu können kann man einfach die Bewertungsfunktion ändern, d.h bei Ungleichheit kann man einen positiven Wert zurückgeben, welcher als Strafe interpretiert werden kann und bei Gleichheit 0 oder einen negativen Wert. Es muss allerdings beachtet werden, dass in der Rekursion bei einer distanzbasierten Bewertungfunktion nicht das Maximum, sondern das Minimum ermittelt werden muss.

Ein Beispiel für eine distanzbasierte Bewertungsfunktion:

$w(x,y)={\begin{cases}0&,x\ =\ y\\1&,x\ \not =\ y\\\end{cases}}$

Probleme

Durch die Konstruktion des Algorithmus in der dargestellten Weise ist es nicht möglich die längste übereinstimmende Teilsequenz beider Sequenzen zu finden. Um dieses Problem zu beheben wurden von T.Smith und M. Waterman einige Modifikationen am Needleman-Wunsch-Algorithmus vorgeschlagen, die genau dies ermöglichen. Der modifizierte Algorithmus wird Smith-Waterman-Algorithmus genannt.

Komplexität

Der Needleman-Wunsch-Algorithmus hat eine Laufzeitkomplexität von $O(n*m)$. Diese Laufzeit entsteht, da der Algorithmus jedes der $n*m$ Felder mit einer Bewertung belegen muss. Die Berechnung der Bewertung funktioniert in der Regel in $O(1)$.

Die Speicherplatzkomplexität liegt ebenfalls bei $O(n*m)$, da alle Werte der Feldes abgespeichert werden.

Der Algorithmus ist daher für lange Alignments eher ungeeignet. Wenn man davon ausgeht, dass man Integerwerte speichert und diese in Java 4 Byte Speicher belegen, kommt man bei einem Alignment von $10.000*10.000$ Buchstaben schon auf $100.000.000*4$ Byte $\approx$ 381 MegaByte. Die Alignierung ganzer Genome lässt sich so allerdings noch nicht optimal durchführen. (Ein mittleres Bakteriengenom hat ca. 1 - 4 Millionen Basenpaare. Das Menschliche Genom z.B. hat ca. 3 Milliarden Basenpaare)

Der Hirschberg-Algorithmus hingegen berechnet ein globales Alignment auf linearem Speicherplatz $O(n+m)$.