3. Short Read Mapping の正確さ:ぱらぱらめくる『Deep Sequencing Data Analysis』

  • マッピングの基礎知識
    • マッピングでは読まれたリードの配列がreference配列のある部分の合致する部分を探して、その位置を決める
    • 長ければ長いほど配列の特異性は高く、特異的にマッピングされる(ゲノム上に1箇所だけが合致箇所と想定される)、短ければその逆。また、読まれたリードの全塩基が揃って合致しないといけないかと言えば、配列多様性があればもちろん、そうはならない。
    • この辺りの、「合致の特異度の強さ」「塩基違い・塩基のin/delの許容度」などはマッピングツールによって異なり、また、1つのマッピングツールでもその評価点のつけ方・閾値の取り方を変えられる。
    • 従って、同一の実験結果でも、結果は色々になる。
    • 言い換えると、どのツールを使っても、そのツールの実行オプションをどう取っても、その結果には正しくマッピングされたもの(True positive)、誤ってマッピングされたもの(False positive)、誤ってマッピングされなかったもの(False negative)、正しくマッピングされなかったもの(True positive:読んだ配列に実験エラーが入っているため合致する箇所がない)があり、その内訳は異なる。
    • また、マッピング箇所は1箇所とは限らず、複数箇所にそれなりの可能性でマッピングされることもある
    • このように結果は複数の意味で確率的
  • マッピング結果を理解するための確率的基礎知識
    • リード長が短ければ、ゲノム上には複数の同一配列が実際に存在する
    • リード長が短ければ、ランダムにマッチすることもある
      • 長さLの塩基配列が無関係な配列にk箇所のミスマッチ(L-k箇所のマッチ)をする確率は\frac{L!}{k!(L-k)!} (\frac{1}{4})^{L-k} (\frac{3}{4})^k
      • これを、L=1,2,...,k=0,1,...,Lで絵にしてみる
      • Lのうちの\frac{1}{4}が一致する場合が多く、それより少ない一致も、それより多い一致も起きにくい

L <- 50
X <- matrix(0,L,L+1)
p <- 1/4
q <- 1-p

for(i in 1:L){
	X[i,1:(i+1)] <- dbinom(0:i,i,q)
}
persp(1:L,0:L,X,phi=50,theta=60)
      • では、ヒトゲノム全長(30億塩基対のペア)を考えると、上記の計算の確率が低くても、数ウチャ当たるので、それなりの割合がたまたまマッピングされうる。その確率は、リード全長のうち一部の塩基が違ってもマッピングされたことにする、という条件にすれば、さらに上がる
  • リードの個々のベースコールのクオリティ・スコア
    • Quality scoreという、各塩基のベースコールの実験シグナル的な信頼度を表した値がある。サンガー法時代のPhred quality score相当値にし、さらにその小数点以下…なる実数的な値を整数に丸めたもの(参考:こちら)
    • コールの正しさの確率がわかる。その確率eでコールされた塩基であり、他の3塩基の確率は1-eを3等分したもの(3等分せずにもう少し工夫することもある)
  • マッピングにあたってクオリティ・スコアで重みづけをする
    • クオリティスコアの高さに応じて、リードのコールとreference配列の塩基との不一致が「本当に不一致」である確率が変わるので、それを総合的に値にして、reference配列上のマッピング候補部位に点数をつける
    • PSSM(Position-Specific Scoring Matrix)。reference配列上の塩基がAのときに、あるサンプルのその位置の塩基がAであるか、それ以外かは読んでみないとわからないが、配列多様性(SNPとか)の程度などの情報を使えば、確率が出る。その確率と、ベースコールのクオリティ・スコアがあると、「そのベースコールがGだったとして、reference配列の塩基がAだったときに、サンプルの塩基がA,T,G,Cである事後確率」が計算できる。(リードが(A,T,G,C)だったときそれぞれについて、そこの本当の塩基が(A,T,G,C)である確率なので4x4行列になる。それをPSSMという
  • マッピングのスコア
    • PSSMがあると、「reference配列がこれこれで与えられたときに、サンプルの本当の配列はわかっていないけれど、確率的にはこんなもの(SNPとかがあるけれど、その情報を確率的に組み込んだ)となっていると想定できる。そのうえで、リードシークエンスが、その仮想的な『サンプルの本当の配列』に照らして、どこの位置に当てると、『うまく当たっている確率』はいくつになるか」という値が計算できるので、ゲノム上の位置が点数化される
  • まとめると
    • 基本的には条件(reference配列、サンプルの配列がreference配列とどのような似かより方か、リードの実験的精度はどうか)に照らして、ゲノム上の位置のどこに何点がつくか、と考える
    • それを実装するに当たっては、どれくらい確率的にやるか、とか、ミスマッチの個数を積み上げることにする、とか、実装ごとのやり方が異なってくる