CI09: Rozdiel medzi revíziami

Aktuálna revízia z 19:22, 18. november 2020

Obsah

1 Hľadanie motívov zadefinovaných pravdepodobnostnou maticou
- 1.1 EM algoritmus
- 1.2 Gibbsovo vzorkovanie (Gibbs sampling)
2 Vzorkovanie z pravdepodobnostného modelu vo všeobecnosti
3 Poriadnejšie Gibbsovo vzorkovanie pre motívy

Hľadanie motívov zadefinovaných pravdepodobnostnou maticou

Mame danych n sekvencii $S=(S_{1}\dots S_{n})$ , kazda dlzky m, dlzku motivu L, nulova hypoteza q (frekvencie nukleotidov v genome)
Hladame motiv vo forme pravdepodobnostneho profilu dlzky L a jeho vyskyt v kazdej sekvencii
Nech $W[a,i]$ je pravdepodobnost, ze na pozicii i motivu bude baza a, W cela matica
$o_{i}$ je pozicia vyskytu v sekvencii $S_{i}$ , $O=(o_{1}\dots o_{n})$ su vsetky vyskyty
je jednoduchý súčin, kde pre pozície v oknách použijeme pravdepodobnosti z W, pre pozície mimo okna použijeme q
- $\Pr(S_{i}|W,o_{i})=\prod _{{j=1}}^{{L}}W[S_{i}[j+o_{i}-1],j]\prod _{{j=1}}^{{o_{i}-1}}q[S_{i}[j]]\prod _{{j=o_{i}+L}}^{m}q[S_{i}[j]]$
- $\Pr(S|W,O)=\prod _{{i=1}}^{n}\Pr(S_{i}|W,o_{i})$
Hľadáme W a O, ktoré maximalizujú tuto vierohodnosť Pr(S|W,O)
- Nepozname efektivny algoritmus, ktory by vedel vzdy najst maximum
- Dali by sa skusat vsetky moznosti O, pre dane O je najlepsie W frekvencie z dat
- Naopak ak pozname W, vieme najst najlepsie O
  - v kazdej sekvencii i skusame vsetky pozicie $o_{i}$ a zvolime tu, ktora ma najvyssiu hodnotu $Pr(S_{i}|W,o_{i})$

EM algoritmus

Iterativne zlepsuje W, pricom berie vsetky O vahovane podla ich pravdepodobnosti vzhladom na W z minuleho kola
Videli sme na prednaske, tu je trochu prepisany:

Inicializácia: priraď každej pozícii j v sekvencii $S_{i}$ nejaké skóre $p_{{i,j}}$
Iteruj:
- Spočítaj W zo všetkých možných výskytov v $S_{1},\dots ,S_{k}$ váhovaných podľa $p_{{i,j}}$
- Prepočítaj všetky skóre $p_{{i,j}}$ tak, aby zodpovedali pomerom pravdepodobností výskytu W na pozícii j v $S_{i}$ , t.j. $p_{{i,j}}$ je umerne $Pr(S_{i}|W,o_{i}=j)$ , pricom hodnoty normalizujeme tak, aby sucet v riadku bol 1

Gibbsovo vzorkovanie (Gibbs sampling)

Inicializácia: Vezmi náhodné pozície výskytov O
Iteruj:
- Spočítaj W z výskytov O
- Vyber náhodne jednu sekvenciu $S_{i}$
- Pre každú možnú pozíciu j v $S_{i}$ spočítaj skóre $p_{{i,j}}$ (ako v EM) výskytu W na tejto pozícii
- Zvoľ $o_{i}$ náhodne s váhovaním podľa $s_{{i,j}}$

Takto dostavame postupnost vzoriek $O^{{(0)}},O^{{(1)}},...$ .
Za sebou iduce vzorky sa podobaju (lisia sa len v jednej zlozke $o_{i}$ ) nie su teda nezavisle
Pre kazdu vzorku $O^{{(t)}}$ najdeme najlepsie $W^{{(t)}}$ a spocitame vierohodnost $\Pr(S|W^{{(t)}},O^{{(t)}})$ . Nakoniec vyberieme O a W, kde bola vierohodnost najvyssia.
Tento algoritmus (s malymi obmenami) bol pouzity v clanku Lawrence, Charles E., et al. (1993) "Detecting subtle sequence signals: a Gibbs sampling strategy for multiple alignment." Science.
- V clanku v kazdej iteracii maticu W rataju zo vsetkych sekvencii okrem $S_{i}$
- Obcas robia krok, kde nahodne skusaju posunut vsetky vyskyty o jedna dolava alebo doprava
- Tento algoritmus nie je uplne matematicky korektne Gibbsovo vzorkovanie (nema ani poradne zadefinovane rozdelenie, z ktoreho vzorkuje). Na spodku stranky pre informaciu uvadzame algoritmus Gibbsovho vzorkovanie pre hladanie motivov z ineho clanku.

Vzorkovanie z pravdepodobnostného modelu vo všeobecnosti

Majme pravdepodobnostny model, kde D su nejake pozorovane data a X nezname nahodne premenne (napr pre nas D su sekvencie S a X su vyskyty O, pripadne aj matica W)
mozeme hladat X pre ktore je vierohodnost Pr(D|X) najvyssia
alebo mozeme nahodne vzorkovat rozne X z Pr(X|D)

Pouzitie vzoriek

spomedzi ziskanych vzoriek zvolime tu, pre ktoru je vierohodnost Pr(D|X) najvacsia (iny pristup k maximalizovaniu vierohodnosti)
ale vzorky nam daju aj informaciu o tom, aka je velka neurcitost v odhade X
- mozeme odhadovat stredne hodnoty a odchylky roznych velicin
- napr. pri hladani motivov mozeme sledovat ako casto je ktora pozicia vyskytom motivu

generovat nezavisle vzorky z Pr(X|D) moze byt tazke
metoda Markov chain Monte Carlo (MCMC) generuje postupnost zavislych vzoriek $X^{{(0)}},X^{{(1)}},\dots$ , konverguje v limite k cielovej distribucii Pr(X|D)
Gibbsovo vzorkovanie je specialnym pripadom MCMC

Markovove reťazce

Markovov reťazec je postupnosť náhodných premenných $X^{{(0)}},X^{{(1)}},\dots ,$ taká, že $\Pr(X^{{(t)}}|X^{{(0)}},\dots ,X^{{(t-1)}})=\Pr(X^{{(t)}}|X^{{(t-1)}})$ , t.j. hodnota v čase $t$ závisí len od hodnoty v čase $t-1$ a nie ďalších predchádzajúcich hodnôt.
Nás budú zaujímať homogénne Markovove reťazce, v ktorých $\Pr(X^{{(t)}}|X^{{(t-1)}})$ nezávisí od $t$ .
Tiez nas zaujimaju len retazce v ktorych nahodne premenne nadobudaju hodnoty z konecnej mnoziny (mozne hodnoty nazyvame stavy)
- Napriklad stavy A,C,G,T
- V Gibbsovom vzorkovani pre motivy je stav konfiguracia premennych O (t.j. mame (m-L+1)^n stavov)
  - Vzorka v kroku t zavisi od vzorky v kroku t-1 (a lisi sa len v hodnote jedneho o_i)

Matica

Pravdepodobnosti prechodu medzi stavmi za jeden krok mozeme vyjadrit maticou pravdepodobnosti P, ktorej prvok oznacuje pravdepodobnost prechodu zo stavu x do stavu y
- Sucet kazdeho riadku je 1, cisla nezaporne
Ako $p_{{x,y}}^{t}$ budeme oznacovat $\Pr(X^{{(t)}}=y|X^{{(0)}}=x)$ , tieto hodnoty dostaneme umocnenim matice P na t

Stacionarne rozdelenie

Rozdelenie $\pi$ na mnozine stavov sa nazyva stacionarne pre Markovov retazec $P$ , ak pre kazde j plati $\sum _{{i}}\pi (i)p_{{i,j}}=\pi (j)\,$ (alebo v maticovej notacii $\pi P=\pi$ )
Ak matica P splna urcite podmienky (je ergodicka), existuje pre nu prave jedno stacionarne rozdelenie $\pi$ . Navyse pre kazde x a y plati $\lim _{{t\to \infty }}p_{{x,y}}^{{t}}=\pi (y)\,$

Priklady Markovovskych retazcov v bioinformatike

V HMM stavy tvoria Markovov retazec
Ine varianty: nekonecne stavove priestory (zlozitejsia teoria), spojity cas (videli sme pri evolucnych modeloch), retazce vyssieho radu, kde urcujeme $\Pr(X_{t}|X_{{t-r}},\dots ,X_{{t-1}})$ a pod.
Pouzitie v bioinformatike: charakterizacia nahodnych sekvencii (nulova hypoteza), pre DNA sa pouzivaju rady az do 5, lepsie ako nezavisle premenne

Ergodické Markovove reťazce

Vravime ze matica je ergodicka, ak $P^{t}$ pre nejake t>0 ma vsetky polozky nenulove
Priklady neergodickych matic

1 0          0.5 0.5          0 1             0.5 0.5
0 1          0   1            1 0             1   0
nesuvisla    slabo suvisla    periodicka      ergodicka

V HMM stavy tvoria Markovov retazec; hladanie genov ergodicky stavovy priestor, profilove HMM nie

Markov chain Monte Carlo MCMC

Chceme generovať náhodné vzorky z nejakeho cieloveho rozdelenia $\pi$ , ale toto rozdelenie je prilis zlozite.
Zostavime ergodicky Markovov retazec, ktoreho stacionarne rozdelenie je rozdelenie $\pi$ , tak aby sme efektivne vedeli vzorkovat $X^{{(t)}}$ ak vieme $X^{{(t-1)}}$ .
Ak zacneme z lubovolneho bodu $X^{{(0)}}$ , po urcitom case t rozdelenie $X^{{(t)}}$ priblizne $\pi$
Ale za sebou iduce vzorky nie su nezavisle!
Vieme vsak odhadovat ocakavane hodnoty roznych velicin ${\frac {1}{t}}\sum _{{i=1}}^{t}f(X^{{(t)}})$ konverguje k $E_{\pi }[f(X)]$

Gibbsovo vzorkovanie

Cielove rozdelenie $\pi (X)$ je cez vektory dlzky n $X=(x_{1},...x_{n})$
V kazdom kroku vzorkujeme jednu zlozku vektora $x_{i}$ z podmienenej pravdepodobnosti $\Pr(x_{i}|x_{1},\dots ,x_{{i-1}},x_{{i+1}},\dots x_{n})$
Ostatne hodnoty nechame rovnake ako v predchadzajucom kroku
Hodnotu $i$ zvolime nahodne alebo periodicky striedame $i=1,2,\dots ,n$

Dôkaz správnosti Gibbsovho vzorkovania

Pozor! Gibbsovo vzorkovanie nie je vzdy ergodicke, ak niektore kombinacie hodnot maju nulovu pravdepodobnost!
Treba dokazat, ze ak je ergodicky, tak ma ako stacionarnu distribuciu nase zvolene $\pi$
Definicia: Vravime, ze matice P a rozdelenie $\pi$ splnaju detailed balance, ak pre kazde stavy (dva vektory hodnot) x a y mame $\pi (x)p_{{x,y}}=\pi (y)p_{{y,x}}$
Lema: ak pre nejaky retazec P a nejake rozdelenie plati detailed balance, je stacionarna distribucia pre P
- Dokaz: $\sum _{x}\pi (x)p_{{x,y}}=\sum _{x}\pi (y)p_{{y,x}}=\pi (y)\sum _{x}p_{{y,x}}=\pi (y)$
Lema: pre retazec Gibbsovo vzorkovania plati detailed balance vzhladom na cielove rozdelnie
- Dokaz: uvazujme dva za sebou iduce vektory hodnot x a y, ktore sa lisia v i-tej suradnici. Nech $x_{{-i}}$ su hodnoty vsetkych ostatnych premennych okrem $x_{i}$
- $\pi (x)p_{{x,y}}=\pi (x)\Pr(y_{i}|x_{{-i}})=\Pr(x_{{-i}})\Pr(x_{i}|x_{{-i}})\Pr(y_{i}|x_{{-i}})=\pi (y)\Pr(x_{i}|x_{{-i}})=\pi (y)\Pr(x_{i}|y_{{-i}})=\pi (y)p_{{y,x}}$

Poriadnejšie Gibbsovo vzorkovanie pre motívy

Uvedene pre zaujimavost - podla clanku Siddharthan R, Siggia ED, van Nimwegen E (December 2005). "PhyloGibbs: a Gibbs sampling motif finder that incorporates phylogeny". PLoS Comput. Biol. 1 (7): e67. doi:10.1371/journal.pcbi.0010067. PMID 16477324.

Pravdepodobnostny model

Rozsirime model, aby aj O a W boli nahodne premenne, takze mame rozdelenie Pr(S,W,O)
- Potom chceme vzorkovat z Pr(O|S) (marginalizujeme cez vsetky hodnoty W)
Vygeneruje sa nahodne matica pravdepodobnosti W (napr z roznomernej distribucie cez vsetky matice)
V kazdej sekvencii i sa zvoli okno $o_{i}$ dlzky L (rovnomerne z m-L+1 moznosti)
V okne sa generuje sekvencia podla profilu W a mimo okna sa generuje sekvencia z nulovej hypotezy (ako predtym)

Gibbsovo vzorkovanie

Mame dane S, vzorkujeme O () (ak treba, z mozeme zostavit maticu )
- zacni s nahodnymi oknami $O^{{(0)}}$
- v kroku t+1 zvol jednu sekvenciu i a pre vsetky pozicie $o'_{i}$ spocitaj $\Pr(o'_{i}|O_{{-i}}^{{(t)}},S)$ (kde $O_{{-i}}=o_{1}\dots o_{{i-1}}o_{{i+1}}\dots o_{n}$ , t.j. všetky pozície výskytov okrem i-tej).
- nahodne zvol jedno $o'_{i}$ umerne k tymto pravdepodobnostiam
- $O^{{(t+1)}}$ dostaneme z $O^{{(t)}}$ vymenou pozicie v sekvencii i za prave zvolenu
- opakuj vela krat
Konverguje k cielovemu rozdeleniu $\Pr(O|S)$ , ale vzorky nie su nezavisle
Dalsie mozne kroky vo vzorkovani: posun vsetky okna o konstantu vlavo alebo vpravo
Dalsie moznosti rozsirenia modelu/algoritmu: pridaj rozdelenie cez L a nahodne zvacsuj/zmensuj L, dovol vynechat motiv v niektorych sekvenciach, hladaj viac motivov naraz,...

Ako spocitat $\Pr(o_{i}|O_{{-i}},S)$ ?

nezaujimaju nas normalizacne konstanty, lahko znormalizujeme scitanim cez vsetky $o'_{i}$
$\Pr(o_{i}|O_{{-i}},S)=\Pr(O|S)/\Pr(O_{{-i}}|S)$ , ale menovatel konstanta
$\Pr(O|S)=\Pr(S|O)\Pr(O)/\Pr(S)$ , kde $\Pr(S)=\sum _{{O'}}\Pr(S|O')\Pr(O')$
Menovatel nas nezaujima (normalizacna konstanta)
$\Pr(O)$ je tiez konstanta (rovnomerne rozdelenie pozicii okien)
Teda mame $\Pr(o_{i}|O_{{-i}},S)$ je umerne $\Pr(S|O)$
Lahko vieme spocitat $\Pr(S|W,O)$ , potrebujeme "zrusit" W, da sa spocitat vzorec...
Skusame vsetky mozne hodnoty $o'_{i}$ , pocitame pravdepodobnost $\Pr(S|O)$ , vzorkujeme umerne k tomu

Dalsie detaily vypoctu $\Pr(S|O)$ :

Nech $S_{o}$ su len sekvencie v oknach a $S_{{-o}}$ mimo okien. Mame $\Pr(S|O)=\Pr(S_{o}|O)\Pr(S_{{-o}}|O)$
$\Pr(S_{{-o}}|O)$ lahko spocitame (nezavisi od W)
$\Pr(S_{o}|O)=\int \Pr(S_{o}|O,W)\Pr(W)dW$ kde integral ide cez hodnoty, kde $w_{{a,i}}\geq 0$ a $\sum _{a}w_{{a,i}}=1\,$
$\Pr(W)$ je konstanta (rovnomerne rozdelenie; nejde o pravdepodobnost ale hustotu), $\Pr(S_{o}|O,W)=\prod _{{i=1}}^{L}\prod _{a}(w_{{a,i}})^{{n_{{a,i}}}}$ , kde $n_{{a,i}}$ je pocet vyskytov bazy a na pozicii i v oknach $o_{1}\dots o_{n}$
$\Pr(S_{o}|O)=\prod _{{i=1}}^{L}3!/(n+3)!\prod _{a}n_{{a,i}}!$ (bez dokazu)