Celera Assembler

Der Celera Assembler, ein Genom-Assembler, wurde ursprünglich von dem Unternehmen Celera entwickelt und wird nun als Open-Source-Projekt weitergeführt. Er wird dazu genutzt, aus vielen kurzen genomischen Fragmenten, die durch eine Whole-Genome-Shotgun-Sequenzierung gewonnen wurden, eine lange zusammenhängende Sequenz zu rekonstruieren.

Dabei werden mit einem Seed-and-Extend-Ansatz Überlappungen zwischen den bereinigten Fragmenten gesucht. Daraus wird der Overlap-Graph konstruiert. Eine Spanning-Forest-Heuristik setzt die Fragmente zu Unitigs zusammen. Diese wiederum werden vom Greedy-Path-Merging-Algorithmus mit Hilfe der Mate-Pair-Informationen und dem Contig-Mate-Graph zu Scaffolds arrangiert. Die Lücken zwischen den Contigs können eventuell mit bisher ungenutzten Fragmenten und Mate-Pairs heuristisch gefüllt werden. Am Ende wird die Konsensussequenz durch ein Multialignment aller Fragmente über den gefundenen Scaffolds berechnet.

Theoretische Betrachtungen

Ein Assembly lässt sich durch folgendes Problem beschreiben.

Es sei eine Menge $F={f_{2},f_{2},\dotsb ,f_{n}}$ von Strings gegeben. Gesucht ist ein String S mit folgenden Eigenschaften:

1. Jedes $f_{i}$ ist Substring von S, $\forall f_{i}\in F:f_{i}\in S$
2. S ist minimal, $\forall S'$ mit 1. gilt $|S|\leq |S'|$

Dieses Problem ist als Shortest Common Superstring bekannt. Es ist NP-vollständig.

Im biologischen Kontext kommt verschärfend hinzu, dass die Sequenzen fehlerhaft sein können und stets zwei mögliche Orientierungen eines Fragmentes zu berücksichtigen sind. Außerdem ist es aus biologischer Sicht nicht immer sinnvoll eine minimale Sequenz zu ermitteln, da Fragmente aus repeathaltigen Regionen zu überkomprimierten Ergebnissen führen können. Das und die kombinatorische Komplexität bewirken, dass das hier vorgestellte Verfahren stark heuristisch ist.

Verfahren

Der Celera Assembler implementiert den Whole-Genome-Shotgun-Ansatz. Der Assembler lässt sich nach [MSD+00] in mehrere, aufeinander aufbauende Stufen einteilen, welche alle heuristisch sind. Das Herzstück des Assemblers bildet der Greedy-Path-Merging-Algorithmus.

Eingabedaten

Die Gewinnung der Sequenzfragmente geschieht mit Hilfe der „Double-Barrel“-Shotgun-Sequenzierung. Dabei wird ein Klon von beiden Seiten ansequenziert und man erhält ein sogenanntes Mate-Pair von Reads. Dabei verwendete Klone haben typische Längen von 2 kb, 10 kb, 50 kb und 150 kb und die sequenzierten Reads eine Durchschnittslänge von 550 Basen. Die ursprünglichen Reads werden daraufhin auf ein Intervall mit Basen 98-prozentiger Sicherheit gekürzt. Die den gelesenen Basen zugeordneten Qualitätswerte werden im weiteren Prozess nicht berücksichtigt. In den Reads enthaltene Teile von Klonierungs- oder Sequenzierungsvektoren werden ebenfalls radikal entfernt. Die verbleibenden Stücke hoher Qualität sind die Eingabe des Assemblers. Diese werden im Folgenden als Fragmente oder Reads bezeichnet.

Screener

Die Sequenzfragmente werden zunächst nach bekannten Repeats durchsucht. Fragmente, welche Repeatmuster enthalten, werden maskiert und später gesondert betrachtet. Die bereinigten Fragmente bilden die Eingabe für den nächsten Schritt.

Overlapper

Der Overlapper sucht Überlappungen zwischen Fragmenten. Eingabe dieser Phase sind die beim Screening bereinigten Fragmente. Ausgegeben werden alle signifikanten Überlappungen zwischen den Reads. Für jedes Fragment $f_{i}$ ist also zu bestimmen, wie es mit irgendeinem anderen Fragment $f_{j}$ oder dessen reversen Komplement ${\overline {f_{j}}}$ überlappt. Zwei beliebige Fragmente $f_{i}$ und $f_{j}$ können auf verschiedene Weise überlappen.

Sie überlappen an den Enden, eines ist im anderen enthalten oder sie überlappen nicht. Alignment durch dynamische Programmierung ist bei dieser großen Anzahl Fragmente zu langsam. Beim humanen Genom wären rund 27 Mio. Reads paarweise zu vergleichen. Hier sind außerdem nur Überlappungen mit einem hohen Alignmentscore, das heißt mit mehr als 96 % Identität, relevant. Daher wird ein BLAST-ähnlicher "Seed-and-Extend"-Ansatz verfolgt.

Alle Fragmente werden zunächst konkateniert, sei also $H=f_{1}f_{2}\dotsb f_{n}$. Anschließend sucht man für alle Fragmente $f_{i}$ exakte Matches der k-mere von $f_{i}$ auf der Sequenz H, die sogenannten Seeds. Der Parameter k wird hier zwischen 18 und 22 gewählt. Nun versucht man die Seeds durch Banded Alignment zu erweitern. Die Laufzeit ist linear und nur wirklich gute Überlappungen zwischen $f_{i}$ und dem entsprechenden $f_{j}$ werden gefunden.

Wurde eine Überlappung zwischen $f_{i}$ und $f_{j}$ gefunden, kann das drei Ursachen haben. Die Fragmente überlappen tatsächlich auf der Originalsequenz, die Enden der Fragmente stammen aus sich wiederholender Sequenz (Repeat) oder die Sequenzgleichheit ist rein zufällig. Letzteres wird dadurch vermieden, dass nur Überlappungen einer Länge von mindestens 40 Basenpaare betrachtet werden. Die durch Repeats erzeugten Overlaps können vermieden werden, indem man vorher die Fragmente penibel nach bekannten Repeats screent und maskiert oder indem man k fest wählt, für jedes k-mer die Häufigkeit bestimmt und extrem oft auftretende k-mere ignoriert. Die berechneten Überlappungen zwischen den Fragmenten werden in einem Graphen repräsentiert.

Der Overlap Graph wird wie folgt definiert: für jedes Fragment $f_{i}$ gibt es zwei Knoten, einen Startknoten $s_{i}$, einen Endknoten $e_{i}$ und eine gerichtete Kante $(s_{i},e_{i})$, um die Orientierung des Reads darzustellen. Jede Überlappung wird durch eine ungerichtete Kante dargestellt. Die Overlap-Kanten inzidieren mit Knoten an den Enden, an denen die entsprechenden Fragmente überlappen.

Unitigger

Der Unitigger nutzt den Overlap Graph, um die Fragmente zunächst anzuordnen, das heißt allen Knoten Koordinaten zuzuweisen (Layout). Als Eingabe dienen die Überlappungen zwischen Fragmenten aus der Overlap-Phase. Eine einfache, auch hier verwendete Heuristik ist die des spannenden Waldes. Zunächst werden alle Read-Kanten markiert. Die Overlap-Kanten sortiert man aufsteigend nach Länge und fügt solange die kürzesten Kanten hinzu, bis in jeder Zusammenhangskomponente ein spannender Baum entstanden ist. Overlap-Kanten, die einen Kreis schließen würden, werden ignoriert.

Aus dieser Teilmenge der Kanten ergibt sich eine relative Anordnung der Reads dieser Zusammenhangskomponente. Solch ein vermeintliches Multialignment heißt Contig.

Durch Repeats erzeugte falsche Overlap-Kanten können zum falschen Zusammensetzen eines Contig (Misassembly) führen.

Ein Spannbaum, der die falsche Kante enthält, würde eine falsche Anordnung der Reads erzeugen. Lässt man die falschen Kanten außer Acht ergibt sich natürlich die richtige Anordnung. Es gibt aber in diesem Fall kein Layout, das mit allen gefundenen Überlappungskanten konsistent ist. Daher definiert man: Ein Unitig (unique contig) ist ein Contig, das mit allen Kanten des Overlap Graphen konsistent ist. Reads, die nicht zu einem Unitig gehören, werden im Folgenden nicht verwendet. Nun kann es vorkommen, dass ein längeres Stück Sequenz wiederholt auftritt. Das führt dazu, dass ein Contig mit allen Überlappungskanten des Overlap Graphen konsistent ist, dessen Fragmente aber nicht aus einer einzigartigen Region der Ursequenz stammen. Man definiert weiterhin: Ein U-Unitig (unique unitig) ist ein Unitig, das einzigartig auf der Ursequenz ist, also nicht teilweise oder vollständig in einer sich wiederholenden Region liegt. U-Unitigs werden mit Hilfe statistischer Betrachtungen identifiziert. Es wird angenommen, die Ankunft der Fragmente wäre Poisson-verteilt. Sei R die Anzahl Reads, G die Länge der Ursequenz und u die Länge eines Unitigs. Dann ist R/G die durchschnittliche Anzahl Reads pro Base und $\lambda ={\frac {u\cdot R}{G}}$ die erwartete Anzahl Reads pro Unitig (Erwartungswert). Die Wahrscheinlichkeit, dass ein Unitig k Fragmente enthält, ist mit der Poisson-Verteilung ${\frac {\lambda ^{k}}{k!}}e^{-\lambda }$. Entstand das Unitig aus den Reads zweier Repeats, verdoppelt sich der Erwartungswert und die Wahrscheinlichkeit ist ${\frac {2\lambda ^{k}}{k!}}e^{-(2\lambda )}$.

Literatur

• "Algorithms in Bioinformatics II", Kapitel 17, Prof. Dr. Daniel Huson, Universität Tübingen, 2004 (engl.) (PDF-Datei; 1,49 MB)
• "Der Celera Assembler", Seminararbeit, 2005 (PDF-Datei; 890 KB)

Weblinks

• Projektseite auf SourceForge