本发明涉及高通量测序数据的处理领域，具体而言，涉及一种高通量测序数据的处理方法、处理装置、存储介质及处理器。
背景技术：
目前，对基因测序来检测突变的方法有很多。其中通过多重扩增对特定目标区域进行特异性扩增，并对产物进行高通量测序是一种高效、经济的常用方法。然而，高通量测序过程会产生大量的序列信息，因而，如何快速准确地处理这些测序数据信息，便成为一个亟待解决的技术问题。
尽管现有技术中有许多高通量测序数据的处理分析方法，但这些方法存在处理结果准确性低的缺陷。因此，仍需要对现有的测序数据的处理方法进行改进。
技术实现要素：
本发明的主要目的在于提供一种高通量测序数据的处理方法、处理装置、存储介质及处理器，以解决现有的处理结果中存在假阳性位点多的问题。
为了实现上述目的，根据本发明的一个方面，提供了一种高通量测序数据的处理方法，该处理方法包括：获取二级测序序列，二级测序序列为高通量测序数据中能够被目的片段扩增引物识别，且去除对应的扩增引物后的测序序列；比对二级测序序列与参考基因组序列，获得初级变异结果；以及利用已知突变数据中的突变数据修正初级变异结果，得到处理结果。
进一步地，获取二级测序序列的步骤包括：过滤下机的高通量测序数据中低质量的测序数据，获得一级测序序列，低质量的测序数据是指q20低于80％或n碱基比例大于10％的测序序列；利用目的片段的扩增引物识别一级测序序列，得到识别序列；以及去除识别序列中对应的扩增引物，得到二级测序序列。
进一步地，比对二级测序序列与参考基因组序列，得到初级变异结果的步骤包括：根据目的片段的扩增引物的位置信息，从参考基因组序列上截取对应的目的片段的参考比对序列；将二级测序序列与参考比对序列进行比对，得到初级变异结果。
进一步地，在将二级测序序列与参考比对序列进行比对之后，以及得到初级变异结果之前，处理方法还包括：将二级测序序列与参考比对序列进行比对，得到比对序列；根据扩增引物的位置信息，判断比对序列中是否存在异常序列，异常序列指比对质量低于第一阈值的序列或者与参考比对序列的信息不一致的序列；若存在，则从比对序列中过滤掉异常序列，并统计剩余序列的每个位置与参考比对序列的异同，得到初级变异结果。
进一步地，利用已知突变数据中的突变数据修正初级变异结果，得到处理结果的步骤包括：从已知突变数据中筛选出目的片段对应区域内的已知突变数据，得到已知突变的局部序列；从初级变异结果中筛选同时也出存在于已知突变数据中的变异位点，形成初级变异结果的局部序列；将初级变异结果的局部序列与已知突变的局部序列进行比对，得到处理结果。
进一步地，将初级变异结果的局部序列与已知突变的局部序列进行比对，得到处理结果的步骤包括：将初级变异结果的局部序列与已知突变的局部序列进行比对，得到二级变异结果；对二级变异结果进行修正，得到处理结果；其中，对二级变异结果进行修正的步骤包括：判断二级变异结果中是否存在邻近的突变位点，如存在，则判断邻近的突变位点的变异频率是否存在显著差异并且是否有支持序列，若不存在显著差异并且有支持序列，则将相邻的突变位点合并，从而得到处理结果。
进一步地，利用目的片段的扩增引物识别一级测序序列，得到识别序列的步骤包括：步骤a，循环目的片段的所有扩增引物，从每条扩增引物的5’端开始，按照长度l截取特异性序列并记录每对扩增引物的特异性序列的数量、对应的特异性序列以及特异性序列后剩余的引物序列的长度；步骤b，变化长度l，重复步骤a，得到所有扩增引物的不同数量的特异性序列的集合，选择特异性序列数量最多的集合所对应的长度l以及对应的特异性序列集合进行后续分析；步骤c，循环处理一级测序序列中的每条序列，截取每条序列的前25～35bp的序列，从5`端开始，按照特异性序列数量最多的集合所对应的长度l去截取序列，得到测序截取序列集合；步骤d，查找特异性序列数量最多的集合中的特异性序列所对应的扩增引物在测序截取序列集合中出现次数最多的扩增引物及对应次数，并当对应次数的最大值大于设定的第二阈值时，即认为该条一级测序序列是由出现次数最多的扩增引物扩增得到，则将该条一级测序序列记为识别序列。
进一步地，去除识别序列中对应的扩增引物，得到二级测序序列的步骤包括：按照出现次数最多的扩增引物的特异性序列在识别序列中最后出现的位置以及最后出现的位置之后剩余的扩增引物序列的长度，去除识别序列中的扩增引物，得到二级测序序列。
进一步地，从已知突变数据中筛选出目的片段对应区域内的已知突变数据，得到已知突变的局部序列的步骤包括：筛选已知突变数据中目的片段对应区域，形成已知突变区域；记录已知突变区域中每个已知突变的起始位置和终止位置，并沿起始位置和终止位置分别向两端延伸，然后记录延伸后的起始位置和终止位置，且当延伸后的起始位置和终止位置位于目的片段对应区域内时，延伸后的起始位置和延伸后的终止位置对应的序列即为已知突变的局部序列；当延伸后的起始位置和/或终止位置超出目的片段对应区域，则将目的片段对应区域的边界作为已知突变的局部序列的起始位置和/或终止位置。
进一步地，从初级变异结果中筛选同时也出存在于已知突变数据中的变异位点，形成初级变异结果的局部序列的步骤包括：从初级变异结果中筛选同时也出存在于已知突变数据中的变异位点，记录每个变异位点的起始位置和终止位置，并沿起始位置和终止位置分别向两端延伸，延伸至已知突变的局部序列对应的位置，即为初级变异结果的局部序列。
进一步地，将初级变异结果的局部序列与已知突变的局部序列进行比对，得到二级变异结果的步骤包括：查找每个目的片段对应的已知突变的局部序列中是否存在初级变异结果；若存在一个初级变异结果，则根据变异结果的起始位置和终止位置，并沿起始位置和终止位置分别向两端延伸，形成一个样本突变局部序列；若存在多个初级变异结果，则判断多个初级变异结果之间的变异频率是否存在显著差异；若全部都存在显著差异，则根据每个初级变异结果的起始位置和终止位置，并沿起始位置和终止位置分别向两端延伸，形成各自的样本突变局部序列；若存在无显著差异的初级变异结果，则初步判定多个初级变异结果为连锁，并将多个初级变异结果合并形成同一样本突变局部序列，而多个初级变异结果中存在显著差异的剩余初级变异结果，则单独生成各自的样本突变局部序列；判断各样本突变局部序列与已知突变的局部序列是否相同，若相同，则将初级变异结果校准为已知突变结果；若不同，则不作校准；将已校准为已知突变结果的突变位点与剩余的未作校准的突变位点合并，得到二级变异结果。
进一步地，对二级变异结果进行修正的步骤包括判断多个初级变异结果为连锁是否存在假阳性的步骤；其中，判断多个初级变异结果为连锁是否存在假阳性的步骤包括：提取同时覆盖多个变异结果的序列，并统计支持同时覆盖多个变异结果的二级测序序列的比例；若支持同时覆盖多个变异结果的二级测序序列的比例与支持多个变异结果中每个变异结果的序列的比例不存在显著差异，则确认多个初级变异结果为连锁出现，并以连锁突变的方式重新计算突变频率，当重新计算后的突变频率满足第三阈值时，得到修正突变结果；若支持同时覆盖多个变异结果的二级测序序列的比例与支持多个变异结果中每个变异结果的序列的比例存在显著差异，则确认多个初级变异结果为连锁存在假阳性，并将合并的多个变异结果进行拆分后重新计算突变频率，当重新计算后的突变频率满足第三阈值时，得到修正突变结果；将修正突变结果与未修正的突变结果合并，得到处理结果。
为了实现上述目的，根据本发明的一个方面，提供了一种高通量测序数据的处理装置，该处理装置包括：二级测序序列获取单元，用于获取二级测序序列，二级测序序列为高通量测序数据中能够被目的片段扩增引物识别，且去除了对应的扩增引物后的测序序列；初级变异结果获取单元，用于比对二级测序序列与参考基因组序列，获得初级变异结果；以及修正单元，用于利用已知突变数据中的突变数据修正初级变异结果，得到处理结果。
进一步地，二级测序序列获取单元包括：过滤模块，用于过滤下机的高通量测序数据中低质量的测序数据，获得一级测序序列，低质量的测序数据是指q20低于80％或n碱基比例大于10％的测序序列；识别模块，用于利用目的片段的扩增引物识别一级测序序列，得到识别序列；以及去除模块，用于去除识别序列中对应的扩增引物，得到二级测序序列。
进一步地，初级变异结果获取单元包括：截取模块，用于根据目的片段的扩增引物的位置信息，从参考基因组序列上截取对应的目的片段的参考比对序列；以及第一比对模块，用于将二级测序序列与参考比对序列进行比对，得到初级变异结果。
进一步地，第一比对模块在将二级测序序列与参考比对序列进行比对之后，以及得到初级变异结果之前，还包括：第一比对子模块，用于将二级测序序列与参考比对序列进行比对，得到比对序列；判断子模块，用于根据扩增引物的位置信息，判断比对序列中是否存在异常序列，异常序列指比对质量低于第一阈值的序列或者与参考比对序列信息不一致的序列；以及过滤子模块，用于判断子模块的判断结果为存在时，从比对序列中过滤掉异常序列，并统计剩余序列的每个位置与参考比对序列的异同，得到初级变异结果。
进一步地，修正单元包括：已知突变的局部序列模块，用于从已知突变数据中筛选出目的片段对应区域内的已知突变数据，得到已知突变的局部序列；初级变异结果的局部序列模块，用于从初级变异结果中筛选出同时也存在于已知突变数据中的变异位点，形成初级变异结果的局部序列；以及第二比对模块，用于将初级变异结果的局部序列与已知突变的局部序列进行比对，得到处理结果。
进一步地，第二比对模块包括：第二比对子模块，用于将初级变异结果的局部序列与已知突变的局部序列进行比对，得到二级变异结果；修正子模块，用于对二级变异结果进行修正，得到处理结果；其中，修正子模块对二级变异结果执行如下进行修正步骤：判断二级变异结果中是否存在邻近的突变位点，如存在，则判断邻近的突变位点的变异频率是否存在显著差异和/或是否有支持序列，若不存在显著差异和/或有支持序列，则将相邻的突变位点合并，从而得到处理结果。
进一步地，识别模块包括：扩增引物特异性序列第一子模块，用于循环目的片段的所有扩增引物，从每条扩增引物的5’端开始，按照长度l截取特异性序列并记录每对扩增引物的特异性序列的数量、对应的特异性序列以及特异性序列后剩余的引物序列的长度；扩增引物特异性序列第二子模块，用于变化长度l，重复执行扩增引物特异性序列第一子模块的步骤，得到所有扩增引物的不同数量的特异性序列的集合，选择特异性序列数量最多的集合所对应的长度l以及对应的特异性序列集合进行后续分析；测序序列截取子模块，用于循环处理一级测序序列中的每条序列，截取每条序列的前25～35bp的序列，从5`端开始，按照特异性序列数量最多的集合所对应的长度l去截取序列，得到测序截取序列集合；查找子模块，用于查找特异性序列数量最多的集合中的特异性序列所对应的扩增引物在测序截取序列集合中出现次数最多的扩增引物及对应次数，并在对应次数的最大值大于设定的第二阈值时，即认为该条一级测序序列是由出现次数最多的扩增引物扩增得到，则将该条一级测序序列记为识别序列。
进一步地，去除模块包括：去除子模块，用于按照出现次数最多的扩增引物的特异性序列在识别序列中最后出现的位置以及最后出现的位置之后剩余的扩增引物序列的长度，去除识别序列中的扩增引物，得到二级测序序列。
进一步地，已知突变的局部序列模块包括：第一筛选子模块，用于筛选已知突变数据中目的片段对应区域，形成已知突变区域；第一记录模块，用于记录已知突变区域中每个已知突变的起始位置和终止位置，并沿起始位置和终止位置分别向两端延伸，然后记录延伸后的起始位置和终止位置；第一已知突变序列生成模块，用于当延伸后的起始位置和终止位置位于目的片段对应区域内时，将延伸后的起始位置和延伸后的终止位置对应的序列记为已知突变的局部序列；以及第二已知突变序列生成模块，用于当延伸后的起始位置和终止位置超出目的片段对应区域，则将目的片段对应区域的边界作为已知突变的局部序列的起始位置和/或终止位置。
进一步地，初级变异结果的局部序列模块包括：第二筛选子模块，用于从初级变异结果中筛选同时也出存在于已知突变数据中的变异位点；第二记录子模块，用于记录每个变异位点的起始位置和终止位置，并沿起始位置和终止位置分别向两端延伸，延伸至已知突变的局部序列对应的位置，即为初级变异结果的局部序列。
进一步地，第二比对子模块包括：第一查找子元件，用于查找每个目的片段对应的已知突变的局部序列中是否存在初级变异结果；第一样本突变局部序列生成元件，用于当第一查找子元件的查找结果为存在一个初级变异结果时，则根据变异结果的起始位置和终止位置，并沿起始位置和终止位置分别向两端延伸，形成一个样本突变局部序列；第二样本突变局部序列生成元件，用于当第一查找子元件的查找结果为存在多个初级变异结果，则判断多个初级变异结果之间的变异频率是否存在显著差异；若全部都存在显著差异，则根据每个初级变异结果的起始位置和终止位置，并沿起始位置和终止位置分别向两端延伸，则分别形成各自的样本突变局部序列；若存在无显著差异的初级变异结果，则初步判定多个初级变异结果为连锁，并将多个初级变异结果合并形成同一样本突变局部序列，而多个初级变异结果中存在显著差异的剩余初级变异结果，则单独生成各自的样本突变局部序列；校准子元件，用于判断各样本突变局部序列与已知突变的局部序列是否相同，若相同，则将初级变异结果校准为已知突变结果；若不同，则不作校准；第一合并子元件，用于将已校准为已知突变结果的突变位点与剩余的未作校准的突变位点合并，得到二级变异结果。
进一步地，修正子模块包括连锁假阳性判断子元件，连锁假阳性判断子元件包括：提取统计子元件，用于提取同时覆盖多个变异结果的序列，并统计支持同时覆盖多个变异结果的二级测序序列的比例；连锁确认子元件，用于当支持同时覆盖多个变异结果的二级测序序列的比例与支持多个变异结果中每个变异结果的序列的比例不存在显著差异，则确认多个初级变异结果为连锁出现，并以连锁突变的方式重新计算突变频率，当重新计算后的突变频率满足第三阈值时，得到修正突变结果；假阳性确认子元件，用于当支持同时覆盖多个变异结果的二级测序序列的比例与支持多个变异结果中每个变异结果的序列的比例存在显著差异，则确认多个初级变异结果为连锁存在假阳性，并将合并的多个变异结果进行拆分后重新计算突变频率，当重新计算后的突变频率满足第三阈值时，得到修正突变结果；以及第二合并子元件，用于将修正突变结果与未修正的突变结果合并，得到处理结果。
根据本发明的另一方面，提供了一种存储介质，存储介质包括存储的程序，其中，在程序运行时控制存储介质所在的设备执行上述任一种处理方法。
根据本发明的另一方面，提供了一种处理器，处理器用于运行程序，其中，程序运行时执行上述任一种处理方法。
应用本发明的技术方案，通过根据已知的引物信息，从高通量测序得到的原始数据中将每条序列中的引物部分去除，减少了扩增产物交叠区域内存在引物突变而造成的假阳性处理结果。此外，通过进行引物的识别还可以将高通量测序数据中一些错误扩增的序列去除，不仅有助于提高后续分析的准确性，降低假阳性结果，而且有助于降低整体数据量，提高后续分析步骤的效率。
附图说明
构成本申请的一部分的说明书附图用来提供对本发明的进一步理解，本发明的示意性实施例及其说明用于解释本发明，并不构成对本发明的不当限定。在附图中：
图1示出了本申请的一种优选实施例中的高通量测序数据的流程示意图；
图2示出了本申请的一种优选实施例中的突变校准示意图；以及
图3示出了本申请的一种优选实施例中的连锁合并示意图。
具体实施方式
需要说明的是，在不冲突的情况下，本申请中的实施例及实施例中的特征可以相互组合。下面将结合实施例来详细说明本发明。
下面详细描述本申请的实施例，具体描述的实施例是示例性的，旨在用于解释本申请，而不能理解为对本申请的限制。本申请的描述中，“一级”、“二级”、“第一”、“第二”等为描述方便，不涉及重要性方面的描述。实施例中未注明具体技术条件者，按照相关领域内文献描述的技术条件或者相关产品说明书进行。所用试剂或仪器没有特殊说明者，均为可以通过市购获得的常规产品。
如背景技术所提到的，现有技术中在处理多重扩增特定目标区域的高通量测序数据时，处理结果中往往存在假阳性位点多的缺陷，为了改善这一状况，在本申请一种典型的实施方式中，提供了一种高通量测序数据的处理方法，该处理方法包括：获取二级测序序列，二级测序序列为高通量测序数据中能够被目的片段扩增引物识别，且去除对应的扩增引物后的测序序列；比对二级测序序列与参考基因组序列，得到初级变异结果；利用已知突变数据中的突变数据修正初级变异结果，得到处理结果。
由于多重扩增方法对目标区域进行特异性捕获时，不同扩增引物的扩增区域之间可能存在交叠，因此扩增引物的存在可能会对交叠区域的突变检测造成干扰。基于这一原因，本申请的上述高通量测序数据的处理方法，通过根据已知的引物信息，从高通量测序得到的原始数据中将每条序列中的引物部分去除，减少了扩增产物交叠区域内存在引物突变而造成的假阳性处理结果。此外，通过进行引物的识别还可以将高通量测序数据中一些错误扩增的序列去除，不仅有助于提高后续分析的准确性，而且有助于降低整体数据量，提高后续分析步骤的效率。
上述处理方法中，获取二级测序序列的步骤中，除了利用已知的扩增引物信息对测序序列进行识别和去除之前，还包括高通量测序数据的常规预处理步骤，比如去除低质量测序序列的步骤。在本申请一种优选的实施例中，如图1所示的本申请的高通量测序数据的处理方法的详细流程图中，上述获取二级测序序列的步骤包括：过滤下机的高通量测序数据中低质量的测序数据，获得一级测序序列，过滤高通量测序数据中低质量的测序数据，获得一级测序序列，低质量的测序数据是指q20低于80％或n碱基比例大于10％的测序序列；利用目的片段的扩增引物识别一级测序序列，得到识别序列；去除识别序列中对应的扩增引物，得到二级测序序列。
此处要说明的是，在本申请中，下机的高通量测序数据是指从测序仪获得的fastq或者bam格式的数据。
上述获取二级测序序列的步骤中，根据测序质量、碱基识别情况对原始测序数据进行过滤筛选，避免测序过程中的低质量数据对后续数据分析造成的干扰，提高后续分析结果的准确性。
本申请的上述处理方法中，比对的步骤采用本领域常规的比对步骤即可实现。为了进一步提高比对速度以及数据处理效率，在本申请一种优选的实施例中，比对二级测序序列与参考基因组序列，得到初级变异结果的步骤包括：根据目的片段的扩增引物的位置信息，从参考基因组序列上截取对应的目的片段的参考比对序列；将二级测序序列与参考比对序列进行比对，得到初级变异结果。
由于多重扩增所得的序列大多数应是目的区域的片段，因此，本申请在进行数据比对时，选择根据引物的扩增区域来截取参考序列。这样不仅能节约计算资源，还能大大加快比对速度。具体地，比对方式为全局比对，其具体算法为：
(1)参数设定：对比对过程中的碱基相同、碱基错配、碱基插入缺失、碱基插入缺失延伸的分值进行定义；
(2)打分矩阵初始化：
a.将参考序列的每一个碱基作为打分矩阵的一列，首列留空；
b.将测序序列的每一个碱基作为打分矩阵的一行，首行留空；
c.打分矩阵填充：将打分矩阵从左向右、自上而下按照以下规则进行填充：
d.每个空位分别计算由左侧、上方、左上方延伸所得得分。其中，来自左上的情况需要判断当前位置的测序碱基与参考碱基是否相同。相同则加上碱基相同的分值，不同则加上碱基错配的分值；对于来自左侧或上方的情况，需判断前一步是否也为插入缺失。如果是则加上插入缺失延伸的分值，否则加上插入缺失的分值。
(3)将3种情况计算得到的最高分作为此处的比对分值，并记录最高分的路径来源。
(4)最优路径回溯：从打分矩阵的右下角进行回溯，根据每个位点的路径来源，得到比对结果，选取最佳比对结果。
在本申请一种优选的实施例中，在将二级测序序列与参考比对序列进行比对之后，以及得到初级变异结果之前，上述处理方法还包括：将二级测序序列与参考比对序列进行比对，得到比对序列；根据扩增引物的位置信息，判断比对序列中是否存在异常序列，异常序列指比对质量低于阈值的序列或者与参考比对序列的信息不一致的序列；若存在，则从比对序列中过滤掉异常序列，并统计剩余序列的每个位置与参考比对序列的异同，得到初级变异结果。
由于参考序列被截取，因此部分非目的扩增序列可能会被强行比对到截取后的参考序列上，这容易对后续突变检测造成干扰。此外，每对引物所扩增的序列应为对应扩增区域内的序列，因此其比对位置应和扩增目的区域的位置基本一致。基于这两点，上述优选实施例中，序列比对后对序列比对结果进行初步过滤，有助于提高后续分析结果的准确性。
上述异常序列中比对质量低于第一阈值的序列，是指根据比对算法得到比对质量值，设定阈值，低于阈值的比对质量值即为比对质量过低的序列。根据实践经验，第一阈值通常为5。上述与参考比对序列的信息不一致的序列，是指根据比对算法，与参考比对序列完全比对不上的二级序列即为信息不一致的序列。
在现有的高通量测序数据的变异检测处理方法中，突变的具体信息一般是以比对结果为准。这种方法的缺陷是如果某种突变在比对位置附近存在多种比对可能时，输出的比对结果和已有的数据库信息可能不一致，导致无法与现有数据库进行后续关联。此外，对于相对复杂的突变，比对过程可能会将其拆分成数个相对较小的突变以得到最佳比对得分，与真实变化不符。因而，最终的分析结果也不准确。
为了改进上述情况，在本申请一种优选的实施例中，利用已知突变数据中的突变数据修正初级变异结果，得到处理结果的步骤包括：从已知突变数据中筛选出目的片段对应区域内的已知突变数据，得到已知突变的局部序列；从初级变异结果中筛选同时也出存在于已知突变数据中的变异位点，形成初级变异结果的局部序列；将初级变异结果的局部序列与已知突变的局部序列进行比对，得到处理结果。
通过选择根据已知突变的局部序列来对初级变异结果进行修正校准，进而使得处理得到的最终结果更准确。具体地，需要预先生成扩增区域对应的已知突变的局部序列，可以根据扩增引物的具体信息从已知突变数据中进行截取获得。而将初级变异结果中的序列在已知突变区域内的，形成初级变异结果的局部序列，而不在已知突变区域内的不做修正校准。
虽然根据已知突变数据能够对一些可能出现的复杂突变进行了校准，但是在实际检测中，可能还会出现一些并未记录在已知突变数据库中的复杂突变。因此，在得到突变结果后，判断其中是否存在可能连锁出现的多个突变，并对其进行修正，以进一步提高分析处理结果的准确性。在本申请一种优选的实施例中，如图1和图2所示，将初级变异结果的局部序列与已知突变的局部序列进行比对，得到处理结果的步骤包括：将初级变异结果的局部序列与已知突变的局部序列进行比对，得到二级变异结果；对二级变异结果进行修正，得到处理结果；其中，对二级变异结果进行修正的步骤包括：判断二级变异结果中是否存在邻近的突变位点，如存在，则判断邻近的突变位点的变异频率是否存在显著差异并且是否有支持序列，若不存在显著差异并且有支持序列，则将相邻的突变位点合并，从而得到处理结果。
上述处理方法中，利用已知扩增引物识别测序序列的基本原理是使用每条扩增引物的特异性序列来作为对应引物的特异性标记。当某对引物的特异性序列在测序序列的前25～35bp中多次出现时，可以认为该序列是由该对应的扩增引物扩增得到。在识别出对应的扩增引物后，即可根据扩增引物的长度去除对应的扩增引物。根据上述原理，可以设计不用的具体算法来实现上述引物识别与去除的功能。
在本申请一种优选的实施例中，利用目的片段的扩增引物识别一级测序序列，得到识别序列的步骤包括：步骤a，循环目的片段的所有扩增引物，从每条扩增引物的5’端开始，按照长度l截取特异性序列并记录每对扩增引物的特异性序列的数量、对应的特异性序列以及特异性序列后剩余的引物序列的长度；步骤b，变化长度l，重复步骤a，得到所有扩增引物的不同数量的特异性序列的集合，选择特异性序列数量最多的集合所对应的长度l以及对应的特异性序列集合进行后续分析；步骤c，循环处理一级测序序列中的每条序列，截取每条序列的前25～35bp的序列，从5`端开始，按照特异性序列数量最多的集合所对应的长度l去截取序列，得到测序截取序列集合；步骤d，查找特异性序列数量最多的集合中的特异性序列所对应的扩增引物在测序截取序列集合中出现次数最多的扩增引物及对应次数，并在次数的最大值大于设定的第二阈值(通常为3)时，即认为该条一级测序序列是由出现次数最多的扩增引物扩增得到，则将该条一级测序序列记为识别序列。
本申请的这种扩增引物识别和去除的方法，能够有效地将扩增引物突变带来的干扰以及建库和/或测序过程中产生的错误扩增序列进行去除，一方面提高处理结果的准确性，另一方面能够降低整体数据量，提高后续分析步骤的效率。
需要说明的是，上述优选的实施例中，特异性序列是指所有扩增引物中不存在与之相同的序列，其是唯一的。即使不同扩增引物会截取出序列相同的序列，而这样的序列不能被称为特异性序列。
在本申请一种优选的实施例中，去除识别序列中对应的扩增引物，得到二级测序序列的步骤包括：按照出现次数最多的扩增引物的特异性序列在识别序列中最后出现的位置以及最后出现的位置之后剩余的扩增引物序列的长度，去除识别序列中的扩增引物，得到二级测序序列。
具体地，某对扩增引物可以截取的特异性序列的具体数量，一方面跟扩增引物的具体长度有关，另一方面跟所设定的截取的长度有关。扩增引物的具体长度越长，所能够截取的特异性序列的数量就越多。设定的截取长度越长，所能够截取的特异性序列的数量就越少。此处去除识别序列中的扩增引物的具体长度并非与构建文库时所用的扩增引物的全长完全一致。由于测序错误或其他未知原因导致的错误，在每条二级测序序列中，能够被某对扩增引物的特异性序列识别的最后位置之前的序列，以及其后剩余长度的序列，无论其实际是否与扩增引物序列相同，均视为能够被扩增引物所识别的引物序列，因而，都需要去除。也就是说，能够被某对扩增引物的特异性序列识别的最后位置之后剩余长度的序列，有可能与该扩增引物实际剩余的序列并不完全相同，但相同长度的碱基序列也需要去除。而本申请经过多次试验结果验证，这种去除对最终处理结果并无影响。
在本申请一种优选的实施例中，如图2所示，从已知突变数据中筛选出目的片段对应区域内的已知突变数据，得到已知突变的局部序列的步骤包括：筛选已知突变数据中目的片段对应区域，形成已知突变区域；记录已知突变区域中每个已知突变的起始位置和终止位置，并沿起始位置和终止位置分别向两端延伸(优选各延伸10～15bp)，然后记录延伸后的起始位置和终止位置，且当延伸后的起始位置和终止位置位于目的片段对应区域内时，延伸后的起始位置和延伸后的终止位置对应的序列即为已知突变的局部序列；当延伸后的起始位置和/或终止位置超出目的片段对应区域，则将目的片段对应区域的的边界作为已知突变的局部序列的起始位置和/或终止位置。
如图2中所示，根据已知突变数据库中记载的突变信息(包括突变起始位置及变异碱基类型)，通过从起始位置向前延伸以及同时从终止位置向后延伸，分别生成局部参考序列及对应的局部变异序列。而如果前后延伸后区域超出了与扩增引物的目的片段所对应的已知突变区域，目的片段对应区域的边界作为已知突变的局部序列的起始位置和/或终止位置。此处举例说明：假设在目的区域内有一个点突变，左右各延伸10bp，如果延伸的终止位置或者起始位置在目的片段范围内，则该突变的局部序列长度为21bp。假设该突变向右延伸后，终止位置超出目的片段5bp，则这5bp会被扔掉，则该突变的局部序列长度为16bp。通过形成局部序列进行比对，能够确保比对结果在目的片段区域内，从而使比对结果也更准确，避免了局部多种比对结果对突变检测造成的干扰，为后续关联各种数据库提供了方便。
需要注意的是，该生成已知突变的局部序列的步骤进行一次即可。每次分析时，如果目的片段区域没有变化，则不需每次都重新生成已知突变的局部序列。而且，生成已知突变的局部序列的时机不限，只要在比对步骤之前形成即可。
由于并非所有的初级变异结果中的变异位点均在已知突变数据中存在，为了将变异结果与已知变异数据进行比对，并且避免多个初级变异位点比对到不同已知突变数据的位置，可以将同时存在在已知突变数据中及初级变异结果中的变异位点的局部序列单独与已知突变的局部序列进行比对。在本申请一种优选的实施例中，从初级变异结果中筛选同时也出存在于已知突变数据中的变异位点，形成初级变异结果的局部序列的步骤包括：从初级变异结果中筛选同时也出存在于已知突变数据中的变异位点，记录每个变异位点的起始位置和终止位置，并沿起始位置和终止位置分别向两端延伸，延伸至已知突变的局部序列对应的位置，即为初级变异结果的局部序列。
在本申请一种优选的实施例中，将初级变异结果的局部序列与已知突变的局部序列进行比对，得到二级变异结果的步骤包括：查找每每个目的片段对应的已知突变的局部序列中是否存在初级变异结果，若存在一个初级变异结果，则根据变异结果的起始位置和终止位置，并沿起始位置和终止位置分别向两端延伸，形成一个样本突变局部序列；若存在多个初级变异结果，则判断多个初级变异结果之间的变异频率是否存在显著差异；若全部都存在显著差异，则根据每个初级变异结果的起始位置和终止位置，并沿起始位置和终止位置分别向两端延伸，形成分别形成各自的样本突变局部序列；若存在无显著差异的初级变异结果，则初步判定多个初级变异结果为连锁，并将多个初级变异结果合并形成同一样本突变局部序列，而多个初级变异结果中存在显著差异的剩余初级变异结果，则单独生成各自的样本突变局部序列；判断各样本突变局部序列与已知突变的局部序列是否相同，若相同，则将初级变异结果校准为已知突变结果；若不同，则不作校准；将已校准为已知突变结果的突变位点与剩余的未作校准的突变位点合并，得到二级变异结果。
具体地，上述优选的实施例中，如图2所示，初级变异结果中，chra的突变类型与已知突变数据库中的突变类型相同，都是atcg缺失，但初级变异结果中记录的起始位置为x4和y4，而已知突变数据中的起始位置为x1和y1，因而，根据已知突变数据库的突变信息，将该突变位点的初级变异结果校准为已知突变结果。而图2中，由于chrd的突变信息在已知突变数据库中并不存在，因而不进行校准。
在本申请一种优选的实施例中，如图3所示，对二级变异结果进行修正的步骤包括判断多个初级变异结果为连锁是否存在假阳性的步骤：其中，判断多个初级变异结果为连锁是否存在假阳性的步骤包括：提取同时覆盖多个变异结果的序列，并统计支持同时覆盖多个变异结果的二级测序序列的比例；若支持同时覆盖多个变异结果的二级测序序列的比例与支持多个变异结果中每个变异结果的序列的比例不存在显著差异，则确认多个初级变异结果为连锁出现，并以连锁突变的方式重新计算突变频率，当重新计算后的突变频率满足第三阈值时，得到修正突变结果；若支持同时覆盖多个变异结果的二级测序序列的比例与支持多个变异结果中每个变异结果的序列的比例存在显著差异，则确认多个初级变异结果为连锁存在假阳性，并将合并的多个变异结果进行拆分后重新计算突变频率，当重新计算后的突变频率满足第三阈值(通常为2％)时，得到修正突变结果；将修正突变结果与未修正的突变结果合并，得到处理结果。
将二级变异结果中多个单独的存在真实连锁关系的突变经过上述步骤修正为连锁突变，而将以连锁假阳性的形式合并的突变拆分成单独的突变，使得突变结果更准确。
在本申请另一种典型的实施方式中，提供了一种高通量测序数据的处理装置，该处理装置包括：二级测序序列获取单元，用于获取二级测序序列，二级测序序列为高通量测序数据中能够被目的片段扩增引物识别，且去除了对应的扩增引物后的测序序列；初级变异结果获取单元，用于比对二级测序序列与参考基因组序列，得到初级变异结果；修正单元，用于利用已知突变数据中的突变数据修正初级变异结果，得到处理结果。
本申请的上述处理装置，通过执行二级测序序列获取单元获取能够被目的片段扩增引物识别，且去除了对应的扩增引物后的测序序列，然后执行初级变异结果获取单元将获取到的二级测序序列与参考基因组序列进行比对，得到的初级变异结果经过修正单元执行修正步骤后，得到的处理结果中变异结果更准确。
在本申请一种优选的实施例中，二级测序序列获取单元包括：过滤模块，用于过滤下机的高通量测序数据中低质量的测序数据，获得一级测序序列，过滤高通量测序数据中低质量的测序数据，获得一级测序序列，低质量的测序数据是指q20低于80％或n碱基比例大于10％的测序序列；识别模块，用于利用目的片段的扩增引物识别一级测序序列，得到识别序列；去除模块，用于去除识别序列中对应的扩增引物，得到二级测序序列。
上述二级测序序列获取单元，根据测序质量、碱基识别情况对原始测序数据进行过滤筛选，避免测序过程中的低质量数据对后续数据分析造成的干扰，提高后续分析结果的准确性。
在本申请一种优选的实施例中，初级变异结果获取单元包括：截取模块，用于根据目的片段的扩增引物的位置信息，从参考基因组序列上截取对应的目的片段的参考比对序列；第一比对模块，用于将二级测序序列与参考比对序列进行比对，得到初级变异结果。
由于多重扩增所得的序列大多数应是目的区域的片段，因此，上述初级变异结果获取单元在进行数据比对时，选择根据引物的扩增区域来截取参考序列，这样不仅能节约计算资源，还能大大加快比对速度。
在本申请一种优选的实施例中，第一比对模块在将二级测序序列与参考比对序列进行比对之后，以及得到初级变异结果之前，还包括：第一比对子模块，用于将二级测序序列与参考比对序列进行比对，得到比对序列；判断子模块，用于根据扩增引物的位置信息，判断比对序列中是否存在异常序列，异常序列指比对质量低于第一阈值的序列或者参考比对序列与扩增引物的信息不一致的序列；过滤子模块，用于判断子模块的判断结果为存在时，从比对序列中过滤掉异常序列，并统计剩余序列的每个位置与参考比对序列的异同，得到初级变异结果。
由于截取模块是从参考序列中截取目的片段扩增序列，因此部分非目的扩增序列可能会被强行比对到截取后的参考序列上，这容易对后续突变检测造成干扰。此外，每对引物所扩增的序列应为对应扩增区域内的序列，因此其比对位置应和扩增目的区域的位置基本一致。基于这两点，上述优选实施例中，设置判断子模块和过滤子模块，分别对序列比对结果进行异常判断和初步过滤，有助于提高后续分析结果的准确性。
在本申请一种优选的实施例中，修正单元包括已知突变的局部序列模块，用于从已知突变数据中筛选出目的片段对应区域内的已知突变数据，得到已知突变的局部序列，初级变异结果的局部序列模块，用于从初级变异结果中筛选出同时也存在于已知突变数据中的变异位点，形成初级变异结果的局部序列；第二比对模块，用于将初级变异结果的局部序列与已知突变的局部序列进行比对，得到处理结果。
通过修正单元根据已知突变的局部序列来对初级变异结果的局部序列进行修正，进而使得处理得到的变异结果更准确。
在本申请一种优选的实施例中，第二比对模块包括：第二比对子模块，用于将初级变异结果的局部序列与已知突变的局部序列进行比对，得到二级变异结果；修正子模块，用于对二级变异结果进行修正，得到处理结果；其中，修正子模块对二级变异结果执行如下进行修正步骤：判断二级变异结果中是否存在邻近的突变位点，如存在，则判断邻近的突变位点的变异频率是否存在显著差异和/或是否有支持序列，若不存在显著差异和/或有支持序列，则将相邻的突变位点合并，从而得到处理结果。
虽然上述第二比对模块，根据已知突变数据能够对一些可能出现的复杂突变进行了修正校准，但是在实际检测中，可能还会出现一些并未记录在已知突变数据库中的复杂突变。因此，在得到突变结果后，通过执行上述修正子模块判断其中是否存在可能连锁出现的多个突变，并对其进行修正，以进一步提高分析处理结果的准确性。
在本申请一种优选的实施例中，识别模块包括：扩增引物特异性序列第一子模块，用于循环目的片段的所有扩增引物，从每条扩增引物的5’端开始，按照长度l截取特异性序列并记录每对扩增引物的特异性序列的数量、对应的特异性序列以及特异性序列后剩余的引物序列的长度；扩增引物特异性序列第二子模块，用于变化长度l，重复执行扩增引物特异性序列第一子模块的步骤，得到所有扩增引物的不同数量的特异性序列的集合，选择特异性序列数量最多的集合所对应的长度l以及对应的特异性序列集合进行后续分析；测序序列截取子模块，用于循环处理一级测序序列中的每条序列，截取每条序列的前25～35bp的序列，从5`端开始，按照特异性序列数量最多的集合所对应的长度l去截取序列，得到测序截取序列集合；查找子模块，用于查找特异性序列数量最多的集合中的特异性序列所对应的扩增引物在测序截取序列集合中出现次数最多的扩增引物及对应次数，并在次数的最大值大于设定的第二阈值(通常为3)时，即认为该条一级测序序列是由出现次数最多的扩增引物扩增得到，则将该条一级测序序列记为识别序列。
在本申请一种优选的实施例中，去除模块包括：去除子模块，用于按照出现次数最多的扩增引物的特异性序列在识别序列中最后出现的位置以及最后出现的位置之后剩余的扩增引物序列的长度，去除识别序列中的扩增引物，得到二级测序序列。
在本申请一种优选的实施例中，已知突变的局部序列模块包括：第一筛选子模块，用于筛选已知突变数据中目的片段对应区域，形成已知突变区域；第一记录模块，用于记录已知突变区域中每个已知突变的起始位置和终止位置，并沿起始位置和终止位置分别向两端延伸，然后记录延伸后的起始位置和终止位置，且第一已知突变序列生成模块，用于当延伸后的起始位置和终止位置位于目的片段对应区域内时，将延伸后的起始位置和延伸后的终止位置对应的序列记为已知突变的局部序列；第二已知突变序列生成模块，用于当延伸后的起始位置和终止位置超出目的片段对应区域，则将目的片段对应区域的边界作为已知突变的局部序列的起始位置和/或终止位置。
在另一种优选的实施例中，上述初级变异结果的局部序列模块包括：第二筛选子模块，用于从初级变异结果中筛选同时也出存在于已知突变数据中的变异位点；第二记录子模块，用于记录每个变异位点的起始位置和终止位置，并沿起始位置和终止位置分别向两端延伸，延伸至已知突变的局域序列对应的位置，即为初级变异结果的局部序列。
在本申请一种优选的实施例中，第二比对子模块包括：第一查找子元件，用于查找每个目的片段对应的已知突变的局部序列中是否存在初级变异结果，第一样本突变局部序列生成元件，用于当第一查找子元件的查找结果为存在一个初级变异结果时，则根据变异结果的起始位置和终止位置，并沿起始位置和终止位置分别向两端延伸，形成一个样本突变局部序列；第二样本突变局部序列生成元件，用于当第一查找子元件的查找结果为存在多个初级变异结果，则判断多个初级变异结果之间的变异频率是否存在显著差异；若全部都存在显著差异，则根据每个初级变异结果的起始位置和终止位置，并沿起始位置和终止位置分别向两端延伸(优选各延伸10～15bp)，形成各自的样本突变局部序列；若存在无显著差异的初级变异结果，则初步判定多个初级变异结果为连锁，并将多个初级变异结果合并形成同一样本突变局部序列，而多个初级变异结果中存在显著差异的剩余初级变异结果，则单独生成各自的样本突变局部序列；校准子元件，用于判断各样本突变局部序列与已知突变的局部序列是否相同，若相同，则将初级变异结果校准为已知突变结果；若不同，则不作校准；第一合并子元件，用于将已校准为已知突变结果的突变位点与剩余的未作校准的突变位点合并，得到二级变异结果。
在本申请一种优选的实施例中，修正子模块包括：连锁假阳性判断子元件，连锁假阳性判断子元件包括：提取统计子元件，用于提取同时覆盖多个变异结果的序列，并统计支持同时覆盖多个变异结果的二级测序序列的比例；连锁确认子元件，用于当支持同时覆盖多个变异结果的二级测序序列的比例与支持多个变异结果中每个变异结果的序列的比例不存在显著差异，则确认多个初级变异结果为连锁出现，并以连锁突变的方式重新计算突变频率，当重新计算后的突变频率满足第三阈值时，得到修正突变结果；假阳性确认子元件，用于当支持同时覆盖多个变异结果的二级测序序列的比例与支持多个变异结果中每个变异结果的序列的比例存在显著差异，则确认多个初级变异结果为连锁存在假阳性，并将合并的多个变异结果进行拆分后重新计算突变频率，当重新计算后的突变频率满足第三阈值(通常为2％)时，得到修正突变结果；以及第二合并子元件，用于将修正突变结果与未修正的突变结果合并，得到处理结果。
在本申请第三种典型的实施方式中，还提供了一种存储介质，该存储介质包括存储的程序，其中，在程序运行时控制存储介质所在的设备执行上述高通量测序数据的处理方法。
在本申请第四种典型的实施方式中，还提供了一种处理器，该处理器用于运行程序，其中，程序运行时执行上述高通量测序数据的处理方法。
通过以上的实施方式的描述可知，所描述的装置仅仅是示意性的，例如所述单元的划分，可以为一种逻辑功能划分，实际实现时可以有另外的划分方式，例如多个单元或组件可以结合或者可以集成到另一个系统，或一些特征可以忽略，或不执行。另一点，所显示或讨论的相互之间的耦合或直接耦合或通信连接可以是通过一些接口，单元或模块的间接耦合或通信连接，可以是电性或其它的形式。
作为分离部件说明的单元可以是或者也可以不是物理上分开的，作为单元显示的部件可以是或者也可以不是物理单元，即可以位于一个地方，或者也可以分布到多个单元上。可以根据实际的需要选择其中的部分或者全部单元来实现本实施例方案的目的。
另外，在本发明各个实施例中的各功能单元可以集成在一个处理单元中，也可以是各个单元单独物理存在，也可以两个或两个以上单元集成在一个单元中。上述集成的单元既可以采用硬件的形式实现，也可以采用软件功能单元的形式实现。
集成的单元如果以软件功能单元的形式实现并作为独立的产品销售或使用时，可以存储在一个计算机可读取存储介质中。
本领域的技术人员可以清楚地了解到本申请可借助软件加必需的通用硬件平台的方式来实现。基于这样的理解，本发明的技术方案本质上或者说对现有技术做出贡献的部分或者该技术方案的全部或部分可以以软件产品的形式体现出来，该计算机软件产品存储在一个存储介质中，包括若干指令用以使得一台计算机设备(可为个人计算机、服务器或者网络设备等)执行本发明各个实施例所述方法的全部或部分步骤。而前述的存储介质包括：u盘、只读存储器(rom，read-onlymemory)、随机存取存储器(ram，randomaccessmemory)、移动硬盘、磁碟或者光盘等各种可以存储程序代码的介质。
本说明书中的各个实施例均采用递进的方式描述，各个实施例之间相同相似的部分互相参见即可，每个实施例重点说明的都是与其他实施例的不同之处。尤其，对于系统实施例而言，由于其基本相似于方法实施例，所以描述的比较简单，相关之处参见方法实施例的部分说明即可。
下面将结合具体的实施例来进一步说明本申请的有益效果。
实施例1：检测人类基因组靶向测序产物中的突变信息
运用本申请的方法和装置，对216例由靶向测序得到的33个基因测序数据进行分析。其中包括170份癌症患者样本(全血、ffpe、胸腔积液、新鲜组织)、32例健康的无偿献血者的全血样本、14例质控品(室间质评样本、horizon标准品)。采用ionpgm测序平台进行测序，得到216例样本bam测序文件。
其中1例癌症患者(ffpe)、1例horizon标准品、1例健康人通过本发明检测的具体结果与torrentsuite具体检测结果比较如下表1和表2。其中，表1显示的是样本1(癌症患者，ffpe)利用本申请的方法进行检测的结果。
表1：
*：yes(～2％)：yes表示与实际突变结果相符；括号中的百分数表示实际突变频率。表明低频突变采用本申请的方法都能检测到，且检测结果准确，与实际突变结果相符。
表2显示的是样本1(癌症患者，ffpe)采用torrentsuite进行检测的结果。
表2：
附：yes的涵义同表1。synonymous：表示同义突变。由于目前临床上同义突变对用药的指导意义不明确，因而尽管本申请的方法也能检测出来，但在表1中未明确列出此类结果。
从样本1采用本申请的方法与现有方法相比较的表1和表2的结果可以看出，本申请的方法不仅能够将相近的连锁突变以连锁的形式显示突变结果，使得突变结果显示更准确，而且在灵敏度方面也高于现有方法，能够准确检测5％以下的低频突变。
样本2(horizon标准品)，采用本申请的方法进行检测的结果见下表3。
表3：
样本2(horizon标准品)，采用torrentsuite进行检测的结果见下表4。
表4：
从表3和表4的结果比较可以看出，样本2采用本申请及现有技术的处理方法进行处理后，所得到的结果不仅能够准确显示连锁的突变，而且能够排除假阳性突变。
样本3(健康人，全血)，采用本申请的方法进行检测的结果见下表5。
表5：
样本3(健康人，全血)，采用torrentsuite的方法进行检测的结果见下表6。
表6：
从表5和表6的结果比较可以看出，样本3采用本申请及现有技术的处理方法进行处理后，所得到的结果能够排除假阳性突变。
统计216例样本hotspot位点结果，本申请方法的检测结果与torrentsuite检测结果比较如下表7。
表7：
从以上的描述中，可以看出，与目前常规的多重扩增测序的突变检测方法相比，现有技术中使用全基因组作为参考序列进行比对运行速度较慢。此外，对于插入突变和缺失突变，比对位置的不同可能导致和数据库中的结果存在差异，进而导致无法与现有数据库中的生物学意义直接关联。而本申请的测序数据的处理方法和处理装置的优势在于：
1)高效性。通过采用根据特异性序列识别引物的步骤，利用特异性序列进行引物识别而非进行比对，可以快速高效的识别出测序序列对应的引物，极大的节省了计算资源。
2)准确性。本申请使用已知突变的对突变结果进行校准，并且在校准时使用突变的局部序列，能够有效地避免局部多种比对结果对突变检测造成的干扰，为后续关联各种数据库提供了方便。
此外，需要说明的是，本申请可用于众多通用或专用的计算系统环境或配置中。例如：个人计算机、服务器计算机、手持设备或便携式设备、平板型设备、多处理器系统、基于微处理器的系统、置顶盒、可编程的消费电子设备、网络pc、小型计算机、大型计算机、包括以上任何系统或设备的分布式计算环境等等。
显然，本领域的技术人员应该明白，上述的本发明的各模块或各步骤可以用通用的计算装置来实现，它们可以集中在单个的计算装置上，或者分布在多个计算装置所组成的网络上，可选地，它们可以用计算装置可执行的程序代码来实现，从而，可以将它们存储在存储装置中由计算装置来执行，或者将它们分别制作成各个集成电路模块，或者将它们中的多个模块或步骤制作成单个集成电路模块来实现。这样，本发明不限制于任何特定的硬件和软件结合。
以上所述仅为本发明的优选实施例而已，并不用于限制本发明，对于本领域的技术人员来说，本发明可以有各种更改和变化。凡在本发明的精神和原则之内，所作的任何修改、等同替换、改进等，均应包含在本发明的保护范围之内。