一种基于群体基因组的序列比对方法技术

一种基于群体基因组的序列比对方法，具体涉及一种基于群体基因组的人体DNA序列比对方法，为了解决以单一参考基因组作为基准的序列比对方法存在精度偏差，导致测序数据变异检测失败或得到错误变异位点的问题，首先构建基于种子序列的群体基因组索引，提取测序序列的种子序列，利用索引获取种子序列在基因组的位置；其次，根据位置利用稀疏动态规划方法得到种子序列在基因组上最优和次优路径集合；最终，根据最优和次优路径集合得到测序序列在基因组上对应的序列，采用KSW2算法比对测序序列与对应的序列，得到比对得分，根据比对得分得到结果比对文件。属于序列比对领域。到结果比对文件。属于序列比对领域。到结果比对文件。属于序列比对领域。

【技术实现步骤摘要】
一种基于群体基因组的序列比对方法


[0001]本专利技术涉及一种序列比对方法，具体涉及一种基于群体基因组的人体DNA序列比对方法，属于序列比对领域。

技术介绍

[0002]基因组测序数据分析已经成为了生命科学、精准健康和个性化医疗等领域的重要研究方向和发展重点。国际上已开展了包括国际千人基因组计划(1KG)、英国万人基因组计划(UK10K)、美国精准医学计划(All of Us)、新加坡万人基因组计划(SG10K)和中国十万人基因组计划在内的多项大规模人类基因组测序计划。基因组计划产生的测序数据规模将达到PB
‑
ZB级，基因组数据分析需求呈爆炸性增长。
[0003]基因组序列比对是基因组序列分析的一个重要基础步骤，同时也是消耗大量计算时间与资源的瓶颈步骤，在数据量到达PB量级的时候问题愈专利技术显。现有第二代测序技术不能直接生成完整基因组，因此“基因组序列比对”的目的是将大量无序的测序片段在参考基因组上定位，从而获取每个位点的变异信息的过程，而现有常用的以单一参考基因组作为基准的序列比对方法存在精度系统性偏差的问题。基因组序列比对上游承接测序数据的输入，下游接续变异分析与疾病关联分析，所以序列比对的质量对于下游分析有极大影响，不准确的序列比对会导致变异检测失败或给出错误的变异位点，直接影响下游分析的正确性与可信性。

技术实现思路

[0004]本专利技术为了解决以单一参考基因组作为基准的序列比对方法存在精度偏差，导致测序数据变异检测失败或得到错误变异位点的问题，进而提出了一种基于群体基因组的序列比对方法。
[0005]本专利技术采取的技术方案是：
[0006]它包括以下步骤：
[0007]S1、构建基于种子序列的群体基因组索引，基于种子序列的群体基因组索引包括基于种子序列的参考基因组的图模型索引和基于种子序列的变异数据的局部变异序列索引；
[0008]S2、提取测序序列的种子序列，利用局部变异序列索引获取种子序列在变异数据上的变异序列，得到变异序列在参考基因组上的位置，即种子序列在参考基因组上的位置；
[0009]S3、根据种子序列在参考基因组上的位置，利用稀疏动态规划方法在参考基因组上进行路径搜索，得到种子序列在参考基因组上的最优路径和次优路径集合；
[0010]S4、根据最优路径和次优路径集合得到S2的测序序列在参考基因组上对应的基因组序列，采用KSW2算法将测序序列与对应的基因组序列进行比对，得到比对得分，根据比对得分得到结果比对文件；
[0011]S5、获取待测序序列，执行S1
‑
S4对待测序序列进行比对，得到结果比对文件。
[0012]进一步的，S1中构建基于Minimizer的群体基因组索引，基于Minimizer的群体基因组索引包括基于参考基因组的图模型索引和基于变异数据的局部变异序列索引，具体过程为：
[0013]S11、获取人类参考基因组数据和对应的包含已知不同类型变异的变异数据，并对两种数据进行去除冗余数据等规范化操作预处理；
[0014]S12、构建参考基因组的de Bruijn图模型索引表示，得到基于参考基因组的de Bruijn图模型索引；
[0015]S13、构建变异数据的局部变异序列索引表示，得到变异数据的局部变异序列索引；
[0016]S14、根据de Bruijn图模型索引的唯一路径和局部变异序列索引构建基于Minimizer的群体基因组索引表示。
[0017]进一步的，S2中提取测序序列的种子序列，利用局部变异序列索引获取种子序列在变异数据上的变异序列，得到变异序列在参考基因组上的位置，即种子序列在参考基因组上的位置，具体过程为：
[0018]针对测序序列中每一对双端测序序列，对每端序列执行如下操作：
[0019]1)、定义滑窗，并提取种子序列Minimizer：
[0020]利用滑窗将一端序列分为多个固定长度的区间，每个区间内均包括多个k
‑
mer序列，定义一端序列上的相邻滑窗之间相互覆盖，且相邻滑窗的起始位置之间相差1bp；
[0021]从头开始遍历一端序列，根据定义的滑窗大小选取每个滑窗内ASC2值或哈希值最小的k
‑
mer序列作为种子序列Minimizer，得到一端序列上的所有种子序列Minimizer；
[0022]2)、获取种子序列Minimizer在参考基因组上的位置：
[0023]如果某一对双端测序序列方向相对的两端序列上找到的种子序列Minimizer的总数量小于预设阈值，预设阈值为测序序列所有种子序列Minimizer总数的1/3，则利用局部变异序列索引查找某一对双端测序序列的种子序列Minimizer在变异数据上对应的变异序列，若在变异数据上找到对应的变异序列，即表示所述变异序列与种子序列Minimizer一致，记录变异序列在参考基因组上对应位置的四元组{Unipath
start
，Unipath
end
，Read
start
，Read
end
}的集合，其中，Unipath
start
表示变异序列在参考基因组的de Bruijn图模型的唯一路径上的起始位置，Unipath
end
表示变异序列在参考基因组的de Bruijn图模型的唯一路径上的终止位置，Read
start
表示变异序列在测序序列上的起始位置，Read
end
表示变异序列在测序序列上的终止位置；如果在变异序列上未找到对应的变异序列，根据已得到的种子序列Minimizer在参考基因组位置上对应的四元组集合执行S3；
[0024]如果某一对双端测序序列方向相对的两端序列上找到的种子序列Minimizer的总数量大于等于预设阈值，得到种子序列Minimizer在参考基因组上对应的四元组集合，执行S3。
[0025]进一步的，S3中根据种子序列在参考基因组上的位置，利用稀疏动态规划方法在参考基因组上进行路径搜索，得到种子序列在参考基因组上的最优路径和次优路径集合，具体过程为：
[0026]S31、判断双端测序序列中每端测序序列的方向：
[0027]根据种子序列Minimizer在参考基因组上对应的四元组集合，定义双端测序序列
的paired
‑
end reads为L
1f
、L
1r
或L
2f
、L
2r
，L
1f
表示paired
‑
end reads的read1的正向序列覆盖长度，L
1r
表示paired
‑
end reads的read1为反向互补序列覆盖长度，L
2f
表示paired
‑
end reads的read2为正向序列覆盖长度，本文档来自技高网...

【技术保护点】

【技术特征摘要】
1.一种基于群体基因组的序列比对方法，其特征在于：它包括以下步骤：S1、构建基于种子序列的群体基因组索引，基于种子序列的群体基因组索引包括基于种子序列的参考基因组的图模型索引和基于种子序列的变异数据的局部变异序列索引；S2、提取测序序列的种子序列，利用局部变异序列索引获取种子序列在变异数据上的变异序列，得到变异序列在参考基因组上的位置，即种子序列在参考基因组上的位置；S3、根据种子序列在参考基因组上的位置，利用稀疏动态规划方法在参考基因组上进行路径搜索，得到种子序列在参考基因组上的最优路径和次优路径集合；S4、根据最优路径和次优路径集合得到S2的测序序列在参考基因组上对应的基因组序列，采用KSW2算法将测序序列与对应的基因组序列进行比对，得到比对得分，根据比对得分得到结果比对文件；S5、获取待测序序列，执行S1
‑
S4对待测序序列进行比对，得到结果比对文件。2.根据权利要求1中所述的一种基于群体基因组的序列比对方法，其特征在于：S1中构建基于Minimizer的群体基因组索引，基于Minimizer的群体基因组索引包括基于参考基因组的图模型索引和基于变异数据的局部变异序列索引，具体过程为：S11、获取人类参考基因组数据和对应的包含已知不同类型变异的变异数据，并对两种数据进行去除冗余数据等规范化操作预处理；S12、构建参考基因组的de Bruijn图模型索引表示，得到基于参考基因组的de Bruijn图模型索引；S13、构建变异数据的局部变异序列索引表示，得到变异数据的局部变异序列索引；S14、根据de Bruijn图模型索引的唯一路径和局部变异序列索引构建基于Minimizer的群体基因组索引表示。3.根据权利要求2中所述的一种基于群体基因组的序列比对方法，其特征在于：S2中提取测序序列的种子序列，利用局部变异序列索引获取种子序列在变异数据上的变异序列，得到变异序列在参考基因组上的位置，即种子序列在参考基因组上的位置，具体过程为：针对测序序列中每一对双端测序序列，对每端序列执行如下操作：1)、定义滑窗，并提取种子序列Minimizer：利用滑窗将一端序列分为多个固定长度的区间，每个区间内均包括多个k
‑
mer序列，定义一端序列上的相邻滑窗之间相互覆盖，且相邻滑窗的起始位置之间相差1bp；从头开始遍历一端序列，根据定义的滑窗大小选取每个滑窗内ASC2值或哈希值最小的k
‑
mer序列作为种子序列Minimizer，得到一端序列上的所有种子序列Minimizer；2)、获取种子序列Minimizer在参考基因组上的位置：如果某一对双端测序序列方向相对的两端序列上找到的种子序列Minimizer的总数量小于预设阈值，预设阈值为测序序列所有种子序列Minimizer总数的1/3，则利用局部变异序列索引查找某一对双端测序序列的种子序列Minimizer在变异数据上对应的变异序列，若在变异数据上找到对应的变异序列，即表示所述变异序列与种子序列Minimizer一致，记录变异序列在参考基因组上对应位置的四元组{Unipath
start
，Unipath
end
，Read
start
，Read
end
}的集合，其中，Unipath
start
表示变异序列在参考基因组的de Bruijn图模型的唯一路径上的起始位置，Unipath
end
表示变异序列在参考基因组的de Bruijn图模型的唯一路径上的终止位置，Read
start
表示变异序列在测序序列上的起始位置，Read
end
表示变异序列在测
序序列上的终止位置；如果在变异序列上未找到对应的变异序列，根据已得到的种子序列Minimizer在参考基因组位置上对应的四元组集合执行S3；如果某一对双端测序序列方向相对的两端序列上找到的种子序列Minimizer的总数量大于等于预设阈值，得到种子序列Minimizer在参考基因组上对应的四元组集合，执行S3。4.根据权利要求3中所述的一种基于群体基因组的序列比对方法，其特征在于：S3中根据种子序列在参考基因组上的位置，利用稀疏动态规划方法在参考基因组上进行路径搜索，得到种子序列在参考基因组上的最优路径和次优路径集合，具体过程为：S31、判断双端测序序列中每端测序序列的方向：根据种子序列Minimizer在参考基因组上对应的四元组集合，定义双端测序序列的paired
‑
end reads为L
1f
、L
1r
或L
2f
、L
2r
，L
1f
表示paired
‑
end reads的read1的正向序列覆盖长度，L
1r
表示paired
‑
end reads的read1为反向互补序列覆盖长度，L
2f
表示paired
‑
end reads的read2为正向序列覆盖长度，L
2r
表示paired
‑
end reads的read2为反向互补序列覆盖长度；若L
1f
+L
2r
＞L
1r
+L
2f
+偏移值，偏移值为测序序列长度的1/10，则认为正确的方向为read1为正向序列，来源于正链，read2为反向互补序列，来源于反链，返回对应方向包含的所有Minimizer的集合；若L
1r
+L
2f
＞L
1f
+L
2r
+偏移值，则认为正确的方向为read2为正向序列，来源于正链，read1为反向互补序列，来源于反链，返回对应方向包含的所有Minimizer的集合；否则，认为无法区分正确方向，返回两个方向对应的包含所有Minimizer的集合；如果两端测序序列的总覆盖长度小于指定阈值，指定阈值为两端测序序列总长度的1/3，则利用局部变异序列进行Minimizer比对，得到两端测序序列的方向，并返回对应方向包含的所有Minimizer集合；S32、根据每端测序序列的方向对其种子序列Minimizer进行延伸与合并：a、延伸种子序列Minimizer：获取当前种子序列Minimizer在测序序列上的起始位置和终止位置，以及当前种子序列Minimizer在deBruijn图模型中唯一路径上的起始位置和终止位置，将种子序列Minimize在测序序列和唯一路径上对应的位置均进行向前和向后的延伸，直至遇到测序序列上某碱基匹配错误或到达测序序列的最后一个碱基位置或到达唯一路径的最后一个碱基位置，记录种子序列Minimizer在测序序列和唯一路径上的前后延伸位置和延伸后在测序序列上的覆盖长度，按照上述方法，遍历测序序列上的所有的种子序列Minimizer；在遍历时，如果当前种子序列Minimizer在测序序列上的终止位置小于前一个种子序列Minimizer的向后延伸位置，说明当前种子序列Minimizer被包含在前一个种子序列Minimizer的延伸覆盖区域内，无需进行延伸；否则，执行上述的延伸；每完成一个种子序列Minimizer的延伸就得到一个候选种子序列，记录候选种子序列对应的唯一路径编号ID
Unipath
、在测序序列上的偏移位置、在参考基因组上的位置、在测序序列上的覆盖长度、候选种子序列编号ID
...

【专利技术属性】
技术研发人员：国宏哲，李高阳，刘博，王亚东，
申请(专利权)人：哈尔滨工业大学，
类型：发明
国别省市：