九州大学大学院 情報学専攻特別講義 （２） 配列アラインメント

Presentation on theme: "九州大学大学院 情報学専攻特別講義 （２） 配列アラインメント"— Presentation transcript:

1 九州大学大学院 情報学専攻特別講義 （２） 配列アラインメント
九州大学大学院　情報学専攻特別講義 （２） 配列アラインメント 阿久津　達也 京都大学　化学研究所 バイオインフォマティクスセンター

2 講義目的、参考書 講義目的 バイオインフォマティクスにおける先端的なアルゴリズムについて理解する
計算時間および解の最適性もしくは近似精度に理論的保証のあるアルゴリズムを主対象とする 参考書 阿久津達也：バイオインフォマティクスの数理とアルゴリズム、共立出版、2007 Anthony: Discerete Mathematics on Neural Networks, SIAM, 2001　（第８回のみ） その他は講義ノートにおいて参考文献を記載

3 講義内容 バイオインフォマティクス概論（資料なし） 配列アラインメント 配列解析 RNA二次構造予測 タンパク質立体構造の比較と予測
固定パラメータアルゴリズムと部分k木 グラフの比較と列挙 ニューラルネットワークの離散モデル ブーリアンネットワークの解析と制御 講義の進展状況によっては内容に変更の可能性あり

4 配列アラインメント

5 配列検索 バイオインフォマティクスにおける基本原理 配列検索の利用法 配列が似ていれば機能も似ている ただし、例外はある
実験を行い機能未知の配列が見つかった データベース中で類似の配列を検索 機能既知の類似の配列が見つかれば、その配列と似た機能を持つと推定

6 配列アラインメント バイオインフォマティクスの最重要技術の一つ ２個もしくは３個以上の配列の類似性の判定に利用
文字間の最適な対応関係を求める（最適化問題） 配列長を同じにするように、ギャップ記号（挿入、欠失に対応）を挿入 ２個の配列に対するアラインメント：　　　　ペアワイズ・アラインメント ３個以上の配列に対するアラインメント：　マルチプル・アラインメント

7 … ペアワイズ・アラインメント ２本の配列に対するアラインメント
大域アラインメント： 配列全体にわたるアラインメント 列ごとにスコアが定義され、各列のスコアの和が最大となる最適アラインメントを計算 入力配列 ACGT ATCCT アラインメント … A C G T ー A ー C G T A C ー G T A T C C T A T C C T －6 1 －1 スコア スコアの定義 同じ文字： １ 違う文字： －１ ギャップ： －１

8 大域アラインメントと格子状グラフ C A G T 入力文字列から格子状グラフを構成 アラインメントと左上から右下へのパスが一対一対応
（縦横の辺の重みはすべてギャップスコア、斜めの辺は対応する文字間のスコア） アラインメントと左上から右下へのパスが一対一対応 最長経路＝最適アラインメント C A G T -1 1 ー 最適アラインメント 非最適アラインメント

9 動的計画法による最適アラインメントの計算
アラインメントの個数：指数関数のオーダー 動的計画法を用いれば O(mn) 時間 D[i,j] は始点(0,0)から(i,j)までの最適パスのスコア アラインメントの復元（トレースバック） D[m,n] から再帰式で＝となっている頂点を逆にたどる D[i,j] D[i-1,j] D[i-1,j-1] D[i,j-1] -d w(s[i],t[j]) s,t: 入力配列 　s[i]: 配列 s の i 番目の文字　 m=|s|, n=|t| w(x,y): 文字 x,y 間のスコア　 -d: ギャップ記号のスコア

10 局所アラインメント

11 局所アラインメント というアラインメントを計算 配列の一部のみ共通部分があることが多い ⇒共通部分のみのアラインメント 問題の定義
　　⇒共通部分のみのアラインメント 例えば、AATGCATT　と　GATCG　の場合、 　　　　　 　　　A T G C 　　　　　 　　　A T － 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)時間

12 局所アラインメントに対する動的計画法 大域アラインメントに対する動的計画法を少し修正するだけでOK

13 局所アラインメント・アルゴリズムの正当性
証明のアイデア 始点と終点を表す２個の頂点を格子状グラフに追加 始点から終点へのパスと局所アラインメントが１対１対応

14 配列検索の実用的プログラム

15 配列検索の実用プログラム（１） O(mn): m は数百だが、n は数ＧＢにもなる ⇒実用的アルゴリズムの開発
　⇒実用的アルゴリズムの開発 FASTA: 短い配列（アミノ酸の場合、1,2文字、DNAの場合、4-6文字）の完全一致をもとに対角線を検索し、さらにそれを両側に伸長し、最後にＤＰを利用。 BLAST: 固定長（アミノ酸では3, DNAでは11）の全ての類似単語のリストを生成し、ある閾値以上の単語ペアを探し、それをもとに両側に伸長させる。ギャップは入らない。伸長の際に統計的有意性を利用。

16 配列検索の実用プログラム　（２） FASTA: 短い配列（アミノ酸の場合、1,2文字、DNAの場合、4-6文字）の完全一致をもとに対角線を検索し、さらにそれを両側に伸長し、最後にＤＰを利用。 BLAST: 固定長（アミノ酸では3, DNAでは11）の全ての類似単語のリストを生成し、ある閾値以上の単語ペアを探し、それをもとに両側に伸長。

17 配列検索の実用プログラム（３） SSEARCH: 局所アラインメント＋アフィンギャップコスト（Smith-Waterman-Gotohアルゴリズム）をそのまま実行 PSI-BLAST: ギャップを扱えるように拡張したBLASTを繰り返し実行。「BLASTで見つかった配列からプロファイルを作り、それをもとに検索」という作業を繰り返す。

18 マルチプル・アラインメント

19 マルチプル・アラインメント：定式化 S(mi) = -∑cia log pia （cia= i列におけるaの出現回数,
３本以上の配列が与えられた時、長さが同じで、かつ、スコアが最適となるように各配列にギャップを挿入したもの スコアづけ　（全体スコアは基本的に各列のスコアの和:∑S(mi)） 最小エントロピースコア S(mi) = -∑cia log pia　　　（cia= i列におけるaの出現回数, pia = i列におけるaの生起確率） SPスコア(Sum-of-Pairs) S(mi)=∑k＜lw(mk[i],ml[i]) 　　　　　　　　　　　　　　　　　　（ mk[i]= アラインメント後のi列, k行目の文字）

20 多次元DPによるマルチプル・アラインメント
N個の配列に対するマルチプル・アラインメント N次元DPによりO(2NnN)時間 （各配列の長さはO(n)を仮定） 例：N=3 一般の N に対しては NP困難 (i,j,k) (i-1,j,k) (i,j-1,k) (i,j-1,k-1)

21 マルチプル・アラインメントの実用的計算手法
プログレッシブ・アラインメント CLUSTAL-W（広く利用されているソフト）などで採用 逐次改善法との組み合わせが、より有効 逐次改善法 シミュレーテッドアニーリング 遺伝的アルゴリズム HMMによるアラインメント 分枝限定法 １０配列程度なら最適解が計算可能

22 近似アルゴリズム NP困難問題への対処法 近似アルゴリズム 固定パラメータアルゴリズム 指数時間アルゴリズム（O(an)で底aを小さくする）
平均的に高速なアルゴリズム ヒューリスティックなアルゴリズム 最適解との比率の最悪の場合の上限を理論的に保証 最適解がわからないのにどうやって保証するか？ 最適解の下限を理論的に見積もる（最小化問題の場合） 近似解の上限を理論的に見積もる 「近似解の上限／最適解の下限」が比率になる

23 近似アルゴリズム： センタースター法 アルゴリズム ∑k≠iS(sk,si)が最小となる sk を計算
スコアに三角不等式を仮定し、最小化問題として定義 SPスコアを使用 アイデア： 中心となる配列を定め、それと各配列とのアラインメント結果をまとめる アルゴリズム ∑k≠iS(sk,si)が最小となる sk を計算 S(sk,si)をもとにギャップを適切に挿入し、マルチプルアラインメントA を構成 図中の xi はsi を表す

24 センタースター法の解析 定理 A のSPスコア ≦ (2-2/N)・最適解の SPスコア 証明 N・Starのスコア ≦2・最適解のスコア
　 (図2でx1に（=s1）接続しない辺のスコアの合計はスター(N-2)個分以下) 図中の xi はsi を表す

25 アラインメントに関するまとめ ペアワイズ・アラインメント マルチプル・アラインメント 補足
動的計画法で O(n2) 時間、線形領域も可能 局所アラインメント、線形ギャップスコアでも同様 マルチプル・アラインメント NP困難だが、距離で定義した場合、2近似が可能 補足 ペアワイズ・アラインメントは O(n2/log n) 時間で可能 [Crochemore et al.: Proc. SODA 2002]。でも、SETHのもとで、　 O(n2-ε）時間は不可能。 様々なギャップスコアや疎行列の場合の動的計画法についても多くの研究 [Galil, Park: Theoret. Comp. Sci. 1992] マルチプル・アラインメントの近似率は 2-K/N まで改善（K は 任意の定数） [Bafna et al.: Theoret. Comp. Sci. 1997] N に関係なく 2 より良くできるかは研究課題

26 編集距離計算量の下限

27 編集距離 入力： 出力： x を y に変換するのに必要な最小の編集操作数 （ペアワイズの全域配列アラインメントと本質的に等価）
編集操作: (A) 文字の置換 (B) 1文字挿入 (C) 1文字削除 ED(x,y): 文字列 x, y 間の編集距離 x=ATGAT A T G C ー ED(x,y)=3 y=ACTCT 挿入 置換 削除 定理：編集距離はSETHのもとでO(n2-ε)時間で計算不可 [Backurs & Indyk, STOC 2015] PT(x,y): 文字列 x と、 y の連続部分文字列との最小編集距離 （近似文字列マッチングの誤差と本質的に等価）

28 SAT(充足可能性問題) 用語 ＳＡＴ k-SAT: すべての節が k 個以下のリテラルの論理和
リテラル： 変数もしくは変数の否定(NOT) 節： リテラルの論理和(OR) 充足可能： 節(集合)を真とする変数への 0-1 割り当てが存在 ＳＡＴ 入力： 節集合　S 問題： S は充足可能(すべての節を真とする割り当てが存在する)か? k-SAT: すべての節が k 個以下のリテラルの論理和

29 SETH（Strong Exponential Time Hypothesis）
SATに関する結果（乱拓アルゴリズムを含む） 3-SAT: O(20.388n) [Hertli, SCICOMP, 2014] 4-SAT: O(20.555n) [Hertli, SCICOMP, 2014] k-SAT： 2(1-1/O(k))n時間 [Paturi et al., JACM, 2005] 一般の場合： 2(1-1/O(log(m/n)))n 時間 （m：節の個数）[Calabro et al., Proc. CCC, 2006] SETH n 変数のSATは、どの定数 ε>0 に関しても、2(1-ε)n時間では解けない（という仮説） 直交ベクトル問題 入力：　d次元0-1ベクトルの集合A, B　(|A|=|B|=N) 出力: x・y=0（内積が0）となる (x,y)∊A×B はあるか？ SETHが成立⇒この問題は O(N2-εd) 時間では解けない

30 SATから直交ベクトル問題への還元 解析 SATの変数を XA={x1,…, xn/2}と XB={xn/2+1,…,xn}に２分割
XA（resp., XB）へのすべての0-1割り当てについて この割り当てで充足される節を０、充足されない節を１として、それらを並べた0-1ベクトルを作成し、 そのベクトルの集合をA（resp., B）とする (|A|=|B|=N=2n/2) 解析

31 証明の方針 直交ベクトル問題を編集距離問題に還元 直交ベクトルが存在 ⇔ P1 が P2 中で閾値以下でマッチ
直交ベクトルが存在 ⇔ 編集距離が閾値以下 A, B の要素（d次元0-1ベクトル）を文字列に変換 P1 = Aから変換した文字列をつなげたもの P2 = Bから変換した文字列をつなげたもの＋α 直交ベクトルが存在 ⇔ P1 が P2　中で閾値以下でマッチ 直交ベクトル（から変換した文字列）でちょうど重なるように、ずらしてマッチ

32 証明の概略（１）： 0-1の文字列化 Aのベクトル中の0,1 ⇒ CG1(0), CG1(1), g
Bのベクトル中の0,1 ⇒ CG2(0), CG2(1) 1-1 間の距離は大きく（3l0)、他のビット間の距離は小さい（l0 ） g と、B の0,1 間の距離はほんの少し大きい（ l0 +1）

33 証明の概略（２）： ベクトルの文字列化 a∊A, b∊B に対し、VG1(a), VG2(b) を以下のように構成
証明の概略（２）：　ベクトルの文字列化 a∊A, b∊B に対し、VG1(a), VG2(b) を以下のように構成 a,bが直交すれば、ED(VG1(a),VG2(b))=l+2l2+dl0 :=Es 直交しなければ、 ED(VG1(a),VG2(b))=Es+d :=Eu 　　（直交すれば CG どうしがマッチ、それ以外は g…g とマッチ）

34 証明の概略（３）：ベクトル集合の文字列化 A, B に対し、P1, P2 を下図のように構成 （|A|≦|B|）
直交する (a,b) があれば、PT(P1,P2)≦(|A|-1) Eu + Es なければ、 PT(P1,P2) = |A| Eu ⇒ P’1=3 |P’2| P13 |P’2| , P’2= P2 として ED(P’1,P’2) を計算

35 編集距離計算下限に関するまとめ 編集距離計算、ペアワイズアラインメント 補足：多くの結果（下限）が前後して出現
SETHが成立 ⇒ （任意のε>0に対し）O(n2-ε） 時間では不可能 補足：多くの結果（下限）が前後して出現 局所アラインメント: O(n2-ε) [Abboud et al., ICALP 2014] Dynamic Time Warping Distance: O(n2-ε ) [Abboud et al., FOCS 2015] Longest Common Subsequence (for k strings): O(nk-ε ) [上と同じ] 部分木同型性： O(n2-ε） [Abboud et al., SODA 2016] 文脈自由文法解析, RNA二次構造予測： O(nω-ε ) [Abboud et al., FOCS 2015] クリーク問題を還元（ωは行列乗算に関する係数） 確率文脈自由文法解析： O(n3-ε ) [Saha, FOCS 2015] （all pairs）最短経路問題を還元 順序木の編集距離計算: O(n3-ε ) [Bringmann et al., SODA 2018] （all pairs）最短経路問題、および、クリーク問題を還元 SETH を（妥当な仮定のもとで）証明する、もしくは、否定するのは重要な研究課題