第2回 2019/4/15
- 講師：後藤 修
本講義の概要
- 配列アラインメントとは
- 動的計画法を用いた最良アラインメント計算アルゴリズム
- 大域法と局所法
- アルゴリズムの実装
- 線形メモリーアルゴリズム
- 高速アルゴリズムゲノム解析
- 多重配列アラインメント
- 確率的アラインメント、ペアHMM、プロファイルHMM
配列アラインメントとは
定義
複数\((N\geq 2)\) の文字列 \(\mathbf{s}_1,\mathbf{s}_2,\ldots,\mathbf{s}_N\) が与えられたとき、それぞれの文字列の間で対応する要素（文字、残基）を、適宜空白を入れながら、縦一列に整列させること、あるいは整列させられたもの
目的
- 単独配列（遺伝子、タンパク質）の未知の機能・構造を類似配列の既知情報をもとに予測（2配列アラインメント） keywords: ホモロジー検索、ホモロジーモデリング、遺伝子予測
- 同一ファミリー、オーソログ間で共通する機能・構造単位を予測（多重アラインメント） keywords: プロファイル、ドメイン、モチーフ、系統フットプリント、進化トレース
- 進化系統樹の作成
種類と方法
基本
| 配列数 | 大域的・半大域的 | 局所的 |
|---|---|---|
| 2配列(厳密) | NWG法 | SWG法 |
| 2配列(高速) | PH Cgaln FASTS BLAST |
|
| 多重 | 累進法 反復改善法 整合法 |
期待値最大化Gibbsサンプリング |
応用
| スレディング | スプライスアラインメント |
|---|---|
| 高次構造と配列とのアラインメント | 真核生物ゲノム配列とcDNAやアミノ酸配列とのアラインメント |
動的計画法を用いた最良アラインメント計算アルゴリズム
動的計画法による文字列比較
- 最長部分配列(Longest Common Subsequence)
- 最小編集距離(Minimum Edit Distance)（配列 \(a\) から配列 \(b\) に変換する置換・欠失・挿入の数を最小化）
- 最大一致(Maximal Match) = LCS-Gap Penalty（狭義のNeedleman-Wunsch法）
- 大域的・半大域的(Global, Semi-global)アラインメント（広義のNeedleman-Wunsch法）
- 局所的(Local)アラインメント（Smith-Waterman法）
Minimum Edit Distance
$$
D_{i,j} = \min
\begin{cases}
D_{i-1, j-1}+d\left(a_{i}, b_{j}\right), \\
D_{i-1, j}+1, \\
D_{i, j-1}+1
\end{cases}\\
d\left(a_{i}, b_{j}\right)=a_{i}==b_{j} ? 0 : 1
$$
Longest Common Subsequence
$$
C_{i,j} = \min
\begin{cases}
C_{i-1, j-1}+d\left(a_{i}, b_{j}\right), \\
C_{i-1, j}+1, \\
C_{i, j-1}+1, \\
0
\end{cases}\\
d\left(a_{i}, b_{j}\right)=a_{i}==b_{j} ? 0 : 1
$$
Pairwise Sequence Alignment \(O(N^2)\) with Affine Gap Penalty (\(g(k)=-uk-v\))
$$
\begin{aligned}
H_{i,j} &= \max\left\{H_{i, j}^1,H_{i, j}^2,H_{i, j}^3\right\}\\
H_{i, j}^1 &= H_{i-1, j-1}^1 + s(a_i,b_j)\\
H_{i, j}^2 &= \max\left\{H_{i-1,j}+g(1), H_{i-1, j}^2-u\right\}\\
H_{i, j}^3 &= \max\left\{H_{i,j-1}+g(1), H_{i, j-1}^3-u\right\}
\end{aligned}
$$
大域法と局所法
局所的最良アラインメント法
Smith-Waterman法 = NW法 + 最大部分配列問題（偏った組成の領域を探す）
最大部分列問題
以下の sum を考える。
Kadane-Gries \(O(n)\) Algorithm:
F[0] = H[0] = 0; // 初期化
for (i=1; i<=N; ++i){
H[i] = max{H[i-1] + S[i], 0}; // 0より小さくなったらその次から再スタート
F[i] = max{F[i-1], H[i]}; // 最大値を記録
}Smith-Waterman-Gotoh Algorithm (\(g(k)=-uk-v\))
Recursion Eq.
$$
\begin{aligned}
P_{i,j} &= \max\left\{H_{i-1,j}-v,P_{i-1, j}\right\}-u\\
Q_{i,j} &= \max\left\{H_{i,j-1}-v,Q_{i, j-1}\right\}-u\\
H_{i,j} &= \max\left\{H_{i-1,j-1}+s(a_i,b_j),P_{i,j},Q_{i,j},0\right\}-u\\
\end{aligned}
$$
Initial Condition
$$
\begin{aligned}
& H_{i,0}=H_{0,j}=0;\ P_{i,0} = P_{0,j} = Q_{i,0} = Q_{0,j} = -\infty\\
& u>0; v>0;\\
& s(a_i,b_j)\geq0\ \mathrm{if}\ a_i == b_i\\
& s(a_i,b_j)\leq0\ \mathrm{otherwise}
\end{aligned}
$$
$$\mathrm{Score} = H_{i^*-1,j^*-1} = \underset{0\leq i\leq m, 1\leq j\leq n}{\max}\left\{H_{i,j}\right\}$$
編集距離アルゴリズムの実装
漸化式(recurrent method)
for (i=0; i<m; ++i){
for (j=0; j<n; ++j){
d[i,j] = min{d[i-1,j-1] + s(a[i],b[j]),
d[i-1,j],
d[i,j-1]}
}
}再帰式(recursive call)
double D(i,j){
if(i==0) return(j);
if(j==0) return(i);
return (min{D(i-1,j-1) + s(a[i],b[i]),
D(i-1,j),
D(i, j-1)}
)
}上の式の場合、\(D(i,j)(i\neq0)\) を計算する時に、再帰式の中で \(D(i-1,j-1)\) を \(3\) 回計算することになるため同じ値の計算が指数関数的に増えてしまい、無駄が多くなっている。
進行方向
for (i=0; i<m; ++i){
for (j=0; j<n; ++j){の方向は、かなり実装が楽。
一方で、斜め方向の場合並列処理が可能になるので、大規模なDP matrixを考えた時は効率が良い。
トレースバック法
DP行列が求められても、どのパスが最良のアラインメントに対応するのかがわからない。そこで、DP行列から最良のアラインメントを求める。
- スコア法
- ビットマップ法
- 連結リスト法
- 線形メモリー法
線形メモリー法
上図のように終点からたどる方法もあるが、その場合2次元配列がメモリとして必要になってしまう。そこで、ここでは一次元配列のみを利用してトレースバックを行う手法を紹介する。
アルゴリズムは簡単で、逆向きからも同じようにDP行列を作成し、最も値が大きいものを繋げば良い、というだけである。
A problem with Affine Gap penalties
実はアフィンギャップの場合はギャップコストが異なる関係でうまくいかない。
この場合は、「ギャップ開始か伸長を記録する」ことで、後で補正すれば問題なく同じアルゴリズムが使える。
配列アラインメントのまとめ
- バイオインフォマティクスにおいては最も実用的な計算機利用のひとつ
- 対象・目的により様々なバリエーション
- 歴史は古いが現在でも盛んに研究されている
- タンパク質・機能RNAの立体構造予測
- タンパク質・機能RNAの活性部位予測
- ゲノム配列上の遺伝子予測
- ゲノム配列上の制御領域予測
- 配列より高次な構造（木、グラフ）への拡張