超快速检测人类基因组单碱基突变和微插入缺失的方法

摘要

本发明提出一种超快速检测人类基因组单碱基突变和微插入缺失的方法，是一种能从人基因组DNA测序结果中快速地检测出单碱基突变、微插入缺失的可行方法。本发明通过把人参考基因组序列科学有效地切分为小的子参考序列块，把人重测序中的几乎全部步骤(包括分析时间较长的步骤)都切分为计算复杂度大大降低的子任务块，而各子任务块之间相互不影响，最后把从各子参考序列块中得到的遗传多态性信息进行去冗余、校正，然后过滤，从而得到原人重测序流程中需要获取的遗传多态性信息。本发明解决了人重测序生物信息分析时间过长的问题，开创了一种新的分析模式。

说明书

技术领域

本发明涉及生物医学数据分析领域，具体涉及一种超快速检测人类基因组 单碱基突变和微插入缺失的方法。

背景技术

基因是遗传的物质基础，生命体一切生命活动(生老病死)，基因都起着重 要的基础作用。随着人类基因组计划、国际单体型计划、千人基因组计划的顺 利完成，特别是高通量测序技术(High Through-put Sequencing)的迅猛发展、 日臻成熟，越来越多的遗传多态性位点(主要是单碱基突变和微插入缺失)与 疾病健康的关系被陆续发现。目前，人们已经可以通过比较低廉的价格来对人 类的基因组进行测序、分析，获得该个体的这些多态性位点在全基因组水平上 的基因型。有了这些遗传信息，可以辅助医生们对可疑病患的确诊(某些疾病 上)，告诉医生哪个药物对该病患的疗效最好，甚至告诉人们自己对哪些疾病具 有较高的罹患风险。这些都是个体化医疗的明天，而廉价、快速地获取上诉多 态性的信息则成为了个体化医疗的关键和基础所在。

为了达到上述目的，人们构建了相关的人重测序分析的方法，该系统一般 要求对人的全基因组达到起码30×的覆盖深度，也就是大概90G碱基的输入数 据量(一般以FASTQ格式存储)。在计算机分析资源充足的情况下，该流程从 输入FASTQ文件到最终得到突变(或说遗传多态性)结果，起码需要7天的时 间，在不考虑基因测序的时间下，这样较长的数据分析时间已经远远超出了目 前大多数的医学检测服务周期，也与个体化医疗中应有的“快速”的特性不符。

人们为了缩短人重测序的分析周期，有的从最耗费时间的序列比对入手， 通过运用基于GPU的CUDA架构来进行序列比对，然而这样做却无法检测到 微插入与缺失，而且时效提高也没有太多，最多是让序列比对这个步骤比之前 的快了3倍。况且其他耗费时间的步骤没有得到提升，故对整体流程分析时间 的缩短比较有限。

另外一种思路是，对于主要耗费时间的关键步骤，运用基于分布式的存储 架构，把一个复杂的任务切分成多份来分析，从而大幅降低该步骤的分析时间。 这个方法是可以实现分析周期的大大缩短，但是该方法限制较大，对计算资源 的要求较高。主要表现在：需要搭建分布式的存储架构，而且每一个小份任务 对内存的要求几乎与整一个任务的无异，故需要大量的计算机内存资源。

发明内容

为了克服上述现有技术的不足，本发明提出一种超快速检测人类基因组单 碱基突变和微插入缺失的方法，本方法通过对人基因组参考序列科学巧妙地切 分，对人重测序分析的整体目标，从串行或部分并行地在大小为3G的人基因 组参考序列中寻找基因组多态性事件，改为从参考序列的各个可合并小块中， 并行化地进行分析，最后通过去除冗余和失真信息，从而获得全基因组水平的 多态性信息(主要是单碱基突变和微插入缺失的信息)。

为了解决上述问题，本发明的技术方案为：

一种超快速检测人类基因组单碱基突变和微插入缺失的方法，包括：

获取人全基因组序列的测序结果，对人类基因组DNA样品进行扩增、建 库和测序处理，得到人的全基因组测序结果；

对得到的人全基因组测序结果进行来料质控，以确保测序得到的序列，其 测序质量符合进行信息分析的要求；

对来料质控后的人的全基因组测序结果进行过滤接头、低质量处理，得到 已过滤序列，对已过滤序列进行来料质控，以确定这些可用的已过滤序列从质 与量上符合信息分析的标准；

将人参考基因组切分成若干小块的参考基因组块，且对每一小块分别建立 适用于比对软件的二进制索引库；

将上述待比对序列，序列比对到切分好的某参考基因组块，对其进行人重 测序生物信息分析，得到测序序列在该参考基因组块中的比对结果；

对该比对结果进行排序、重比对、比对质量与比对得分值校正，得到修正 后的比对结果；

对修正后的比对结果进行去重，得到能够用于检测遗传多态性的比对结果 文件；

基于该结果文件，检测该参考基因组块中的某个小块的遗传多态性SNP和 Indel；然后分别对该参考基因组块中其他小块进行同样的分析；不同参考基因 组块的分析能够并行处理，且相互不影响；上述所有分析都完成后，得到若干 个分别对应于参考基因组块的遗传多态性结果；

对上述遗传多态性的结果进行合并，得到全基因组的遗传多态性结果；对 得到全基因组的遗传多态性结果进行校正、过滤，并对过滤后的多态性结果进 行结果质控，以确定该结果是正常无误的人的遗传多态性结果。

本发明要解决的技术问题，是提供一种能从人基因组DNA测序结果中快 速地检测出单碱基突变、微插入缺失的可行方法。该方法，通过把人参考基因 组序列科学有效地切分为小的子参考序列块，把人重测序中的几乎全部步骤(包 括分析时间较长的步骤)都切分为计算复杂度大大降低的子任务块。而各子任 务块之间相互不影响，最终通过一定规则和算法，把从各子参考序列块中得到 的遗传多态性信息进行去冗余、校正，然后过滤，从而得到原人重测序流程中 需要获取的遗传多态性信息。解决了人重测序生物信息分析时间过长的问题， 开创了一种新的分析模式。

优选的，所述将人参考基因组切分成若干小块的参考基因组块的具体过程 为：是根据人参考基因组各染色体长度、基因组中的N序列区的分布、基因分 布，对参考基因组序列进行切分，得到若干块各自独立又互有联系的原参考基 因组的小块。继而分别对这些小块，根据后续所用到的比对软件建立二进制索 引库。

优选的，所述对得到的人全基因组测序结果进行来料质控，是对不同测序 长度的测序结果分别以下述指标进行质控，以保障测序结果的序列符合人重测 序流程分析对于测序数据的最低质量标准，其指标包括：

Q20、Q30的百分率，测序错误率，AT-GC分离率，每个位置在整体序列 中被测出为N的比例。

优选的，所述对来料质控后的人的全基因组测序结果进行过滤接头、低质 量处理，得到待比对的序列的具体方式为：是将测序结果中含有的、由测序建 库过程中引入的接头序列去除，将平均测序质量低的小部分序列予以去除。判 断测序质量是否偏低，一般的阈值为平均测序错误率>3％，可根据具体分析项 目来灵活变动，建议阈值起码为5％。

优选的，所述对已过滤序列进行来料质控的过程具体是：

对经过过滤得到的可用数据量占原数据量的比率进行统计及质控，用于判 定此次测序结果的可用率是否达到测序服务的标准；

若无法达到，且数据量经过大致推算后，符合至少30×覆盖度的指标，则 仍能够用该数据进行后续分析处理；否则，则需加测或重测，以达到对数据量 的最低要求。

下面所说的步骤，均为把上述所说的经过滤后得到的待比对序列，并行化 地对于切分好的各参考基因组小块进行生物信息分析。

本发明提供的超快速检测人类基因组多态性和微缺失的方法，在人重测序 生物信息分析步骤中，包含、但却不硬性包含以下子步骤：

对经过滤处理后得到的待比对的序列，采用并行的方式对切分好的参考基 因组块进行分析，其具体过程为：

(1)序列比对，把测序结果比对到某参考基因组的分块中；该步骤中，还 需要根据具体测序类型、测序序列长度、测序建库策略及多态性检测要求设置 相适宜的比对参数；

(2)排序，将上述步骤(1)得到的比对结果按照染色体和比对坐标排序， 其中比对坐标排序按从小到大排序；

(3)重比对、比对质量与比对得分值校正，首先找出比对区域中多态性位 点出现多的区域，对该区域的比对结果进行优化，对比对结果由于高频SNP、 Indel的存在而比对得分低的情况予以校正；

(4)去重，由同一个DNA分子通过PCR扩增而来的其他分子，相对于原 模板分子而言，都是冗余的数据，故在进行信息分析的时候需要予以区分并标 记，在多态性检测的时候与其他非冗余的数据予以区分对待；

(5)SNP检测，基于上述步骤(4)生成的比对结果，对该基因组块区域 进行SNP检测，得到该参考基因组块上的SNP检测结果；

(6)Indel检测，基于上述步骤(4)生成的比对结果，对该基因组块区域 进行Indel检测，得到该参考基因组块上的Indel检测结果。

对于上述步骤，可以根据具体分析要求和数据特点来适当修改或增删，并 不仅仅局限于上述所示的步骤。待比对序列分别对应于切分好的各参考基因组 块的生物信息分析步骤完成。

优选的，上述遗传多态性的结果进行合并，得到全基因组的遗传多态性结 果，是通过相邻参考序列块之间的切分特点、序列比对的特点，把从所有分块 中得到的多态性结果进行去冗余和校正；并根据位点覆盖深度、附近序列覆盖 率、位点平均测序质量、位点附近序列比对的均一性指标对多态性结果进行过 滤，最终得到全基因组的可信的多态性结果。

优选的，在多态性结果中还进行注释与统计，是对上述得到的最终多态性 结果进行注释、统计。这样就可以知道这些多态性位点位于哪个功能元件，或 哪些基因、编码哪个氨基酸，这些多态性(或一说突变)是否与某些疾病或药 物作用位点有关系等。

上述对统计结果进行结果质控是对得到的多态性结果以及统计结果的以下 指标进行质控，以确保经过本方法得到的多态性结果从基因组学和遗传性、生 物信息学角度来看，是基本合理的：

1.与质谱的21个位点的基因分型结果相对比，以确保在整个测序、分析的 过程中无样品污染与样品颠倒。

2.检测出的SNP占公共数据库中的比率较高，以确定SNP检测的准确性。

3.检测出的Indel占公共数据库中的比率较高，以确定Indel检测的准确性。

4.短Indel较长Indel多。

5.CCDS区内，长度为3的倍数的Indel较与之长度相近的Indel多。

与现有技术相比，本发明的有益效果为：

1.通过利用更多计算机处理器，能极大地缩短人重测序的分析时间，而对 计算机内存要求却比现有人重测序分析的低。

2.具有较好的稳点性，因为每个子任务所要花费的时间大大缩短，且所需 的内存较小，即使某几个子任务出错了，也会能快重新运行并完成，不会对整 体分析时间造成太大影响。

3.具有很高的扩展性，在计算资源充足的情况下，可以灵活根据对分析时 间的要求，设置较大的切分块份数，从而达到更高的分析效能、更短的分析时 间。

附图说明

图1为本发明实施例提供的超快速检测人类基因组多态性和微缺失的方法 的总流程图。

图2为本发明实施例提供的超快速检测人类基因组多态性和微缺失的方法 流程中，人重测序生物信息分析步骤的流程图。

图3为基于本发明超快速检测人类基因组多态性和微缺失的方法的实施例 流程图。

图4为本发明超快速检测人类基因组多态性和微缺失的方法的实施例中， 所检测出的SNP示意图。

图5为本发明超快速检测人类基因组多态性和微缺失的方法的实施例中， 所检测出的SNP的统计结果(含部分结果质控项)示意图。

图6为本发明超快速检测人类基因组多态性和微缺失的方法的实施例中， 所检测出的SNP的注释结果示意图。

图7为本发明超快速检测人类基因组多态性和微缺失的方法的实施例中， 所检测出的Indel结果示意图。

图8为本发明超快速检测人类基因组多态性和微缺失的方法的实施例中， 所检测出的Indel的统计结果(含部分结果质控项)示意图。

图9为本发明超快速检测人类基因组多态性和微缺失的方法的实施例中， 所检测出的全部Indel的长度分布统计图。

图10为本发明超快速检测人类基因组多态性和微缺失的方法的实施例中， 所检测出的CCDS区内的Indel的长度分布统计图。

图11为本发明超快速检测人类基因组多态性和微缺失的方法的实施例中， 所检测出的Indel的注释结果示意图。

图12为本发明超快速检测人类基因组多态性和微缺失的方法的实施例与 原人重测序流程，对同样的人全基因组数据，在不同CPU数下面的具体分析时 间图表。

图13为本发明超快速检测人类基因组多态性和微缺失的方法的实施例与 原人重测序流程，对同样的人全基因组数据，在不同CPU数下面的分析时间对 比图。

具体实施方式

下面结合附图对本发明做进一步的描述，但本发明的实施方式并不限于此。

本发明提供了一种超快速检测人类基因组单碱基突变和微插入缺失的方 法，该方法的具体实现过程如图1所示，包括：

获取人全基因组序列的测序结果：对人类基因组DNA样品进行扩增、建 库和测序处理，得到人的全基因组测序结果。对测序结果进行来料质控，以确 保测序得到的序列，其测序质量符合进行后续信息分析的要求。

对来料质控后的人的全基因组测序结果进行过滤接头、低质量处理，得到 已过滤的序列，对已过滤的序列进行来料质控，以确定这些可用的待比对序列 从质与量上符合信息分析的标准。将人参考基因组进行科学有效地切分成若干 小块，且对每一小块分别建立适用于比对软件的二进制索引库。

将上述待比对序列，比对到切分好的某参考基因组块，得到测序序列在该 参考基因组块中的比对结果。对该比对结果进行排序、重比对、比对质量与比 对得分值校正，得到修正后的比对结果。

对该结果进行去重(仅对于测序序列分为多个子文件输出的情况)，得到能 够用于检测遗传多态性的比对结果文件。基于该结果文件，检测该参考基因组 小块中的遗传多态性SNP和Indel。上述从序列比对到检测SNP、Indel，可视 为对之前切分好的参考基因组块中的某个小块的遗传多态性检测，然后分别对 该参考基因组块其他所有小块进行同样的分析。不同参考基因组块之间的分析 可以并行处理，且相互不影响。待上述所有分析都完成后，

将会有若干个分别对应于参考基因组各个小块的遗传多态性(SNP、Indel) 的结果。对这些结果进行合并，得到全基因组的遗传多态性结果。对这些结果 进行校正、过滤，并对统计结果进行结果质控，以确定该结果是正常无误的人 的遗传多态性结果。

在本实施例中，在对参考基因组进行切分的步骤中，根据对人参考基因组 各染色体长度、基因组中的N序列区的分布、基因分布，对参考基因组序列进 行科学有效地切分，得到若干块各自独立又互有联系的原参考基因组的小块。 然后分别对这些小块建立适用于比对软件的二进制索引库。

对参考基因组的切分，可根据项目综合考虑分析时间、分析资源等条件来 确定，一般推荐为大于30份。无论需要切分的份数是多少，我们在切分的时候， 统一遵循以下规则：

(1)各小块的碱基数量基本一致，起码需要在同一个数量级。

(2)若切分的份数较多(如大于50份)，优先根据染色体来切分，即优 先切分成每条染色体各一份，共24(M染色体和Y染色体归并为一 份)份，然后再对较长的染色体进一步切分。

(3)若要对某染色体进行切分，断点尽量选在基因组的N序列区。

若切分断点与其他规则相左，无法选在N序列区，则选在目前数据库无功 能注释的基因间区，且相邻切分区域头尾之间有一定的重合，重合大小一般为 测序平均插入片段3倍方差的大小。

在本实施例中，在对测序结果进行来料质控，是对不同测序长度的测序结 果分别以下述指标进行质控，以保障测序结果的序列符合人重测序流程分析对 于测序数据的最低质量标准：

1.Q20、Q30的百分率。

2.测序错误率。

3.AT-GC分离率。

4.每个位置在整体序列中被测出为N的比例。

本发明提供的超快速检测人类基因组多态性和微缺失的方法的实施例中， 在对来料质控后的人的全基因组测序结果进行过滤接头、低质量处理时，通过 将测序结果中含有的、由测序建库过程中引入的接头序列去除，将平均测序质 量偏低的小部分序列予以去除。

在本实施例中，在对已过滤的序列进行来料质控中，对经过上述过滤步骤 后得到的可用数据量占原数据量的比率进行统计，及质控，以判定此次测序结 果的可用率是否达到测序服务的标准。若无法达到，而数据量经过大致推算后， 仍符合至少30×覆盖度的指标，则仍能够用该数据进行后续分析处理。否则， 则需加测或甚至重测，以达到对数据量的最低要求。

下面所说的步骤，均为把上述所说的经过滤后得到的待比对序列，并行化 地对于切分好的各参考基因组小块的进行生物信息分析。

在本实施例中，在人重测序生物信息分析步骤中，包含、但却不硬性包含 以下子步骤：

1.序列比对步骤。把测序结果比对到某参考基因组的分块中。该步骤中， 还需要根据具体测序类型、测序序列长度、测序建库策略、多态性检测要求等 因素设置相适宜的比对参数。

2.排序步骤。将上述步骤(1)得到的比对结果，按照染色体和比对坐标 排序(从小到大)。该步骤只是为了符合后续分析步骤对输入的比对结果的要 求而做。

3.重比对、比对质量与比对得分值校正步骤。该步骤主要是，首先找出比 对区域中多态性位点出现较多的区域，通过目前人们已经发现的、公共数据库 收录的高频SNP、Indel结果，辅助来对该区域的比对结果进行优化，对比对结 果由于高频SNP、Indel的存在而比对得分偏低的情况予以一定校正，从而让后 续的变异检测软件更准确地对SNP、Indel进行检测。

4.去重步骤。由同一个DNA分子通过PCR扩增而来的其他分子，相对于 原模板分子而言，都是冗余的数据，故本实施例在进行信息分析的时候需要予 以区分并标记，在多态性检测的时候与其他非冗余的数据予以区分对待。

5.SNP检测。基于上述步骤4生成的比对结果，对该基因组块区域进行 SNP检测，得到该参考基因组小块上的SNP检测结果。

6.Indel检测。基于上述步骤4生成的比对结果，对该基因组块区域进行 Indel检测，得到该参考基因组小块上的Indel检测结果。

对于上述步骤，可以根据具体分析要求和数据特点适当修改或增删，并不 仅仅局限于上述所示的步骤。待比对序列分别由切分好的各参考基因组小块的 生物信息分析步骤完成。

在本实施例中，在校正多态性结果的步骤中，通过相邻参考序列小块之间 的切分特点、序列比对的特点，合理有效地把从所有分块中得到的多态性结果 进行去冗余和合理校正。并根据位点覆盖深度、附近序列覆盖率、位点平均测 序质量、位点附近序列比对的均一性等指标对多态性结果进行过滤，最终得到 全基因组的可信的多态性结果。

在本实施例中，还包括对多态性结果进行注释与统计，是对上述得到的最 终多态性结果进行注释、统计。这样就可以知道这些多态性位点都位于哪个功 能元件，或哪些基因、编码哪个氨基酸，这些多态性(或一说突变)是否与某 些疾病或药物作用位点有关系等。

在本实施例中，在对统计结果进行结果质控中，对由上述步骤得到的多态 性结果以及统计结果的以下指标进行质控，以确保经过本方法得到的多态性结 果从基因组学和遗传性、生物信息学角度来看，是基本合理的：

1.与质谱的21个位点的基因分型结果相对比，以确保在整个测序、分析的 过程中无样品污染与样品颠倒。

2.检测出的SNP占公共数据库中的比率较高，以确定SNP检测的准确性。

3.检测出的Indel占公共数据库中的比率较高，以确定Indel检测的准确性。

4.短Indel较长Indel多。

5.CCDS区内，长度为3的倍数的Indel较与之长度相近的Indel多。

如图3所示，一种超快速检测人类基因组多态性和微缺失的方法，包括以下 步骤：

步骤301，测序，以获得人DNA的全基因组测序结果。本实施例中，测序方 法可以采用高通量测序技术，如Illumina Hiseq 2000测序介绍，该技术是基于边 合成边测序技术(SBS,Sequencing By Synthesis)的下一代测序技术，能在一个 分析周期中产出足以进行人重测序分析要求的大量测序数据，且数据质量较高。

本实施例中，采用的是某测序服务项目的数据，该项目是人重测序分析项目， 数据量约为30×的全基因组数据。若需要利用该数据进行测试，请先征得项目人 与数据持有人(即本项目的服务对象)的同意。

步骤302，对不同测序长度的测序结果分别对以下指标进行质控，以保障测 序结果的序列符合人重测序流程分析时对于测序数据的最低质量标准，其指标包 括：(1)Q20、Q30的百分率，(2)测序错误率，(3)AT-GC分离率，(4)每个 位置在整体序列中被测出为N的比例。

本实施例中，运用fastQC软件对上述质控指标进行计算统计并给出质控结 果。

步骤303，通过将测序结果中含有的、由测序建库过程中引入的接头序列去 除，将平均测序质量偏低的小部分序列予以去除。本实施例使用的是SOAPnuke 软件(http://soap.genomics.org.cn)，该软件还可以通过预估数据量，来按照分析 所需的目标数据量把多余的剔除掉。

步骤304，通过对人参考基因组各染色体长度、基因组中的N序列区的分布、 基因分布、适宜后续正确合并相邻块之间的多态性结果等因素的综合评价，对参 考基因组序列进行科学有效地切分，得到若干份各自独立又互有联系的原参考基 因组的小块。然后分别对这些小块，根据后续所用到的比对软件(如：bwa、bowtie) 建立二进制索引库文件。切分参考基因组的软件可以将人类参考基因组有效、正 确地切分成所需的任意份数。

步骤305，对经过上述过滤步骤后得到的可用数据量占原数据量的比率进行 统计及质控，以判定此次测序结果的可用率是否达到一般测序服务的标准。该统 计项也是对测序质量的另一个质控指标。若无法达到，而数据量经过大致推算后， 仍符合至少30×覆盖度的指标，则仍可用该数据进行后续分析处理。否则，则需 加测或甚至重测，以达到对数据量的最低要求。下面实施例中所描述的步骤，均 为把上述所说的经过过滤后得到的待分析序列，并行化地对切分好的各参考基因 组小块的进行生物信息分析。

步骤306，对上述步骤得到的测序结果序列，分别比对到切分好的参考基因 组小块中，不同小块之间的比对是并行进行且相互不影响的。一般用目前使用较 为广泛的bwa比对软件。其他基于Burrows-Wheeler Transfer及后缀数组算法的短 序列比对软件(如：bowtie等)亦可，但需要在切分参考基因组的步骤，即步骤 304中，建立该比对软件指定的二进制索引库文件。步骤306中，还需要根据具体 测序类型、测序序列长度、测序建库策略、多态性检测要求等因素设置相适宜的 比对参数。本实施例所用的是bwa软件及samtools做格式转换，输出的是BAM格 式的文件。

步骤307，比对结果，按照染色体和比对坐标排序(从小到大)。该步骤只是 为了符合后续分析步骤对输入的比对结果的要求而做，若比对结果为BAM格式 的文件，可用samtools软件的calmd或sort命令即可实现。本实施例所用的是 samtools软件的calmd命令，输出BAM格式文件。

步骤308，对上述比对结果进行重比对以及比对质量值的校正。详细方法是， 首先找出比对结果中多态性位点出现较多的区域，通过目前人们已经发现的、公 共数据库收录的高频SNP、Indel结果，辅助来对该区域的比对结果进行优化，对 比对结果由于高频SNP、Indel的存在而比对得分偏低的情况予以一定校正，从而 让后续的变异检测软件更准确地对SNP、Indel进行检测。可用的分析软件有GATK (Genome Analysis Toolkit)，推荐版本为1.6-13以上。所用到的公共数据有： dbSNP、hapmap、OMIM、mills、千人基因组等。本实施例所用的软件是GATK  v1.6-13，公共数据库是dbSNP141、hapmap、OMIM、mills、千人基因组第1-3期 的数据。同样是输出BAM格式文件。

步骤309，对上述比对结果进行去重。由同一个DNA分子通过PCR扩增而来 的其他分子，相对于原模板分子而言，都是冗余的多余数据，在进行信息分析的 时候需要予以区分并标记，在多态性检测的时候与其他非冗余的数据予以区分对 待。

步骤310，基于上述步骤生成的比对结果，对该基因组小块区域进行SNP检 测。可用的检测软件有：GATK、samtools、MAQ snp、SOAPsnp等。不同检测 软件对SNP的检测效能和检测时间等各有特点和优势，权衡各个因素后，本实施 例所用GATK软件。本步骤得到该参考基因组小块中的潜在SNP结果，为VCF格 式。

步骤311，基于上述步骤310生成的比对结果，对该基因组小块区域进行Indel 检测。可用的检测软件有：GATK，samtools的mpileup命令、DIndel等。不同检 测软件对SNP的检测效能和检测时间等各有特点和优势，权衡各个因素后，本实 施例所用的是GATK软件。本步骤得到该参考基因组小块中的潜在Indel结果，为 VCF格式。

本实施例中，待比对序列分析分别由对切分好的各参考基因组小块进行的生 物信息分析步骤完成。

步骤312，通过相邻参考序列小块之间的切分特点、序列比对的特点，合理 有效地把从所有分块中得到的多态性结果进行去冗余和合理校正。并根据位点覆 盖深度、附近序列覆盖率、位点平均测序质量、位点附近序列比对的均一性等指 标对多态性结果进行过滤，最终得到全基因组的可信的多态性结果。过滤用的是 跟多态性检测软件一致的软件，为GATK中的子工具bcftools，以保持分析的连续 性，输出的文件依然为VCF格式。

步骤313，对上述得到的最终多态性结果进行注释、统计。则可知道这些多 态性位点都位于哪个功能元件或哪些基因、编码哪个氨基酸，这些多态性(或一 说突变)是否与某些疾病或药物作用位点有关系等。注释可用的软件有：annovar、 Polyphen、SIFT等，本实施例运用的是annovar软件，输出的文件为csv文件，可 用office的excel软件打开，方便生物学和医学研究者在windows或MAC平台查看 分析结果。

步骤314，对上述步骤得到的多态性结果以及统计结果的以下指标进行质控， 以确保经过本方法得到的多态性结果从基因组学和遗传性、生物信息学角度来 看，是基本合理的，其质控主要体现在：

(1)与质谱的21个位点的基因分型结果相对比，以确保在整个测序、分析 的过程中无样品污染与样品颠倒。

(2)检测出的SNP占公共数据库中的比率较高，以确定SNP检测的准确性。

(3)检测出的Indel占公共数据库中的比率较高，以确定Indel检测的准确性。

(4)短Indel较长Indel多。

(5)CCDS区内，长度为3的倍数的Indel较与之长度相近的Indel多。

到此本实施例的所有分析内容完成，具体检测结果以及相关的统计结果见图 4至图11。针对同样的人全基因组重测序的下机数据(约30×)，在不同的计算资 源(CPU数)下，本实施例整个分析流程与传统人重测序流程的使用时间的对比， 显示出在CPU数较多的情况下，分析时间缩短了3倍以上。具体结果见图12、图 13。

以上所述的本发明的实施方式，并不构成对本发明保护范围的限定。任何 在本发明的精神原则之内所作出的修改、等同替换和改进等，均应包含在本发 明的权利要求保护范围之内。