590 likes | 701 Views
V5 – Analyse von Genomsequenzen. - Genom-Assemblierung finde identische k -Tupel - Genom-Alignment Suche nach MUMs (maximal unique matches) andere wichtige Bereiche, für die wir heute keine Zeit haben - Gene identifizieren Hidden Markov Modelle - Transkriptionsfaktorbindestellen
E N D
V5 – Analyse von Genomsequenzen - Genom-Assemblierung finde identische k-Tupel - Genom-Alignment Suche nach MUMs (maximal unique matches) andere wichtige Bereiche, für die wir heute keine Zeit haben - Gene identifizieren Hidden Markov Modelle - Transkriptionsfaktorbindestellen Position Specific Scoring Matrices (PSSM) - finde Repeat-Sequenzen Suche nach bekannten Repeat-Motiven Suche auf Suffix-Baum Softwarewerkzeuge
Whole Genome Shotgun Assemblierung Es gibt 2 Strategien für die Sequenzierung von Genomen: clone-by-clone Methode whole-genome shotgun Methode (Celera, Gene Myers). Die Shotgun Sequenzierung wurde bereits 1977 von F. Sanger et al. eingeführt und ist seither eine Standardmethode für die Sequenzierung von Genen. Umstritten war jedoch, ob man sie auch für komplette Genome verwenden kann. ED Green, Nat Rev Genet 2, 573 (2001) Softwarewerkzeuge
Arachne Programm • von Serafin Batzoglou (MIT, Doktorarbeit 2000) • konstruiere Graph G für Überlappungen zwischen Paaren von reads aus • Shotgun-Daten • prozessiere G um Supercontigs von gemappten reads zu erhalten. Wichtige Variation der whole-genome shotgun Sequenzierung: sequenziere reads jeweils von beiden Enden eines Klons. Da die Inserts nach ihrer Größe ausgewählt werden, ist damit der ungefähre Abstand zwischen dem Paar von reads bekannt. Man nennt diese earmuff(Ohrenwärmer) Verbindungen. Batzoglou et al. Genome Res 12, 177 (2002) Softwarewerkzeuge
Arachne: erzeuge Überlappungsgraphen Liste von reads R = (r1, ..., rN) , N ist die Anzahl der reads. Jeder read ribesitzt eine Länge li< 1000. Wenn beide reads von den Endpunkten desselben Klons stammen (earmuff link), besitzt rieine Verknüpfung zu einem anderen read rjin einer festen Distanz dij. Erstes Ziel: erzeuge Graphen G der Überlappungen (Kanten) zwischen Paaren an reads (Knoten) dies ergibt die Paare an reads in R, die aligniert werden müssen. Da R sehr lang sein kann, sind N2 alignments nicht praktikabel. erstelle Tabelle für das Vorkommen vonk-Tupel (Strings der Länge k) in den reads, zähle die Anzahl von k-Tupel Treffern für jedes Paar an reads. Führe dann paarweise Alignments zwischen den Paaren an reads durch, die mehr als cutoff gemeinsame k-mere besitzen. Batzoglou PhD thesis (2002) Softwarewerkzeuge
Arachne: Tabelle für Vorkommen von k-meren Ermittle die Anzahl an k-Tupel Treffern in der Vorwärts- und Rückwärts-Richtung zwischen jedem Paar von reads in R. (1) Ermittle alle Triplets (r,t,v) r = Nummer des reads in R t = Index eines k-mers, das in r vorkommt v = Richtung des Auftretens (vorwärts oder rückwärts) (2) sortiere die Menge der Paare nach den k-mer Indices t (3) verwende eine sortierte Liste um eine Tabelle T von Quadrubletts (ri, rj, f, v) zu erstellen, wobei riund ridie reads sind, die mindestens einen gemeinsamen k-mer enthalten, v die Richtung angiebt, und f die Anzahl an gemeinsamen k-mers zwischen ri und rjin Richtung v. Batzoglou PhD thesis (2002) Softwarewerkzeuge
Arachne: Tabelle für Vorkommen von k-mers Hier: k = 3 Batzoglou PhD thesis (2002) Softwarewerkzeuge
Arachne: Tabelle für Vorkommen von k-mers • Wenn ein k-Tupel „zu oft“ auftritt gehört er wahrscheinlich zu einer • Repeat-Sequenz. • Man sollte diese nicht für die Detektion von Überlappungen verwenden. • Implementierung • finde k-Tupel(r,t,v) und sortieren sie in 64 Dateien entsprechen den ersten • drei Nukleotiden jedes k-mers. • Für i=1,64 • lade Datei in den Speicher, sortiere nach t, speichere sortierte Datei ab. • end • lade 64 sortierte Dateien nacheinander in den Speicher, • fülle Tabelle T nacheinander auf. • In der Praxis ist k = 8 bis 24. Batzoglou PhD thesis (2002) Softwarewerkzeuge
Arachne: paarweise read-Alignments Führe paarweise Alignments zwischen den Reads durch, die mehr als Cutoff gemeinsame k-mers besitzen. Sobald man zu häufige k-mers ausschließt (mehr als ein zweiter Cutoff), ist sichergestellt, daß nur O(N) viele paarweise Sequenzalignments durchgeführt werden müssen. Nur eine kleine Anzahl an Basen-Austauschen und Indels ist in einer überlappenden Region zweier alignierter reads erlaubt. Output des Alignment-Algorithmus: für die reads ri, rj gibt es Quadrubletts (b1, b2, e1, e2) für jede detektierte Überlappungsregion mit den Anfangspositionen b1, b2 und Endpositionen e1,e2. Falls eine signifikante Überlappungsregion vorliegt, wird (ri, rj, b1, b2, e1, e2) eine Kante im Überlappungsgraphen G. Batzoglou PhD thesis (2002) Softwarewerkzeuge
Kombination teilweiser Alignments 3 teilweise Alignments der Länge k=6 zwischen einem Paar von reads werden zu einem einzigen vollen Alignment der Länge k=19 kombiniert. Die vertikalen Linien verbinden übereinstimmenden Basen, wogegen x Mismatche sind. Dies ist eine oft auftretende Situation, in der ein ausgedehnter k-mer Treffer ein volles Alignment von zwei reads ist. Batzoglou et al. Genome Res 12, 177 (2002) Softwarewerkzeuge
Repeats erzeugen Mehrdeutigkeit Ohne das Auftreten von Sequen-zierungsfehlern und Repeats wäre es einfach, alle entdeckbaren paarweise Abstände von reads zu finden und den Graph G zu konstruieren. Da es Repeats jedoch sehr häufig auftreten, bedeutet eine Verbindung zwischen zwei reads in G nicht ohne weiteres eine wahre Überlappung. Eine „Repeat-Verbindung“ ist eine Verbindung in G zwischen zwei reads, die aus verschiedenen Regionen des Genoms stammen und in der repetitiven Sequenz überein-stimmen. Batzoglou PhD thesis (2002) Softwarewerkzeuge
Sequence contigs unerläßlich für die Assemblierung ist die ausreichende Überdeckung (mehrfache Sequenzierung = coverage) derselben Genomregionen Batzoglou PhD thesis (2002) Softwarewerkzeuge
Verbinden von Contigs Durch die Löschung von k-mers hoher Frequenz wird einiges an Repetition im Genom vor der Erzeugung von G effizient maskiert. Zur Erkennung von repetitiven Verbindung dienen weitere heuristische Algorithmen, die hier nicht diskutiert werden sollen. Sequenz-Contigs werden gebildet indem Paare von reads verbunden werden, die eindeutig verbunden werden können. Tatsächlich ist die Situation viel schwieriger als hier gezeigt, da Repeats häufig nicht zu 100% zwischen Kopien konserviert sind. Batzoglou PhD thesis (2002) Softwarewerkzeuge
Benutze Überlapp-Paarungen um die reads zu verbinden Arachne sucht nach 2 Plasmiden mit gleicher Insert-Länge, deren Sequenzen an beiden Enden überlappen paired pairs. (A) A paired pair of overlaps. The top two reads are end sequences from one insert, and the bottom two reads are end sequences from another. The two overlaps must not imply too large a discrepancy between the insert lengths. (B) Initially, the top two pairs of reads are merged. Then the third pair of reads is merged in, based on having an overlap with one of the top two left reads, an overlap with one of the top two right reads, and consistent insert lengths. The bottom pair is similarly merged. Unten: eine Menge von paired pairs werden zu contigs zusammengefasst und eine Konsensussequenz erzeugt. Batzoglou et al. Genome Res 12, 177 (2002) Softwarewerkzeuge
Detection of repeat contigs Some of the identified contigs are repeat contigs in which nearly identical sequence from distinct regions are collapsed together. Detection by (a) repeat contigs usually have an unusually high depth of coverage. (b) they will typically have conflicting links to other contigs. Contig R is linked to contigs A and B to the right. The distances estimated between R and A and R and B are such A and B cannot be positioned without substantial overlap between them. If there is no corresponding detected overlap between A and B then R is probably a repeat linking to two unique regions to the right. After marking repeat contigs, the remaining contigs should represent the correctly assembled sequence. Batzoglou et al. Genome Res 12, 177 (2002) Softwarewerkzeuge
Contig assembly If (a,b) and (a,c) overlap, then (b,c) are expected to overlap. Moreover, one can calculate that shift(b,c)=shift(a,c)-shift(a,b). A repeat boundary is detected toward the right of read a, if there is no overlap (b,c), nor any path of reads x1, ..., xksuch that (b,x1), (x1,x2) ..., (xk,c) are all overlaps, and shift(b,x1) + ... + shift(xk,c) shift(a,c) – shift(a,b). Batzoglou et al. Genome Res 12, 177 (2002) Softwarewerkzeuge
Consistency of forward-reverse links • The distance d(A,B) (length of gap or negated length of overlap) between two linked contigs A and B can be estimated using the forward-reverse linked reads between them. • The distance d(B,C) between two contigs B,C that are linked to the same contig A can be estimated from their respective distances to the linked contig. Batzoglou et al. Genome Res 12, 177 (2002) Softwarewerkzeuge
Filling gaps in supercontigs • Contigs A and B are connected by a path p of contigs X1,..., Xk. The distance dp(A,B) between A and B (along the path p) is the length of the sequence in the path that does not overlap A and B. • Contigs Y1 and Y2 share forward-reverse links with the supercontig S. These links position them in the vicinity of the gap between A and B. Therefore, Y1 and Y2 will be used as possible stepping points in the path closing the gap from A to B. Batzoglou et al. Genome Res 12, 177 (2002) Softwarewerkzeuge
Contig Coverage and Read Usage Batzoglou et al. Genome Res 12, 177 (2002) Softwarewerkzeuge
Characterization of Contigs and Supercontigs Batzoglou et al. Genome Res 12, 177 (2002) Softwarewerkzeuge
Base Pair Accuracy base quality x*10 means that (on average) one sequencing error occurs in 10-x bases. Batzoglou et al. Genome Res 12, 177 (2002) Softwarewerkzeuge
Misassemblies Batzoglou et al. Genome Res 12, 177 (2002) Softwarewerkzeuge
Computational Performance Batzoglou et al. Genome Res 12, 177 (2002) Softwarewerkzeuge
Contig Coverage and Read Usage Batzoglou et al. Genome Res 12, 177 (2002) Softwarewerkzeuge
Comparison of different assemblers • you should look out for: • - smallest number of contigs + misassembled contigs • highest possible coverage by contigs • lowest possible coverage by misassembled contigs Pevzner, Tang, Waterman PNAS 98, 9748 (2001) Softwarewerkzeuge
There is no error-free assembler to date Comparative analysis of EULER, PHRAP, CAP, and TIGR assemblers (NM sequencing project). Every box corresponds to a contig in NM assembly produced by these programs with colored boxes corresponding to assembly errors. Boxes in the IDEAL assembly correspond to islands in the read coverage. Boxes of the same color show misassembled contigs. Repeats with similarity higher than 95% are indicated by numbered boxes at the solid line showing the genome. To check the accuracy of the assembled contigs, we fit each assembled contig into the genomic sequence. Inability to fit a contig into the genomic sequence indicates that the contig is misassembled. For example, PHRAP misassembles 17 contigs in the NM sequencing project, each contig containing from two to four fragments from different parts of the genome. „Biologists "pay" for these errors at the time-consuming finishing step“. Pevzner, Tang, Waterman PNAS 98, 9748 (2001) Softwarewerkzeuge
What comes next? Finishing the genome Usually, the assembly of shotgun data is finished with a number of contigs with some remaining gaps. Also, within each contig there are some regions of high error rate. The goal of the finishing phase is then to get a single continuous contig with low error rate. „Finishers“ apply ad hoc rules to decide where additional data is necessary. This experimental data may then be generated in experiments using different chemistry or higher coverage. Autofinish (phrap group) is a program to help humans with deciding which new reads to get. Softwarewerkzeuge
Whole Genome Alignment (WGA) Nachdem die genomische DNA-Sequenz eng verwandter Organismen verfügbar wird, ist die erste Frage, wie das Alignment beider Genome aussieht. Globale Genom-Alignments machen nur für eng verwandte Organismen Sinn. Im anderen Fall muß man erst die genomischen Rearrangements betrachten. Dann kann man die systenischen Regionen (Regionen, in denen Gen-Reihenfolge des nächsten gemeinsamen Vorfahrens in beiden Spezies konserviert blieb) betrachten und lokale Genom-Alignments dieser Regionen produzieren. Softwarewerkzeuge
Vergleich von Maus und Mensch auf Genomebene aus dem Paper des Mouse Genome Sequencing Consortiums „Initial sequencing and comparative analysis of the mouse genome“, Nature 420, 520-562 (5.12.2002). Excellent paper! Well readable! Wichtigste Ergebnisse: * das Mausgenom ist etwa 14% kürzer als das menschliche Genom. Die unterschiedliche Länge liegt wohl an der höheren Deletionsrate in Maus. * über 90% des Maus- und Menschen-Genoms kann in entsprechende Regionen mit konservierter Syntenie eingeteilt werden * auf dem Nukleotid-Level kann etwa 40% des menschlichen Genoms mit dem Maus-Genom aligniert werden (diese am stärksten orthologen Sequenzen blieben wohl in beiden Linien vom gemeinsamen Vorfahren erhalten). Der Rest wurde wohl in einem oder beiden Genomen gelöscht. * die neutrale Substitutionsrate beträgt etwa 0.5 Nucleotid-Substitutionen pro Position seit der Divergenz der beiden Spezien. Etwa doppelt so viele Austausche haben in Maus gegenüber Mensch stattgefunden. Softwarewerkzeuge
Vergleich von Maus und Mensch auf Genomebene Key findings: * der Anteil kurzer (50-100 bp) Segmente in den Säugetier-Genomen, der reinigender Selektion unterliegt, ist etwa 5%, d.h. wesentlich höher als der Anteil der Protein-kodierenden Regionen Genome enthalten viele zusätzliche Eigenschaften wie UTRs (untranslated regions), regulatorische Elemente, nicht-Protein-kodierende Gene, chromosomale Strukturelemente, die unter Selektion für die biologische Funktion stehen. * die Evolution von Säugetier-Genomen verläuft ungleichmäßig. Es gibt deutliche Unterschiede an Divergenz je nach Genomposition. * Sowohl Maus wie Mensch-Genom enthalten etwa 30.000 Gene, die für Proteine kodieren. Der Anteil an Mausgenen mit einem eindeutigen Orthologen im menschlichen Genom ist etwa 80%. Der Anteil der Mausgene ohne ein homologes Gen im menschlichen Genom ist < 1%. Softwarewerkzeuge
Konservierung von Syntenie zwischen Mensch und Maus Ein typisches 510-kb Segment des Maus-Chromosoms 12, das mit einem 600-kb Stück des menschlichen Chromosom 14 verwandt ist. Blaue Linien: reziprok eindeutige Treffer in beiden Genomen. Rote Markierungen kennzeichnen die Länge der passenden Regionen. Die Abstände zwischen diesen „Landmarks“ sind im Maus-Genom kleiner als im Mensch, was mit der 14% kürzeren Gesamtlänge des Genoms übereinstimmt. The mouse genome. Nature420, 520 - 562 Softwarewerkzeuge
Entsprechung syntenischer Regionen 342 Segmente und 217 Blöcke >300 kb mit konservierter Syntenie im Mensch sind im Maus-Genom markiert. Jede Farbe entspricht einem bestimmten menschlichen Chromosom. The mouse genome. Nature420, 520 - 562 Softwarewerkzeuge
Sensitivität Im globalen Mensch:Maus Alignment sind mehr als eine Millionen Regionen stärker als 70% konserviert (auf 100-bp Level) – diese Regionen decken > 200 Million bp ab. Nur 62% von ihnen werden von (lokalen) BLAT-Treffern abgedeckt. Dies bedeutet, daß man 38% der konservierten Abschnitte nur durch das globale Alignment finden kann! Idee: lokales Alignment soll als Anker-Verfahren für anschliessendes globales Alignment dienen. Dadurch hofft man, viele zusätzliche konservierte Regionen ausserhalb der Anker-Regionen zu finden. Couronne, ..., Dubchak, Genome Res. 13, 73 (2003) Softwarewerkzeuge
hohe Sensitivität von globalen Alignments Beispiel: das globale Alignment der mouse finished sequence NT_002570 gegen die Region, die mit BLAT-Ankern gefunden wurde, zeigt konservierte kodierende und nicht-kodierende Elemente, die mit BLAT nicht gefunden wurden. Couronne, ..., Dubchak, Genome Res. 13, 73 (2003) Softwarewerkzeuge
Zusätzliche Informationen aus globalem WGA • Unterschiede in Repeat-Merkmalen • Duplikationen (große Fragmente, chromosomal) • Tandem-Repeats • Große Insertionen und Deletionen • Translokationen von einem Teil des Genoms zu einem anderen • Single Nucleotide Polymorphism Softwarewerkzeuge
Methods for WGA: iterative pairwise global alignment These Methods follow a general strategy of iteratively merging two multiple alignments of two disjoint subsets of sequences into a single multiple alignment of the union of those subsets. Construct a hash table on either the query string, or the database string (or both) for all possible substrings of a pre-specified size (say l) Find exactly matching substrings of length l using this hash table (seeds). In the second phase, these seeds are extended in both directions, and combined if possible, in order to find better alignments. If the global pairwise alignment of two genomic DNA sequences S1 and S2 is computed by standard dynamic programming algorithms (which requires O( | S1|∙| S2| time, where |S| is the length of sequence S) such iterative methods cannot be used in practice to align DNA sequences of entire genomes due to time and memory limitations. examples are: FASTA, BLAST, MegaBLAST, BL2SEQ, Wu-blast, flash,PipMaker (BLASTZ), and PatternHunter Softwarewerkzeuge
Methods for WGA: anchor-based global multiple alignment These methods try to identify substrings of the sequences under consideration that they are likely parts of a global alignment. (As mentioned, these substrings can be obtained from local alignments). These substrings form „anchors“ in the sequences to be aligned. These methods first align the anchors and subsequently close the gaps (align the substrings between the anchors). Anchor-based alignment methods are well suited for aligning very long sequences. MUMmer is a very successful implementation of this strategy for aligning two genome sequences. Softwarewerkzeuge
Was ist MUMmer? • A.L. Delcher et al. 1999, 2002 Nucleic Acids Res. • http://www.tigr.org/tigr-scripts/CMR2/webmum/mumplot • Nimm an, dass zwei Sequenzen eng verwandt sind (sehr ähnlich) • MUMmer kann zwei bakterielle Genome in weniger als 1 Minute alignieren • nutzt Suffix-Bäume um Maximal Unique Matches zu finden • Definition eines Maximal Unique Matches (MUM): • Eine Subsequenz, die in beiden Sequenzen genau einmal ohne Abweichungen vorkommt und in keine Richtung verlängert werden kann. • Grundidee: ein MUM ausreichender Länge wird sicher Teil eines globalen Alignments sein. A maximal unique matching subsequence (MUM) of 39 nt (shown in uppercase) shared by Genome A and Genome B. Any extension of the MUM will result in a mismatch. By definition, an MUM does not occur anywhere else in either genome. Delcher et al. Nucleic Acids Res 27, 2369 (1999) Softwarewerkzeuge
MUMmer: wichtige Schritte • Erkenne MUMs (Länge wird vom Benutzer festgelegt) ACTGATTACGTGAACTGGATCCA ACTCTAGGTGAAGTGATCCA ACTGATTACGTGAACTGGATCCA ACTCTAGGTGAAGTGATCCA 10 1 20 ACTGATTACGTGAACTGGATCCA ACTC--TAGGTGAAGTG-ATCCA 1 10 20 Softwarewerkzeuge
Definition von MUMmers • Für zwei Strings S1 und S2 und einen Parameter l • Der Substring u ist eine MUM Sequenz wenn gilt: • |u| > l • u kommt genau einmal in S1 und genau einmal in S2 (Eindeutigkeit) • Für jeden Buchstaben a kommt weder ua noch au sowohl in S1 als auch in S2 vor (Maximalität) Softwarewerkzeuge
Wie findet man MUMs? • Naiver Ansatz • Vergleiche alle Teilsequenzen von A mit allen Teilsequenzen von B. Dies dauert O(nn) • verwende Suffix-Bäume als Datenstruktur • ein naiver Ansatz, einen Suffix-Baum zu konstruieren hat eine quadratische Komplexität in der Rechenzeit und dem Speicherplatz • durch klevere Benutzung von Pointern gibt es lineare Algorithmen in Rechenzeit und Speicherplatz wie den Algorithmus von McCreight Softwarewerkzeuge
Suffix-Bäume • Suffix-Bäume sind seit über 20 Jahren wohl etabliert. • Einige ihrer Eigenschaften: • ein “Suffix” beginnt an jeder Position I der Sequenz und reicht bis zu ihrem Ende. • Eine Sequenz der Länge N hat N Suffices. • Es gibt N Blätter. • Jeder interne Knoten hat mindest zwei Kinder. • 2 Kanten aus dem selben Knoten können nicht mit dem selben Buchstaben beginnen. • Am Ende wird $ angefügt CACATAG$ Softwarewerkzeuge
Konstruktion eines Suffix-Baums CACATAG$ Suffixes: 1. CACATAG$ C A C A T A G $ 1 Softwarewerkzeuge
Konstruktion eines Suffix-Baums CACATAG$ A Suffixes: 1. CACATAG$ 2. ACATAG$ C A C C A A T T A A G G $ $ 2 1 Softwarewerkzeuge
Konstruktion eines Suffix-Baums CACATAG$ A Suffixes: 1. CACATAG$ 2. ACATAG$ 3. CATAG$ C A C C A A T T T A A A G G G $ $ $ 2 3 1 Softwarewerkzeuge
Konstruktion eines Suffix-Baums CACATAG$ A Suffixes: 1. CACATAG$ 2. ACATAG$ 3. CATAG$ 4. ATAG$ C T A G $ 4 A C C A A T T T A A A G G G $ $ $ 2 3 1 Softwarewerkzeuge
Konstruktion eines Suffix-Baums CACATAG$ A Suffixes: 1. CACATAG$ 2. ACATAG$ 3. CATAG$ 4. ATAG$ 5. TAG$ C T A G $ 4 A C C T A A T T A T A A A G G G G $ $ $ $ 2 3 1 5 Softwarewerkzeuge
Konstruktion eines Suffix-Baums CACATAG$ A Suffixes: 1. CACATAG$ 2. ACATAG$ 3. CATAG$ 4. ATAG$ 5. TAG$ 6. AG$ C T A G $ 4 A C C T A A T G T A T $ A A A G 6 G G G $ $ $ $ 2 3 1 5 Softwarewerkzeuge
Konstruktion eines Suffix-Baums CACATAG$ G $ 7 A Suffixes: 1. CACATAG$ 2. ACATAG$ 3. CATAG$ 4. ATAG$ 5. TAG$ 6. AG$ 7. G$ C T A G $ 4 A C C T A A T G T A T $ A A A G 6 G G G $ $ $ $ 2 3 1 5 Softwarewerkzeuge
Konstruktion eines Suffix-Baums CACATAG$ $ G $ 8 7 A Suffixes: 1. CACATAG$ 2. ACATAG$ 3. CATAG$ 4. ATAG$ 5. TAG$ 6. AG$ 7. G$ 8. $ C T A G $ 4 A C C T A A T G T A T $ A A A G 6 G G G $ $ $ $ 2 3 1 5 Softwarewerkzeuge
Speicherplatz sparen CACATAG$ [8, 1] [7, 2] [2, 1] Suffix-Bäume können sehr groß werden, da ihre Größe mit der des Genoms vergleichbar ist. Es ist daher wichtig, Speicherplatz effizient zu nützen. [5, 4] [1, 2] [7, 2] [3, 6] [5, 4] [3, 6] [5, 4] Softwarewerkzeuge