生命情報学

生命情報学 （５）

隠れマルコフモデル

阿久津 達也

京都大学 化学研究所

内容



配列モチーフ



最尤推定、ベイズ推定、MAP推定



隠れマルコフモデル(HMM)



Viterbiアルゴリズム



EMアルゴリズム



Baum-Welchアルゴリズム



前向きアルゴリズム、後向きアルゴリズム



プロファイルHMM

モチーフ発見



配列モチーフ

: 同じ機能を持つ遺伝子配列などに見

られる共通の文字列パターン

正規表現など文法表現を用いるもの

例： ロイシンジッパーモチーフ L-x(6)-L-x(6)-L-x(6)-L ジンクフィンガーモチーフ C-x(2,4)-C-x(3)-[LIVMFYWC]-x(8)-H-x(3,5)-H 人間にとってわかりやすいが表現力が弱い

確率的な表現法を用いるもの

重み行列（プロファイル） HMM (隠れマルコフモデル) 人間にとってわかりにくいが 一般に表現力は高い A C G T 3.8 1.5 -1.5 0.2 -3.5 1.3 -2.9 -4.1 -0.3 3.7 2.3 -4.6 3.1 -1.3 1.2 4.2 A C G G A C score= 3.8 + 1.3 + 4.2 + 3.1

モチーフの例

• ジンクフィンガーモチーフ

C-x(2,4)-C-x(3)-[LIVMFYWC]-x(8)-H-x(3,5)-H

• ロイシンジッパーモチーフ

局所マルチプルアラインメント



複数配列と長さ L が与えられた時、スコア最大と

なるように各配列から長さ L の部分列を抽出



モチーフ発見

などに有用

A T C G G A T A T C C G A T T T C G G A A Sequence 1 Sequence 2 Sequence 3

相対エントロピースコアのもとでの

局所マルチプルアラインメント



相対エントロピースコアの定義



f

_j

(a): （モチーフ領域の）j列目におけるaの出現頻度



p(a): aの出現頻度（事前確率）



L: モチーフ領域の長さ

∑ ∑

=

L

₌

j

j

a

j

a

p

a

f

a

f

score

₁

)

(

)

(

log

)

(



実用的アルゴリズム



Gibbsサンプリング

, EMアルゴリズム

Gibbs サンプリング

1. 各配列 x

_j

からランダム

に部分配列 t

_j

を選ぶ

2. 1個の配列 x

_i

をランダム

に選ぶ

３. x

_i

の部分列 t

_i

’ を

に比例する確率で選ぶ

4. t

_i

をt

_i

’ でおきかえる

5. ステップ2-4を十分な回

数だけ繰り返す

( t_i[j]: 部分列t_i のj列目の文字 ) Step 1 x1 x3 t1 t3 x2 t2 Step 2-3 x2 t2’ Probabilistic Search x1 x3 t1 t3 x2 t2’ Step 4

])

[

'

(

])

[

'

(

1

p

t

j

j

t

f

i i j L j

∏

=

最尤推定



P(D|

θ)

（

尤度

）



モデルパラメータ

θ のもとでのデータ D の出現確率



最尤法



P(D|

θ) を最大化する θ を選ぶ



例



コインを5回投げて、表が3回出た後、裏が2回出た



p(表)=a, p(裏)=1-a とすると、P(D|

θ)=a

3

(1-a)

2 

a=3/5の時、 P(D|

θ) は最大

ベイズ推定と

MAP推定



ベイズ推定

：尤度とモデル（パラメータ）の事前確率

から、ベイズの定理により、事後確率を推定



最大事後確率（MAP）推定

 P(D|θ)P(θ) を最大化する θ を計算  P(θ) が一様分布なら最尤推定と同じ

が連続値の時）

（

ただし、

_P

_D

_P

_D

θ'

_P

θ'

θ

D

P

θ

P

θ

D

P

D

θ

P

)

(

)

|

(

)

(

)

(

)

(

)

|

(

)

|

(

θ'

∫

=

=

不正サイコロのベイズ推定



公正サイコロと不正サイコロ

 公正： P(i|公正)=1/6  不正： P(6|不正)=1/2, P(i|不正)=1/10 for i≠6  P(公正)=0.99, P(不正)=0.01 

６が３回続けて出た場合の事後確率

21

.

0

)

99

.

0

(

)

(

)

01

.

0

(

)

5

.

0

(

)

01

.

0

(

)

5

.

0

(

)

666

(

)

(

)

|

666

(

)

666

|

(

3 6 1 3 3

=

+

=

=

P

P

P

P

不正

不正

不正

隠れマルコフモデル

(

HMM)



HMM

≒

有限オートマトン＋確率



定義

 出力記号集合Σ  状態集合 S={1,2,…n}  遷移確率(k→ｌ)

a

_kl  出力確率

e

_k

(b)

 （開始状態= 終了状態= 0) 0.4 0.6 0.3 0.7 0.5 0.5

1

2

3

A: 0.2 B: 0.8 A: 0.7 B: 0.3 A: 0.1 B: 0.9

HMMにおける基本アルゴリズ

ム



Viterbiアルゴリズ

ム

 出力記号列から 状態列を推定  構文解析 

Baum-Welchアル

ゴリズム

(

EMアルゴリズム

）

 出力記号列から パラメータを推定  学習 0.4 0.6 0.3 0.7 0.5 0.5 1 2 3 A: 0.2 B: 0.8 A: 0.7 B: 0.3 A: 0.1 B: 0.9 BABBBAB 2312131 1 2 3 BABBBAB ABBABBAAB BBBABBABAB BAABBBBA 0.4 0.6 0.3 0.7 0.5 0.5 1 2 3 A: 0.2 B: 0.8 A: 0.7 B: 0.3 A: 0.1 B: 0.9

時々いかさまをするカジノ

 サイコロの出目だけが観測可能、どちらのサイコロを振ってい るかは観測不可能  サイコロの出目から、どちらのサイコロを振っているかを推定  6,2,6,6,3,6,6,6, 4,6,5,3,6,6,1,2 →不正サイコロ  6,1,5,3,2,4,6,3, 2,2,5,4,1,6,3,4 →公正サイコロ  6,6,3,6,5,6,6,1, 5,4,2,3,6,1,5,2 →途中で公正サイコロに交換 1: 1/6 2: 1/6 3: 1/6 4: 1/6 5: 1/6 6: 1/6 1: 1/10 2: 1/10 3: 1/10 4: 1/10 5: 1/10 6: 1/2 0.05 0.1 0.95 0.9 公正サイコロ 不正サイコロ

Viterbiアルゴリズム(1)

 観測列（出力配列データ） x=x₁…x_Lと状態列π=π₁…π_Lが与えら

れた時、その同時確率は

P(x,π)=a_{0 π1}Πe_πi (x_i)a_πiπi+1 但し、π_L+1=0

 x が与えられた時、最も尤もらしい状態列は π*=argmax_π P(x,π)  例：どちらのサイコロがいつ使われたかを推定 0.05 公 不 公 不 公 不 0.05 0.1 0.95 0.9 0.95 0.05 0.1 0.9 0.95 0.9 0.1 x2=3 0 ₀ 公 不 0.5 0.5 x1=4 x3=2 6 1 6 1 6 1 3 2 1 3 2 1 π ( ,π) ( , ) 0.5 0.95 0.95 max P

_x

_x

_x

= P

_x

_x

_x

公公公 = ⋅ ⋅ ⋅ ⋅ ⋅

Viterbiアルゴリズム(2)



x から、

π

*

=argmax

_π

P(x

,π)

を計算



そのためには x

₁

…x

_i

を出力し、状態 k に至る確率

最大の状態列の確率

_v

_k

_(i)

を計算



v

_k

(i)は以下の式に基づき

動的計画法

で計算

)

)

(

(

)

(

)

1

(

max

1

v

a

x

e

v

k kl k i l l

i

+

=

+

i











∏

₋

=

=

_e

_x

_a

v

j

_j

_j

i

j

j

i

k

π 1

π

1

π

π

)

(

_max

(

)

Viterbiアルゴリズム(3)

0.05 公 不 公 不 公 不 0.05 0.1 0.95 0.9 公 不 0.95 0.05 0.1 0.9 0.95 0.9 0.1

4

6

1

i

_-1

i

i

+1

}

)

(

1

.

0

)

1

(

),

(

95

.

0

)

1

(

max{

)

1

(

i

_e

_v

i

_e

_v

i

v

公

+

=

公

⋅

⋅

公 公

⋅

⋅

不

EM(Expectation Maximization)アルゴリズム



「

欠けているデータ

」のある場合の

最尤推定

のため

の一般的アルゴリズム



最大化は困難であるので、

反復により尤度を単調

増加

させる（

θ

t

より

θ

t+1

を計算）



HMMの場合、「

欠けているデータ

」は

状態列

の最大化

目標：

パラメータ集合

欠けているデータ、

観測データ、

∑

=

y

)

(

log

)

(

log

:

:

:

x,y|

θ

P

x|

θ

P

θ

y

x

EMアルゴリズムの導出

とすれば尤度は増大 よって、 、 ロピーで常に正なので 最後の項は相対エント とおくと、 右辺第１項を についての和をとり、 をかけて 両辺に ) | ( max arg ) | ( ) | ( ) | ( log ) | ( log ) , | ( ) , | ( log ) , | ( ) | ( ) | ( ) | ( log ) | ( log ) | ( ) , | ( log ) , | ( ) | , ( log ) , | ( ) | ( log ) , | ( ) , | ( log ) | , ( log ) | ( log 1 t θ t t t t t y t t t t t t t t y y t t θ θ Q θ θ θ Q θ θ Q θ x P θ x P θ x y P θ x y P θ x y P θ θ Q θ θ Q θ x P θ x P θ θ Q θ x y P θ x y P θ y x P θ x y P θ x P y θ x y P θ x y P θ y x P θ x P = − ≥ − + − = − − = − = +

∑

∑

∑

EMアルゴリズムの一般形

1.

初期パラメータ

Θ

0

を決定。t=0とする

2.

Q(

θ|θ

t

)=∑P(y|x, θ

t

) log P(x,y|

θ) を計算

3.

Q(

θ|θ

t

)を最大化する

θ

*

を計算し、

θ

t+1

₌

θ

*

とする。t=t+1とする

前向きアルゴリズム



配列 x の生成確率

P(x)=

∑P(x,π)

を計算



Viterbiアルゴリズム

と類似



f

_k

(i)=P(x

₁

…x

_i

,π

_i

=k)

をDPにより計算

a

a

x

e

k k k kl k k i l l k

L

f

x

P

i

f

i

f

f

f

0 0

)

(

)

(

)

1

(

)

(

)

(

0

)

0

(

,

1

)

0

(

∑

∑

=

−

=

=

=

1

2

3

1

2

3

a

11

a

21

a

31 ) ) 1 ( ) 1 ( ) 1 ( ( ) ( ) ( 31 3 21 2 11 1 1 1

a

a

a

x

e

i f i f i f i f _i − + − + − =

後向きアルゴリズム

∑

∑

=

+

=

=

+ k l l l k kl l i k k k k

b

x

e

a

b

x

e

a

b

a

b

x

P

i

i

L

)

1

(

)

(

)

(

)

1

(

)

(

)

(

)

(

1 0 1 0 

b

_k

(i)=

P(x

_i+1

…x

_L

|π

_i

=k)

をDPにより計算



P

(π

_i

=k|x) =

f

_k

(i)b

_k

(i)/P(x)

1

2

3

1

2

3

a

11

a

12

a

13 )) 1 ( ) ( ) 1 ( ) ( ) 1 ( ) ( ( ) ( 3 1 3 13 2 1 2 12 1 1 1 11 1 + + + + + = + + + i i i i

b

x

e

a

b

x

e

a

b

x

e

a

b

i i i

Viterbi と前向きアルゴリズムの比較



Viterbiアルゴリズム



Forwardアルゴリズム

i i i i n i x

e

a

θ

x

P

_{π π π ,} 1 1

)

|

π

,

(

₋ =

∏

=

{

(

,

π

|

)

}

max

π

P

x

θ

{

(

,

π

|

)

}

π

θ

x

P

∑

1

2

a

12

a

21

A

e

1,A

B

e

1,B

A

e

2,A

C

e

2,C Π for “BACA”= { 1121, 1122, 1221,1222 }

a

11

a

22

HMMに対するEMアルゴリズム

（

Baum-Welchアルゴリズム）

{ }

∑

∑

∑

∑

= + = + = b i j k j k j j i j l j i l kl j k j j x i i f P b k i i f P kl j j

b

x

E

b

x

e

a

x

A

| 1 ) ( ) ( ) ( 1 ) ( ) 1 ( ) ( ) ( ) ( 1 値 から現れる回数の期待 が状態 文字 番目の配列 値 が使われる回数の期待 k b b k j kl

E

x

a

A

kl j : ) ( :

:

∑

∑

=

=

' ' '

(

'

)

)

(

)

(

ˆ

ˆ

b k k k l kl kl kl

_b

b

b

E

E

e

A

A

a

パラメータの更新式

Baum-WelchのEMによる解釈

[

]

∑

∑

∑

∑∑

∑∑

∑

∑

∑∑

∑

∑

∏∏

∏∏

= = = + =       + × = = = = = = = = = = = = ' ' 0 1 1 1 0 1 0 1 ) ( ) , ( 1 , ) ' ( ) ( ) ( log log ) ( log ) ( log ) ( ) ( log ) , ( ) , | ( ) | ( ) | , ( log ) , | ( ) | ( ) ( ) | , ( l kl kl kl b k k k i i i _{i i} M k M l kl kl k M k b k M k M k M l kl kl b k k t π t t π t M k M l π kl π b M k b k

A

A

a

E

E

e

p

q

q

p

a

A

e

E

a

A

e

E

a

e

b b b b b π b π b θ x π P θ θ Q θ π x P θ x π P θ θ Q A E b θ π x P kl k の時、最大より、 は ここで より、 および

配列アラインメント



２個もしくは３個以上の配列の類似性の判定に利用

 ２個の場合：ペアワイズアラインメント  ３個以上の場合：マルチプルアラインメント 

文字間の最適な対応関係

を求める（

最適化問題

）



配列長が同じになるよう、

ギャップ記号

を挿入

HBA_HUMAN VGAHAGEY HBB_HUMAN VNVDEV MYG_PHYCA VEADVAGH GLB5_PETMA VYSTYETA LGB2_LUPLU FNANIPKH GLB1_GLYDI IAGADNGAGV HBA_HUMAN V G A - - H A G E Y HBB_HUMAN V - - - - N V D E V MYG_PHYCA V E A - - D V A G H GLB5_PETMA V Y S - - T Y E T A LGB2_LUPLU F N A - - N I P K H GLB1_GLYDI I A G A D N G A G V

プロファイル

HMM (1)



配列をアラインメントするためのHMM



タンパク質配列分類やドメイン予測

などに有用

 例：ドメインの種類ごとにHMMを作る  PFAM(http://pfam.wustl.edu/) 

一致状態(M)、欠失状態(D)、挿入状態(I)

を持つ

M M M I I I I D D D BEGIN END

プロファイル

HMM (2)

M M M I I I I D D D BEGIN END M M ・ M A A － A － G － G A － A － G A A － － A A － A G － C C － C C － － ・ ・ こうもり ラット ネコ ハエ ヤギ マルチプル アラインメント プロファイル HMM

プロファイル

HMM (3)



各配列ファミリーごとに HMM を作成



スコア最大のHMMのファミリーに属すると予測

HMM1 HMM2 class 1 class 2

Known Seq. (Training Data)

EM EM New Seq. (Test Data) Score= 19.5 Score= 15.8 win ! Viterbi Viterbi

まとめ



配列モチーフ



局所マルチプルアラインメント



Gibbsサンプリング



HMMによる配列解析



最尤推定、ベイズ推定、MAP推定



隠れマルコフモデル

(HMM)



Viterbiアルゴリズム



Baum-Welchアルゴリズム

 EMアルゴリズムに基づく  前向きアルゴリズム、後向きアルゴリズム 

プロファイルHMM