置換スコア行列 | アラインメントから置換スコア行列を作る方法

塩基配列に着目したとき、ある塩基が他の塩基に置換する事象の起こりやすさを、行列にまとめたものを置換スコア行列という。このスコアは定数として与えても良く、比較的に保存されているマルチプルアラインメントから確率を計算して与えてもよい。前者の場合、例えば、ある塩基が他の塩基に置換しなかったときのスコアを α、他の塩基に置換したときのスコアを β とおくと、置換スコア行列は次のように書くことができる。

\[ \mathbf{S} = \begin{pmatrix} \alpha & \beta & \beta & \beta \\ \beta & \alpha & \beta & \beta \\ \beta & \beta & \alpha & \beta \\ \beta & \beta & \beta & \alpha \\ \end{pmatrix} \]

このとき α および β は何らかの基準を用いて与える必要がある。blastn のデフォルトオプションを確認すると reward が 2 で、penalty が -3 となっている。それぞれが α = 2 および β = -3 に対応している。

置換速度行列

このように定数を使った置換スコア行列は多様にあり、何らかのマルチプルアラインメントから各塩基の置換速度を計算して、その速度を行列の要素とすることもできる。この場合、この行列を置換速度行列などと呼ばれたりする。置換速度行列には Junks-Cantor (JC) モデル、

\[ \mathbf{S} = \begin{pmatrix} -\alpha & \frac{\alpha}{3} & \frac{\alpha}{3} & \frac{\alpha}{3} \\ \frac{\alpha}{3} & -\alpha & \frac{\alpha}{3} & \frac{\alpha}{3} \\ \frac{\alpha}{3} & \frac{\alpha}{3} & -\alpha & \frac{\alpha}{3} \\ \frac{\alpha}{3} & \frac{\alpha}{3} & \frac{\alpha}{3} & -\alpha \\ \end{pmatrix} \]

Kimura モデル

\[ \mathbf{S} = \begin{pmatrix} -2\beta -\alpha & \beta & \beta & \alpha \\ \beta & -2\beta -\alpha & \alpha & \beta \\ \beta & \alpha & -2\beta -\alpha & \beta \\ \alpha & \beta & \beta & -2\beta -\alpha \\ \end{pmatrix} \]

などがある。

対数オッズスコア

置換スコア行列の各要素を対数オッズスコアで埋める方法もある。この対数オッズスコアはアラインメントから計算する。アライメントの取り方は、目的に応じて近縁種の塩基配列やアミノ酸配列を用いたり、あるいは遠縁種のそれを用いたりする。アライメントから置換スコアを求めるために 2 つの確率モデルを仮定する。それぞれをランダムモデル（R）と一致モデル（M）という。これを説明するために、（塩基配列などを表す）文字列 x と文字列 y を考え、両者からなるペアワイズアライメントを L とする。

ランダムモデル（R）において、文字 a が出現する確率を q_a と仮定する。このとき、2 つの文字列 x と y を考えたとき、文字列 x に a が出現する確率は q_a であり、文字列 y に文字 b が出現する確率は q_b である。これに対して、一致モデルにおいて、文字列 x の文字 a の位置に文字列 y の文字 b がアラインメントされる確率を p_ab と表す。このとき、文字 a から文字 b に置換する確率は次のように求められる。

\[ \frac{p_{ab}}{q_{a}q_{b}} \]

この確率を対数化したものが対数オッズスコアである。すなわち、文字 a から文字 b に置換するときのスコア s(a, b) は次のように計算される。

\[s(a,b) = \log \frac{p_{ab}}{q_{a}q_{b}} \]

この s(a, b) の値で置換スコア行列を埋めることができる。この方法で作られる置換スコア行列には PAM スコア行列や BLOSUM スコア行列などがある。

PAM スコア行列

PAM スコア行列は、近縁でかつよく似たタンパク質のアミノ酸配列からスコアを計算して作られたものである。手順としては、まずこれらのタンパク質のアミノ酸配列から系統樹を作成し祖先配列を計算し、続けて系統樹の枝を辿って置換回数を集計し、文字 a が文字 b に置換する変異確率 M_ab を計算する。ただし、マルチプルアラインメントから計算される（一致モデルから計算される） a と b の対応が生じた確率を q(a,b) と表し、偶然に（ランダムモデルから計算される） a が生じた確率を p(a) と表す。

\[ M_{ab} = \frac{q(a, b)}{p(a)} \]

次に、100 個のアミノ酸あたりに 1 個が置換するように、M_ab を調整する。この M_ab を a から b への変異スコアとして表した行列 (M)_ab が PAM スコア行列の基本となる。

また、行列 (M)_ab は 100 個のアミノ酸あたりに 1 個が置換するスコア（確率）なので、これを単位時間あたりの進化確率と見なすことができる。そのため、行列 (M)_ab を n 乗することで、遠縁の進化を表すことができるようになる。このように、一般化された PAM 行列の各成分 s(a, b) は最終的には次のように計算できる。

\[ s(a,b) = \log\frac{q(a,b)}{p(a)p(b)} = \log\frac{\left(M^{N}\right)_{ab}}{p(b)} \]

よく使われる n には 30, 70, 250 などがある。それぞれを PAM30、PAM70、PAM250 などのようによぶ。

BLOSUM スコア行列

BLOSUM はマルチプルアラインメントデータベースである BLOCKS に登録された配列から算出されるスコア行列である。まず、BLOCKS にあるマルチプルアラインメントを配列一致率でクラスタリングし、一致率が閾値以上ならば一つのクラスタにまとめる操作を行う。次に、各クラスタ間の置換数を集計し、確率 q(a, b) を推定する。これを利用してスコア行列を計算する。

\[ s(a,b) = \log\frac{q(a,b)}{p(a)p(b)} = \log\frac{\left(M^{N}\right)_{ab}}{p(b)} \]

クラスタリングする際の配列一致率の閾値として 45%、62%、80% などが用いられる。それぞれの閾値から計算されるスコア行列を BLOSUM45、BLOSUM62、BLOSUM80 などとよぶ。この数値は配列一致率であるので、この値が大きいほど、近縁の進化を表していることになる。