跳到内容
基因组学研究中心

使用DRAGEN Bio-IT平台从PopGen数据集中准确高效地调用小变量和大变量

Andrew Gross, Sorina Maciuca, Anthony Cox, Duke Tran，邱云江，黄卓毅，Jennifer Del Giudice

分享本文

高通量测序已越来越多地应用于临床实践和群体遗传学(PopGen)研究，将挑战从获得技术转移到以高效而准确的方式从生成的数据中提取最多的信息。凭借其全面的产品，该公司Illumina DRAGEN(基因组学动态读取分析)生物it平台能够以高精度和高速调用小型和大型变体，从而能够在一个平台内高效地交付有意义的见解。

纽约基因组中心发布了1000个基因组计划的高覆盖率全基因组测序数据1允许我们将这些不同的呼叫者置于一个有代表性的队列中。这使得我们不仅可以从许多不同的样本中观察变量调用，而且还可以深入到覆盖数据不均匀或违反变量调用者假设的区域。

在这里，我们在1000个基因组计划的数据上部署了DRAGEN平台，并且:

  1. 识别小的和大的变体，聚合和公开可用。2
  2. 通过例子，展示了DRAGEN为调查更大的变体提供的各种见解。
  3. 利用DRAGEN在队列样本上的高精度和过滤能力，以标记潜在的伪产物或违反孟德尔假设的小变异。标记的变体也可以公开使用。3.

I)在1000基因组计划数据集上使用的DRAGEN速度和变体调用器

DRAGEN平台具有多个连接使用的管道，以便从队列抽样中高精度地调用小变量和大变量。表1描述了变体调用可用的DRAGEN管道。

表1。DRAGEN管道用于1000个基因组项目分析

CYP2D6 *4、5和SMN 1/26、7目前使用v3.5.7b对1000基因组计划数据集进行的分析中不包括管道，但在DRAGEN v.3.7 (CYP2D6)和DRAGEN v.3.8 (SMN1/2)中可用。

对于1000基因组计划数据集分析，所有使用的管道都通过一个DRAGEN命令进行编排，并为每个样本部署在一个端到端运行中。为了提供DRAGEN高处理速度的示例，表2和图1显示了使用DRAGEN v3.5.7b在Illumina Connected Analytics上执行分析的2504个样本的平均运行时间8云计算平台利用f1.4 4xlarge AWS实例。

表2。ICA上DRAGEN v3.5.7b的平均运行时间(f1.4 4xlarge实例)
(n=2504，参考:hg38 Alt敏感)

图1所示。Sample DRAGEN (v3.5.7b)定时比较

DRAGEN使用不同的基因组引用、更大的AWS实例(例如，f1.16xlarge)或在本地DRAGEN服务器上运行，可以加快分析运行时间。运行时的改进也可以通过新的DRAGEN软件版本获得。9

II)对用于1000基因组计划数据的大变体调用的DRAGEN管道的额外见解

结构变异(SV)是在个体的基因组DNA中发现的相对于参考序列的大基因组变异。虽然比小核苷酸变异(SNVs)频率低，但这些大变异可以破坏基因的功能，它们在疾病中的影响已被记录在案。10

为了识别所有不同类型的大变异，DRAGEN平台集成了多种工具，用于串联收集对个人基因组结构的洞察。DRAGEN采用互补的方法，使用读取深度、分割和配对读取信息，以及针对已知基因组位点的定制算法。

在本文中，我们将展示如何将DRAGEN- sv调用器与其他DRAGEN管道结合使用，以调用sv的子类型，如拷贝数变化(Copy Number Variations, cnv)。我们还描述了DRAGEN如何部署靶向方法来调用难以调用区域的变体(例如，使用Expansion Hunter用于STRs，使用CYP2D6调用用于细胞色素P450 2D6酶编码基因中的变体，使用SMN1和2调用用于在脊髓性肌萎缩中发挥作用的生存运动神经元基因中的变体)。

DRAGEN SV呼叫

DRAGEN SV调用程序派生自开源的Manta SV调用程序。11像Manta一样，它利用分裂和配对读取的证据来发现和评分结构变异，包括删除、插入、串联复制和断裂。按照惯例，结构变异的报告尺寸要低至50bp。对于200碱基或更大的删除和串联重复，SV调用者额外断言变体及其侧翼区域之间的深度变化与变体类型一致，进一步提高了这些变体类型的精度。

在1000个基因组计划数据集(n=2504)上，DRAGEN-SV管道能够唯一地调用36961删除(1kb-20kb)和7534唯一串联重复(1kb-20kb)分布在不同的频率阈值上，如图2所示。

图2:DRAGEN-SV管道调用的删除和串联重复在不同频率阈值上的分布。

德拉根CNV呼叫

DRAGEN CNV调用者使用基于深度的方法将基因组分割为相邻种系拷贝数的区域。这是在种系假设下进行的，因此不适合检测镶嵌变异(尽管高纯度的镶嵌CNV可能由CNV调用者的输出VCF文件报告)。一个经验法则是，对于一个30倍的基因组，这个调用者在10kb或更大的变异上有很好的表现，但这大致随覆盖深度而变化。一个非常大的基因组CNV可以在CNV调用者中表示，也可以反映在其源染色体上的倍性变化。此外，我们还可以在DRAGEN-CNV调用器中看到一个CNV，如果DRAGEN-SV调用器除了基于深度的支持之外还支持拆分或配对读取，则还可以看到它。 

在1000基因组计划数据集上，DRAGEN CNV管道能够调用2891个大小大于10kb的唯一CNV，其频率分布如图3所示。

在样本中聚合cnv，我们看到每个样本中大约有150个大于10kb的调用。与小变异一样，我们注意到许多cnv非常常见，在1%的频率阈值下约有25个变异，每个样本有6个单例变异。

图3:DRAGEN-CNV管道调用的大于1kb的cnv在不同频率阈值下的分布
DRAGEN-CNV调用器的一个特性是基因组被分割成命令行上指定的覆盖箱约定，这在每次运行中都是一致的。在这种情况下，我们可以简单地将多个覆盖文件(DRAGEN输出以*.target.counts.gc_normalized结尾)合并到一个表结构中，以便在总体规模上有效地存储和查询覆盖更改(图4)。
图4:一个常见缺失的队列可视化。

热图中的垂直切片表示覆盖箱，每个切片的颜色强度表示该给定箱在1000个基因组队列中覆盖样本深度的分布。突出显示的区域表示公共删除的坐标。

DRAGEN ExpansionHunter

str /Repeat扩展是一类特殊的插入变体，其中插入序列包含重复母题的附加副本。由于SV调用的方法(依赖于插入序列的从头组装)，这仍然是一个很难调用的变体类。尽管如此，DRAGEN平台可以利用变量结构的先验知识，使用序列图(ExpansionHunter)准确地调用重复展开。12、13、14DRAGEN团队从遗传学文献中策划了临床相关的重复扩增，并将大多数已知的重复扩增报告作为该变体调用包的一部分。

我们在1000基因组计划的样本中生成了一些STRs的大小分布(图5)。

图5:1000基因组计划样本中CSTB、AFF2、ATXN10和JPH3基因中STRs的等位基因频率。

DRAGEN Ploidy Caller

在许多病理中，可以观察到从基线倍性的变化，DRAGEN倍性调用者报告了所有常染色体和性染色体的平均倍性水平。该调用者能够报告种系全染色体倍性改变或镶嵌变异，镶嵌纯度低至15%左右。这些变化是通过它们的DNA拷贝数来观察的，可能被认为是狭义上的结构变异，也可能不被认为是结构变异(许多人将它们称为结构变化)，但在这个讨论中包括它们，因为它们可以在实际解释的意义上与sv集中在一起 

我们可以使用1000 Genomes的数据来检查不同样本阵列上倍性调用的分布(图6)，使我们能够更好地理解整个染色体倍性变化的检测极限，而这在通常由不到几个样本组成的传统基准数据集中是不可行的。对这个样本池的染色体覆盖率进行评估，可以让我们看到绝大多数染色体位于预期覆盖率的90%到110%之间，而少数异常值的覆盖率有所增加。 

图6:1000基因组计划数据集样本中倍性调用的分布。

DRAGEN ROH来电者

DRAGEN也支持称为大纯合度(ROH)。该调用者对小变异基因型数据进行操作，并突出显示一些纯合变异比杂合变异富集的区域。大的ROH呼叫可以用来强调单亲同二体或亲缘关系的存在。由于ROH在人群中的自然变异，1000基因组计划数据集对于在典型人群中设置ROH区块的基线期望非常有指导意义，并可用于设置启发式阈值，将临床病例标记为异常，以便进一步调查(图7)。

图7:在1000基因组计划数据集中，大型ROH区块的数量分布以及大型ROH区块调用中snv的比例。

III)利用DRAGEN平台组合多种数据类型:不平衡易位示例

1000基因组计划数据集的广度使我们能够观察到标准基准样本无法获得的罕见变异类型的多样性。一个例子是NA20533中出现的不平衡易位。这些类型的大规模基因组重排非常罕见，当在临床样本中观察到时，几乎总是致病的。在这个例子中，我们怀疑不平衡易位是一个体细胞系的人为因素。尽管如此，我们可以将其作为临床相关变异类型的一个例子。

查看倍性调用者输出，我们看到了一个值得进一步研究的异常染色体:

倍性估计，，常染色体中位数覆盖率，36.47
…
倍性估计，13中位数/常染色体中位数，1.24
…

在这个例子中，我们使用DRAGEN的target.counts。通过gc-corrected '文件来可视化整个基因组的覆盖率(图8)。将其与染色体覆盖率的背景分布进行比较，我们可以看到这是一个明显的异常值:

图8:NA20533覆盖的全基因组可视化。

数据被表示为热图，图像中的切片表示深度箱的分布(从.target.counts。Gc_normalized文件)在100kb基因组间隔。

在上面的图8中，我们可以看到在13号染色体上有一个非常大的重复，在17号染色体上有一个很大的末端缺失，这表明存在不平衡易位。转到CNV调用(表3)，我们可以看到删除表示在一个LOSS调用中，而复制被分割到5个单独的调用中。这种分裂是由于存在常见的拷贝数变化或噪声将变体调用分裂成片段，这对于如此大的cnv是常见的。

表3:CNV的子集。NA20533的VCF文件对应于与不平衡易位相邻的CNVs。

回到我们的覆盖数据，我们可以查看这些CNV的边缘，看看我们是否信任CNV提供的断点，或者我们是否想要改进它们(图9)。
图9:NA20533在易位断点处的覆盖可视化。

点表示来自`.target.counts的规范化、分箱化覆盖数据。gc_normalized”文件。

看到这些覆盖数据，让我们有信心，我们确实对这种不平衡易位的断点有一个很好的估计。最后，由于我们知道这些易位会导致基因组重排以及拷贝数变化，我们可以在CNV定义的断点附近查询结构变体调用文件' .sv.vcf.gz '，通过配对和分割读取数据寻找该变体的证据(表4):

表4:NA20533易位断点的VCF记录。

这个例子展示了DRAGEN管道如何促进在许多尺度上观察罕见的大型变体，并允许我们从读取深度和分割读取的角度观察基因组。通过CNV调用者和相关的覆盖数据，我们得到了非常可靠的调用和对变异的清晰解释，而通过结构变异调用者，我们得到了断点分辨率和对这种特定基因组重排的机制形成的理解。

IV)小变量调用的准确性

为了统一跨样本的变体表示并进行队列分析，我们使用了 的gVCF基因型组件Illumina DRAGEN (Dynamic Read Analysis for GENomics)生物信息技术平台， 3.6.3版本。

gVCF基因型从队列中的每个样本中提取小变量调用者输出，并在所有队列成员中对任何样本中看到的所有变体进行基因型。对于不具有特定变体的样本，可以根据小变量调用者输出中的深度信息估计其纯合子参考信度，但gVCF基因型器不会试图根据来自种群的信息调整基因型。输出是在标准多样VCR15格式，可用于下游队列分析。

我们以队列调用集的形式发布其输出，该队列调用集包括来自1000个基因组计划的2504个样本，所有样本的变体基因分型和注释的群体频率。用于构建此调用集的每个样本小变量调用在以前的版本中可用。16

通过与使用广泛采用的GATK最佳实践工作流独立生成的调用集进行比较，我们展示了结果数据集的高质量17用于种系小变异。

在1000个基因组计划数据上使用DRAGEN PopGen管道进行可伸缩队列调用

DRAGEN队列发布由每个染色体的多样本VCF组成，包含2504个不相关样本。在整个基因组中，该数据集总共包含1.51亿个位点，其中1.38亿个snp和1800万个indel。值得注意的是，与其他变体调用器不同，DRAGEN输出的候选等位基因有一些证据，但不会以足够高的可信度被调用。因此，输出中的一小部分位点具有等位基因计数(AC) 0，其余的一些位点包含AC=0的等位基因，称为等位基因。这些等位基因的潜在价值如图13所示。如果需要，可以使用以下命令删除它们:

bcftools view -a {in.vcf.gz} | bcftools filter -e 'ALT="。"' -Oz -o {out.vcf.gz}

排除AC=0的记录，有1.27亿个位点，1.14亿个snp, 1500万个indel和800万个多等位变异。相比之下，类似的GATK调用集包含1.2亿个位点，1.08亿个snp, 1200万个indel和900万个多等位变异。图10显示了这些变体在等位基因频谱中的分布，揭示了DRAGEN调用了明显更罕见的变体——在DRAGEN输出中有1300多万个频率低于5%的变体。

图10:等位基因频谱中被称为变异的分布。

注意，这个图使用了完整的调用集——在VQSR之后，GATK将包含更少的变量过滤。

DRAGEN队列 呼叫精度

我们使用各种指标对DRAGEN和GATK获得的总体调用集的准确性进行了基准测试:对不同真值集的错误率，偏离Hardy-Weinberg均衡和三人中的孟德尔违反。

错误率

首先，我们计算了具有真值变异的特征良好的样本中的假阳性/假阴性计数，该样本由NIST瓶中基因组(GIAB)发布18)财团。这个真实样本，被称为NA12878，是最初1000个基因组计划队列的一部分。从多样本VCF中提取代表NA12878的列，并使用来自NIST数据集3.3.2版本的高置信区域将变量与真值集进行比较。图11显示了来自DRAGEN调用(“DRAGEN- gg”)和来自GATK输出的2个调用集的结果错误计数:联合基因分型后的所有变体(“GATK- jg”)和仅通过变体质量评分重新校准的变体(“GATK- vqsr”)。DRAGEN-GG在snp和INDELs中具有最低的假阳性和假阴性数量。

图11 a,b:针对(a) SNP和(b) INDEL的3个不同变体调用管道的错误计数:DRAGEN-GG, GATK-JG, GATK-VQSR。
一个公平的问题是，这个样本有多具有代表性，对于其余人群的准确性有真值变异。在缺乏更多真值数据的情况下，这个问题很难回答，但我们可以看看队列中NA12878错误的普遍程度，以了解其他样本中发生了什么。图12显示，即使在VQSR过滤之后，DRAGEN假阳性(图11a中的FP)通常也不如GATK假阳性常见。
图12:每个管道的NA12878 SNP FP，所有2504个样本的等位基因频率。
类似地，对于假阴性，我们通过在其他样本中寻找它们的存在，检查了NA12878中遗漏的真变异(图11a中的FN状态)在总体水平上仍未被检测到的数量。在其余队列中，只有15%的DRAGEN FNs和21%的GATK FNs未被检测到，这表明我们正在接近单一真实样本所能显示的总体敏感性的极限。对于在队列中检测到的fn，我们可以计算它们的等位基因频率。图13表明DRAGEN发现的额外变异往往具有罕见的频率。特别令人感兴趣的是等位基因数为0的那些，它们有一些支持证据，但没有足够的信心被称为。注意，此分析仅在snp中执行，以避免因INDEL站点的真值集和队列调用集之间的不同变体表示而造成混淆。
图13:每个管道的NA12878个SNP FN变体，所有2504个样本的等位基因频率。VQSR不包括在内，因为它的FNs是未经过滤的GATK的超集。

偏离Hardy-Weinberg平衡

接下来，我们测量了Hardy-Weinberg平衡的偏差，以比较 每个 调用集中显示伪像迹象的站点的比例。使用GATK VariantAnnotator计算每个位点的ExcessHet度量(过剩杂合度)，代表Hardy-Weinberg平衡检验的phredscale p值。更高的值意味着更高的技术工件的机会。图14显示DRAGEN数据集具有最高比例的不违反Hardy-Weinberg平衡的站点，突出显示了其调用的高精度。我们标记为overesshet > 28.69的过滤站点，其p值与Hardy-Weinberg期望相差超过3个标准差。这导致不到0.05%的DRAGEN站点被标记。我们在DRAGEN队列调用附带的站点专用VCF文件中提供这些信息，其中包括ExcessHet指标，以及其他常用计算的站点注释。如果用户喜欢将所有信息放在一个文件中，他们可以使用以下命令合并DRAGEN输出和站点专用VCF:

bcftools注解-a {sites.vcf.gz} -c INFO，+FILTER -Oz -o {output.vcf.gz} {dragen_calls.vcf.gz}

图14:根据overesshet度量被调用站点的分布。

当样本不相关时，该值越高，技术伪影的可能性就越大。

孟德尔的错误

最后，我们计算了样本NA20891、NA20882、NA20900中隐藏的三个一组的孟德尔错误率。与 truth 集合相比，家庭关系的违反是一个更广泛地评估准确性的有用度量标准，因为它们不局限于基因组高置信区域内的变体。表5显示了孟德尔误差在三个位点中至少一个成员变异位点总数上的数量。DRAGEN-GG 和 GATK-VQSR的性能相似，以较少的调用为代价，在 上略微优于GATK-VQSR。

表5所示。每个管道的孟德尔错误率

致谢:感谢Shyamal Mehtalia、Egor Dolzhenko、Christopher Saunders、Heidi Norton和Rami Mehio对本文测试的参与和支持，感谢他们的审阅和对数据的访问。

外部链接

https://github.com/Illumina/manta

https://github.com/Illumina/canvas

https://github.com/Illumina/ExpansionHunter

https://github.com/Illumina/SMNCopyNumberCaller

https://github.com/Illumina/gvcfgenotyper

https://github.com/Illumina/Cyrius

DRAGEN对AWS上1000个基因组数据集的再分析

有两种方式访问数据:通过AWS CLI或AWS管理控制台:

数据集的链接

小变种2504个样品:

https://s3.console.aws.amazon.com/s3/buckets/1000genomes-dragen?prefix=data/dragen-3.5.7b/hg38_altaware_nohla-cnv-anchored/gVCF-genotyper-3.6.3-2/hg38_alt_aware_nohla/2504samples/

小变种NA12878结果:

https://s3.console.aws.amazon.com/s3/buckets/1000genomes-dragen?prefix=data/dragen-3.5.7b/hg38_altaware_nohla-cnv-anchored/gVCF-genotyper-3.6.3-2/hg38_alt_aware_nohla/NA12878_annotation/

CNV/SV 2504样品频率:

https://s3.console.aws.amazon.com/s3/buckets/1000genomes-dragen?region=us-west-2&prefix=data/dragen-3.5.7b/hg38_altaware_nohla-cnv-anchored/cnv-sv-frequency/2504samples/

AWS的这篇文章解释了这些数据集:

https://aws.amazon.com/blogs/industries/dragen-reanalysis-of-the-1000-genomes-dataset-now-available-on-the-registry-of-open-data/

参考文献
  1. 玛尔塔Byrska-Bishop;et al。扩大的1000个基因组计划队列的高覆盖率全基因组测序，包括602个三人组。bioRxiv2021.02.06.430068
  2. 小变种2504个样品。
  3. 结果是NA12878的小变种。
  4. 陈晓霞，沈芳，等。。Cyrius:精确的CYP2D6基因分型使用全基因组测序数据。药物基因组学J(2021)。
  5. Cyrius文章
  6. 陈X，桑奇-胡安，A.， French, C.E.等。基于基因组测序数据的脊髓性肌萎缩症诊断和载体筛选。麝猫地中海22, 945-953(2020)。
  7. SMA的文章
  8. Illumina Connected Analytics
  9. 德拉根重新分析了1000个基因组数据集，现在可以在注册开放数据。图3
  10. Medhat Mahmoud等人。结构变体的召唤:长和短的召唤。基因组生物学(2019)20:246 doi: 10.1186/s13059-019-1828-7
  11. 陈晓宇等。Manta:用于生殖系和癌症测序应用的结构变异和索引的快速检测。生物信息学，2016年4月15日;32(8):1220-2。Doi: 10.1093 /生物信息学/ btv710。Epub 2015年12月8日PMID: 26647377。
  12. Egor Dolzhenko等人。从无pcr全基因组序列数据中检测长重复扩增。Genome Res, 2017年11月;27日(11):1895 - 1903。doi: 10.1101 / gr.225672.117。Epub 2017 9月8日。
  13. Egor Dolzhenko等人。ExpansionHunter:一种基于序列图的工具，用于分析短串联重复区域的变异。生物信息学，第35卷，第22期，2019年11月15日，第4754-4756页。doi: 10.1093 /生物信息学/ btz431
  14. 评论者/扩展猎人文章
  15. VCF (Variant Call Format) 4.2版本规范
  16. 德拉根重新分析了1000个基因组数据集，现在可以在注册开放数据。
  17. 数据收集1000个基因组
  18. 瓶子里的基因组
前一篇文章
下一篇文章