变异探测算法深度解析

一、somatic与germline区别

image

胚系突变(Germline Mutation)和体细胞突变(Somatic Mutation)在WGS、WES、Gene Panel检测时常常遇到,二者最大的区别是胚系突变本质上是遗传性,可以遗传给后代,而体细胞突变可能对****细胞造成损害,癌变或细胞死亡,****突变不能遗传给后代

Germline与Somatic 生信检测方法

Germline突变频率主要集中在50%(杂合) , 100%(纯合),一般是比对错误或者测序错误等原因产出噪音,干扰变异探测。我们可以使用统计学分析或者机器学习等方法,对潜在变异位点进行区分,得到突变位点。

Somatic一般使用对照样本进行分析, 例如,取肿瘤细胞作为Tumor样本,取血液中的白细胞或者癌旁正常细胞作为对照样本来分析,Tumor中特有的突变为Somatic mutation 。

二、变异检测方法

image

变异检测方法:

  1. 测序原始数据经过清洗后的fastq文件, 比对到参考基因组上,获得BAM文件。

  2. 根据比对质量,过滤低比对质量的数据。

  3. 提取每个位点的等位基因。潜在SNV位点,出现某些等位基因与参考等位基因不同, 如图,第5、23和28位。这些等位基因可能是SNVs或测序错误

  4. 对于每一个候选SNV,可以应用不同的统计方法来确定该基因座中是否存在SNV

  5. 将探测到的snv存储在VCF文件中

image

当前探测SNV的统计方法:

  1. 通过计算等位基因(counting alleles)

  2. 通过二项分布( binomial distribution)

  3. 通过泊松二项分布(Poisson-binomial distribution)

  4. 通过贝叶斯方法(Bayesian)

1. 基于等位基因的变异探测方法

基本方法是基于对等位基因的计数。

设D0是D中可信度高的基的子集。通常我们过滤碱基的质量分数≥20。因此,D0 = {bi∈D | qi≥20}。然后,在D0的所有底数中计算每个等位基因出现的次数。

  • 如果D0参考等位基因的比例低于θlow(一般,20%),它被称为纯合非参序列等位基因(homozygous non-reference allele) ;

  • 如果D0参考等位基因的比例高于θhigh(一般,80%),它被称为纯合子参考序列等位基因( homozygous reference allele);

  • 否则,它被称为杂合基因型

image

有三个非参考序列等位基因至少出现一次。

对于位置5,T出现的时间少于20%。我们预测该基因型为AA(非参考序列等位基因纯合子)。对于位置23,参考基A出现的频率超过80%。我们预测该基因型为AA(参考序列等位基因纯合子)。

对于28号位点,75%的reads包含参考基t。我们预测该位点的基因型为GT(杂合子位点)。

这种方法被用于许多商业软件程序,包括Roches GSMapper, CLC基因组工作台和DNSTAR Lasergene。测序深度高(> 20×)时效果较好。

然而,这种方法缺陷也比较明显:

  1. 简单的质量过滤可能会导致信息丢失。

  2. 这种方法不能提供不确定性的度量。

  3. 这种方法可能会低估杂合基因型

  4. 不能给出p值

2. 基于二项分布的变异探测方法

令D = {b1,…, bn}是覆盖特定位置的一组碱基。

设随机变量X为n个碱基中,突变的个数。Prn (X = k)为观察D中k个突变的概率。

假设D中有K个非参考序列碱基,假设n个碱基的排序误差是独立的。当序列误差概率p已知时(比如p = 0.01), X服从二项分布。然后,我们有

image

三条read覆盖位置 j,位置 j 的三个碱基为D = {A, G, A}。
image

注意,有两个基是非参考序列碱基。观察两个非参考变量的p值为

image

该方法虽然确定了不确定性概率,但没有利用每个基的质量分数

3. 基于泊松分布的变异探测方法

二项分布假定对同一位置上的每个碱基的测序错误率是相同的。然而,不同碱基的测序错误率实际上是不同的。每个碱基的测序错误率可以通过PHRED质量评分来估计。

设随机变量X为非参考序列碱基个数,碱基总数为n 。表示P rn(X = k) D = {b1中k个变量的概率,…,在零模型下为bn}。

我们将二项分布推广到一个泊松二项分布,其中不同的基的序列误差概率是不同的。然后,我们有

image

同样,实际计算如下

image
image

4. 基于贝叶斯的变异探测方法

D代表观测数据(即,特定位点的碱基), G代表位点的基因型 ,有10个可能的基因型:AA, CC, GG, TT, AC, AG, AT, CG, CT, GT 。D = {b1,…, bd}和G基因型A1A2。我们的目标是计算Pr(G|D) ,然后,我们目标是使Pr(G|D)最大化的基因型G。

根据贝叶斯

image

Pr(D|G)是后验概率,Pr(G)先验概率 。

** Pr(D|G) 后验概率计算方法**

由于碱基来自不同的read,所以read base是独立的

image

假设G=A1A2, Pr(bi|G) 通过下边公式计算

image

ei为错误概率,根据碱基质量分数得出的。(参见本公众号前期文章, 《测序数据质控报告分析》文章,有详细介绍)

image

Pr(G)先验概率计算方法

G有10种可能的基因型。先验概率Pr(G)受其为纯合子参考型、杂合子或纯合子非参考基因型的身份的影响。

设r为参考序列碱基,s为替代等位基因。

通常设置

Homozygous SNP rate = altHOM = 0.0005

Heterozygous SNP rate = altHET = 0.001

(例如, r=G and s=A.)

image

许多方法利用额外的生物信息来提高Pr(G)的估计。例如,我们可以使用已知的数据dbSNP来得出 。

计算步骤

Pr(b1=A|AG)=1/2(Pr(b1=A|A)+Pr(b1=A|G))=1/2((1-10-2)+10-2/3)=0.49667

Pr(b2=G|AG)=1/2(Pr(b2=G|A)+Pr(b2=G|G))=1/2(10-1/3+(1-10-1))=0.466667

Pr(b3=A|AG)=1/2(Pr(b3=A|A)+Pr(b3=A|G))=1/2((1-10-5)+10-5/3)=0.499997

Pr(D|AG) = 0.496670.4666670.499997 = 0.115888

Pr(AG|D) = Pr(D|AG)*Pr(AG)=0.000116

因此,我们预测该基因型为AG。

image

基因组分析 微信公众号推出 《50篇文章深入理解NGS》系列文章, 第三篇文章 《变异探测算法深度解析》,争取每周更新一篇高质量生信干货帖子。

请点击 关注微信公众号 ,**转发 **给同学和同事,您的认可,是对我最大的支持 ,任何问题,后台可以留言。

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

推荐阅读更多精彩内容