読者です 読者をやめる 読者になる 読者になる

アメリエフのブログ

Amelieff Staff Blog

VCF format③

どうもakbです。
VCF formatの連載三日目です。
今日は、FORMATの"PL"について説明しようと思います。

まず、前回の復習ですが、FORMATの"GT"はGenotypeの略で
0:reference 1:alternate(sample1)のように数字に意味を持たせており、
0/0,0/1,1/1はそれぞれREFのホモ、REFとALTのヘテロ、ALTのホモを示します。

まず図をご覧ください。
(例)

--------------------------------------------
REF = A ALT = G GT:PL 0/1:191,0,255
--------------------------------------------
"PL"に対応する値は"191,0,255"のようにカンマ区切りで示されています。
なぜ三つあるかというと"REF = A ALT = G"から考えられるGenotypeごとに phred-scaledの尤度計算をしているからです。

つまり、REF = A ALT = G の場合

AA AG GG
191 0 255

この数字は小さいほど信頼性が高くなりますので、
AG = 0が採用され、0/1、REFとALTのヘテロとなるわけです。

では、REF=T  ALT=A,GのようにALTの候補が2つある場合です。

--------------------------------------------
REF=T  ALT=A,G   GT:PL 1/1:178,82,76,104,0,98
--------------------------------------------
この場合、"REF=T  ALT=A,G"から考えられるGenotypeは

TT TA AA TG AG GG (←この順番になります)
178 82 76 104 0 98

AG = 0が採用され、REFとALTのヘテロ、つまりA/Gのように見えます。
しかし、ここは注意が必要です。
なぜならGTは1/1を示しており、ALTのホモを表しているからです。
この場合、ALT同士の組み合わせになっているため
"1/1"という表記になっていますが、ALTのヘテロです。

(2014年11月14日訂正)


三回にわたってVCF formatについて連載してきました。
今回、ご紹介できなかったID等は、下記URLをご参照下さい。

http://samtools.github.io/hts-specs/VCFv4.1.pdf

vcf-variant-call-format-version-41

---------------------------------------------------
(2014年11月14日追記)
【訂正】
上記のGTの説明につきまして、誤りがあるとご指摘をいただきました。
ありがとうございます。
当ブログをご覧のみなさまには、ご迷惑をおかけしまして申し訳ありません。

GTの数字の意味は、以下の通りになります。

0: REFアリル
1: ALTの1つめのアリル
2: ALTの2つめのアリル
:
(以下同)


そのため、GTが1/1ということは、ここではA/Aのホモということになります。

GATKのガイドによりますと、PLは、GTで示したジェノタイプの確率を1.0(PLでは0と示す)とし、それをもとに他のジェノタイプの尤度を相対的に算出しているため、上記の例のようなケース(GTとPLでジェノタイプが食い違うケース)は実際にはないと考えられます。

一例ですが、下記のような組み合わせがあると考えられます。

【ジェノタイプがAAの場合】
--------------------------------------------------------
REF=T  ALT=A,G   GT:PL 1/1:178,82,0,104,76,98
--------------------------------------------------------

【ジェノタイプがAGの場合】
--------------------------------------------------------
REF=T  ALT=A,G   GT:PL 1/2:178,82,76,104,0,98
--------------------------------------------------------