DAR Approximate string matching Casus: biological sequence alignment
1
Text search Approx string matching dynamic programming, edit distance example application: Google search
Text indexing inverted list q-grams
Casus: protein strings similarity, specific for biological domein application of q-gram indexing: BLAST
2 2
Biological sequence alignment
3 3
Biological sequence alignment Ik heb een sequentie S die bij soort A eigenschap X bepaalt. Kan ik in de genomen van andere soorten sequences vinden die grote gelijkenis vertonen met S? Ik heb een aantal patiënten met aandoening Y. Kan ik in hun genoom een sequentie vinden die zij gemeenschappelijk hebben en die afwijkt in vergelijking met niet-patiënten?
4
DNA sequencing: mijlpalen
1953 Crick & Watson ontdekken moleculaire
structuur van het DNA (dubbele helix) 1975 Sanger ontdekt sequencing techniques 2000-2003 human genome sequenced Anno nu: onderzoek naar individuele genetische verschillen
5 5
Sequence basics
Model: een DNA sequence is een string over het
alphabet {A,C,G,T} DNA segmenten worden gekopieerd naar een sequentie van messenger RNA om een eiwit te produceren (U in RNA = T in DNA) triplet structure definieert vertaling naar sequentie van aminozuren (codons)
6 6
Sequence basics
7 7
Sequence basics Een gen is een DNA-fragment dat een proteine codeert (versimpeld) Een proteine = sequentie van aminozuren (20 letter alphabet) ATGACCAGGATCTTTAAGTGA start codon = ATG stop codon = TGA codons: ATG-ACC-AGG-ATC-TTT-AAG-TGA methionine-threonine-arginine-… MTR… 6 manieren om een DNA-string af te lezen
8 8
Biological sequence databases
SWISSPROT, GENBANK, … Zowel proteines als ACGT-sequenties Generieke vorm query:
find protein sequences similar to MKYMTVTDLNNAGATV
9 9
Biological sequence databases
SWISSPROT example entries (FASTA format) >gi|1171675|sp|P42268|NDD_BPR70 NUCLEAR DISRUPTION PROTEIN MKYMTVTDLNNAGATVIGTIKGGEWFLGTPHKDILSKPGFYFLVSEFDGSCV SARFYVGNQRSKQGFSAVLSHIRQRRSQLARTIANNNMAYTVFYLPASKM KPLTTGFGKGQLALAFTRNHHSEYQTLEEMNRMLADNFKFVLQAY >gi|123527|sp|P05228|HRP2_PLAFA HISTIDINE-RICH PROTEIN PRECURSOR (CLONE PFHRP-III) MVSFSKNKVLSAAVFASVLLLDNNNSEFNNNLFSKNAKGLNSNKRLLHESQA HAGDAHHAHHVADAHHAHHVADAHHAHHAANAHHAANAHHAANAHHAANA HHAANAHHAANAHHAANAHHAANAHHAANAHHAANAHHAANAHHAANAHH AANAHHAADANHGFHFNLHDNNSHTLHHAKANACFDDSHHDDAHHDGAHH DDAHHDGAHHDDAHHDGAHHDGAHHDGAHHNATTHHLH
10 10
Protein based data similarity Twee fragmenten DNA zijn homoloog als ze
overeenkomsten vertonen die gebaseerd zijn op gemeenschappelijk afstamming We zoeken een concept sequence similarity dat ‘gapping’ ondersteunt en rekening houdt met ‘letter distance’ GSAQVKGHGKKVA G+ +VK+HGKKV GNPKVKAHGKKVL
HV---D--DMPNAL ++ +L QLQVTGVVVTDATL
11 11
Similarity
Doel 1: een probabilistisch model opstellen ter
onderbouwing van notie van sequence similarity Doel 2: algoritmen opstellen om optimale alignments te vinden
12 12
Scoring model
Twee fragmenten DNA zijn homoloog als ze
overeenkomsten vertonen die gebaseerd zijn op gemeenschappelijk afstamming Men heeft collecties van paren strings verzameld die corresponderen met homologieën Op grond van deze collecties kunnen we tellingen doen betreffende de frequenties van matching van paren aminozuren Door (bio)chemische oorzaken zijn niet alle matches even waarschijnlijk
13 13
Scoring model
String x bestaat uit een reeks symbolen xi Symbool xi heeft een waarschijnlijkheid qxi Volgens random model R (er is geen sprake van een
homologie) hebben twee strings x en y een waarschijnlijkheid
14 14
Scoring model Volgens het Match model komen paren strings voor met een gezamenlijke waarschijnlijkheid, uitgaande van gemeenschappelijke afstamming van een symbool c
De odds-ratio (voor wel/geen homologie) is dan
15 15
Scoring model Om mathematisch-technische redenen werken we liever met log-odds
De log-odds ratio voor strings x en y is
16 16
BLOSUM matrices
Op basis van tellingen zijn matrices ter representatie van s(a,b) opgesteld; hieronder een fragment van de BLOSUM50-matrix D, E en K charged; V, I en L hydrofoob Variatie t.b.v. sensitiviteit (50, 62, 80)
17 17
Scoring model: gaps Penalty voor gap met lengte g optie 1
optie 2 (d = gap open, e = gap extension)
18 18
Alignment algoritmen
Aantal mogelijke alignments voor sequences met
lengte m en n is erg groot Dynamic programming: druk de optimale alignment van twee sequences uit in de optimale oplossingen voor deelsequences
19 19
Needleman-Wunsch (global)
Globale match Eenvoudige gap penalty We zoeken de optimale match van
x1 .. xm en y1 .. yn We hebben optimale matches van x1 .. xm-1 en y1 .. yn-1 x1 .. xm-1 en y1 .. yn x1 .. xm en y1 .. yn-1
20 20
Needleman-Wunsch
De optimale match van x1 .. xm en y1 .. yn verkrijg je door de beste optie te kiezen uit: x1 .. xm-1 en y1 .. yn-1 ,
waarbij xm gepaard wordt aan yn x1 .. xm-1 en y1 .. yn : xm gepaard aan gap x1 .. xm en y1 .. yn-1 : yn gepaard aan gap
21 21
Needleman-Wunsch
De keuze uit deze drie mogelijkheden wordt bepaald aan de hand van de scorefunctie F
22 22
Needleman-Wunsch De keuze uit deze drie mogelijkheden wordt bepaald aan de hand van de scorefunctie F In de matrix:
F (i-1, j-1)
F (i , j-1)
F (i-1, j)
F (i , j)
Via trace-back pijlen geven we de herkomst van de oplossing weer 23 23
Needleman-Wunsch Voorbeeld Match HEAGAWGHEE en PAWHEAE d = 8 ---
H
E
A
0
P A W
24 24
Needleman-Wunsch Voorbeeld Match HEAGAWGHEE en PAWHEAE d = 8 --
H
E
A
--
0
-8
-16
-24
P
-8
A
-16
W
-24
25 25
Needleman-Wunsch Voorbeeld Match HEAGAWGHEE en PAWHEAE d = 8 --
Needleman-Wunsch Alignment is klaar als je het veld rechtsonderin bereikt. Deze bevat de globale score van de match Time complexity: O(mn)
28 28
Smith-Waterman (local) Voorbeeld Globale match: score = +1 HEAGAWGHE-E --P-AW-HEAE Lokale match: score = +21 HEA HEA \
Lokale matches zijn interessanter
29 29
Smith-Waterman (local) Vraag Hoe pas je Needleman-Wunsch aan zodanig dat je lokale matches vindt? Antwoord 1. Zorg dat je op elk punt in de matrix ‘met een schone lei’ begint: zet negatieve scores direct op 0 2. Je mag overal stoppen
30 30
Smith-Waterman (local)
De keuze uit de mogelijkheden wordt bepaald aan de hand van de scorefunctie F
31 31
Smith-Waterman Voorbeeld Match HEAGAWGHEE en PAWHEAE d = 8 --
H
E
A
--
0
0
0
0
P
0
A
0
W
0
32 32
Smith-Waterman Voorbeeld Match HEAGAWGHEE en PAWHEAE d = 8
--
H
E
A
--
0
0
0
0
P
0
0
0
0
A
0
0
0
W
0
0
33 33
Smith-Waterman Voorbeeld Match HEAGAWGHEE en PAWHEAE d = 8
--
H
E
A
--
0
0
0
0
P
0
0
0
0
A
0
0
0
5
W
0
0
0
0
34 34
Smith-Waterman (local) Alignment kan op elk punt in de matrix eindigen
Als twee strings een 3-gram delen, noemen we dat
een hit HEAGAWGHEE en PAWHEAE
Idee van filtering:
Zoek via 3-gram hits ‘kansrijke posities’ Pas op zo’n kansrijke positie een lokaal expansiealgoritme toe
37 37
Blast
Als twee strings een q-gram delen, noemen we dat
een hit HEAGAWGHEE en PAWHEAE
Het hebben van één hit geeft veel false positives; Hoe kunnen we spelen met selectiviteit en
sensitiviteit?
38 38
Blast
Extensie van hit-notie Twee q-grams zijn similar als hun onderlinge matchscore boven een bepaalde drempel komt (default = 11) HEAGAWGHEE en PAWHEAE Effecten: verhoging sensitivity (recall) verlaging selectivity (precision)
39 39
BLAST2 (PBLAST) Eis voor selectie is Two-hit diagonal principe: paar hits met dezelfde onderlinge afstand in query en db-string query = CWYWRWYYC dbstr = RRWYWAWYYRR query = CWYWRWYYC dbstr = RRWYWABCWYYRR
OK!
FOUT!
Wat gebeurt er met sensitivity/selectivity?
40 40
Blast: expansie
Eenvoudige ungapped expansie: ga uit van match op hits schuif naar rechts totdat match-score een val van
8 punten heeft gemaakt ten opzichte van optimum idem naar links
41 41
Blast: expansie
Gapped expansie: pas lokaal een variant van Smith-Waterman
(dynamic programming) toe
42 42
43
PBLAST: an RDB approach
44 44
PBLAST: an RDB approach
Oefening: Druk BLAST-filtering uit in Relationele algebra (of SQL)