一种基因序列比对方法及系统技术方案

本发明专利技术公开一种基因序列比对方法及系统，该方法包括：将参考基因组序列和查询基因组序列存储于分布式存储系统；在Spark异构分布式计算平台框架下：按行偏移量分割参考基因组序列，并进行预处理，得到多个预处理参考数据集；采用后缀数组算法为每一预处理参考数据集建立索引，并合并所有建立索引后的预处理参考数据集，得到参考序列索引文件；采用种子扩展算法将查询基因组序列中的各个片段与参考序列索引文件进行CUDA细粒度序列比对，确定每个片段在参考序列索引文件中的位置信息；合并所有片段在参考序列索引文件中的位置信息得到基因序列比对结果。本发明专利技术提高了大规模序列比对算法计算的速度和精度。

【技术实现步骤摘要】
一种基因序列比对方法及系统
本专利技术涉及基因比对
，特别是涉及一种基于异构云平台的基因序列比对方法及系统。
技术介绍
自第一代基因组测序技术专利技术以来，基因组测序技术已经经过几次技术迭代。第一代基因组测序技术又被称为Sanger测序法。其特点是测序序列长，长度可达到1kbp(kilobasespair)，及准确率高，但是由于其价格昂贵及通量低的缺点，逐渐被第二代基因组测序技术(NextGenerationSequencing,NGS)取代。各国学者针对二代测序数据的特点，设计了各类序列比对、序列拼接及序列分类算法。在这些序列比对算法研究中，大部分都是基于Seed-and-Extend框架，如后缀数组，BWT等。随着二代测序技术的广泛应用，计算分子生物学家获得了宝贵的研究数据。由于第二代测序序列的长度较短，二代测序数据在特定的基因组学问题上的效果仍不理想。第三代测序技术使用单分子测序技术进行测序，其拥有测序序列长的优点，其测序序列平均可长达10kbp。由于，第三代测序技术产生的序列特征与第一二代测序数据显著不同，且拥有序列读长长的特点。因此全基因组测序可以发现更多的遗传变异信息，第三代测序的序列比对软件国外有：BLASR、Minimap2等。BLASR是第一个专门为三代测序而设计的序列比对方法。分别采用FM-index和后缀数组作为索引结构，通过索引查找测序片段在参考基因组上的最长公共前缀作为种子，随后通过串链的方法定位候选区域。但由于该索引的设计是面向精确比对的，因此在查找带有误配的序列时，计算性能较低。Minimap2使用种子索引(seed)、种子链接(chain)再进行序列比对(align)的策略。Minimap2是目前比对第三代序列最快的软件之一，是三代测序数据分析中主要对标的软件，但其序列比对精度不高。针对第三代测序技术，中国学者也开发了一些达到国际先进水平的序列比对工具，例如由哈工大生物信息研究院开发的LAMSA以及rHAT[20]等，rHAT算法在小规模基因组上性能不错，但是rHAT在大规模基因上速度较慢，有很大提升空间。综上，虽然国内外研究学者针对第三代测序序列比对提出了不同的算法，算法实现方法、实现细节等有所不同。但均存在因为比对规模大，导致索引建立时内存不足而比对结果精度不高，或者由于序列比对时间长而采用非精确匹配导致的匹配误差大等问题。因此，如何根据三代测序的数据特性，设计高效的序列比对算法，加速序列比对的过程，减少比对过程的时间消耗是三代测序序列比对面临的主要问题。三代数据在处理重复区域、复杂区域的序列比对，检测结构体变异，和重构转录本等方面有得天独厚的优势。第三代测序技术的reads长度较二代测序长，每个reads平均长度可达到10kbp。面向二代测序数据的算法已经不能满足三代测序数据分析的要求。再者，由于数据量大，计算复杂度高，一些应用程序的运行时间非常长，如SNP检测和GWAS比对也可能需要几天甚至几个月才能完成一个数据集的处理。而如何对这些测序得到的大量基因组数据进行处理分析并设计出高性能的处理算法成为亟需解决的问题。面向三代数据的序列比对算法近年来不断发展，但是无论在速度还是准确性方面，都有很大的提升空间。
技术实现思路
本专利技术的目的是提供一种基因序列比对方法及系统，以提高大规模序列比对算法计算的速度和精度。为实现上述目的，本专利技术提供了一种基因序列比对方法，该方法包括：将参考基因组序列和查询基因组序列存储于分布式存储系统；在Spark异构分布式计算平台框架下：按行偏移量分割所述参考基因组序列，并进行预处理，得到多个预处理参考数据集；采用后缀数组算法为每一所述预处理参考数据集建立索引，并合并所有建立索引后的预处理参考数据集，得到参考序列索引文件；采用种子扩展算法将所述查询基因组序列中的各个片段与所述参考序列索引文件进行CUDA细粒度序列比对，确定每个所述片段在所述参考序列索引文件中的位置信息；合并所有所述片段在所述参考序列索引文件中的位置信息得到基因序列比对结果。可选的，所述按行偏移量分割所述参考基因组序列，并进行预处理，得到多个预处理参考数据集，具体包括：将参考基因组按行偏移量利用zipWithIndex函数为每个字符串ri分配唯一标号，并按所述唯一标号将所述参考基因组切割成多个片段序列；利用sortByKey函数算子对多个所述片段序列按照所述片段序列在所述参考基因组中的位置进行排序，得到RDD1数据集；其中，所述片段序列在所述参考基因组中的位置为sortByKey函数算子中的key；并行地，利用RDD.sortByValue函数算子对多个所述片段序列按照所述片段序列在所述参考基因组中的行位置进行排序，得到RDD2数据集；采用RDD.joinByKey函数算子对所述RDD1数据集和所述RDD2数据集进行连接操作，得到连接数据集；在所述连接数据集中确定各个所述片段在参考基因组的相对位置，以所述相对位置为values对所述连接数据集中的各个所述片段进行排序，得到预处理初始数据集；利用Action算子触发持久化将所述预处理初始数据集存储到分布式存储系统中，得到预处理参考数据集。可选的，所述采用后缀数组算法为每一所述预处理参考数据集建立索引，并合并所有建立索引后的预处理参考数据集，得到参考序列索引文件，具体包括：将所述预处理参考数据集中的每一片段序列以预设偏移量分割成多个种子，所述种子的长度为k；将每一所述种子按A、G、C、T分为四类，并分别建立索引数组；对每一所述索引数组用<key,value>按value的值进行排序，得到参考序列索引文件。可选的，所述采用种子扩展算法将所述查询基因组序列中的各个片段与所述参考序列索引文件进行CUDA细粒度序列比对，确定每个所述片段在所述参考序列索引文件中的位置信息，具体包括：将所述查询基因组序列中的每个所述片段按最大允许误配k切分成查询片段种子；按照所述查询片段种子的首尾字母在所述参考序列索引文件匹配位置，作为候选位置；在所述候选位置上扩展所述查询片段种子，并将所述候选位置前后的附近序列与扩展后的查询片段种子进行近似匹配，判断所述候选位置是否为真实匹配位置，若是，则确定所述查询片段与参考序列片段相匹配，并将相匹配的所述查询片段与所述参考序列片段列举出；所述附近序列为按照所述查询片段种子的扩展长度确定的参考基因组序列；所述扩展长度小于或等于最大允许误配k。可选的，所述在所述候选位置上扩展所述查询片段种子，并将所述候选位置前后的附近序列与扩展后的查询片段种子进行近似匹配，判断所述候选位置是否为真实匹配位置，具体包括：构建(n+1)×(m+1)的得分矩阵H；所述得分矩阵H的元素为比对结果，n为查询片段的长度，m为参考序列片段的长度；得分矩阵中的元素值根据所述元素值的左面邻接点、上面邻接点和左上角面邻接点元素的值计算；即当矩阵H(i-1,j-1)、H本文档来自技高网...

【技术保护点】
1.一种基因序列比对方法，其特征在于，所述方法包括：/n将参考基因组序列和查询基因组序列存储于分布式存储系统；在Spark异构分布式计算平台框架下：/n按行偏移量分割所述参考基因组序列，并进行预处理，得到多个预处理参考数据集；/n采用后缀数组算法为每一所述预处理参考数据集建立索引，并合并所有建立索引后的预处理参考数据集，得到参考序列索引文件；/n采用种子扩展算法将所述查询基因组序列中的各个片段与所述参考序列索引文件进行CUDA细粒度序列比对，确定每个所述片段在所述参考序列索引文件中的位置信息；/n合并所有所述片段在所述参考序列索引文件中的位置信息得到基因序列比对结果。/n

【技术特征摘要】
1.一种基因序列比对方法，其特征在于，所述方法包括：
将参考基因组序列和查询基因组序列存储于分布式存储系统；在Spark异构分布式计算平台框架下：
按行偏移量分割所述参考基因组序列，并进行预处理，得到多个预处理参考数据集；
采用后缀数组算法为每一所述预处理参考数据集建立索引，并合并所有建立索引后的预处理参考数据集，得到参考序列索引文件；
采用种子扩展算法将所述查询基因组序列中的各个片段与所述参考序列索引文件进行CUDA细粒度序列比对，确定每个所述片段在所述参考序列索引文件中的位置信息；
合并所有所述片段在所述参考序列索引文件中的位置信息得到基因序列比对结果。


2.根据权利要求1所述的基因序列比对方法，其特征在于，所述按行偏移量分割所述参考基因组序列，并进行预处理，得到多个预处理参考数据集，具体包括：
将参考基因组按行偏移量利用zipWithIndex函数为每个字符串ri分配唯一标号，并按所述唯一标号将所述参考基因组切割成多个片段序列；
利用sortByKey函数算子对多个所述片段序列按照所述片段序列在所述参考基因组中的位置进行排序，得到RDD1数据集；其中，所述片段序列在所述参考基因组中的位置为sortByKey函数算子中的key；
并行地，利用RDD.sortByValue函数算子对多个所述片段序列按照所述片段序列在所述参考基因组中的行位置进行排序，得到RDD2数据集；
采用RDD.joinByKey函数算子对所述RDD1数据集和所述RDD2数据集进行连接操作，得到连接数据集；
在所述连接数据集中确定各个所述片段在参考基因组的相对位置，以所述相对位置为values对所述连接数据集中的各个所述片段进行排序，得到预处理初始数据集；
利用Action算子触发持久化将所述预处理初始数据集存储到分布式存储系统中，得到预处理参考数据集。


3.根据权利要求1所述的基因序列比对方法，其特征在于，所述采用后缀数组算法为每一所述预处理参考数据集建立索引，并合并所有建立索引后的预处理参考数据集，得到参考序列索引文件，具体包括：
将所述预处理参考数据集中的每一片段序列以预设偏移量分割成多个种子，所述种子的长度为k；
将每一所述种子按A、G、C、T分为四类，并分别建立索引数组；
对每一所述索引数组用<key,value>按value的值进行排序，得到参考序列索引文件。


4.根据权利要求1所述的基因序列比对方法，其特征在于，所述采用种子扩展算法将所述查询基因组序列中的各个片段与所述参考序列索引文件进行CUDA细粒度序列比对，确定每个所述片段在所述参考序列索引文件中的位置信息，具体包括：
将所述查询基因组序列中的每个所述片段按最大允许误配k切分成查询片段种子；
按照所述查询片段种子的首尾字母在所述参考序列索引文件匹配位置，作为候选位置；
在所述候选位置上扩展所述查询片段种子，并将所述候选位置前后的附近序列与扩展后的查询片段种子进行近似匹配，判断所述候选位置是否为真实匹配位置，若是，则确定所述查询片段与参考序列片段相匹配，并将相匹配的所述查询片段与所述参考序列片段列举出；所述附近序列为按照所述查询片段种子的扩展长度确定的参考基因组序列；所述扩展长度小于或等于最大允许误配k。


5.根据权利要求4所述的基因序列比对方法，其特征在于，所述在所述候选位置上扩展所述查询片段种子，并将所述候选位置前后的附近序列与扩展后的查询片段种子进行近似匹配，判断所述候选位置是否为真实匹配位置，具体包括：
构建(n+1)×(m+1)的得分矩阵H；所述得分矩阵H的元素为比对结果，n为查询片段的长度，m为参考序列片段的长度；
得分矩阵中的元素值根据所述元素值的左面邻接点、上面邻接点和左上角面邻接点元素的值计算；即当矩阵H(i-1,j-1)、H(i,j-1)和H(i-1,j)值计算结束后，H(i,j)值才能计算；在记录矩阵元素H(i,j)值的同时，记录下得到H(i,j)的路径，建立路径矩阵T；
利用计算得到的第g-1条反对角线和第g-2条反对角线上的元素值，执行并行计算得到第...

【专利技术属性】
技术研发人员：郑芳，赵良，田芳，倪福川，汪毅，姚雅鹃，姚娟，章程，
申请(专利权)人：华中农业大学，
类型：发明
国别省市：湖北;42