P

LAN WYKŁADU

1. Modelowanie matematyczne

2. Zało˙zenia i własno´sci modelowania opartego o ukryte modele Markowa

3. Architektura i parametryzacja modelowych układów HMM

4. Zastosowania, zalety i wady modeli HMM

5. Modelowanie rodziny sekwencji biopolimerów

Bioinformatyka 2 (BT172)

Wykład 7

Ukryte modele Markowa

Krzysztof Murzyn

{5,12}.XII.2005

Z

AŁO ˙ZENIA

HMM

obserwacje mog ˛a by´c dowolnego typu (tj. niekoniecznie liczby, ale

równie˙z dane znakowe (np. sekwencje aminokwasowe, nukleotydowe,

etc.)

obserwowana sekwencja zdarze´n jest uporz ˛adkowana

obserwacje w sekwencji posiadaj ˛a

własno´s´c Markowa

: tj. kolejna

obserwacja w sekwencji nie zale˙zy od jej poprzedzaj ˛acej

M

ODELOWANIE MATEMATYCZNE

definicja

: u˙zycie j˛ezyka matematyki (

funkcja

,

rozkład

,

zmienna

, ...) do

opisania zachowania jakiego´s układu

przykłady

: analiza kontroli metabolicznej (MCA), symulacje dynamiki

molekularnej (MD), ocena uliniowie´n sekwencji (rozkład Gumbela),

klasyfikacja białek (ukryte modele Markowa, HMM), etc.

Model matematyczny opisuje dany układ za pomoc ˛a

zmiennych

.

Warto´sci zmiennych mog ˛a nale˙ze´c do ró˙znych zbiorów, tj. by´c

rozmaitego typu: liczb rzeczywistych, całkowitych, warto´sci

logicznych,

ci ˛agów znakowych

i tym podobnych.

Wła´sciwy model to grupa

funkcji

wi ˛a˙z ˛acych ze sob ˛a ró˙zne

zmienne

za pomoc ˛a okre´slonych

parametrów

i w ten sposób opisuj ˛acych

powi ˛azania mi˛edzy wielko´sciami w układzie.

A

RCHITEKTURA

HMM

A

RCHITEKTURA SYSTEMU

V

EIL

Upstream

Start Codon

Exon

Stop Codon

Downstream

3’ Splice Site

Intron

5’ Splice Site

5’ Poly−A Site

M

ODEL

ATG

a

t

g

Upstream

Exon

M

ODEL EGZONU

a

c

g

t

a

c

g

t

a

c

g

t

a

g

g

a

16 Backedges

3’ Splice Site

Start Codon

5’ Splice Site

Downstream

Pocz ˛atkowy etap modelowania procesu

stochastycznego (np. reprezentowanego przez

sekwencj˛e reszt aminokwasowych lub

nukleotydowych) w oparciu o ukryte modele

Markowa obejmuje ustalenie architektury całego

modelu, który mo˙ze obejmowa´c kilka oddzielnych

FSM o zdefiniowanej topologii (tzw.

modularno´s´c

modelu).

T

HE EXON AND STOP CODON MODELS IN

VEIL

This model can be entered in two ways: either just after outputting

a

start codon

, or upon leaving

the 3’ splice site

model, which

follows the intron model. The three central columns of states

correspond to the three codon positions. Each of these 12 states is

labeled with the base that it can output. The

system outputs bases

three at a time

,

looping back

(

16 possibilities

: a[acgt], c[acgt], ...)

after each codon. Note that the paths corresponding to a

stop codon

(TAA, TAG, and TGA) all force the system to exit from the model

(four states at lower right of figure). Alternatively, the system can

exit through the 5’ splice site

, in which case an intron must follow

the exon. The two

blank states

on either end of the model can

output

any base

; these “absorbing states” allow the model to align itself to

the proper reading frame, as splice junctions need not respect

codon boundaries.

Własno´s´c Markowa w genomowym DNA?

Odpowiednia

architektura poszczególnych modeli pozwala rozwi ˛aza´c ten

problem.

W

ŁASNO ´SCI

HMM

HMM modeluje proces stochastyczny, którego pewne własno´sci nie s ˛a

znane/jawne – innymi słowy s ˛a:

ukryte

kolejne obserwacje s ˛a reprezentowane przez ukryte modele Markowa

(

HMM

), ka˙zdy

HMM

jest automatem sko´nczonym (ang. finite state

machine), na który składa si˛e sko´nczona liczba stanów i przej´s´c mi˛edzy

nimi

topologi˛e okre´slonego FSM okre´sla si˛e mianem

architektury

HMM

ka˙zdemu stanowi w

HMM

przypisuje si˛e prawdopodobie´nstwa emisji

warto´sci zmiennych losowych z zadanego, sko´nczonego zbioru (ang.

emission) i prawdopodobie´nstwa przej´s´c do innych stanów (ang.

transition)

P

ARAMETRYZACJA

HMM

po ustaleniu architektury systemu konieczne jest wyznaczenie dla

ka˙zdego stanu w HMM wielko´sci

prawdopodobie´nstw emisji

i

przej´s´c

do innych stanów

parametry wyznaczane s ˛a w procedurze zwanej

trenowaniem

modelu,

która przebiega o odpowiednio przygotowany zbiór danych (por.

techniki nauczania maszynowego); trenowanie prowadzone jest zwykle

w oparciu o

algorytm Bauma-Welsha

trenowanie modelu pozwala dopasowa´c wielko´sci parametrów systemu

w taki sposób aby sekwencjom obserwacji ze zbioru ucz ˛acego

przypisywane były wysokie warto´sci prawdopodobie´nstw

wielko´s´c zbioru ucz ˛acego

: wprost proporcjonalna do liczby parametrów układu,

odwrotnie proporcjonalna do długo´sci sekwencji obserwacji w zbiorze ucz ˛acym:

wiele krótkich sekwencji vs. mniej długich

lokalne maksimum

: wyznaczone warto´sci parametrów mog ˛a by´c suboptymalne;

rozwi ˛azanie

: trenowanie modelu rozpoczyna´c od przypisania ’sensownych’

pocz ˛atkowych warto´sci parametrów

przetrenowanie

: ang. over-fitting, kiedy układ ´swietnie modeluje dane ze zbioru

ucz ˛acego, a kiepsko radzi sobie z modelowaniem danych nie uj˛etych w procesie

trenowania; wierne

odtwarzanie

danych/wzorców vs. ich

generalizowanie

P

ROFILE

HMM

W przypadku, kiedy na

architektur˛e całego modelu

składaj ˛a si˛e powtarzaj ˛ace

si˛e elementy (tj. modele

Markowa o takiej samej

liczbie i rodzajów stanów

poł ˛aczonych w taki sam

sposób) mówimy o

profilu

ukrytych modeli Markowa

(ang. profile HMMs).

HMMER P

LAN

7

4

D

M1

S

N

B

M2

M3

M4

E

C

T

J

I1

I2

I3

1

D2

D3

D

SAM

Start

End

M1

I1

D1

D2

D3

I2

I3

I4

M2

M3

V

ITERBI VS

. forward

Viterbi

: najbardziej prawdopodobna anotacja sekwencji (egzon/intron):

forward

: prawdopodobie´nstwo, ˙ze najlepsza anotacja jest poprawna:

(stosunek

do

sumy prawdopodobie´nstw dla 6 alternatywnych anotacji z

G

(splice) i 8 z

A

w stanie

5

;

analizowany model mo˙ze wygenerowa´c dan ˛a sekwencj˛e nukleotydow ˛a na ł ˛acznie 14

alternatywnych sposobów korzystaj ˛ac z 14 ró˙znych szlaków stanów)

W

YKORZYSTANIE

HMM

wyznaczanie prawdopodobie´nstwa okre´slonej

serii obserwacji

w oparciu

o przyj˛ety model badanego procesu (

algorytm forward

, ’do-przodu’)

– przy modelowaniu rodziny sekwencji białek, algorytm forward wykorzystywany

jest do oceny podobie´nstwa sekwencji kwerendy z modelem (identyfikacja

sekwencji homologicznych); suma prawdopodobie´nstw wygenerowania danej

sekwencji aminokwasowej na wszystkie mo˙zliwe sposoby (szlaki stanów) przez

okre´slony profil HMM

dekodowanie: przypisanie ka˙zdej z kolejnych obserwacji w sekwencji

najbardziej prawdopodobnego

stanu w modelu (

algorytm Viterbiego

)

– przy modelowaniu rodziny białek, algorytm Viterbiego wykorzystywany jest do

dodawania nowej sekwencji do istniej ˛acego

MSA

(ka˙zdej z kolejnych reszt w nowej

sekwncji zostaje przypisany

najbardziej prawdopodobny

stan modelu HMM)

– przy przewidywaniu struktury drugorz ˛adowej białek, dla danej sekwencji mo˙zna

ustali´c poło˙zenie poszczególnych elementów (

-helisy,

-arkusz, etc.)

– przy przewidywaniu struktury genu org. eukariotycznego, pewne odcinki sekwencji

zostan ˛a opisane jako egzony, inne jako introny, etc.

O

BLICZENIA

...

sekwencja

CTTCATGTGAAAGCAGACGTAAGTCA

szlak stanów

EEEEEEEEEEEEEEEEEE

5

IIIIIII

Prawdopodobie´nstwo wygenerowania przez model

!

okre´slonej sekwencji

"

w oparciu o jeden z mo˙zliwych szlaków stanów

#

:

$&%

"'(#*)!,+.-

$&%

#.)/!,+

0

$&%

"1)

#2'!,+

3547698

:2;

prawdopodobie´nstwo obrania szlaku

6

(tj. iloczyn prawdopodobie´nstw przej´s´c mi˛edzy

kolejnymi stanami) przy generowaniu przez HMM zadanej sekwencji

CTTCA..TCA

(18

stanów

E

, jeden

5

i 7

I

)

354<6=8

:2;?>@ACBDEBFAG9HIJDEBFA/@KD@ABLDEBFAGMKDNBFA/@O>QPFAPSRTD@UB9VXW

354ZY[8

6]\^:2;

prawdopodobie´nstwo wygenerowania przez HMM (tj. ilczonyn prawdopodobie´nstw

emisji symboli nukleotydów z kolejnych stanów w szlaku) zadanej sekwencji

CTTCA..TCA

:

3_4`Y[8

6\a:b;>cB[AedfgHihjDEBFAGfLDEB[ACklKDEBFA/@KDEBFACkDEBFA/@mDNBFACkn>@ACkF@NRD@NB9V]Ho

ostatecznie:

piq

$&%

"'(#*)!,+.-

rQsLt?u`v?v

C

O JEST UKRYTE W

HMM

analizowany HMM generuje dwa ci ˛agi informacji

w

sekwencje nukleotydów (tj.

CTTCA..TCA

)

w

sekwencj˛e stanów (tj.

EEE..5..II

poniewa˙z dana jest wył ˛acznie sekwencja nukleotydowa, któr ˛a HMM ma

wygenerowa´c, sposób w jaki to zrobi (tj. szlak stanów: przej´s´c/emisji) jest

nieznany/ukryty

M

ETODY MODELOWANIA ZBIORU SEKWENCJI

Problem

: dysponuj ˛ac zbiorem spokrewnionych ewolucyjnie sekwencji,

zidentyfikowa´c ich cechy wspólne w celu stworzenia reprezentatywnego

modelu umo˙zliwiaj ˛acego klasyfikowanie innych (nowych, nieznanych)

sekwencji jako członków okre´slonej rodziny

Rozwi ˛azania

:

x

uliniowienia wielosekwencyjne (

MSA

, ang. Multiple Sequence

Alignment

x

wyra˙zenia regularne (

RE

, ang. Regular Expressions)

x

sekwencje konsensusowe

x

ukryte modele Markowa (

HMM

, ang. Hidden Markov Models)

x

pozycyjnie zró˙znicowane macierze warto´sciuj ˛ace (

PSSM

, ang.

Position Specific Substitution Matrices, tzw. profile podstawie´n)

Z

ALETY

HMM

modularno´s´c

: zło˙zony problem mo˙ze by´c rozbity na autonomiczne modele

HMM

przezroczysto´s´c

: samo-opisuj ˛aca si˛e architektura układu modelowego

y

łatwo´s´c interpretacji układu i wyników modelowania (por. sieci

neuronowe)

y

w poł ˛aczeniu z modularno´sci ˛a łatwo´s´c wbudowywania w układ

modelowy (architektura) specyficznych dla danego problemu

informacji

W

ADY

HMM

własno´s´c Markowa

: rozwi ˛azania w odpowiedniej architekturze układu

modelowego

wymagania obliczeniowe

: trenowanie, metoda programowania

dynamicznego (Viterbi, forward)

trenowanie

: problem lokalnego maksimum, over-fitting

P

RZYKŁAD

:

SEKWENCJE O TEJ SAMEJ DŁUGO ´SCI

,

´ZRÓDŁOWE

ULINIOWIENIE BEZ PRZERW

Rozwa˙zmy przedstawienie przykładowego

MSA

jako liniowego ci ˛agu stanów

odpowiadaj ˛acych kolejnym kolumnom uliniowienia...

A: 1.0

C: 0.0

T: 0.0

G: 0.0

A: 0.0

C: 0.75

T: 0.0

G: 0.25

A: 0.75

C: 0.25

T: 0.0

G: 0.0

A: 0.25

C: 0.0

T: 0.75

G: 0.0

A: 0.0

C: 0.0

T: 1.0

G: 0.0

A: 0.25

C: 0.75

T: 0.0

G: 0.0

1.0

1.0

1.0

1.0

1.0

ACATTC

ACCTTC

ACATTC

AGAATA

Porównanie wykorzystania modeli rodziny sekwencji reprezentowanych

przez przykładowe 5 sekwencji nukleotydowych

ACA---ATG

TCAACTATC

ACAC--AGC

AGA---ATC

ACCG--ATC

z

sprawdza´c jako´s´c uliniowienia nowych sekwencji do

istniej ˛acego

MSA

z

wykorzysta´c wyra˙zenie regularne:

[AT][CG][AC][ACG-][-C][-T]A[GT][CG]

a

mo˙ze

[AT][CG][AC][ACGT-]{1,3}A[GT][CG]

,

ale rozwa˙zmy porównanie dwóch nowych

hipotetycznych sekwencji:

TGCT--AGG

i

ACAC--ATC

z

dysponuj ˛ac odpowiednim

PSSM

, ocena podobie´nstwa obu ostatnio

porównywanych sekwencji b˛edzie si˛e istotnie ró˙zni´c:

PROFILE

podstawie´n mogłyby tu by´c rzeczywi´scie pomocne (np. PSI-BLAST)

z

sekwencja konsensusowa:

ACA---ATC

? tak, ale zarówno sposób jej

wyznaczenia jak i pó´zniejsze wykorzystanie w krytyczny sposób zale˙z ˛a

od przyj˛etej miary dystansu mi˛edzy sekwencjami (kryteriów oceny

podobie´nstwa)

z

hmm, a mo˙ze

HMM

? a precyzyjniej:

PROFILE

HMM

s

P

RZYKŁAD

:

UWZGL ˛

EDNIAMY DELECJE

...

A: 1.0

C: 0.0

T: 0.0

G: 0.0

A: 0.0

C: 0.8

T: 0.0

G: 0.2

A: 0.75

C: 0.25

T: 0.0

G: 0.0

A: 0.75

C: 0.0

T: 0.25

G: 0.0

A: 0.0

C: 0.0

T: 1.0

G: 0.0

A: 0.2

C: 0.8

T: 0.0

G: 0.0

A: 0.0

C: 0.5

T: 0.0

G: 0.5

e

1.0 − e

1.0

0.8

1.0

1.0

0.2

1.0

0.6

0.4

ACATT−

−−−

C

ACCTT−

−−−

C

ACATT−

−−−

C

AGAATG

CGC

A

AC

−−

TG

CGC

C

{

tu w przypadku insercji:

|~}€

i‚

a prawdopodobie´nstwo przej´scia wynosi

ƒ

‚

(insercje

zachodz ˛a w dwóch spo´sród 5 sekwencji)

{

delecje mogłyby by´c modelowane przez dodanie przej´s´c typu:

„2…†‡„2…‰ˆ‹Š

gdzie

Œ

}cjŽUŽŽ‘

; jednak wprowadzenie dodatkowego stanu nie emituj ˛acego ˙zadnego symbolu

(ang. silent states) podnosi transparentno´s´c modelu i obni˙za liczb˛e jego parametrów; por.

’*“•”?–˜—™’*“›šœgNš

vs.

žZŸ

ƒE 

, gdzie

’

jest długo´sci ˛a profilu, np. przy

’

}

”U¡¡

: mamy 390 vs.

4950 parametrów

P

RZYKŁAD

:

UWZGL ˛

EDNIAMY INSERCJE

A: 1.0

C: 0.0

T: 0.0

G: 0.0

A: 0.0

C: 0.75

T: 0.0

G: 0.25

A: 0.75

C: 0.25

T: 0.0

G: 0.0

A: 0.25

C: 0.0

T: 0.75

G: 0.0

A: 0.0

C: 0.0

T: 1.0

G: 0.0

A: 0.25

C: 0.75

T: 0.0

G: 0.0

A: 0.0

C: 0.5

T: 0.0

G: 0.5

e

1.0 − e

0.25

1.0

1.0

1.0

1.0

0.75

ACATT−

−−−

C

ACCTT−

−−−

C

ACATT−

−−−

C

AGAATG

CGC

A

¢

insercje s ˛a tymi odcinkami sekwencji, które nie przystaj ˛a do istniej ˛acego profilu

¢¤£

: prawdopodobie´nstwo wydłu˙zania przerwy (emisji reszty na pozycji odpowiadaj ˛acej

przerwom w innych sekwencjach); tu:

£~¥€¦

§¨

(12 mo˙zliwo´sci insercji w 4 sekwencjach z

czego 3 ’wykorzystane’)

¢

uwzgl˛ednienie jednego dodatkowego stanu w celu opisywania insercji jest rozwi ˛azaniem

niewystarczaj ˛acym (

£

nie zale˙zy od długo´sci insercji), por.

afiniczne punktowanie przerw

:

– wiele stanów insercji; zró˙znicowane wielko´sci prawdopodobie´nstw przej´s´c

©2ª«­¬®ª‰¯

§

vs.

¬®ª«­¬°ª/¯

§

K

LASYFIKACJA SEKWENCJI

M

ODELOWANIE

PRZYPADKOWYCH
SEKWENCJI

M

start

end

w przypadku sekwencji

aminokwasowych,

prawdopodobie´nstwa emisji

odpowiadaj ˛a składowi

sekwencji np. w bazie

SWISSPROT, np.:

A

: 0.0764,

E

: 0.0647,

W

: 0.0121

...

±

poniewa˙z układ HMM jest modelem probabilistycznym,

mo˙ze on generowa´c (niemal) dowolne sekwencje reszt

±

wytrenowanie okre´slonego układu HMM dla danej rodziny

sekwencji powoduje jednak, ˙ze pewne sekwencje b˛ed ˛a

generowane z wy˙zszym prawdopodobie´nstwem

±

ocena podobie´nstwa pary sekwencji

²

; dla zadanej

sekwencji kwerendy wyznacza si˛e prawdopodobie´nstwo jej

wygenerowania przez model odpowiedniej rodziny

³5´¶µ^·U¸]¹

º~»¼»¼½

oraz prawdopodobie´nstwo wygenerowania

takiej samej sekwencji przez model przypadkowy

³5´¶µ^·U¸]¹

¾À¿FÁ/Á½

operuj ˛acy na cz˛esto´sciach wyst˛epowania

poszczególnych reszt w typowych białkach/genomach etc.:

²ÄÃ

Á/ÅÆ

³5´™µa·U¸¹

ºL»Ç»¼½

³5´¶µ^·¸¹

¾È¿Á/Á½

pozytywna ocena

podobie´nstwa par sekwencji (model

rodziny sekwencji vs. model przypadkowych

sekwencji): prawdopodobny homolog

negatywna ocena

podobie´nstwa sekwencji: niskie

prawdopodobie´nstwo, ˙ze dany model rodziny sekwncji

mógł tak ˛a sekwencj˛e wygenerowa´c