Bio, C#, Root

03 Der Needleman-Wunsch-Algorithmus

  1. Einführung in die Bioinformatik
  2. Scoring von Alignments
  3. Der Needleman-Wunsch-Algorithmus

Die F-Matrix (abstrakt)

Der Needleman-Wunsch-Algorithmus geht bei der Berechnung des idealen Alignements in zwei Schritten vor, um exponentiell skalierende Rechenzeit zu vermeiden. Zunächst wird die sogenannte F-Matrix berechnet, aus welcher anschließend das ideale Alignment zusammengesetzt wird.

Zur Veranschaulichung alignen wir die Sequenzen GACFC und GCFHC.

Jedes Feld in der Matrix (jede Zelle der Tabelle) beinhaltet einen Score. Diese Punktzahl ist jene, die auf dem Weg zur Berechnung des Feldes gesammelt wurde. Zur Berechnung einer Zelle sind immer die Felder links oben, oben und links notwendig:

Unter den drei Feldern wird jenes gewählt, durch das sich die Höchste Punktzahl erreichen lässt:

  • Oben, wenn eine Lücke in der zweiten Sequenz (Spalten der Tabelle) die meisten Punkte bringt. Punktzahl = oberes Feld + gap penalty = -8 + -8 = -16
  • Oben-links, wenn das entsprechende Aminosäure-Paar die meisten Punkte bringt. Punktzahl = oben-links + Score(G-G) = 0 + 8 = 8
  • Links, wenn eine Lücke in der ersten Sequenz (Zeilen der Tabelle) die meisten Punkte bringt. Punktzahl = linkes Feld + gap Penalty = -8 + -8 = -16

Die höchste Punktzahl ist 8, also wird das Feld oben-links gewählt.

Wir gehen nun ein Feld nach rechts und berechnen wieder die möglichen Punktzahlen:

Die blauen Striche zeigen den Rechenweg, der gegangen wurde, um die Punktzahl eines Feldes zu berechnen. Horizontale und Vertikale Striche verwenden die gap Penalty. Beim diagonalen Rechenweg wird dagegen die Punktzahl des Pärchens verwendet.

Das Paar GC hat einen Score von -3. (siehe BLOSUM50)

  • -16 + -8 = -24
  • -8 + -3 = -11
  • 8 + -8 = 0

Die höchste Punktzahl lässt sich also durch eine Lücke in der ersten Sequenz erreichen.

Auf diese Weise wird die Tabelle Zeile für Zeile gefüllt:

Nur ein Rechenweg hat zur höchstmöglichen Punktzahl geführt. Er beinhaltet sowohl Lücken in der ersten, als auch Lücken in der zweiten Sequenz.

Auswählen des Alignments (abstrakt)

Der zum letzten Feld führende Rechenweg muss nun aufgespürt werden. Dazu wird immer das Nachbarfeld ausgewählt, welches für die Berechnung verwendet wurde. (Die blaue Linie kennt der Algorithmus natürlich nicht.)

Ähnlich wie bei der Berechnung der Tabelle, wird das Nachbarfeld gewählt, durch das die (höchste) Punktzahl errechnet wurde. In diesem Fall das links-oben-Feld. Das entsprechende Aminosäuren-Paar wird dem fertigen Alignment (im Bild unten rechts) hinzugefügt. (Merke: diagonal: Match, links: Lücke im Ersten, links: Lücke im Zweiten)

Das gewählte Feld ist als nächstes dran:

Hier ergibt sich die höchste Punktzahl durch eine Lücke in der ersten Sequenz, also geht’s mit dem linken Feld weiter.

Schließlich gelangt man wieder beim ersten Feld an und hat das korrekte Alignment zusammengesetzt:

Du kannst das Programm auch herunterladen und die Berechnung schrittweise nachverfolgen:

https://www.dropbox.com/s/tadha32wpub5p4c/Needleman-Wunsch%20Visualisierung.zip?dl=0 (falls es nicht starten sollte, musst du wahrscheinlich das .NET Framework 4 Client Profile nachinstallieren – der Installer ist auch dabei)

Wenn du im Programm bist, kannst du auswählen, dass die Berechnung nur schrittweise durchgeführt wird (None, FMatrix, Alignment, Both). Mit dem “calculate”-Knopf wird gestartet und mit “next” geht’s einen Rechenschritt weiter. Wenn du mit der Maus über eine der Zahlen fährst, erscheint die Berechnung des Feldes.

Die F-Matrix (Code)

Der Code für den Needleman-Wunsch Algorithmus ist etwas ausführlicher, als die vorherigen Beispiele, daher habe ich ihn hier auf drei Methoden aufgeteilt:

  • AlignSequenccesAsync -> wird aufgerufen um die Sequenzen zu berechnen
  • CalculateFMatrix -> berechnet aus zwei Sequenzen die F-Matrix
  • SelectAlignment -> wählt aus einer F-Matrix das Alignment aus

Beginnen wir zunächst mit dem Berechnen der F-Matrix:

Oben haben wir gesehen, wie die F-Matrix manuell berechnet werden kann. Bei Sequenzen mit mehreren Tausend Aminosäuren ist das natürlich nicht mehr machbar.

In C# verwenden wir für die F-Matrix ein zweidimensionales int-Array, welches in Zeile 35 initialisiert wird.

Zunächst werden in Zeile 37-41 die jeweils erste Zeile bzw. Spalte der Tabelle gefüllt. Die erste Zelle (gap vs. gap) bekommt hier ausnahmsweise keine gap-penalty – gap-gap (am Anfang des Alignments) macht schließlich auch keinen Sinn. Für die restlichen Zellen der Spalte bzw. Zeile, wird für die Punktzahl zur Berechnung nur die gap-penalty verwendet, schließlich treten auf dem Weg zur Zelle keine echten Paare auf.
Für die übrigen Zellen wird der im vorherigen Kapitel beschriebene Rechenweg durchgeführt: In den Zeilen 47-49 werden die möglichen Punktzahlen aus den umliegenden Feldern berechnet. Die höchste wird anschließend ausgewählt und in der Tabelle (f) gespeichert (Zeile 55).

Auswählen des Alignments (Code)

Die F-Matrix enthält nun alle Informationen, die zur Auswahl des idealen Alignments benötigt werden. Beginnend mit der letzten Zelle, wird nun jeweils die Differenz zu den drei anderen Zellen berechnet (scoreDiag, scoreUp, scoreLeft in Zeile 65-68). Anschließend wird geprüft, ob die jeweilige Nachbarzelle für die Berechnung verwendet wurde. Für die oben, bzw. links liegende Zelle muss die Differenz zur aktuellen Zelle der gap-penalty entsprechen, während die Differenz zur oben-links-Zelle dem Score für das aktuelle Paar entspricht.

Das fertige Alignment (alignmentA, alignmentB) wird mit dem entsprechenden Paar verlängert. Die Zählvariablen i und j werden nun verringert, sodass das gewählte Nachbarfeld als nächstes dran ist. (Für diagonal beide verringern; nach oben j verringern; nach links i verringern) (Zeile 69-87)

Sobald eine der beiden Variablen i, oder j bei 0 angekommen ist (man ist am Ende einer Sequenz angelangt), muss der Rest der verbleibenden Sequenz mit gaps alignt werden. (Zeile 89-100)

alignmentA und alignmentB enthalten nun die berechneten Sequenzen, sind allerdings noch in umgekehrter Reihenfolge. Das Umkehren der Reihenfolge wird in der AlignSequencesAsync-methode durchgeführt. Dazu wird die Methode Reverse-String verwendet, auf die ich jetzt nicht weiter eingehe.

AlignSequencesAsync verknüpft CalculateFMatrix und SelectAlignment:

Anfangs (Zeile 15+16) sind die berechneten Sequenzen leer (“”). Die eigentliche Berechnung wird auf einem neuen Thread durchgeführt. Das wird so gemacht, damit das Programm auch bei längeren Berechnung noch ‘responsive’ bleibt und es nicht so aussieht, als sei es abgestürzt.

In den Zeilen 24 und 25 wird das berechnete Alignment umgekehrt – SelectAlignment geht schließlich rückwärts durch die Tabelle.

Nachdem alignmentTask fertig ist (Zeile 28) wird das fertige Alignment als string[] zurückgegeben (Zeile 29).

So viel zur Theorie. Das Beispielprojekt zum Selbermachen folgt in Kürze.

3 comments

  1. joe

    Dear Kuchenzeit
    Excuse me, but I understood that
    Could you please provide code to Learning ?

    1. thecake

      Hey,
      you can take a look at the code in this project: https://www.dropbox.com/s/tadha32wpub5p4c/Needleman-Wunsch%20Visualisierung.zip?dl=0

      It’s the step-by-step Needleman-Wunsch algorithm that I used to create the pictures above.

      cheers

  2. plusnoir

    Hey! Für meinen Blog (http://pressbit.wordpress.com/) schreibe ich gerade einen Beitrag über eben diesen Algorithmus. Bin bei meinen Recherchen über deinen Eintrag gestolpert. Cooler Artikel!

Leave a Reply

Your email address will not be published. Required fields are marked *