九州大学大学院 情報学専攻特別講義 (4) RNA二次構造予測 阿久津 達也 京都大学 化学研究所 バイオインフォマティクスセンター
講義内容 バイオインフォマティクス概論(資料なし) 配列アラインメント 配列解析 RNA二次構造予測 タンパク質立体構造の比較と予測 固定パラメータアルゴリズムと部分k木 グラフの比較と列挙 ニューラルネットワークの離散モデル ブーリアンネットワークの解析と制御 講義の進展状況によっては内容に変更の可能性あり
RNA二次構造予測
RNA二次構造予測 RNA二次構造予測(基本版) 塩基対 入力: RNA配列 a=a[1]…a[n] 出力: 以下を満たし、スコア 二次構造 U G C 塩基対 二次構造 二次構造でない 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} (i,j) ,(h,k) ∈M となる i ≤h ≤j ≤k が存在しない 塩基対: B={{a,u},{g,c}} スコア関数(最も単純なもの) μ(a[i],a[j])= -1 if {a[i],a[j]} ∈B μ(a[i],a[j])= 0 otherwise スコアが最小でないものも二次構造とよび、最小のものを最適二次構造とよぶこともある {g,u} も塩基対に含まれる場合がある
RNA二次構造の二種類の表現
RNA二次構造の例 RNA配列
Nussinovアルゴリズム
予測アルゴリズム(Nussinovアルゴリズム) 入力配列: a=a[1]…a[n] 動的計画法 初期化 メインループ 最適解(= -塩基対の個数)
アルゴリズムの説明 メインループ
Nussinovアルゴリズムの解析 定理: 上記アルゴリズムは O(n3) 時間で最適解を計算 略証: テーブル W(i,j) のサイズはO(n2)。1個のテーブル要素の計算にO(n)時間(最後の行)。
RNA二次構造予測と 確率文脈自由文法
RNA二次構造予測と確率文脈自由文法 (1) 確率文脈自由文法(SCFG): 導出確率が最大となる構文解析木を計算 ⇒ 確率の代わりにスコアを用いる 文法表現としては 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 の行列乗算にかかる計算時間 高速行列乗算に基づくRNA二次構造予測 基本的に 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) M(n)はn×n行列の乗算の計算量
Valiantアルゴリズムの概略(4) メインルーチン: 下図のとおり 時間計算量:
まとめ RNA二次構造予測 補足 動的計画法により O(n3)時間 Valiant アルゴリズムなどの利用により少しだけ改善 Valiantアルゴリズムは、RNAアラインメント・構造同時予測問題や結合RNA二次構造予測問題にも適用可 [Zakov et al.: Alg. Mol. Biol. 2011] 擬似ノットなしRNA二次構造予測の(クリークによる)下限は Ω(nω-ε) 時間[Abboud et al.: FOCS 2015]。現状では上限 とのギャップ。 ⇒ 肯定的に解決された(O(n2.8…) time) [Bringmann et al.: FOCS 2016]
RNA二次構造予測の O(n2.8…polylog(n))時間アルゴリズム 原論文: [Bringmann et al.: FOCS 2016]
背景 Valiant アルゴリズムにより二次構造は計算可能 ただし、行列乗算において、(×,+)を(+,max)に置き換える必要がある ⇒ 高速乗算が利用できない 一般の行列に対して、(+,max)の行列乗算がO(n3-ε)時間で計算できるかは未解決 一方、行列の要素の絶対値が M以下であれば、 時間で計算できることは既知(#1) [Alon et al.:JCSS, 1997] RNA二次構造予測では行列に bounded difference (BD) という制約(隣接要素の差異が小さい)が存在(#2) #1と#2の性質を利用し、さらに乱拓アルゴリズムを巧みに用いることにより、BD行列(+,max)乗算の 時間を達成 さらに、脱乱拓可が可能( 時間) は を表す ω は通常((×,+))の高速行列乗算の計算量のべき (ω≦2.3729)
アルゴリズムの概略 Phase 1 Phase 2: 以下をρ回繰り返す Phase 3 (C=A・B を計算) (i,k,j)がΔの整数倍の要素についてのみ計算を行い、C=A・B の近似 D を計算 ⇒ すべての要素の誤差が 4ΔW 以内 Phase 2: 以下をρ回繰り返す ランダムに(ir,jr)を選び、絶対値が48ΔW以下の要素からなるAr, Br を作成し、 Ar・Br の結果をもとに D を更新 ( 時間行列乗算を利用) ⇒ 正しく計算されない要素の個数は Phase 3 正しく計算されなかった(i,k,j)(ただし、Δの倍数)に対し、 それを含む近傍(Δ×Δ×Δ)について再計算
Phase 1 I(i’) := {i ∊ [1..n] | i’-Δ < i ≦ i’} 補題1: |Ci,j – Di,j| ≦ 4ΔW n=8, Δ=2
Phase 2: (+,max)行列乗算に関する性質 A’: 行列 Ai,k の各列 k の要素に同じ数 Fk を加算 B’: 行列 Bk,j の各行 k の要素から同じ数 Fk を減算 すると、A’・B’=A・B が成立 Fk をうまく定義することにより、(重要な部分について)要素の絶対値の小さな行列を得る ⇒ 時間行列乗算を適用
Phase 2:アルゴリズム 基本的に が成立するので、 が、min{Ai,k+Bk,j} を与える k について成立すれば より、Ci,j を正しく計算 (少し余裕をみて)どの程度の個数の (i,k,j) についてρ回のうち1回も が成立しないかを見積もる
Phase 2:解析 補題2 (i,k,jr), (ir,k,jr), (ir,k,j) が weakly relevant であれば、 strongly relevant: k は mink を与える weakly relevant: mink となる k の候補の数の上限を与える 補題2 (i,k,jr), (ir,k,jr), (ir,k,j) が weakly relevant であれば、 (i,k,j) は weakly uncovered でない 補題3 weakly relevant で、かつ、weakly covered でない (i,k,j) の個数は、(非常に高い確率で) 証明の方針:補題2をもとに、(kを固定した場合の)weakly relevant な辺 (i,j)からなる二部 グラフにおける長さ4の閉路の個数に関する結果を利用
Phase 3:アルゴリズム 補題3から、 個の (i,k,j) について、再計算すれば十分 ⇒ しかし、そのような (i,k,j) の探索にO(n3)時間かかる そこで、まずΔの倍数の (i,k,j) についてチェックし、該当しそうなものについて、その近傍を再計算 ⇒ 少し条件を変える必要 近傍の再計算 (i,k,j)=(4,6,8) Δ=2
Phase 3:正当性 補題4: 各(i,k,j)について以下が成立 定義2: (i,k,j)∊ I(i’)×I(k’)×I(j’)について 補題4: 各(i,k,j)について以下が成立 strongly relevant ⇒ approximately relevant approximately relevant ⇒ weakly relevant strongly uncovered ⇒ approximately uncovered approximately uncovered ⇒ weakly uncovered よって、Phase 3により、 ・strongly relevant はすべて対応(ただし、既にstrongly uncoveredでないのは対応済み) ・approximately relevant の個数は weakly relevant の個数以下 ・strongly relevant で strongly uncovered であるものはすべて対応 ・(approximately relevant のうち) approximately uncovered の個数は weakly uncovered の個数以下 strongly relevant なものがもれなく計算されていればOK ⇒ Phase 3までで、すべて計算される
時間計算量の解析 Phase 1: Phase 2: Phase 3: + よって、全体で ここで、 とすると、 を得る ここで、 とすると、 を得る W が O(1) であれば、この値は 最悪、O(n3)時間で済むので、平均計算時間も同様 定理 塩基対スコアの和に基づくRNA二次構造予測問題は (高確率および平均で) 時間で解ける
補題3:証明(1) 補題5:G=(U⋃V,E) を |U|=|V|=n, |E|=m である二部グラフとする。m≧2n1.5 であれば、長さ4の閉路数は m4/(32n4) 個以上 補題3の証明 z=c(n2/ρ)ln(n) とする。 Gk を E(Gk)={(i,j) | (i,k,j)は weakly relevant} により初期化。その後、weakly uncovered でなくなったら、その辺を削除。 (i,j) を含む長さ 4 の閉路 i→jr→ir→j→i が存在し、(ir,jr) が選択されれば、(i,j) は削除される z 個以上の閉路に含まれる (i,j) が削除されない確率は よって、c>4 とすれば高確率で、z 個以上の閉路に含まれる (i,j)は、すべて Gk から削除されている
補題3:証明(2) 証明(続き) mk を(最終的な) Gk の辺数とする。 (#)が成立するので、mk≧2n1.5 を満たす k のみを考える そのような各 k においては、1個の辺について z 以下の閉路しかない(それ以上であれば Gk から削除されているはず) よって、 Gk における閉路の総数は mkz 個以下。一方、補題5を考えると (mk/n)4/32 < mkz が成立。よって、以下が成立 よって、これをすべての k について足し合わせ、それに (#) を加えて、以下を得る 補題3 weakly relevant で、かつ、weakly covered でない (i,k,j) の個数は、(非常に高い確率で)
まとめ 行列の隣接値の差が定数以内であれば、(min,+)による行列乗算が(高確率および平均で) 時間で解ける ⇒ Valiant アルゴリズムとの組み合わせでRNA二次構造予測も同様の計算時間 さらに、これらの結果は、Phase 2, Phase 3の改良により 時間まで改良されている また、その脱乱拓化もでき、 時間の決定性アルゴリズムも開発されている クリークによる下限Ω(nω-ε) [Abboud et al.: FOCS 2015]とはギャップが存在。また、通常の行列乗算がO(n2)に近づいたとしても、計算量は と議論されている