情報生命科学特別講義III
(11) RNA二次構造予測
阿久津 達也
京都大学 化学研究所
バイオインフォマティクスセンター
講義予定















第1回: 文字列マッチング
第2回: 文字列データ構造
第3回: たたみ込みとハッシュに基づくマッチング
第4回: 近似文字列マッチング
第5回: 配列アラインメント
第6回: 配列解析
第7回: 進化系統樹推定
第8回: 木構造の比較:順序木
第9回: 木構造の比較:無順序木
第10回: 文法圧縮
第11回: RNA二次構造予測
第12回: タンパク質立体構造の予測と比較
第13回: 固定パラメータアルゴリズムと部分k木
第14回: グラフの比較と列挙
第15回: まとめ
RNA二次構造予測
RNA二次構造予測
RNA二次構造予測(基本版)
入力: RNA配列 a=a[1]…a[n]
出力: 以下を満たし、スコア
Σ(i,j)∈M μ(a[i],a[j]) が最小となる塩基対
の集合 M={(i,j)|1≤i+1<j≤n,{a[i],a[j]}∈B}

塩基対
A
U
G
C
二次構造
(a[i],a[j]) ,(a[h],a[k]) ∈M となる i ≤h ≤j ≤k
が存在しない
A G A G C U


塩基対: B={{a,u},{g,c}}
スコア関数(最も単純なもの)


二次構造でない
μ(a[i],a[j])= -1 if {a[i],a[j]} ∈B
μ(a[i],a[j])= 0 otherwise
A G A G C U
スコアが最小でないものも二次構造とよび、最小
のものを最適二次構造とよぶこともある
{g,u} も塩基対に含まれる場合がある
RNA二次構造の二種類の表現
RNA二次構造の例
RNA配列
Nussinovアルゴリズム
予測アルゴリズム(Nussinovアルゴリズム)


入力配列: a=a[1]…a[n]
動的計画法
初期化
W (i, i)  0, W (i, i  1)  0
W (i  1, j )


W (i, j  1)

W (i, j )  min
W (i  1, j  1)   (a[i], a[ j ])

 min { W (i, k )  W (k  1, j ) }
i k  j 1
メインループ
最適解(= -塩基対の個数)
W (1, n)
 1, x, y {{A, U}, {G, C}}の時
それ以外の時
 0,
 ( x, y)  
W (i, j )  部分配列a[i]a[i  1]...a[ j 1]a[ j ] の最適解の値
アルゴリズムの説明
W (i  1, j )


W (i, j  1)

W (i, j )  min
W (i  1, j  1)   (a[i], a[ j ])

 min { W (i, k )  W (k  1, j ) }
i k  j 1
メインループ
W (i, j )  部分配列a[i]a[i  1]...a[ j 1]a[ j ] の最適解の値
Nussinovアルゴリズムの解析
定理: 上記アルゴリズムは O(n3) 時間で最適解を計算
略証: テーブル W(i,j) のサイズはO(n2)。1個のテーブル要素の計
算にO(n)時間(最後の行)。
W (i  1, j )


W (i, j  1)

W (i, j )  min
W (i  1, j  1)   (a[i], a[ j ])

 min { W (i, k )  W (k  1, j ) }
i k  j 1
RNA二次構造予測と
確率文脈自由文法
RNA二次構造予測と確率文脈自由文法 (1)

確率文脈自由文法(SCFG): 導出確率が最大となる構文解
析木を計算 ⇒ 確率の代わりにスコアを用いる
rule
X→ε
X→a
X→u
X→g
X→c
X→YZ
X→ a Y u
X→ u Y a
X→ g Y c
X→ c Y g
score
0
0
0
0
0
0
1
1
1
1
score for X
0
0
0
0
0
score(X)+score(Y)
score(Y)+1
score(Y)+1
score(Y)+1
score(Y)+1
文法表現としては X→aYu, X→XY などではなく、S→aSu, S→SS などが正式
RNAの場合はスコアは1ではなく、-1に置き換えることが必要
RNA二次構造予測と確率文脈自由文法 (2)


スコア最大(≒確率最大)の構文解析木 ⇔ 最適二次構造
実際、NussinovアルゴリズムはCYKアルゴリズム(文脈自由
文法の構文解析アルゴリズム)に類似
Valiantアルゴリズムの利用
[Akutsu: J. Comb. Opt. 1999] [Zakov et al.: Alg. Mol. Biol. 2011]
二次構造予測の計算量の改良

文脈自由文法の構文解析



高速行列乗算に基づくRNA二次構造予測





高速行列乗算に基づく Valiant アルゴリズムを用いれば O(nω) 時間
O(nω) は n×n の行列乗算にかかる計算時間
基本的に Valiant アルゴリズムを適用
 しかし、行列乗算の基本演算を (+,×) から (max,+) に変える必要
 (max,+) の行列演算は、ほんの少し O(n3) より良くなるだけ
O(n3((log log n)/(log n))1/2)時間 [Akutsu: J. Comb. Opt. 1999]
O(n3(log3(log n))/log n)時間 [Zakov et al.: Alg. Mol. Biol. 2011]
SCFGの内側アルゴリズム[Akutsu: J.Comb.Opt. 1999]、外側アルゴリズム[Zakov
et al.:Alg. Mol. Boil. 2011] は(+,×)演算で済むので O(nω) 時間で可能
Four-Russian アルゴリズムに基づくRNA二次構造予測

O(n3/log n)時間 [Frid, Gusfield: Proc. WABI 2009]
ω は二十数年ぶりに 2.3737 から 2.3736 へ、さらに、2.327 へ改善された [Wiliianms: Proc. STOC 2012]
Valiantアルゴリズムの概略(1)




基本的に分割統治
W(i,j)=Σk W(i,k)×W(k+1,j) の計算(青のベクトルと赤の
ベクトルの乗算)を
高速化
行列乗算を適用する
ため、複数の W(i,j)
の計算をまとめて
実行
左下三角は計算不要
Valiantアルゴリズムの概略(2)
基本戦略
白と黄が計算済みとして、青を計算(青は一部計算済み)
 ピンクの部分の行列積を計算後、青と加算(⇒結果は緑)
 緑と黄色からなる行列を作り、再帰計算により青を計算
(左下の白は不要なのですべて0としてOK)
高速
乗算
再帰計算
Valiantアルゴリズムの概略(3)
時間計算量
T2(n)= T2(n/2)+ 2T3(3n/4)+ T4(n)+ O(n2)
T3(n)= M(n/3)+ T2(2n/3)+ O(n2)
T4(n)= 2M(n/4)+ T2(n/2)+ O(n2)
⇒
T2(n)= 4T2(n/2)+ 4M(n/4)+ O(n2)
log n
T2 (n)  O(n log n)  4M (n / 4)   2( 2 ) k
2
k 0
M(n)はn×n行列の乗算の計算量
 O ( n )
Valiantアルゴリズムの概略(4)
メインルーチン: 下図のとおり
時間計算量: T (n)  2T (n / 2)  T2 (n)  O(n 2 )
log n
T (n)  O (n 2 )  T2 (n)   2  k  2T2 (n)  O (n 2 )
k 0
 O ( n )
二次構造予測の
平均計算時間の改良
[Wexler et al.:J. Comp. Biol. 2007]
平均計算時間の改良: アイデア
最悪の時間計算量の O(n3) からの本質的改良は極
めて困難
 高速行列乗算は実用的でない
 Nussinovアルゴリズムや他の動的計画法アルゴリズ
ムは平均的にも O(n3) 時間かかる
⇒ 平均計算時間の改良 [Wexler et al.:J. Comp. Biol. 2007]
 ブランチループの計算がボトルネックとなっていた



Valiant 型のアルゴリズムでは行列乗算により改良
アイデア: 必要な k のみをチェックすることにより改良
min{ W (i, k )  W (k  1, j ) }
ik  j
詳細なエネルギーモデル(1)


W(i,j): 部分配列 a[i..j] に対する最適解
V(i,j): a[i]とa[j]が塩基対として結合する場合の最適解

W (i, j )  min V (i, j ),W (i  1, j ),W (i, j  1), min{ W (i, k )  W (k  1, j ) }
ik  j
V (i, j )  mineh(i, j ), es(i, j )  V (i  1, j  1),VBI (i, j ),VM (i, j )
VBI (i, j )  min ebi(i, j, i' , j '}  V (i' , j ' )
i i ' j ' j
VM (i, j )  min W (i  1, k )  W (k  1, j  1) b
i  k  j 1

詳細なエネルギーモデル(2)

前述のモデルを j≧i+4 の場合のみを考えて簡略化
W (i, j )  W ' (i, j ) j  i  4
V ' (i, j )  V (i, j )

j  i  4
W ' (i, j )  min V ' (i, j ), minW ' (i, k )  W ' (k  1, j )
ik  j




W ' (i, j )  min V ' (i, j ), min V ' (i, k )  W ' (k  1, j ) 
ik  j


定理: V’(i,j)≧V’(i,k)+W’(k+1,j) がある k (i < k < j)について
成立すれば、すべての j’>j について
V’(i,j)+W’(j+1,j’) ≧ V’(i,k)+W’(k+1,j’) が成立

必要な k のみの計算
定理: V’(i,j)≧V’(i,k)+W’(k+1,j) がある k (i < k < j)について
成立すれば、すべての j’>j について
V’(i,j)+W’(j+1,j’) ≧ V’(i,k)+W’(k+1,j’) が成立
V’(i,j’) の計算においては、V’(i,j)≦V’(i,k)+W’(k+1,j) (i < k < j)
が成立する j について計算すれば良い ( j’>j )
アイデア:塩基対を作ることによりエネルギーが減る場所のみを k の候補とする
アルゴリズム


V’(i,j)<W’(i,j) となる
j のみを以降では k
の候補として採用
ψ(n): 長さ n の配列
に対する L の大きさ
の最大値の期待値
procedureCandidateFold
for i  n to 1 do
L  {};
for j  i to n do
W ' (i, j )  min V ' (i, k )  W ' (k  1, j ) ;
kL
if (V ' (i, j )  W ' (i, j )) th e n
W ' (i, j )  V ' (i, j );
L  L  { j}
定理: CandidateFold は平均的に O(n2 ψ(n)) 時間で動作
妥当な現実的な仮定(polymer-zeta property)のもとで ψ(n) は
定数になることが知られている
⇒ RNA二次構造予測は平均的に O(n2) 時間で実行可能
単純擬似ノットつき
二次構造の予測
[Akutsu: Disc. Appl. Math. 2000]
擬似ノット
擬似ノット
i ≤h ≤j ≤k を満たす塩基対ペア (a[i],a[j]) ,(a[h],a[k]) ∈M
A G A G C U
単純擬似ノット: 下の図に示される擬似ノット(定義は省略)
単純擬似ノットに対する動的計画法 (1)
もとの W(i,j) に加え、3種類のテーブルを用いる
WL(i,j,k): a[i] と a[j] が塩基対をなす場合
WR(i,j,k): a[j] と a[k] が塩基対をなす場合
WM(i,j,k): a[i] , a[j], a[k] のどのペアも塩基対をなさない場合
 WL (i  1, j  1, k )

WL (i, j , k )   (a[i ], a[ j ])  minWM (i  1, j  1, k )
W (i  1, j  1, k )
 R
WL (i, j  1, k  1)

WR (i, j , k )   (a[ j ], a[k ])  minWM (i, j  1, k  1)
W (i, j  1, k  1)
 R
WL (i, j , j )   (a[i ], a[ j ]) for all i  j
WR (i0 , j , k )   (a[ j ], a[ j  1]) for all j
WL (i0  1, j , k )  WL (i0  1, j , k )  WL (i0  1, j , k )
 0 for all otherk  j , k  j  1
WM (i  1, j , k ), WM (i, j  1, k ), WM (i, j , k  1)

WM (i, j , k )  min
WL (i  1, j , k ), WL (i, j  1, k ),

WR (i, j  1, k ), WR (i, j , k  1)

W pseudo (i0 , k0 ) 
min
i0 i  j  k  k0
W pseudo (i, j )


W (i  1, j )


W (i, j  1)
W (i, j )  min
W (i  1, j  1)   (a[i ], a[ j ])

W (i, k )  W (k  1, j )
 min
ik  j
WL (i, j, k ), WM (i, j, k ), WR (i, j, k )
単純擬似ノットに対する動的計画法(2)
WL(i,j,k) 計算の説明
WL (i  1, j  1, k )

WL (i, j , k )   (a[i ], a[ j ])  minWM (i  1, j  1, k )
W (i  1, j  1, k )
 R
より複雑な擬似ノットつき二次構造の予測
木接合文法に基づく方法 [Uemura et al.: Theoret. Comp. Sci. 1999]

O(n4) 時間、 O(n5) 時間(再帰的構造を含む場合)
PKNOTSアルゴリズム [Rivas, Eddy: J. Mol. Biol. 1999]

擬似ノットを組み合わせた構造にも対応、O(n6)時間
平面的擬似ノット:
NP困難 [Akutsu: Disc. Appl. Math. 2000]
まとめ

RNA二次構造予測

動的計画法により O(n3)時間
Valiant アルゴリズムなどの利用により少しだけ改善

Polymer-zeta propertyを仮定すると、平均的にO(n2) 時間


擬似ノットつきRNA二次構造予測

計算量は対象とする擬似ノットの複雑さに依存
補足

Valiantアルゴリズムは、RNAアラインメント・構造同時予測問
題や結合RNA二次構造予測問題にも適用可 [Zakov et al.: Alg. Mol. Biol.
2011]

擬似ノットなしRNA二次構造予測の O(n3-ε) 時間アルゴリズム
の開発は研究課題
ダウンロード

RNA二次構造予測