情報生命科学特別講義

情報生命科学特別講義III
(5)配列アラインメント
阿久津 達也
京都大学 化学研究所
バイオインフォマティクスセンター
講義予定















第1回: 文字列マッチング
第2回: 文字列データ構造
第3回: たたみ込みとハッシュに基づくマッチング
第4回: 近似文字列マッチング
第5回: 配列アラインメント
第6回: 配列解析
第7回: 進化系統樹推定
第8回: 木構造の比較:順序木
第9回: 木構造の比較:無順序木
第10回: 文法圧縮
第11回: RNA二次構造予測
第12回: タンパク質立体構造の予測と比較
第13回: 固定パラメータアルゴリズムと部分k木
第14回: グラフの比較と列挙
第15回: まとめ
配列アラインメントとは?
配列検索

バイオインフォマティク
スにおける基本原理



配列が似ていれば機能
も似ている
機能未知の配列
VLPIKSKLP......
配列検索
配列データベース
ただし、例外はある
配列検索の利用法



実験を行い機能未知の配列
が見つかった
データベース中で類似の配列
を検索
機能既知の類似の配列が見
つかれば、その配列と似た機
能を持つと推定
DFECILTSKLG.....
.
ACILTSTRE......
VLPIKSDLP......
HPFACILPDEL......
類似配列
VLPIKSDLP......
配列アラインメント
• バイオインフォマティクスの
最重要技術の一つ
• 2個もしくは3個以上の配
列の類似性の判定に利用
• 文字間の最適な対応関係
を求める(最適化問題)
• 配列長を同じにするように
、ギャップ記号(挿入、欠失
に対応)を挿入
A L G F G S L Y G
A L G G V S V G
A L G F G
A L G
S L Y G
G V S V
2個の配列に対するアラインメント:
ペアワイズ・アラインメント
3個以上の配列に対するアラインメント: マルチプル・アラインメント
G
ペアワイズ・アラインメント
ペアワイズ・アラインメント
• 2本の配列に対するアラインメント
• 大域アラインメント: 配列全体にわたるアラインメント
• 列ごとにスコアが定義され、各列のスコアの和が最
大となる最適アラインメントを計算
入力配列
ACGT
ATCCT
アラインメント
A C G T ー ー
A ー C G T
A C ー G T
ー A T C C T
A T C C T
A T C C T
スコア
-6
スコアの定義
1
同じ文字: 1
違う文字: -1
-1
ギャップ: -1
…
大域アラインメントと格子状グラフ



入力文字列から格子状グラフを構成
アライメントと左上から右下へのパスが一対一対応
最長経路=最適アラインメント
A
A
1
C
G
-1
-1
T
-1
最適アラインメント
-1
A ー C G T
-1
1
T
-1
C
-1
1
-1
-1
-1
C
-1
1
-1
-1
-1
T
-1
-1
-1
-1
-1
-1
-1
1
-1
-1
-1
A T C C T
非最適アラインメント
A C G T ー ー
ー A T C C T
アラインメントの個数
• 漸化式による解析
• アラインメントの個数=格子状グラフにおけるパス数
m  0 or n  0
 1,
A(m, n)  
 A(m  1, n  1)  A(m, n  1)  A(m  1, n), m  0 and n  0
• 組み合わせ的解析
• ギャップ文字が現れない列の個数を k とすると
• アラインメント全体の長さは m+n-k
• (si,-) という型の列の個数は m-k
• (-,tj) という型の列の個数は n-k
• よって、アラインメントの個数は m+n-k 個の要素をk個、m-k個、n-k個に
分割する組み合わせの数
A(m, n) 
min(m, n )

k 0
( m  n  k )!
k !( m  k )!( n  k )!
• Delannoy数として知られている
入力配列: s, t (|s|=m, |t|=n)
動的計画法による最適アラインメントの計算
• アラインメントの個数:指数関数のオーダー
• 動的計画法を用いれば O(mn) 時間
• D[i,j] は始点(0,0)から(i,j)までの最適パスの長さ
• アラインメントの復元(トレースバック)
• D[m,n] から再帰式で=となっている頂点を逆にたどる
D[i,0]  i  ( d )
D[0, j ]  j  ( d )
i  0,  , m
j  0,  , n
D[i  1, j ]  d


D[i, j ]  max
D[i, j  1]  d
 D[i  1, j  1]  w( s[i ], t[ j ])

D[i-1,j-1]
w(s[i],t[j])
D[i-1,j]
-d
D[i,j-1]
-d
D[i,j]
スコア行列
• 残基間(アミノ酸文字間)の類似性を表す行列
– PAM250, BLOSUM45 など
A
A
R
N
D
C
Q
E
G
H
I
L
K
M
F
P
S
T
W
Y
V
5
-2
-1
-2
-1
-1
R
N
D
C
Q
E
G
H
I
L
K
M
F
P
T
W
Y
V
-2 -1 -2 -1 -1 -1 0 -2 -1 -2 -1 -1 -3 -1 1 0
7 -1 -2 -4 1 0 -3 0 -4 -3 3 -2 -3 -3 -1 -1
-1 7 2 -2 0 0 0 1 -3 -4 0 -2 -4 -2 1 0
-2 2 8 -4 0 2 -1 -1 -4 -4 -1 -4 -5 -1 0 -1
-4 -2 -4 13 -3 -3 -3 -3 -2 -2 -3 -2 -2 -4 -1 -1
1 0 0 -3 7 2 -2 1 -3 -2 2 0 -4 -1 0 -1
-3
-3
-4
-5
-5
-1
-2
-1
-2
-3
-3
-1
0
3
-3
-4
-1
-3
BLOSUM50 スコア行列
(置換行列)の一部分
S
局所アラインメント
局所アラインメント

配列の一部のみ共通部分があることが多い
⇒共通部分のみのアラインメント


というアラインメントを計算
問題の定義



例えば、AATGCATE と GATCG の場合、
ATG C
AT-C
入力: 2個の配列 s, t スコア関数 w(x,y)
出力: Sopt(s[h…k],t[h’…k’]) が最大となる部分文字列の
組(s[h…k],t[h’…k’])に対する最適アラインメント
大域アラインメントを繰り返すとO(m3n3)時間
⇒Smith-WatermanアルゴリズムならO(mn)時間
局所アラインメントに対する動的計画法

大域アラインメントに対する動的計画法を少し修正
するだけでOK
D[i,0]  0
D[0, j ]  0
i  0, , m
j  0,, n
0
 D[i  1, j ]  d

D[i, j ]  max
 D[i, j  1]  d
 D[i  1, j  1]  w( xi , y j )
maxD[i, j ]}
局所アラインメント・アルゴリズムの正当性

証明のアイデア
 始点と終点を表す2個の頂点を格子状グラフに追加
 始点から終点へのパスと局所アラインメントが1対1対応
0
 D[i  1, j ]  d

D[i, j ]  max
 D[i, j  1]  d
 D[i  1, j  1]  w( xi , y j )
maxD[i, j ]}
0
0
0
(一部の辺は
省略)
ギャップコスト
ギャップペナルティ
ギャップ (g =3)

線形コスト




-gd



L Y G
A V G V S D L
g: ギャップ長
d: ギャップペナルティ
この図の例では、コスト= -3d
アフィンギャップコスト

A L G
G
–d – e(g-1)
d: ギャップ開始ペナルティ
e: ギャップ伸張ペナルティ
この図の例では、コスト= -d - 2e
よく利用されるペナルティ (d,e)=(12,2),(11,1)
アフィンギャップコストによるアラインメント
 D[i  1, j ]  d
DX [i, j ]  max
 DX [i  1, j ]  e
 D[i, j  1]  d
DY [i, j ]  max
 DY [i, j  1]  e
DX [i, j ]


D[i, j ]  max
DY [i, j ]
 D[i  1, j  1]  w( s[i ],t[ j ])

三種類の行列を用いる動的
計画法によりO(mn)時間
 Smith-Watermanアルゴリズ
ムとの組み合わせが広く利
用されている
⇒ Smith-Waterman-Gotoh
アルゴリズム

任意ギャップコストによるアラインメント

動的計画法(下式)により、O(n 3 )時間
(ただし、m=O(n)とする)

 D[i  1, j  1]  w( s[i ], t[ j ])

D[i, j ]  max  max D[k , j ]   (i  k )
 k 0,...,i 1
 max D[i, k ]   ( j  k )
 k 0,...j 1
ギャップコストと計算時間の関係




線形: O(n2)時間
アフィン:O(n2)時間
凸: O(n2α(n))時間
任意: O(n3)時間
線形
O(n 2 )時間
2
凸 O(n α(n))時間
α(n)はアッカーマン関数の逆関数
(実用的には定数と同じ)
アフィン
任意
O(n 2 )時間
3
O(n )時間
線形領域アラインメント
線形領域アラインメント



スコアの計算だけならO(m+n)領
域で簡単に可能
トレースバックが難しい
アイデア: 分割統治法を利用
計算時間
C×mn
LinearSpaceAlign(i1 , i2 , j1 , j2 )
1. i2  i1が定数以下なら通常の DPで計算・
出力して終了
C×(mn/2)
2. ih  (i1  i2 ) / 2  i1とする
3. (i1 , j1 )から (ih , j )までのの最適スコア D1[ j ]、
(ih , j )から (i2 , j2 )までの最適スコア D2 [ j ]を
通常のDPで計算
4. D1[ j ]  D2 [ j ]が最大となる jを jhとする
5. LinearSpaceAlign(i1 , ih , j1 , jh )を実行
6. (ih , jh )を出力
7. LinearSpaceAlign(i1 , ih , jh , j2 )を実行
C×(mn/4)
配列検索の実用的アルゴリズム
配列検索の実用的アルゴリズム
• データベース検索: O(mn): m は数百だが、n は数GBにもなる
⇒実用的アルゴリズムの開発
• FASTA: 短い配列の完全一致(アミノ酸の場合、1,2文字、DNAの場合、46文字)をもとに対角線を検索し、さらにそれを両側に伸長し、最後にDPを
利用
• BLAST: 固定長(アミノ酸では3, DNAでは11)の全ての類似単語のリストを
生成し、ある閾値以上の単語ペアを探し、それをもとに両側に伸長させる。
ギャップは入らない。伸長の際に統計的有意性を利用
• SSEARCH: 局所アラインメント(Smith-Waterman-Gotohアルゴリズム)を
そのまま実行
• PSI-BLAST: ギャップを扱えるように拡張したBLASTを繰り返し実行。「
BLASTで見つかった配列からプロファイルを作り、それをもとに検索」という
作業を繰り返す
• PatternHunter: 穴あきシードを用いる(連続した文字ではなく飛び飛びの文
字の完全一致をもとに検索)
マルチプル・アラインメント
マルチプル・アラインメント:定式化

3本以上の配列が与えられた時、長さが同じで、かつ、スコアが最適となるよ
うに各配列にギャップを挿入したもの
HBA_HUMAN
HBB_HUMAN
MYG_PHYCA
GLB5_PETMA
LGB2_LUPLU
GLB1_GLYDI

VGAHAGEY
VNVDEV
VEADVAGH
VYSTYETA
FNANIPKH
IAGADNGAGV
HBA_HUMAN
HBB_HUMAN
MYG_PHYCA
GLB5_PETMA
LGB2_LUPLU
GLB1_GLYDI
V
V
V
V
F
I
G
E
Y
N
A
A
A
S
A
G
A
D
H
N
D
T
N
N
A
V
V
Y
I
G
G
D
A
E
P
A
E
E
G
T
K
G
Y
V
H
A
H
V
スコアづけ (全体スコアは基本的に各列のスコアの和:∑S(mi))
 最小エントロピースコア


(cia= i列におけるaの出現回数,
pia = i列におけるaの生起確率)
S(mi) = -∑cia log pia
SPスコア(Sum-of-Pairs)

S(mi)=∑k<lw(mk[i],ml[i])
( mk[i]= アラインメント後のi列, k行目の文字)
多次元DPによるマルチプル・アラインメント

N個の配列に対するマルチプ
ル・アラインメント
(i,j-1,k-1)
N次元DPによりO(2NnN)時間
(各配列の長さはO(n)を仮定)

例:N=3
 D[i  1, j  1, k  1] 

w( s1[i ], s2 [ j ],s 3 [k ])

 D[i, j  1, k  1]  w(, s2 [ j ],s 3 [k ])

 D[i  1, j , k  1]  w( s1[i ], s 3 [k ])
D[i, j , k ]  max 
 D[i  1, j  1, k ]  w( s1[i ], s2 [ j ],)
 D[i, j , k  1]  w(,, s 3 [k ])

 D[i, j  1, k ]  w(, s2 [ j ],)
 D[i  1, j , k ]  w( s [i ],,)
1


一般の N に対しては NP困難
(i,j-1,k)
(i-1,j,k)
(i,j,k)
マルチプル・アラインメントの実用的計算手法

プログレッシブ・アラインメント







CLUSTAL-W(広く利用されているソフト)などで採用
逐次改善法との組み合わせが、より有効
逐次改善法
シミュレーテッドアニーリング
遺伝的アルゴリズム
HMMによるアラインメント
分枝限定法

10配列程度なら最適解が計算可能
近似アルゴリズム

NP困難問題への対処法






近似アルゴリズム
固定パラメータアルゴリズム
指数時間アルゴリズム(O(an)で底aを小さくする)
平均的に高速なアルゴリズム
ヒューリスティックなアルゴリズム
近似アルゴリズム


最適解との比率の最悪の場合の上限を理論的に保証
最適解がわからないのにどうやって保証するか?



最適解の下限を理論的に見積もる(最小化問題の場合)
近似解の上限を理論的に見積もる
「近似解の上限/最適解の下限」が比率になる
近似アルゴリズム




スコアに三角不等式を仮定
し、最小化問題として定義
SPスコアを使用
アイデア: 中心となる配列を
定め、それと各配列とのアラ
インメント結果をまとめる
アルゴリズム
1. ∑k≠iS(sk,si)が最小とな
る sk を計算
2. S(sk,si)をもとにギャッ
プを適切に挿入し、マ
ルチプルアラインメン
トA を構成
x1
x2
STAR-1
STAR-2
STAR-3
STAR-4
x3
ACGT
A-GT
A-CGT
AACGT
ACG-T
A-GGT
x4
A-CG-T
A--G-T
AACG-T
A--GGT
図中の xi はsi を表す
近似アルゴリズムの解析
定理
 A のSPスコア ≦
(2-2/N)・最適解の
SPスコア
証明
 N・Starのスコア ≦2・最
適解のスコア
(図1でNスターにより、各辺
を2回ずつカバー)
 A のスコア ≦
(N-1)・Starのスコア
(図2でx1に(=s1)接続しない
辺のスコアの合計はスター
(N-2)個分以下)
図1
x1
x2
x3
x4
図2
三角不等式
x1
x2
xk
xi
x3
xj
x4
図中の xi はsi を表す
まとめ

ペアワイズ・アラインメント



動的計画法で O(n2) 時間、線形領域も可能
局所アラインメント、線形ギャップスコアでも同様
マルチプル・アラインメント

NP困難だが、距離で定義した場合、2近似が可能
補足
 ペアワイズ・アラインメントは O(n2/log n) 時間で可能 [Crochemore
et al.: Proc. SODA 2002]


様々なギャップスコアや疎行列の場合の動的計画法について
も多くの研究 [Galil, Park: Theoret. Comp. Sci. 1992]
マルチプル・アラインメントの近似率は 2-K/N まで改善(K は
任意の定数) [Bafna et al.: Theoret. Comp. Sci. 1997]

N に関係なく 2 より良くできるかは研究課題