一种面向大规模基因数据的读段定位方法
一种面向大规模基因数据的读段定位方法
【专利摘要】本发明公开了一种面向大规模基因数据的读段定位方法；属于生物信息分析领域。该方法包括如下步骤：基因读段数据随机分割；数据的负载平衡；读段的空间索引；子读段不跨越剪切位定位；子读段跨越剪切位定位；子读段拼接；读段定位信息统计。本发明采用基于空位种子的读段定位方法，对该方法进行了改进，使之能够处理跨越剪切位读段定位的情况，设计了并行化程序并且在MapReduce框架下实现，提高了读段定位的效率。此外，本发明提出了一种负载平衡的解决方案，对个别节点执行时间较长，降低整个作业执行效率的情况具有良好的处理能力，因此具有较高的使用价值。
【专利说明】一种面向大规模基因数据的读段定位方法
【技术领域】
[0001]本发明属于生物信息分析领域，特别涉及一种面向大规模基因数据的读段定位方法。
【背景技术】
[0002]基因序列匹配(读段定位)是生物信息学的一项重要研究课题，也是基因数据分析的关键环节。基因序列需要在进行读段定位过程之后，才能展开后续一系列的生物信息分析的步骤，如基因表达水平评估、选择性剪切事件识别、聚类等等。因此，基因读段定位已受到生物信息领域众多研究者的广泛关注。
[0003]新一代基因测序技术的快速发展产生了海量的基因读段数据，这给传统的读段定位算法带来了极大的挑战。庞大的基因读段数据使用传统串行算法在单核机器上运行，已经远远不能满足研究者的需求。一般情况下定位几十GB的基因读段数据需要一周左右的时间，甚至更长，才能获得最终的定位结果，从而导致了整个基因数据分析流程过于缓慢。目前的读段长度一般较长，使用传统不跨越剪切位的空位种子算法可能会遗漏一些正确的定位，得到的定位结果往往不够精确。此外，在并行化过程中总会存在某个节点的执行时间较长，使得整个作业等待个别节点，这种情况严重影响了读段定位的效率。目前的研究热点也是针对如何解决上述问题而展开。
[0004]典型的跨越剪切位的读段定位方法包括数据集分割、读段分割、子读段不跨越剪切位定位、子读段跨越剪切位定位和拼接。良好的分割方法和高效的子读段跨越剪切位定位能在一定程度上消除负载不均衡的情况以及降低算法时间复杂度，达到提高读段定位效率的效果。
【发明内容】
[0005]本发明为了解决基因读段定位过程中存在的效率低下的问题以及数据分割中负载不均衡问题，提出了一种面向大规模基因数据的读段定位方法，能有效地提高基因读段定位的效率。
[0006]本发明采用的技术方案如下:
[0007]—种面向大规模基因数据的读段定位方法，包括如下步骤:
[0008]一种面向大规模基因数据的读段定位方法，其特征在于，包括如下步骤:
[0009]步骤1，基因读段数据随机分割:将给定的基因读段数据集随机分割成指定的块数；
[0010]步骤2，数据的负载平衡:利用MapReduce (映射归约)框架，从分割好的数据块中随机抽取少量基因读段尝试执行，检测出可能执行时间较长的数据块，将其与其他数据块重新分割，均衡负载；
[0011]步骤3，读段的空间索引:采用MapReduce框架，在Map过程中将读段分割成指定段数的子读段，对每个子读段生成空位种子模式，将每个模式的哈希值与对应的读段信息记录在哈希表中，构建读段的空间索引；
[0012]步骤4，子读段不跨越剪切位定位:在Map过程中利用已有的读段定位软件将子读段中能够进行连续定位的读段首先定位，记录定位信息；
[0013]步骤5，子读段跨越剪切位定位:经过步骤4定位后，能够连续定位的子读段已经获得了定位信息，在Map过程中根据已经定位的子读段信息，推断不能连续定位的子读段的定位信息，实现跨越剪切位的读段定位；
[0014]步骤6，子读段拼接:经过步骤5所有读段分割成的子读段的定位信息已全部被记录，在Map过程中将所有定位的子读段进行组装，记录所有能够连接成原始读段的定位信息；
[0015]步骤7，读段定位信息统计:在Reduce过程中将所有通过Map过程传递过来的读段定位的信息进行汇总，统计读段定位在参考序列上的整体信息。
[0016]所述步骤2数据的负载平衡具体包括如下步骤=MapReduce框架将整个数据集默认分割为N块，记为(S1, S2,, SJ，从每块中随机抽取一定数量的基因数据，形成相应的子块，执行读段定位算法，检测出潜在的可能执行时间较长的分块；重新分割检测出的数据块，并且将其分别与其他未分割的数据块合并，最后再次分割成指定的块数。
[0017]所述步骤3读段的空间索引具体包括如下步骤:在Map过程中，将对应每个数据集子块中的读段分割为M个子读段，对每个子读段创建种子模式；每一个种子模式必然穷尽了所有可能的含有指定错配个数碱基的方式，将每一个种子模式转换成哈希值，以键值对的形式将哈希值和读段信息存储在哈希表里；在定位过程中就可以快速的从这张表中通过读段的哈希值获取指定读段的信息，这张哈希表就是所有读段的空间索引。
[0018]所述步骤5子读段跨越剪切位定位具体包括如下步骤:在每一个Map中，根据读段的长度，将子读段能够跨越的剪切位点的个数限定在指定的数量以内；根据已有的生物信息，限定内含子最大长度，提高定位效率；根据已经定位的子读段信息，通过延展尝试匹配的方式来确定不能够连续定位的子读段位置信息；如果原始读段的所有子读段都能够进行定位，则保留定位信息；否则舍弃该读段。
[0019]所述步骤6子读段拼接具体包括如下步骤:首先，每个子读段都记录它能够定位到参考序列上的所有可能位置，用一个数组存储所有位置信息；然后，从第一个子读段的第一个定位信息，即从第一个子读段的位置数组中取出第一个位置，与第二个子读段的所有位置信息进行检测，查看是否与其中一个位置具有顺接关系；如果有，则继续检测第二个子读段顺接的那个位置与第三个子读段所有位置信息的关系，直至最后一个子读段；如果没有，则停止检测；接着继续从第一个子读段的第二个位置进行检测，直至检查完第一个子读段的所有位置信息；最后，提取所有能够完整匹配的读段定位信息。
[0020]本发明是专门针对基因读段定位而提出的方法。与现有技术相比，本发明具有以下特征:
[0021](I)本发明针对现有基因数据量大，传统串行基因读段定位算法效率低下这一问题，提出了并行化的读段定位算法，并且利用MapReduce框架加以实现，在不影响读段定位正确性的条件下，提高了基因读段定位的效率；
[0022](2)本发明对基于空位种子的不跨越剪切位的读段定位算法进行改进，使之能够处理跨越剪切位的读段定位情况；此外，加入了生物信息:内含子的长度和剪切位点的数量，缩小了算法的搜索空间；改进读段定位时的算法，降低了算法的时间复杂度，从而进一步提高了基因读段定位的效率。
[0023](3)本发明在读段定位的过程中，专门考虑了负载平衡问题，提出了一种实际的解决方案，有效处理了个别节点执行时间长，降低整体效率的情况，因此具有较高的使用价值。
【专利附图】
【附图说明】
[0024]图1为本发明整体流程图；
[0025]图2为本发明中数据的负载平衡步骤子流程图；
[0026]图3为本发明中子读段不跨越剪切位定位步骤子流程图；
[0027]图4为本发明中子读段跨越剪切位定位步骤子流程图。
【具体实施方式】
[0028]以下结合【专利附图】
【附图说明】本发明的【具体实施方式】。
[0029]如图1所示，本发明公开了一种面向大规模基因数据的读段定位方法，具体步骤如下:
[0030]步骤1，基因读段数据随机分割:将给定的基因读段数据集随机分割成指定的块数；
[0031]步骤2，数据的负载平衡:采用MapReduce框架，从分割好的数据块中随机抽取少量基因读段尝试执行，检测出可能执行时间较长的数据块，将其重新分割，分别添加到其余数据块中，最后将其余的数据块再次进行分割，实现数据的均衡负载；
[0032]步骤3，读段的空间索引:采用MapReduce框架，在Map过程中将读段分割成指定段数的子读段，对每个子读段生成空位种子模式，将每个模式的哈希值与对应的读段信息记录在哈希表中，构建读段的空间索引；
[0033]步骤4，子读段不跨越剪切位定位:在Map过程中利用已有的读段定位软件将子读段中能够进行连续定位的读段首先定位，记录定位信息；
[0034]步骤5，子读段跨越剪切位定位:在Map过程中根据已经定位的子读段信息，推断不能连续定位的子读段的定位信息，实现跨越剪切位的读段定位；
[0035]步骤6，子读段拼接:在Map过程中将所有定位的子读段进行组装，记录所有能够连接成原始读段的定位信息；
[0036]步骤7，读段定位信息统计:在Reduce过程中将所有通过Map过程传递过来的读段定位的信息进行汇总，统计读段定位在参考序列上的整体信息，如定位数量、位置等；
[0037]假设读段定位过程中允许错配数量为2。
[0038]步骤I中，基因读段数据随机分割可以利用MapReduce框架实现。假设随机将读段数据集分割为N块，则可以指定Reduce的数量为N ;而在Map中可以设置一个能够产生1_N之间数的随机产生器，每次遍历一条读段时，产生一个随机数，该随机数对应着Reduce，表示将该条读段分配到相应的Reduce，最后只需要取出Reduce中的读段就完成数据分割了。
[0039] 步骤2中，利用MapReduce框架将整个数据集默认分割为N ±夹，记为(SpS2,...-,:!，从每块(64MB)中随机抽取一定数量的基因数据，形成相应的子块，执行读段定位算法，检测出潜在的可能执行时间较长的分块(假设检测出的是S1, S2)。
[0040]重新分割数据块，平衡负载:将S1, S2各自分割成N块，分别记为{Sn, S12,S1J，{S21, S22,S2J。合并其余N-2 ±夹，再重新分割成N ±夹，记为{S31, S32,...，S3J。将S11, S21, S31合并成Y i，按照同样的合并方法得到N个最终分块，记为{S' 2,...,S/ N}。
[0041]其中从每块中随机抽取一定数量的样本可以通过一个MapReduce任务实现，具体步骤:
[0042]Map 阶段:
[0043](I)根据路径读取读段数据集；
[0044](2) while (抽取的读段数量已满足要求)
[0045]①随机产生一个读段下标；
[0046]②将随机抽取到的读段送入reduce ；
[0047]③删除已抽过的读段；
[0048]④抽取的读段数量加I ;
[0049](3)结束 while 循环。
[0050]Reduce 阶段:
[0051](I)遍历定位数量集合；
[0052]①输出参考序列名称和最终的定位数量；
[0053](2)结束。
[0054]步骤2流程如图2所示，从各个子读段数据块中随机抽取一部分读段，做为测试数据集；在MapReduce框架下以得到的子数据集做为输入，执行读段定位算法，以检测出可以影响整体执行时间的个别节点；以Web浏览方式查看每一个Map节点的执行时间，根据时间的长度推断可能负载较大的节点；将检测出的负载较大节点中的读段数据集重新分割；将重新分割的子数据集添加到其他未分割的数据集中，达到平衡负载的目的。
[0055]步骤3中，在Map过程中，将对应每个数据集子块中的所有读段都分割为N个子读段。
[0056]针对每个子读段创建种子模式，主要思想如下:假设将子读段再分割为连续不交叉的k份，最多允许m个错配(m〈k)，对于每一个低于错配指定值的匹配，在已分的k个部分中，至少有k-m个部分没有错配。这k-m个部分属于k个部分，因此共有Cf"组合。而可以选择的种子模式就应该有01个。每一个种子模式必然穷尽了所有可能的含有指定错配个数碱基的方式。
[0057]将每一个种子模式转换成哈希值，以键值对的形式将哈希值和对应原始读段信息存储在哈希表里。在定位过程中就可以快速的从这张表中通过读段的哈希值获取指定读段的信息，这张哈希表就是所有读段的空间索引。
[0058]步骤4中，在Map过程中利用已有的读段定位软件将子读段中能够进行连续定位的读段首先定位，记录定位信息。本发明中采用的是SeqMap软件中的读段定位算法，将读段分割为4个子读段，错配2个，则每个子读段的种子模式就有6种。构建所有子读段的空间索引，遍历所有参考序列，将能够连续定位的子读段定位到参考序列上。
[0059]步骤4流程如图3所示，先将子读段再次进行分割，本发明中分割为4段，创建6中空位种子模式；然后将所有种子模式转换为哈希值，并且将哈希值与对应的读段信息存储到哈希表中，形成所有子读段的空间索引；依次遍历参考序列中的每一个碱基；判断步骤15中形成的哈希值在空间索引中是否存在；如果空间索引中存在该哈希值，则将对应的读段取出来与参考序列中的碱基进行进一步匹配；如果不存在，则中断此次尝试匹配；继续遍历下个碱基。
[0060]步骤5中，经过步骤4定位后，能够连续定位的子读段已经获得了定位信息；
[0061]在每一个Map中，根据读段的长度，将子读段能够跨越的剪切位点的个数限定在指定的数量以内(假设读段长84bp，可以限定最多允许跨越2个剪切位点)，因为跨越剪切位点数量越多，算法就需要更多的时间去检测更多可能的定位方式；根据已有的生物信息，限定内含子最大长度(比如已知的哺乳动物基因组中只有非常少的内含子长度大于400Kb)，这样可以减少算法的搜索空间，提高定位效率；
[0062]根据已经定位的子读段信息，通过延展尝试匹配的方式来确定不能够连续定位的子读段位置信息。跨越剪切位的读段定位可以分为两种情况:
[0063]第一种情况是未定位子读段的上下游子读段都已经定位，这种情况是比较容易处理的。本发明中分别从上游读段的最后一个碱基和下游读段的第一个碱基开始匹配，分别记录下错配0、1、2的前一个位置，在记录下的位置中查找能够连接的，且错配最少的位置，记为剪切位点。通过这种只需一次遍历查找到剪切位点的方法，可以之间将这一过程的时间复杂度降低至O (N)。
[0064]第二种情况是未定位子读段的上下游子段只有一个已经定位。对于这种情况，本发明中截取未定位子读段末尾一段长度，记为h-mer ;首先通过检测h-mer在上游子读段或下游子读段中的位置，然后再利用第一种情况的处理方式进行读段定位。而在本发明中h-mer的长度取值区别于一般算法中的取值(直接将h-mer长度固定为2)。本发明动态获取h-mer的长度，首先通过已定位的上游或下游子读段进行延展匹配直至出现2个错配的位置，记录此时延展的长度L，h-mer的长度就取为未定位子读段长度减去L。通过这一做法能够尽可能减少未定位子读段尝试匹配次数，并且使得改进算法将时间复杂度O从原来的0(N3)降低至0(N2)。
[0065]进行上述跨越剪切位定位后，如果原始读段的所有子读段都能够进行定位，则保留定位信息；否则舍弃该读段。
[0066]步骤5流程如图4所示，用已有的基于空位种子的读段定位软件将能够连续定位的子读段定位，记录位置信息；根据已经定位的子读段信息，推断未定位的子读段；未定位子读段的可能的第一种情况，它的上游和下游子读段均已定位成功，此时通过延展上下游子读段来确定未定位子读段的位置信息；未定位子读段可能的第二种情况，它的上游或下游子读段之一未定位，这种情况下可以截取未定位子读段末尾一段长度，检测它在上游子读段或下游子读段中的位置，这样就将第二种情况转化为第一种情况了，再通过第一种情况的处理方式就可以完成定位；拼接子读段步骤将所有子读段定位信息进行分析，拼接成完整读段定位信息。
[0067]步骤6中，经过步骤5所有读段分割成的子读段的定位信息已全部被记录；
[0068]首先，每个子读段都记录它能够定位到参考序列上的所有可能位置(可以多余I个)，可以用一个数组存储所有位置信息；[0069] 然后，从第一个子读段的第一个定位信息，即从第一个子读段的位置数组中取出第一个位置，与第二个子读段的所有位置信息进行检测，查看是否与其中一个位置具有顺接关系；如果有，则继续检测第二个子读段顺接的那个位置与第三个子读段所有位置信息的关系，直至最后一个子读段；如果没有，则停止检测；接着继续从第一个子读段的第二个位置进行检测，直至检查完第一个子读段的所有位置信息；
[0070]最后，提取所有能够完整匹配的读段定位信息。
[0071 ] 步骤7中，在Reduce过程中将所有通过Map过程传递过来的读段定位的信息进行汇总，统计读段定位在参考序列上的整体信息，如定位数量、位置等。
[0072]实施例:
[0073]本实施例所使用的云计算平台共有21台服务器(I个master节点，20个slave节点)，每个节点配置相同。节点的处理器为Intel (R)Xeon(R)CPU E5620，主频为2.40GHZ，操作系统为64位Debian Linux操作系统。Hadoop平台版本为hadoop-0.20.2。
[0074]采用的基因数据是从拟南芥菜基因组中抽取10条基因做为参考序列，平均长度约为12000bp。读段数据集主要有两个来源:1.用ES培养液培养拟南芥菜，重复实验3次，获得三次情况下的拟南芥菜RNA，经测序得到3个读段样本集。2.在缺P情况下培养拟南芥菜，同样重复实验3次，获得3个读段样本集。两种情况下获得的读段数据集共25.2GB,我们分别称这6个读段样本集为ES1、ES2、ES3、P1、P2和P3。
[0075]本实施例包括以下部分:
[0076]1.基因读段数据随机分割:
[0077]基因读段数据随机分割可以利用MapReduce框架实现。假设随机将读段数据集分害I]为N ±夹，则可以指定Reduce的数量为N ;而在Map中可以设置一个能够产生1-N之间数的随机产生器，每次遍历一条读段时，产生一个随机数，该随机数对应着Reduce，表示将该条读段分配到相应的Reduce，最后只需要取出Reduce中的读段就完成数据分割了。主要代码如下:
[0078]Map 阶段:
[0079](l)read_data = get_read_dataset(read_path)；
[0080](2) random = Math, rand () *N+1 ；
[0081](3)context, write(random, read_info)；
[0082]Reduce 阶段
[0083](I) for (count = 0 ；count<read_list.length ；count++)
[0084]① read = read_list.get (count)；
[0085]② context, write (task_num, read)；
[0086](2) end.[0087]2.数据的负载平衡:
[0088]利用MapReduce框架将整个数据集默认分割为N ±夹，记为(S1, S2, , Sj，从每块(64MB)中随机抽取一定数量的基因数据，形成相应的子块，执行读段定位算法，检测出潜在的可能执行时间较长的分块。
[0089]主要代码如下:
[0090](I) split_to_N-2 (S1, S2, N_2) ^/S1, S2 是潜在负载较高的块；[0091](2)add(list_S1, list_S2, N_2) ；//将重新分割的块添加到其余N_2块中；
[0092](3) re split (read_N-2_path, N) ；// 将 N_2 块重新划分为 N 块。
[0093]其中从每块中随机抽取一定数量的样本可以通过一个MapReduce任务实现，具体步骤:
[0094]Map 阶段:
[0095](I)read_data = get_readdata(read_path)；
[0096](2)while(read_num< = fix_readNum)
[0097]① random = Math, rand ()*read_totalNum+l ；
[0098]② read = read_data.get (random)；
[0099]③ context, write(fix_reduce, read)；；
[0100]④ read_data.move (random)；
[0101]⑤read_num++;
[0102](3) end whi le.[0103]Reduce 阶段:
[0104](I)for(count = 0 ；count<read_data.length ；count++)；
[0105]O lcontext.write (read_name, map_num)；
[0106](2) end.[0107]3.读段的空间索引:
[0108]本发明中读段的长度是84bp，至多允许读段跨越2个(根据生物信息学的经验在84bp长度的读段中很少会有跨越2个剪切位点以上的读段)剪切位点。因此，在Map过程中，将对应每个数据集子块中的所有读段都分割为4个(适中)子读段。
[0109]针对每个子读段创建种子模式，这种种子模式来自于SeqMap软件，具体可以参考文献:Jiang H, Wong W H.SeqMap:mapping massive amount of oligonucleotides to thegenome[J].Bioinformatics, 2008, 24(20):2395-2396.[0110]主要代码如下:
[0111]Map 阶段:
[0112](I) list_sub_read = split_read (read, f ix_read_num) ；// 将读段分割为指定段数；
[0113](2) list_seed_modes = create_modes (list_sub_read) -J/ 创建种子模式；
[0114](3) for (count = 0 ；count<list_seed_modes.length ；count++)
[0115]①hashcode = hashcode(list_seed_modes.get(count))；
[0116]②hash_map.put(hashcode, read_info)；
[0117]Il将哈希值和对应原始读段信息存储在哈希表里；
[0118](3) end.[0119]4.子读段不跨越剪切位定位:
[0120]在Map过程中利用已有的读段定位软件将子读段中能够进行连续定位的读段首先定位，记录定位信息。本发明中采用的是SeqMap软件中的读段定位算法，将读段分割为4个子读段，错配2个，则每个子读段的种子模式就有6种。构建所有子读段的空间索引，遍历所有参考序列，将能够连续定位的子读段定位到参考序列上。[0121]主要代码如下:
[0122](I) for (int count = O ;count〈ref—list, length ;count++)//遍历参考序列数据
集
[0123]① ref = ref—list, get (count)；
[0124]②for (int 1c = 0 ；loc<ref.length ;loc++)// 遍历参考序列每个减基
[0125]a.seed—mode = create—mode (ref.substring (loc，loc+84))；
[0126]b.1f (hash—map.get (seed—mode) ! = null) match (ref, read)；
[0127]③ end;
[0128](2) end.[0129]其中math()函数代码如下:
[0130](I) for (int count = 0 ；count<read.length ；count++)
[0131]① if( ! read, get (count).equals (ref.get (count))) break ；
[0132]else context, write (read—name，one)；
[0133](2) end.[0134]5.子读段跨越剪切位定位:
[0135]经过步骤4定位后，能够连续定位的子读段已经获得了定位信息；根据已经定位的子读段信息，通过延展尝试匹配的方式来确定不能够连续定位的子读段位置信息。
[0136]主要代码如下:
[0137](I) if (pre—subread.isFixed&&next_subread.1sFixed) // 第一种情况
[0138]① brige—match (read，ref)；
[0139](2) end ；
[0140](3) if ( ! pre—subread.isFixed | | ! next—subread.isFixed))// 第二种情况
[0141]①h—mer = read, substring (read.length-L) ；// 动态获取 h—mer
[0142]②loc = search (h_mer, ref) ；// 查找 h—mer 位置
[0143]③ brige—match (h—mer，ref)；
[0144](4) end.[0145]其中brige_match()函数代码如下:
[0146](I) for (int loc = 0 ；lco<read.length ；loc++)
[0147]① if ( ! read, get (loc).equals (ref.get (pre—subread.length+loc)))
[0148]a.mis—num++ ；
[0149]b.list—loc = add (loc)；
[0150]②if (mis—num>2) break ；
[0151](2) end.[0152](3) for (int loc = read, length ；loc> = 0 ；loc—)
[0153]① index = I ；
[0154]② if ( ! read, get (loc).equals (ref.get (next—subread.loc-1ndex)))
[0155]a.mis—num++ ；
[0156]b.list—loc = add (loc)；
[0157]③if (mis—num>2) break ；[0158](4) end,
[0159](5) spliced_site = find_min_mistake (list_loc).// 查找错配最少的连接点
[0160]6.子读段拼接，提取所有能够完整匹配的读段定位信息。
[0161]7.统计读段定位在参考序列上的整体信息。
[0162]以上对本发明所提供的一种面向大规模基因数据的读段定位方法进行了详细介绍。值得注意的是，具体实现该技术方案的方法和途径有很多，以上所述仅是本发明的优选实施方式，只用于帮助理解本发明的方法及核心思想；同时，对于本领域的一般技术人员，在本发明核心思想的基础上，做出的修改和调整都将视为本发明的保护范围。综上所述，本说明书内容不应理解为对本发明的限制，本发明的保护范围应由所附的权利要求来限定。
【权利要求】
1.一种面向大规模基因数据的读段定位方法，其特征在于，包括如下步骤: 步骤1，基因读段数据随机分割:将给定的基因读段数据集随机分割成指定的块数； 步骤2，数据的负载平衡:利用MapReduce框架，从分割好的数据块中随机抽取少量基因读段尝试执行，检测出可能执行时间较长的数据块，将其与其他数据块重新分割，均衡负载； 步骤3，读段的空间索引:采用MapReduce框架，在Map过程中将读段分割成指定段数的子读段，对每个子读段生成空位种子模式，将每个模式的哈希值与对应的读段信息记录在哈希表中，构建读段的空间索引； 步骤4，子读段不跨越剪切位定位:在Map过程中利用已有的读段定位软件将子读段中能够进行连续定位的读段首先定位，记录定位信息； 步骤5，子读段跨越剪切位定位:经过步骤4定位后，能够连续定位的子读段已经获得了定位信息，在Map过程中根据已经定位的子读段信息，推断不能连续定位的子读段的定位信息，实现跨越剪切位的读段定位； 步骤6，子读段拼接:经过步骤5所有读段分割成的子读段的定位信息已全部被记录，在Map过程中将所有定位的子读段进行组装，记录所有能够连接成原始读段的定位信息； 步骤7，读段定位信息统计:在Reduce过程中将所有通过Map过程传递过来的读段定位的信息进行汇总，统计读段定位在参考序列上的整体信息。
2.根据权利要求1所述的一种面向大规模基因数据的读段定位方法，其特征在于，所述步骤2数据的负载平衡具体包括如下步骤: MapReduce框架将整个数据集默认分割为N ±夹，记为(S1, S2, , Sj，从每块中随机抽取一定数量的基因数据，形成相应的子块，执行读段定位算法，检测出潜在的可能执行时间较长的分块；重新分割检测出的数据块，并且将其分别与其他未分割的数据块合并，最后再次分割成指定的块数。
3.根据权利要求1所述的一种面向大规模基因数据的读段定位方法，其特征在于，所述步骤3读段的空间索引具体包括如下步骤: 在Map过程中，将对应每个数据集子块中的读段分割为M个子读段，对每个子读段创建种子模式；每一个种子模式必然穷尽了所有可能的含有指定错配个数碱基的方式，将每一个种子模式转换成哈希值，以键值对的形式将哈希值和读段信息存储在哈希表里；在定位过程中就可以快速的从这张表中通过读段的哈希值获取指定读段的信息，这张哈希表就是所有读段的空间索引。
4.根据权利要求1所述的一种面向大规模基因数据的读段定位方法，其特征在于，所述步骤5子读段跨越剪切位定位具体包括如下步骤: 在每一个Map中，根据读段的长度，将子读段能够跨越的剪切位点的个数限定在指定的数量以内；根据已有的生物信息，限定内含子最大长度，提高定位效率； 根据已经定位的子读段信息，通过延展尝试匹配的方式来确定不能够连续定位的子读段位置信息；如果原始读段的所有子读段都能够进行定位，则保留定位信息；否则舍弃该读段。
5.根据权利要求1所述的一种面向大规模基因数据的读段定位方法，其特征在于，所述步骤6子读段拼接具体包括如下步骤:首先，每个子读段都记录它能够定位到参考序列上的所有可能位置，用一个数组存储所有位置信息； 然后，从第一个子读段的第一个定位信息，即从第一个子读段的位置数组中取出第一个位置，与第二个子读段的所有位置信息进行检测，查看是否与其中一个位置具有顺接关系；如果有，则继续检测第二个子读段顺接的那个位置与第三个子读段所有位置信息的关系，直至最后一个子读段；如果没有，则停止检测；接着继续从第一个子读段的第二个位置进行检测，直至检查完第一个子读段的所有位置信息； 最后，提 取所有能够完整匹配的读段定位信息。
【文档编号】G06F19/24GK103971031SQ201410185387
【公开日】2014年8月6日 申请日期:2014年5月4日 优先权日:2014年5月4日
【发明者】杨明, 涂金金, 高阳 申请人:南京师范大学
【专利摘要】本发明公开了一种面向大规模基因数据的读段定位方法；属于生物信息分析领域。该方法包括如下步骤：基因读段数据随机分割；数据的负载平衡；读段的空间索引；子读段不跨越剪切位定位；子读段跨越剪切位定位；子读段拼接；读段定位信息统计。本发明采用基于空位种子的读段定位方法，对该方法进行了改进，使之能够处理跨越剪切位读段定位的情况，设计了并行化程序并且在MapReduce框架下实现，提高了读段定位的效率。此外，本发明提出了一种负载平衡的解决方案，对个别节点执行时间较长，降低整个作业执行效率的情况具有良好的处理能力，因此具有较高的使用价值。
【专利说明】一种面向大规模基因数据的读段定位方法
【技术领域】
[0001]本发明属于生物信息分析领域，特别涉及一种面向大规模基因数据的读段定位方法。
【背景技术】
[0002]基因序列匹配(读段定位)是生物信息学的一项重要研究课题，也是基因数据分析的关键环节。基因序列需要在进行读段定位过程之后，才能展开后续一系列的生物信息分析的步骤，如基因表达水平评估、选择性剪切事件识别、聚类等等。因此，基因读段定位已受到生物信息领域众多研究者的广泛关注。
[0003]新一代基因测序技术的快速发展产生了海量的基因读段数据，这给传统的读段定位算法带来了极大的挑战。庞大的基因读段数据使用传统串行算法在单核机器上运行，已经远远不能满足研究者的需求。一般情况下定位几十GB的基因读段数据需要一周左右的时间，甚至更长，才能获得最终的定位结果，从而导致了整个基因数据分析流程过于缓慢。目前的读段长度一般较长，使用传统不跨越剪切位的空位种子算法可能会遗漏一些正确的定位，得到的定位结果往往不够精确。此外，在并行化过程中总会存在某个节点的执行时间较长，使得整个作业等待个别节点，这种情况严重影响了读段定位的效率。目前的研究热点也是针对如何解决上述问题而展开。
[0004]典型的跨越剪切位的读段定位方法包括数据集分割、读段分割、子读段不跨越剪切位定位、子读段跨越剪切位定位和拼接。良好的分割方法和高效的子读段跨越剪切位定位能在一定程度上消除负载不均衡的情况以及降低算法时间复杂度，达到提高读段定位效率的效果。
【发明内容】
[0005]本发明为了解决基因读段定位过程中存在的效率低下的问题以及数据分割中负载不均衡问题，提出了一种面向大规模基因数据的读段定位方法，能有效地提高基因读段定位的效率。
[0006]本发明采用的技术方案如下:
[0007]—种面向大规模基因数据的读段定位方法，包括如下步骤:
[0008]一种面向大规模基因数据的读段定位方法，其特征在于，包括如下步骤:
[0009]步骤1，基因读段数据随机分割:将给定的基因读段数据集随机分割成指定的块数；
[0010]步骤2，数据的负载平衡:利用MapReduce (映射归约)框架，从分割好的数据块中随机抽取少量基因读段尝试执行，检测出可能执行时间较长的数据块，将其与其他数据块重新分割，均衡负载；
[0011]步骤3，读段的空间索引:采用MapReduce框架，在Map过程中将读段分割成指定段数的子读段，对每个子读段生成空位种子模式，将每个模式的哈希值与对应的读段信息记录在哈希表中，构建读段的空间索引；
[0012]步骤4，子读段不跨越剪切位定位:在Map过程中利用已有的读段定位软件将子读段中能够进行连续定位的读段首先定位，记录定位信息；
[0013]步骤5，子读段跨越剪切位定位:经过步骤4定位后，能够连续定位的子读段已经获得了定位信息，在Map过程中根据已经定位的子读段信息，推断不能连续定位的子读段的定位信息，实现跨越剪切位的读段定位；
[0014]步骤6，子读段拼接:经过步骤5所有读段分割成的子读段的定位信息已全部被记录，在Map过程中将所有定位的子读段进行组装，记录所有能够连接成原始读段的定位信息；
[0015]步骤7，读段定位信息统计:在Reduce过程中将所有通过Map过程传递过来的读段定位的信息进行汇总，统计读段定位在参考序列上的整体信息。
[0016]所述步骤2数据的负载平衡具体包括如下步骤=MapReduce框架将整个数据集默认分割为N块，记为(S1, S2,, SJ，从每块中随机抽取一定数量的基因数据，形成相应的子块，执行读段定位算法，检测出潜在的可能执行时间较长的分块；重新分割检测出的数据块，并且将其分别与其他未分割的数据块合并，最后再次分割成指定的块数。
[0017]所述步骤3读段的空间索引具体包括如下步骤:在Map过程中，将对应每个数据集子块中的读段分割为M个子读段，对每个子读段创建种子模式；每一个种子模式必然穷尽了所有可能的含有指定错配个数碱基的方式，将每一个种子模式转换成哈希值，以键值对的形式将哈希值和读段信息存储在哈希表里；在定位过程中就可以快速的从这张表中通过读段的哈希值获取指定读段的信息，这张哈希表就是所有读段的空间索引。
[0018]所述步骤5子读段跨越剪切位定位具体包括如下步骤:在每一个Map中，根据读段的长度，将子读段能够跨越的剪切位点的个数限定在指定的数量以内；根据已有的生物信息，限定内含子最大长度，提高定位效率；根据已经定位的子读段信息，通过延展尝试匹配的方式来确定不能够连续定位的子读段位置信息；如果原始读段的所有子读段都能够进行定位，则保留定位信息；否则舍弃该读段。
[0019]所述步骤6子读段拼接具体包括如下步骤:首先，每个子读段都记录它能够定位到参考序列上的所有可能位置，用一个数组存储所有位置信息；然后，从第一个子读段的第一个定位信息，即从第一个子读段的位置数组中取出第一个位置，与第二个子读段的所有位置信息进行检测，查看是否与其中一个位置具有顺接关系；如果有，则继续检测第二个子读段顺接的那个位置与第三个子读段所有位置信息的关系，直至最后一个子读段；如果没有，则停止检测；接着继续从第一个子读段的第二个位置进行检测，直至检查完第一个子读段的所有位置信息；最后，提取所有能够完整匹配的读段定位信息。
[0020]本发明是专门针对基因读段定位而提出的方法。与现有技术相比，本发明具有以下特征:
[0021](I)本发明针对现有基因数据量大，传统串行基因读段定位算法效率低下这一问题，提出了并行化的读段定位算法，并且利用MapReduce框架加以实现，在不影响读段定位正确性的条件下，提高了基因读段定位的效率；
[0022](2)本发明对基于空位种子的不跨越剪切位的读段定位算法进行改进，使之能够处理跨越剪切位的读段定位情况；此外，加入了生物信息:内含子的长度和剪切位点的数量，缩小了算法的搜索空间；改进读段定位时的算法，降低了算法的时间复杂度，从而进一步提高了基因读段定位的效率。
[0023](3)本发明在读段定位的过程中，专门考虑了负载平衡问题，提出了一种实际的解决方案，有效处理了个别节点执行时间长，降低整体效率的情况，因此具有较高的使用价值。
【专利附图】
【附图说明】
[0024]图1为本发明整体流程图；
[0025]图2为本发明中数据的负载平衡步骤子流程图；
[0026]图3为本发明中子读段不跨越剪切位定位步骤子流程图；
[0027]图4为本发明中子读段跨越剪切位定位步骤子流程图。
【具体实施方式】
[0028]以下结合【专利附图】
【附图说明】本发明的【具体实施方式】。
[0029]如图1所示，本发明公开了一种面向大规模基因数据的读段定位方法，具体步骤如下:
[0030]步骤1，基因读段数据随机分割:将给定的基因读段数据集随机分割成指定的块数；
[0031]步骤2，数据的负载平衡:采用MapReduce框架，从分割好的数据块中随机抽取少量基因读段尝试执行，检测出可能执行时间较长的数据块，将其重新分割，分别添加到其余数据块中，最后将其余的数据块再次进行分割，实现数据的均衡负载；
[0032]步骤3，读段的空间索引:采用MapReduce框架，在Map过程中将读段分割成指定段数的子读段，对每个子读段生成空位种子模式，将每个模式的哈希值与对应的读段信息记录在哈希表中，构建读段的空间索引；
[0033]步骤4，子读段不跨越剪切位定位:在Map过程中利用已有的读段定位软件将子读段中能够进行连续定位的读段首先定位，记录定位信息；
[0034]步骤5，子读段跨越剪切位定位:在Map过程中根据已经定位的子读段信息，推断不能连续定位的子读段的定位信息，实现跨越剪切位的读段定位；
[0035]步骤6，子读段拼接:在Map过程中将所有定位的子读段进行组装，记录所有能够连接成原始读段的定位信息；
[0036]步骤7，读段定位信息统计:在Reduce过程中将所有通过Map过程传递过来的读段定位的信息进行汇总，统计读段定位在参考序列上的整体信息，如定位数量、位置等；
[0037]假设读段定位过程中允许错配数量为2。
[0038]步骤I中，基因读段数据随机分割可以利用MapReduce框架实现。假设随机将读段数据集分割为N块，则可以指定Reduce的数量为N ;而在Map中可以设置一个能够产生1_N之间数的随机产生器，每次遍历一条读段时，产生一个随机数，该随机数对应着Reduce，表示将该条读段分配到相应的Reduce，最后只需要取出Reduce中的读段就完成数据分割了。
[0039] 步骤2中，利用MapReduce框架将整个数据集默认分割为N ±夹，记为(SpS2,...-,:!，从每块(64MB)中随机抽取一定数量的基因数据，形成相应的子块，执行读段定位算法，检测出潜在的可能执行时间较长的分块(假设检测出的是S1, S2)。
[0040]重新分割数据块，平衡负载:将S1, S2各自分割成N块，分别记为{Sn, S12,S1J，{S21, S22,S2J。合并其余N-2 ±夹，再重新分割成N ±夹，记为{S31, S32,...，S3J。将S11, S21, S31合并成Y i，按照同样的合并方法得到N个最终分块，记为{S' 2,...,S/ N}。
[0041]其中从每块中随机抽取一定数量的样本可以通过一个MapReduce任务实现，具体步骤:
[0042]Map 阶段:
[0043](I)根据路径读取读段数据集；
[0044](2) while (抽取的读段数量已满足要求)
[0045]①随机产生一个读段下标；
[0046]②将随机抽取到的读段送入reduce ；
[0047]③删除已抽过的读段；
[0048]④抽取的读段数量加I ;
[0049](3)结束 while 循环。
[0050]Reduce 阶段:
[0051](I)遍历定位数量集合；
[0052]①输出参考序列名称和最终的定位数量；
[0053](2)结束。
[0054]步骤2流程如图2所示，从各个子读段数据块中随机抽取一部分读段，做为测试数据集；在MapReduce框架下以得到的子数据集做为输入，执行读段定位算法，以检测出可以影响整体执行时间的个别节点；以Web浏览方式查看每一个Map节点的执行时间，根据时间的长度推断可能负载较大的节点；将检测出的负载较大节点中的读段数据集重新分割；将重新分割的子数据集添加到其他未分割的数据集中，达到平衡负载的目的。
[0055]步骤3中，在Map过程中，将对应每个数据集子块中的所有读段都分割为N个子读段。
[0056]针对每个子读段创建种子模式，主要思想如下:假设将子读段再分割为连续不交叉的k份，最多允许m个错配(m〈k)，对于每一个低于错配指定值的匹配，在已分的k个部分中，至少有k-m个部分没有错配。这k-m个部分属于k个部分，因此共有Cf"组合。而可以选择的种子模式就应该有01个。每一个种子模式必然穷尽了所有可能的含有指定错配个数碱基的方式。
[0057]将每一个种子模式转换成哈希值，以键值对的形式将哈希值和对应原始读段信息存储在哈希表里。在定位过程中就可以快速的从这张表中通过读段的哈希值获取指定读段的信息，这张哈希表就是所有读段的空间索引。
[0058]步骤4中，在Map过程中利用已有的读段定位软件将子读段中能够进行连续定位的读段首先定位，记录定位信息。本发明中采用的是SeqMap软件中的读段定位算法，将读段分割为4个子读段，错配2个，则每个子读段的种子模式就有6种。构建所有子读段的空间索引，遍历所有参考序列，将能够连续定位的子读段定位到参考序列上。
[0059]步骤4流程如图3所示，先将子读段再次进行分割，本发明中分割为4段，创建6中空位种子模式；然后将所有种子模式转换为哈希值，并且将哈希值与对应的读段信息存储到哈希表中，形成所有子读段的空间索引；依次遍历参考序列中的每一个碱基；判断步骤15中形成的哈希值在空间索引中是否存在；如果空间索引中存在该哈希值，则将对应的读段取出来与参考序列中的碱基进行进一步匹配；如果不存在，则中断此次尝试匹配；继续遍历下个碱基。
[0060]步骤5中，经过步骤4定位后，能够连续定位的子读段已经获得了定位信息；
[0061]在每一个Map中，根据读段的长度，将子读段能够跨越的剪切位点的个数限定在指定的数量以内(假设读段长84bp，可以限定最多允许跨越2个剪切位点)，因为跨越剪切位点数量越多，算法就需要更多的时间去检测更多可能的定位方式；根据已有的生物信息，限定内含子最大长度(比如已知的哺乳动物基因组中只有非常少的内含子长度大于400Kb)，这样可以减少算法的搜索空间，提高定位效率；
[0062]根据已经定位的子读段信息，通过延展尝试匹配的方式来确定不能够连续定位的子读段位置信息。跨越剪切位的读段定位可以分为两种情况:
[0063]第一种情况是未定位子读段的上下游子读段都已经定位，这种情况是比较容易处理的。本发明中分别从上游读段的最后一个碱基和下游读段的第一个碱基开始匹配，分别记录下错配0、1、2的前一个位置，在记录下的位置中查找能够连接的，且错配最少的位置，记为剪切位点。通过这种只需一次遍历查找到剪切位点的方法，可以之间将这一过程的时间复杂度降低至O (N)。
[0064]第二种情况是未定位子读段的上下游子段只有一个已经定位。对于这种情况，本发明中截取未定位子读段末尾一段长度，记为h-mer ;首先通过检测h-mer在上游子读段或下游子读段中的位置，然后再利用第一种情况的处理方式进行读段定位。而在本发明中h-mer的长度取值区别于一般算法中的取值(直接将h-mer长度固定为2)。本发明动态获取h-mer的长度，首先通过已定位的上游或下游子读段进行延展匹配直至出现2个错配的位置，记录此时延展的长度L，h-mer的长度就取为未定位子读段长度减去L。通过这一做法能够尽可能减少未定位子读段尝试匹配次数，并且使得改进算法将时间复杂度O从原来的0(N3)降低至0(N2)。
[0065]进行上述跨越剪切位定位后，如果原始读段的所有子读段都能够进行定位，则保留定位信息；否则舍弃该读段。
[0066]步骤5流程如图4所示，用已有的基于空位种子的读段定位软件将能够连续定位的子读段定位，记录位置信息；根据已经定位的子读段信息，推断未定位的子读段；未定位子读段的可能的第一种情况，它的上游和下游子读段均已定位成功，此时通过延展上下游子读段来确定未定位子读段的位置信息；未定位子读段可能的第二种情况，它的上游或下游子读段之一未定位，这种情况下可以截取未定位子读段末尾一段长度，检测它在上游子读段或下游子读段中的位置，这样就将第二种情况转化为第一种情况了，再通过第一种情况的处理方式就可以完成定位；拼接子读段步骤将所有子读段定位信息进行分析，拼接成完整读段定位信息。
[0067]步骤6中，经过步骤5所有读段分割成的子读段的定位信息已全部被记录；
[0068]首先，每个子读段都记录它能够定位到参考序列上的所有可能位置(可以多余I个)，可以用一个数组存储所有位置信息；[0069] 然后，从第一个子读段的第一个定位信息，即从第一个子读段的位置数组中取出第一个位置，与第二个子读段的所有位置信息进行检测，查看是否与其中一个位置具有顺接关系；如果有，则继续检测第二个子读段顺接的那个位置与第三个子读段所有位置信息的关系，直至最后一个子读段；如果没有，则停止检测；接着继续从第一个子读段的第二个位置进行检测，直至检查完第一个子读段的所有位置信息；
[0070]最后，提取所有能够完整匹配的读段定位信息。
[0071 ] 步骤7中，在Reduce过程中将所有通过Map过程传递过来的读段定位的信息进行汇总，统计读段定位在参考序列上的整体信息，如定位数量、位置等。
[0072]实施例:
[0073]本实施例所使用的云计算平台共有21台服务器(I个master节点，20个slave节点)，每个节点配置相同。节点的处理器为Intel (R)Xeon(R)CPU E5620，主频为2.40GHZ，操作系统为64位Debian Linux操作系统。Hadoop平台版本为hadoop-0.20.2。
[0074]采用的基因数据是从拟南芥菜基因组中抽取10条基因做为参考序列，平均长度约为12000bp。读段数据集主要有两个来源:1.用ES培养液培养拟南芥菜，重复实验3次，获得三次情况下的拟南芥菜RNA，经测序得到3个读段样本集。2.在缺P情况下培养拟南芥菜，同样重复实验3次，获得3个读段样本集。两种情况下获得的读段数据集共25.2GB,我们分别称这6个读段样本集为ES1、ES2、ES3、P1、P2和P3。
[0075]本实施例包括以下部分:
[0076]1.基因读段数据随机分割:
[0077]基因读段数据随机分割可以利用MapReduce框架实现。假设随机将读段数据集分害I]为N ±夹，则可以指定Reduce的数量为N ;而在Map中可以设置一个能够产生1-N之间数的随机产生器，每次遍历一条读段时，产生一个随机数，该随机数对应着Reduce，表示将该条读段分配到相应的Reduce，最后只需要取出Reduce中的读段就完成数据分割了。主要代码如下:
[0078]Map 阶段:
[0079](l)read_data = get_read_dataset(read_path)；
[0080](2) random = Math, rand () *N+1 ；
[0081](3)context, write(random, read_info)；
[0082]Reduce 阶段
[0083](I) for (count = 0 ；count<read_list.length ；count++)
[0084]① read = read_list.get (count)；
[0085]② context, write (task_num, read)；
[0086](2) end.[0087]2.数据的负载平衡:
[0088]利用MapReduce框架将整个数据集默认分割为N ±夹，记为(S1, S2, , Sj，从每块(64MB)中随机抽取一定数量的基因数据，形成相应的子块，执行读段定位算法，检测出潜在的可能执行时间较长的分块。
[0089]主要代码如下:
[0090](I) split_to_N-2 (S1, S2, N_2) ^/S1, S2 是潜在负载较高的块；[0091](2)add(list_S1, list_S2, N_2) ；//将重新分割的块添加到其余N_2块中；
[0092](3) re split (read_N-2_path, N) ；// 将 N_2 块重新划分为 N 块。
[0093]其中从每块中随机抽取一定数量的样本可以通过一个MapReduce任务实现，具体步骤:
[0094]Map 阶段:
[0095](I)read_data = get_readdata(read_path)；
[0096](2)while(read_num< = fix_readNum)
[0097]① random = Math, rand ()*read_totalNum+l ；
[0098]② read = read_data.get (random)；
[0099]③ context, write(fix_reduce, read)；；
[0100]④ read_data.move (random)；
[0101]⑤read_num++;
[0102](3) end whi le.[0103]Reduce 阶段:
[0104](I)for(count = 0 ；count<read_data.length ；count++)；
[0105]O lcontext.write (read_name, map_num)；
[0106](2) end.[0107]3.读段的空间索引:
[0108]本发明中读段的长度是84bp，至多允许读段跨越2个(根据生物信息学的经验在84bp长度的读段中很少会有跨越2个剪切位点以上的读段)剪切位点。因此，在Map过程中，将对应每个数据集子块中的所有读段都分割为4个(适中)子读段。
[0109]针对每个子读段创建种子模式，这种种子模式来自于SeqMap软件，具体可以参考文献:Jiang H, Wong W H.SeqMap:mapping massive amount of oligonucleotides to thegenome[J].Bioinformatics, 2008, 24(20):2395-2396.[0110]主要代码如下:
[0111]Map 阶段:
[0112](I) list_sub_read = split_read (read, f ix_read_num) ；// 将读段分割为指定段数；
[0113](2) list_seed_modes = create_modes (list_sub_read) -J/ 创建种子模式；
[0114](3) for (count = 0 ；count<list_seed_modes.length ；count++)
[0115]①hashcode = hashcode(list_seed_modes.get(count))；
[0116]②hash_map.put(hashcode, read_info)；
[0117]Il将哈希值和对应原始读段信息存储在哈希表里；
[0118](3) end.[0119]4.子读段不跨越剪切位定位:
[0120]在Map过程中利用已有的读段定位软件将子读段中能够进行连续定位的读段首先定位，记录定位信息。本发明中采用的是SeqMap软件中的读段定位算法，将读段分割为4个子读段，错配2个，则每个子读段的种子模式就有6种。构建所有子读段的空间索引，遍历所有参考序列，将能够连续定位的子读段定位到参考序列上。[0121]主要代码如下:
[0122](I) for (int count = O ;count〈ref—list, length ;count++)//遍历参考序列数据
集
[0123]① ref = ref—list, get (count)；
[0124]②for (int 1c = 0 ；loc<ref.length ;loc++)// 遍历参考序列每个减基
[0125]a.seed—mode = create—mode (ref.substring (loc，loc+84))；
[0126]b.1f (hash—map.get (seed—mode) ! = null) match (ref, read)；
[0127]③ end;
[0128](2) end.[0129]其中math()函数代码如下:
[0130](I) for (int count = 0 ；count<read.length ；count++)
[0131]① if( ! read, get (count).equals (ref.get (count))) break ；
[0132]else context, write (read—name，one)；
[0133](2) end.[0134]5.子读段跨越剪切位定位:
[0135]经过步骤4定位后，能够连续定位的子读段已经获得了定位信息；根据已经定位的子读段信息，通过延展尝试匹配的方式来确定不能够连续定位的子读段位置信息。
[0136]主要代码如下:
[0137](I) if (pre—subread.isFixed&&next_subread.1sFixed) // 第一种情况
[0138]① brige—match (read，ref)；
[0139](2) end ；
[0140](3) if ( ! pre—subread.isFixed | | ! next—subread.isFixed))// 第二种情况
[0141]①h—mer = read, substring (read.length-L) ；// 动态获取 h—mer
[0142]②loc = search (h_mer, ref) ；// 查找 h—mer 位置
[0143]③ brige—match (h—mer，ref)；
[0144](4) end.[0145]其中brige_match()函数代码如下:
[0146](I) for (int loc = 0 ；lco<read.length ；loc++)
[0147]① if ( ! read, get (loc).equals (ref.get (pre—subread.length+loc)))
[0148]a.mis—num++ ；
[0149]b.list—loc = add (loc)；
[0150]②if (mis—num>2) break ；
[0151](2) end.[0152](3) for (int loc = read, length ；loc> = 0 ；loc—)
[0153]① index = I ；
[0154]② if ( ! read, get (loc).equals (ref.get (next—subread.loc-1ndex)))
[0155]a.mis—num++ ；
[0156]b.list—loc = add (loc)；
[0157]③if (mis—num>2) break ；[0158](4) end,
[0159](5) spliced_site = find_min_mistake (list_loc).// 查找错配最少的连接点
[0160]6.子读段拼接，提取所有能够完整匹配的读段定位信息。
[0161]7.统计读段定位在参考序列上的整体信息。
[0162]以上对本发明所提供的一种面向大规模基因数据的读段定位方法进行了详细介绍。值得注意的是，具体实现该技术方案的方法和途径有很多，以上所述仅是本发明的优选实施方式，只用于帮助理解本发明的方法及核心思想；同时，对于本领域的一般技术人员，在本发明核心思想的基础上，做出的修改和调整都将视为本发明的保护范围。综上所述，本说明书内容不应理解为对本发明的限制，本发明的保护范围应由所附的权利要求来限定。
【权利要求】
1.一种面向大规模基因数据的读段定位方法，其特征在于，包括如下步骤: 步骤1，基因读段数据随机分割:将给定的基因读段数据集随机分割成指定的块数； 步骤2，数据的负载平衡:利用MapReduce框架，从分割好的数据块中随机抽取少量基因读段尝试执行，检测出可能执行时间较长的数据块，将其与其他数据块重新分割，均衡负载； 步骤3，读段的空间索引:采用MapReduce框架，在Map过程中将读段分割成指定段数的子读段，对每个子读段生成空位种子模式，将每个模式的哈希值与对应的读段信息记录在哈希表中，构建读段的空间索引； 步骤4，子读段不跨越剪切位定位:在Map过程中利用已有的读段定位软件将子读段中能够进行连续定位的读段首先定位，记录定位信息； 步骤5，子读段跨越剪切位定位:经过步骤4定位后，能够连续定位的子读段已经获得了定位信息，在Map过程中根据已经定位的子读段信息，推断不能连续定位的子读段的定位信息，实现跨越剪切位的读段定位； 步骤6，子读段拼接:经过步骤5所有读段分割成的子读段的定位信息已全部被记录，在Map过程中将所有定位的子读段进行组装，记录所有能够连接成原始读段的定位信息； 步骤7，读段定位信息统计:在Reduce过程中将所有通过Map过程传递过来的读段定位的信息进行汇总，统计读段定位在参考序列上的整体信息。
2.根据权利要求1所述的一种面向大规模基因数据的读段定位方法，其特征在于，所述步骤2数据的负载平衡具体包括如下步骤: MapReduce框架将整个数据集默认分割为N ±夹，记为(S1, S2, , Sj，从每块中随机抽取一定数量的基因数据，形成相应的子块，执行读段定位算法，检测出潜在的可能执行时间较长的分块；重新分割检测出的数据块，并且将其分别与其他未分割的数据块合并，最后再次分割成指定的块数。
3.根据权利要求1所述的一种面向大规模基因数据的读段定位方法，其特征在于，所述步骤3读段的空间索引具体包括如下步骤: 在Map过程中，将对应每个数据集子块中的读段分割为M个子读段，对每个子读段创建种子模式；每一个种子模式必然穷尽了所有可能的含有指定错配个数碱基的方式，将每一个种子模式转换成哈希值，以键值对的形式将哈希值和读段信息存储在哈希表里；在定位过程中就可以快速的从这张表中通过读段的哈希值获取指定读段的信息，这张哈希表就是所有读段的空间索引。
4.根据权利要求1所述的一种面向大规模基因数据的读段定位方法，其特征在于，所述步骤5子读段跨越剪切位定位具体包括如下步骤: 在每一个Map中，根据读段的长度，将子读段能够跨越的剪切位点的个数限定在指定的数量以内；根据已有的生物信息，限定内含子最大长度，提高定位效率； 根据已经定位的子读段信息，通过延展尝试匹配的方式来确定不能够连续定位的子读段位置信息；如果原始读段的所有子读段都能够进行定位，则保留定位信息；否则舍弃该读段。
5.根据权利要求1所述的一种面向大规模基因数据的读段定位方法，其特征在于，所述步骤6子读段拼接具体包括如下步骤:首先，每个子读段都记录它能够定位到参考序列上的所有可能位置，用一个数组存储所有位置信息； 然后，从第一个子读段的第一个定位信息，即从第一个子读段的位置数组中取出第一个位置，与第二个子读段的所有位置信息进行检测，查看是否与其中一个位置具有顺接关系；如果有，则继续检测第二个子读段顺接的那个位置与第三个子读段所有位置信息的关系，直至最后一个子读段；如果没有，则停止检测；接着继续从第一个子读段的第二个位置进行检测，直至检查完第一个子读段的所有位置信息； 最后，提 取所有能够完整匹配的读段定位信息。
【文档编号】G06F19/24GK103971031SQ201410185387
【公开日】2014年8月6日 申请日期:2014年5月4日 优先权日:2014年5月4日
【发明者】杨明, 涂金金, 高阳 申请人:南京师范大学
相关技术
网友询问留言
已有0条留言
- 还没有人留言评论。精彩留言会获得点赞！
1