E値とは何か？ ~BLASTの統計学~

 BLASTを用いたアライメント結果はどの程度信頼できるのか？という抽象的な問いに答

える一つの指標がE値です．その正体については良く理解されないまま都合の良い閾値と して扱われていることが多いのですが，E値はつまるところ一般的な統計分野で使われて いる”期待値”と意味は同じですので難しい事はありません．そこでBLASTでのE値の扱わ れ方を今一度クリアに理解するために改めて次の具体例について考えてみましょう．!

!

1. まず一本のクエリ配列を用意する．!

2. 次に10,000本の配列を含むデータベースを用意する．ここでデータベースはクエリ配

列のホモログ配列(=クエリ配列と進化的な類縁関係のある配列)一本と類縁のない配

列9,999本からなるものとする．!

3. そしてクエリ配列とデータベースで総当たりで1対1のアライメントを行なってそれ ぞれスコアを算出する．!

4. クエリ配列と類縁関係のない9,999本の配列の内，ホモログ配列とクエリ配列のアラ イメントスコアに到達した配列の数をカウント．これがE値が見つもっている数．!

5. その数を9,999で割った値はP値に相当．E値が十分小さいときはE値 P値．!

!

こうして文にしてみると何の事はないですね．しかしこのイメージをキチンと持っておく と，配列解析についてより実感の伴った議論ができると思います．例えば，E=0.05は多く の統計解析において十分な有意水準として用いられていますが，配列解析では”ギリギリ” な値である事が分かると思います．類縁関係の有無を議論できる水準は局所アライメント ではE=10^-10，大域アライメントでE=10^-100程度であると言われています．!

!

＜今更聞けないここだけの話＞!

  ところで本題に進む前の念のための確認なのですが，よくE値は1E-10とか1E-100とか いうように出力されますが，1E-10の中にEが出てくるからE値というのでは決してありま せん！”そんなの当たり前だろ”，と思った方は以下読み飛ばして下さって結構です．”まじ で?!”，と驚かれた方は続きをそっと読んでおきましょう．!

!

  まず1E-10のEはE値のEではありません．ネイピア数とかでもありません．e^-10と1E-10 では意味が全く違います．このEは”科学表記法”というもので例えば1E-10は10の-10乗を 意味する表記法です．つまり!

です．文献によってはEの代わりにD等も使われていたりしますが意味は同じです．某

◯ahoo!◯恵袋でEは自然対数の底ですよ〜とかいう回答がありましたが真っ赤な嘘です．


e

¹⁰

=

_eeeeeeeeee¹

1E 10 = 0.0000000001

  さて，BLASTの出力ファイルを見てみると，E値の他にS値やビットスコア等が目に入 ると思います．これらの値は密接に関係していて，以下のようにして互いに変換する事が できます．これらのスコア体系について以下，簡単に確認していきましょう．!

!

 

!

E値とP値の導出 

「アライメント結果からどのようにしてE値を導出するのか」という問題はすなわち「ア ライメントスコアをいかにして確率の世界に引き込むのか」という問題に書き換える事が できます．まずアライメント結果を順に一残基ずつ見ていきましょう．するとそれは”当た り”か”当たりでない”かの二通りの試行の連続として捉える事ができます．それをコインの 表と裏で例えるとちょうどコイン投げの問題に帰着できることになります．つまりコイン を投げて”表”が連続して出る最長の長さの期待値はいくらか？という問題に置き換える事 ができるということです．!

!

1. コイン投げ!

  この問題については都合の良い事にエルデシュ(エルデシュ数のエルデシュです)とレ ニーによって1970年にすでに解かれていて，表が出る確率がpで試行回数がn回の時，表 が連続する長さの最大値をRnとすると，!

となることが知られています．これは直観的には以下のようにして理解する事ができるで しょう．まず連続する表の長さが1になる回数を考えます．するとpの確率で表がでますの でn回コイン投げしたときに表になっている所はだいたいnpあると見積もる事ができま す．あとは同様に表が連続する長さが2の場合はnp²，長さが3はnp³...となって長さiはnpⁱ カ所ある，と見積もることができます．ここで最長の連続は一カ所しか無いと仮定すると，

その時の長さをlとした場合np^l=1と書く事ができますので，それを解くと，!

R

_n

log

¹

p

(n) 1

Fig.23-BLASTスコアの体系(Korf I et al., “BLAST”より引用)

エルデシュの式が示されました．もちろん厳密な証明はきちんとしなくてはいけません が，イメージはこんなもんでしょう．!

!

2. アライメントへの応用!

 さて，このコイン投げの式を使うと長さm，nの二つの配列のアライメントにおいて一 致する最長の長さMは以下のようにモデル化する事ができます．!

実際には，天然配列の近隣残基は完全に独立しているわけではないのでその期待値と分散 は以下のような近似式になることがArratiaやWatermanによって示されています．!

ここでγ=0.577...(オイラーの数)，q=1-pです．ここは頑張って計算したらとこうなったと 捉えておいてください(詳しくはAttatiaら(1986)をどうぞ)．するとE(M)は!

のようにまとめることができます．ここでλ=loge(1/p)，Kは雑多な定数を今一度置き直し たもので，両方ともデータベースと塩基組成に依存する決まる定数であると言えます．最 後に一致する長さMはアライメントスコアSと置き換える事ができるので，!

となります．これがアライメントスコアの期待値の式です．!

np

^l

= 1 p

^l

= 1

n log

_p

p

^l

= log

_p

1

n l = log

¹

p

n

M = log

¹

p

(mn)

E (M ) = log

¹

p

(mn) + log

¹

p

(q) + log

¹

p

(e) 1 2 V ar(M (m, n)) log

¹

p

(e)

6 + 1

12

E (M ) log

_e

(Kmn)

E (S ) = log

_e

(Kmn)

3.ポアソン分布を使ったP値の算出!

  類縁配列のスコアをxとしてSがそれを超える確率を考えます．ここで類縁の無い配列の スコアがxを超える事はめったにない(x>>E(S))と仮定します．するとこの事象はxを平均 としたポアソン分布で近似できますので，Sがxをn回超える確率は，!

で与えられます．従ってSがxを一度も超えない確率は!

となりますので，Sがxを超える確率は!

で求める事ができます．ここで”Sがxを超えない”を”S-E(S)がある定数cを超えない”と言い 換えると，上記の式は以下のように書き直す事ができます．!

x=c+E(S)としてこれを解いていきますと，!

これでP値を求める事ができました．なお，この式変形の途中で!

の関係式を使っています．さらに!

P

n

=

^e _n!^xxn

P

₀

= e

^x

1 P

₀

= 1 e

^x

P (S E(S ) c) = 1 e

^pc

P (S x) = 1 e

^{p(x E(S))}

= 1 e

^px

ln(Kmn)

= 1 e

^{p1 ( x ln(Kmn))}

= 1 e

^p

1

(

^{ln(Kmn)e x}

)

= 1 e

^{Kmne x}

= ln 1 p e = 1

p

P = e

ですので，x>2の範囲では!

がP値の良い近似式となります．さらに!

なるS’を定義すれば上記の式は!

まで簡略化できます．なお，!

のようにS’の底を2に変換すると以後の統計解析に使いやすいという事で，このように定 義されたビットスコアはE値と並んでアライメントの有意性を評価するのによく使われて います．!

!

4. 指数極値分布を使ったE値の算出!

  ちなみにポアソン分布の裏返しである指数分布を使うと今度はE値を計算する事ができ ます．詳細は上記の一連の式の説明とかぶりますので主に式変形だけ追っていきますと，

まずSの分布は以下のように指数分布で近似する事ができます．!

するとSがxを超える確率は以下のような指数分布の定積分で求める事ができます．!

! !

x

lim 1 e

^{e x}

e

^x

P (S x) Kmne

^x

S = S ln Kmn

P (S x) e

^x

Bit score = S ln(2)

f (S ) = e

^S

P (S x) =

x

f (S )

=

x

e

^S

dS

= e

^x

と，このようにしてアライメントスコアxの探索空間KmnにおけるE値を求める事ができま した．面白い事にP値の近似式とE値は同じ式の形をしていますね．!

という関係式がこれまでの議論をまとめると導けますが，実際確かめてみると分かる通 り，xが十分大きければ，すなわちP(x)とE(x)が十分小さければ!

となります．!

!

5. ガンベル極値分布を使ったP値の算出!

  さて，P値の導出方法はもう一つあります．これから説明する方法とポアソン分布を 使った方法を比べる事でアライメントの意味がよく見えてきます．その戦略自体は簡単で!

!

1. まずランダムな配列同士のアライメントスコアを適切な確率分布に当てはめる．! 2. 次にある類縁配列同士のアライメントスコアαを算出する．!

3. そして当てはめた確率分布について-∞からαまで積分することでその類縁配列同士の アライメントの有意性を求める(つまりランダムな配列のスコアが類縁配列同士のス コアを超える確率Pを求める)．!

!

ということでランダム試行の結果を確率分布に当てはめれば良いのですが，実はこの分布 は正規分布には従わない事が知られています．なぜなら配列アライメントとは「与えられ た二つの配列についてその無数の並べ方パターンの中から最適なものを採用し，それ以外 を棄却する作業」だからです．もし全ての並べ方パターンからランダムに値を抽出して分 布を描くならば，それは当然正規分布になるのですが，このように最大値を選んで抽出し た場合はガンベル極値分布に従うことが知られています．


E (x) = Kmne

^S

P (x) = 1 e

^E(x)

P (x) E (x)

ガンベル極値分布Yevは以下の式で与えられます．!

上で述べた通り，P(S<x)はYevについて-∞からxまでを積分する事で求める事ができるので!

となり，P(S≧x)は1からP(S<x)を引いて!

で求められます．ここでポアソン分布の時と同様にxを以下のようにして標準化します．!

ここでuはSの期待値，すなわちu=E(S)です．代入して整理すると，!

P値が算出されました．ポアソン分布の時 と一致している事を確認しましょう．!

! !

Y

_ev

= e

^{x e x}

P (S < x) = e

^{e x}

P (S x) = 1 e

^{e x}

X = (x u)

P (S X ) = 1 e

^{e (x u)}

= 1 e

^{e (x}

loge(Kmn))

= 1 e

^{e x+loge(Kmn)}

= 1 e

^{Kmne x}

Fig.24-正規分布とガンベル極値分布の比較．黒:正 規分布(μ=0,σ²=1)，赤:ガンベル極値分布(μ=0,η=1)

ドキュメント内 春プロ勉強会14.pages (ページ 67-74)