HaplotypeCaller是怎么计算PL的?(VCF解密)

贝聿铭先生设计的卢浮宫

写在前面

PL是HaplotypeCaller等GATK变异检测软件在sample层面给出的注释,一般记录在VCF文件的FORMAT/sample一栏中。它代表了根据某位点变异情况给出的基因型判定不正确的可能性。

计算方法

计算PL的基本公式如下:
PL=−10∗logP(Data|Genotype)
公式右面的P(Data | Genotype)指的是:根据观察到的数据 D 计算出的基因型为 G 的条件概率(想知道P(D|G)是怎么计算的请戳这儿)。将P(D|G)取log值并乘以-10后将其转换为Phred-scale格式即为PL,然后将所有基因型的PL进行归一化,使得最有可能的基因型PL为0。

举个例子

假如某位点的参考等位基因是A,现在观察到一个read在该位点为T,并且我们现在有HaplotypeCaller根据这个read计算出来的各基因型的条件概率P(D|G)(当然如果有多个read,它们的贡献将会相加):

Alleles
Reference: A
Read: T

Conditional probabilities calculated by HC
P(AA | Data) = 0.000001
P(AT | Data) = 0.000100
P(TT | Data) = 0.010000

计算初始PL值:

我们要分别计算基因型为0/0, 0/1 和 1/1时的初始PL值,根据前述的公式,结果如下:

Genotype A/A A/T T/T
Raw PL -10 * log(0.000001) = 60 -10 * log(0.000100) = 40 -10 * log(0.010000) = 20

可以发现P(D|G)最低的基因型在转换为PL后变为最大值。这是在我们的意料之中的,因为PL指的就是这个基因型是不正确的概率。一个基因型的Raw PL值越小,代表越有可能是真实的基因型。

标准化:

在将PL值写入VCF之前,还要对Raw PL进行一个小小的计算:对所有的PL值进行标准化,使得最小的PL值为零。很简单,取所有基因型的PL值和最小PL值之间的差值即可。

Genotype A/A A/T T/T
Normalized PL 60 - 20 = 40 40 - 20 = 20 20 - 20 = 0

我们从中可以发现,最终的PL值和原始的P(G|D)之间的关系:我们取的原始P(G|D)之间依次相差100倍,最终的PL之间最终相差20,也就是100取Phred-scale后的值。这样我们扫一眼VCF里面PL值的大小,就可以方便的比较各个基因型之间可能性的差异了。

PL和GQ之间的联系和区别

GQ的值其实就是PL次小值和最小值之间的差异,由于PL的最小值总是0,所以GQ就是PL的次小值。在上面的例子中,GQ = 20 - 0 = 20。需要注意的是,为了节省计算空间,在GATK中,GQ值最大为99,就算实际计算出的GQ大于99,VCF中也只会记录为99哦。

参考资料

https://software.broadinstitute.org/gatk/documentation/article?id=5913

最后编辑于
©著作权归作者所有,转载或内容合作请联系作者
  • 序言:七十年代末,一起剥皮案震惊了整个滨河市,随后出现的几起案子,更是在滨河造成了极大的恐慌,老刑警刘岩,带你破解...
    沈念sama阅读 217,509评论 6 504
  • 序言:滨河连续发生了三起死亡事件,死亡现场离奇诡异,居然都是意外死亡,警方通过查阅死者的电脑和手机,发现死者居然都...
    沈念sama阅读 92,806评论 3 394
  • 文/潘晓璐 我一进店门,熙熙楼的掌柜王于贵愁眉苦脸地迎上来,“玉大人,你说我怎么就摊上这事。” “怎么了?”我有些...
    开封第一讲书人阅读 163,875评论 0 354
  • 文/不坏的土叔 我叫张陵,是天一观的道长。 经常有香客问我,道长,这世上最难降的妖魔是什么? 我笑而不...
    开封第一讲书人阅读 58,441评论 1 293
  • 正文 为了忘掉前任,我火速办了婚礼,结果婚礼上,老公的妹妹穿的比我还像新娘。我一直安慰自己,他们只是感情好,可当我...
    茶点故事阅读 67,488评论 6 392
  • 文/花漫 我一把揭开白布。 她就那样静静地躺着,像睡着了一般。 火红的嫁衣衬着肌肤如雪。 梳的纹丝不乱的头发上,一...
    开封第一讲书人阅读 51,365评论 1 302
  • 那天,我揣着相机与录音,去河边找鬼。 笑死,一个胖子当着我的面吹牛,可吹牛的内容都是我干的。 我是一名探鬼主播,决...
    沈念sama阅读 40,190评论 3 418
  • 文/苍兰香墨 我猛地睁开眼,长吁一口气:“原来是场噩梦啊……” “哼!你这毒妇竟也来了?” 一声冷哼从身侧响起,我...
    开封第一讲书人阅读 39,062评论 0 276
  • 序言:老挝万荣一对情侣失踪,失踪者是张志新(化名)和其女友刘颖,没想到半个月后,有当地人在树林里发现了一具尸体,经...
    沈念sama阅读 45,500评论 1 314
  • 正文 独居荒郊野岭守林人离奇死亡,尸身上长有42处带血的脓包…… 初始之章·张勋 以下内容为张勋视角 年9月15日...
    茶点故事阅读 37,706评论 3 335
  • 正文 我和宋清朗相恋三年,在试婚纱的时候发现自己被绿了。 大学时的朋友给我发了我未婚夫和他白月光在一起吃饭的照片。...
    茶点故事阅读 39,834评论 1 347
  • 序言:一个原本活蹦乱跳的男人离奇死亡,死状恐怖,灵堂内的尸体忽然破棺而出,到底是诈尸还是另有隐情,我是刑警宁泽,带...
    沈念sama阅读 35,559评论 5 345
  • 正文 年R本政府宣布,位于F岛的核电站,受9级特大地震影响,放射性物质发生泄漏。R本人自食恶果不足惜,却给世界环境...
    茶点故事阅读 41,167评论 3 328
  • 文/蒙蒙 一、第九天 我趴在偏房一处隐蔽的房顶上张望。 院中可真热闹,春花似锦、人声如沸。这庄子的主人今日做“春日...
    开封第一讲书人阅读 31,779评论 0 22
  • 文/苍兰香墨 我抬头看了看天上的太阳。三九已至,却和暖如春,着一层夹袄步出监牢的瞬间,已是汗流浃背。 一阵脚步声响...
    开封第一讲书人阅读 32,912评论 1 269
  • 我被黑心中介骗来泰国打工, 没想到刚下飞机就差点儿被人妖公主榨干…… 1. 我叫王不留,地道东北人。 一个月前我还...
    沈念sama阅读 47,958评论 2 370
  • 正文 我出身青楼,却偏偏与公主长得像,于是被迫代替她去往敌国和亲。 传闻我的和亲对象是个残疾皇子,可洞房花烛夜当晚...
    茶点故事阅读 44,779评论 2 354