CN103093121A

CN103093121A - 双向多步deBruijn图的压缩存储和构造方法

Info

Publication number: CN103093121A
Application number: CN2012105870596A
Authority: CN
Inventors: 孟金涛; 魏彦杰; 成杰峰; 冯圣中
Original assignee: Shenzhen Institute of Advanced Technology of CAS
Current assignee: Shenzhen Institute of Advanced Technology of CAS
Priority date: 2012-12-28
Filing date: 2012-12-28
Publication date: 2013-05-08
Anticipated expiration: 2032-12-28
Also published as: CN103093121B

Abstract

本发明涉及一种双向多步deBruijn图的压缩存储和构造方法，包括压缩存储步骤，和deBruijn图构造步骤。本发明提供的双向多步deBruijn图的压缩存储和构造方法，（1）结合DNA序列互补双螺旋结构的特点，对de Bruijn图结构进行结构优化，使用双向多步de Bruijn图，可以将需要存储的图的节点减半；（2）针对序列组装的原始de Bruijn图规模异常庞大，对内存造成存储压力等问题，使用双向多步de Bruijn图的压缩存储技术，使得存储双向多步de Bruijn图的内存消耗控制在参考序列大小的100倍以内；（3）构造双向多步de Bruijn图，可将DNA序列组装问题分解为边融合子问题，更适合并行计算。

Description

双向多步deBruijn图的压缩存储和构造方法

【技术领域】

本发明涉及基因测序领域，尤其涉及一种双向多步deBruijn图的压缩存储和构造方法。

【背景技术】

新一代基因测序技术所产生的序列片段具有序列短、高覆盖率、额外的双端信息等特点，使得原有传统的序列拼接技术无法使用，加大了序列拼接技术的复杂度和难点。面对海量的数据，高效的序列拼接技术成为处理测序数据的关键。

序列拼接技术问题主要转化为图论的知识进行求解。当前，序列拼接技术问题解决途径主要有两条。一种途径是传统的先重叠后扩展方法，即OLC（Overlap-Layout-Consensus）方法。该方法主要是将序列拼接技术转化为哈密尔顿路径进行求解。OLC的拼接过程主要分为三个步骤：1）Overlap：对获得的所有read顶点进行排序，通过序列比对算法，寻找DNA序列片段间的重叠信息；2）Layout：根据read之间的重叠信息，排列所有DNA序列片段，形成新的链接体结构，该步骤主要是为了寻找经过每个顶点一次且仅一次的路径，即寻找哈密尔顿路径；3）Consensus：根据新的链接体结构中原始质量数据，在链接体中寻找质量最重的路径，从而组合成最终的DNA序列。

基于OLC的序列拼接在技术上容易实现，但是该方法仅局限于比对大于某个阈值的read之间的信息，忽略了多个read之间的相互信息，从而使得该方法在处理重复区域问题上受到了很大的限制；另外，该方法进行序列比对时，无论使用贪心算法或者BWA算法，都需要耗费大量的内存，同时在任意两个read序列之间比对，使得该序列拼接的算法复杂度为O(n²)，其中n为序列片段的个数。因此，从时间和空间上考虑，该算法均难以用于拼接长达百万碱基的基因组样本或者由百万条以上的DNA序列片段组成的测序数据。基于OLC的拼接主要适用于基于Sanger测序原理的拼接。测序样本也主要是基因组较小的生物,例如病毒、真菌。目前，PHRAP、TIGR、CAP3、CELERA、ARACHNE、PHUSION、SSAKE、VCAKE、SHARCGS等序列拼接技术均是在此基础上开发的。

DNA序列拼接的另外一种途径，是基于de Bruijn图的序列拼接方法。该方法主要是将DNA序列拼接问题转化为欧拉路径进行求解；通过构造并简化deBruijn图来实现整个拼接过程。基于de Bruijn图的拼接方法理论上需要找到一条欧拉路径来重构源基因组参考序列，但实际上只能通过将构造好的de Bruijn图中的无分支路径进行收缩最终合并为完整的contig进行输出。目前，Velvet、Soapdenovo、Idba、Abyss、Pasha等序列拼接技术均是在此基础上开发的。

与OLC相比，基于de Bruijn图的序列拼接技术有很多的优点。例如，在过滤错误信息、重复区域发现、解耦以及利用双端信息对contig进行扩展等方面。然而，由于大基因组构造的de Bruijn图异常庞大，以至于单个计算机的内存无法释放。例如，对人类基因组的数据进行拼接，构建的de Bruijn图中有30G的顶点，而存储这样的图通常需要消耗大约500G到1T的内存。同时，由于在图的数据结构中顶点关联的随机性，不能使用数据预取技术，导致最终对图的简化也将耗费几周的时间。在拼接大基因组和宏基因组测序数据时，上述拼接技术的性能依然无法满足需求。

【发明内容】

本发明要解决的技术问题在于提高现有技术序列拼接运行速度、降低单机内存消耗。

为此，本发明提供一种双向多步deBruijn图的压缩存储和构造方法，包括

压缩存储步骤，具体为

S11、读取一个序列s；

S12、将序列s用滑动窗口切割为多个片段t；

S13、对每个片段t，使用核酸编码表进行编码，并表示为一个64位的整数a；

S14、将片段t进行反转，使用对称互补表将反转的片段互补处理，得到互补片段，并再次使用步骤S13中的核酸编码表将互补片段进行编码，并表示为一个64位的整数b；

S15、取整数a和整数b的最大数，作为片段t和互补片段v的k分子的标志数；

S16、重复步骤S11-S15，直至所有序列完成；

和deBruijn图构造步骤，具体为

S21、读取一个序列s；

S22、将序列s用滑动窗口切割为多个片段t，选取一片段t其标志数为cur、并标记其前、后的片段的标志数分别为pre、lat；

S23、若t的编码小于其互补片段编码，则交换pre,lat的值;

S24、在cur的正向位置映射表的相应bit位置1来表示指向pre的边;

S25、在cur的反向位置映射表的相应bit位置1来表示指向lat的边;

S26、重复步骤S22-S25，处理序列s的其他片段t，直至完成序列s的全部片段t，执行步骤S27；

S27、读取一个新的序列s，重复步骤S22-S26；直至处理完所有的序列，执行步骤S28;

S28、完成双向多步de Bruijn图的构造。

进一步地，所述步骤S12、S22中的滑动窗口为长度为k的滑动窗口，其中0<k<32且k为奇数。

进一步地，所述步骤S13中的核酸编码表为{A:00,C:01,G:10,T:11}。

进一步地，所述步骤S14中的对称互补表为{A->T,C->G， G->C,T->A}。

进一步地，所述步骤S14具体为，将片段t的字符串进行反转，使用对称互补表将反转的字符串中每个字符变为其互补字符，得到互补字符的字符串v，并再次使用步骤S13中的核酸编码表将字符串v进行编码，并表示为一个64位的整数b；

进一步地，所述步骤S22中，若片段t没有之前或之后的片段，若片段t没有之前或之后的片段，则对pre或者lat值赋为空或NULL。

进一步地，步骤S24中正向位置映射表为{A:0，C:1，G:2，T:3}，位置查询字符为pre的最后一个字符。

进一步地，步骤S25中反向位置映射表为{A:4，C:5，G:6，T:7}，位置查询字符为lat的第一个字符的互补字符。

本发明提供的双向多步deBruijn图的压缩存储和构造方法，采用压缩存储步骤和deBruijn图构造步骤，（1）结合DNA序列互补双螺旋结构的特点，对deBruijn图结构进行结构优化，使用双向多步de Bruijn图，可以将需要存储的图的节点减半；（2）针对序列组装的原始de Bruijn图规模异常庞大，对内存造成存储压力等问题，使用双向多步de Bruijn图的压缩存储技术，使得存储双向多步de Bruijn图的内存消耗控制在参考序列大小的100倍以内；（3）构造双向多步de Bruijn图，可将DNA序列组装问题分解为边融合子问题，更适合并行计算。本发明的方法提高了现有技术序列拼接运行速度、降低了单机内存消耗。

【附图说明】

图1示出本发明的压缩存储步骤流程图。

图2示出本发明的deBruijn图构造步骤流程图。

【具体实施方式】

下面结合附图和具体实施例对本发明作进一步详细说明。

本发明提供一种双向多步deBruijn图的压缩存储和构造方法，包括

压缩存储步骤，所需原始数据包括第一代，第二代和新一代的测序仪器产生出来的FASTA格式文件，将FASTA文件中的序列逐个切割成k分子并且用二进制编码进行压缩存储为一个64位的长整型k分子的标志数。

如图1所示，具体为

S11、读取一个序列s；其中，序列s取自FASTA格式文件；

S12、将序列s用滑动窗口切割为多个片段t；

S16、重复步骤S11-S15，直至所有序列完成。

通过上述步骤将两个传统的De Brujin图中的kmer，转化为一个64位的k分子的标志数来存储。该步骤可以将其他软件例如velvet、IDBA、SOAPdenovo里的两个压缩kmer存储为一个压缩k分子的标志数，并且在得到k分子的标志数后也可以反过来求出该k分子的长度为k的片段t和它的互补片段v。

和deBruijn图构造步骤，1、使用上述压缩存储步骤中计算k分子的标志数，2、将每个片段以及和它前后相邻的片段的扩展字符作为该k分子和其前后相邻的片段的对应的k分子的边并初始化k分子数据结构的边；3、将初始化后的k分子数据结构以k分子的标志数为关键值存入hash_map。

如图2所示，具体为

S21、读取一个序列s；

S23、若t的编码小于其互补片段编码，则交换pre,lat的值;

S28、完成双向多步de Bruijn图的构造。

实施例1

本发明实施例提供的一种双向多步deBruijn图的压缩存储和构造方法，包括

压缩存储步骤，具体为

S11、读取一个序列s；

S12、将序列s用长度为k的滑动窗口切割为多个片段t；其中0<k<32且k为奇数。

S13、对每个片段t，使用核酸编码表{A:00,C:01,G:10,T:11}进行编码，并表示为一个64位的整数a；

S14、将片段t的字符串进行反转，使用对称互补表{A->T，C->G， G->C,T->A}将反转的字符串中每个字符变为其互补字符，得到互补字符的字符串v，并再次使用步骤S13中的核酸编码表将字符串v进行编码，并表示为一个64位的整数b；

S16、重复步骤S11-S15，直至所有序列完成。

和deBruijn图构造步骤，具体为

S21、读取一个序列s；

S23、若t的编码小于其互补片段编码，则交换pre,lat的值;

S28、完成双向多步de Bruijn图的构造。

实施例2

本发明另一实施例提供的一种双向多步deBruijn图的压缩存储和构造方法，包括

压缩存储步骤，具体为

S11、读取一个序列s；

S16、重复步骤S11-S15，直至所有序列完成。

和deBruijn图构造步骤，1、使用上述压缩存储步骤中计算k分子的标志数，2、将每个片段以及和它前后相邻的片段的扩展字符作为该k分子和其前后相邻的片段的对应的k分子的边并初始化k分子数据结构的边；3.将初始化后的k分子数据结构以k分子的标志数为关键值存入hash_map。具体为

S21、读取一个序列s；

S23、若t的编码小于其互补片段编码，则交换pre,lat的值;

S28、完成双向多步de Bruijn图的构造。

其中所述的k分子数据结构如下：

本发明分别选取四种规模大小不等的DNA模拟数据进行测试。其中，DNA模拟数据分别为Swinepox virus（Swinepox）、Escherichia coli str.K-12substr(E.coli)、Saccharomyces cerevisiae(Yeast)、Caenorhabditis elegans(C.elegans)四种基因组数据，其实验数据来源于基因库文件（NC_003389，NC_000913，NC_001133-NC001148，NC_003279-NC_003284），得到DNA数据信息如表1所示。

表1.DNA模拟数据表

表1可以看出，经过DNA测序仪处理之后，测序数据规模巨大。如果直接构造传统的de Bruijn图造成了系统资源的极大浪费；同时，带来了单机存储的内存压力。由表2中实验数据对比可知，本发明所提出双向多步de Bruijn图的压缩存储和构造方法能节省内存近50%，同时内存消耗控制在源数据规模的100倍以内。

表2.构图所消耗的内存（G）

以上所述本发明的具体实施方式，并不构成对本发明保护范围的限定。任何根据本发明的技术构思所作出的各种其他相应的改变与变形，均应包含在本发明权利要求的保护范围内。

Claims

1.一种双向多步deBruijn图的压缩存储和构造方法，其特征在于，包括压缩存储步骤，具体为

S11、读取一个序列s；

S12、将序列s用滑动窗口切割为多个片段t；

S14、将片段t进行反转，使用对称互补表将反转的片段互补处理，得到互补片段v，并再次使用步骤S13中的核酸编码表将互补片段进行编码，并表示为一个64位的整数b；

S16、重复步骤S11-S15，直至所有序列完成；

和deBruijn图构造步骤，具体为

S21、读取一个序列s；

S23、若t的编码小于其互补片段编码，则交换pre,lat的值;

S28、完成双向多步de Bruijn图的构造。

2.根据权利要求1所述的双向多步deBruijn图的压缩存储和构造方法，其特征在于，所述步骤S12、S22中的滑动窗口为长度为k的滑动窗口，其中0<k<32且k为奇数。

3.根据权利要求1所述的双向多步deBruijn图的压缩存储和构造方法，其特征在于，所述步骤S13中的核酸编码表为{A:00,C:01,G:10,T:11}。

4.根据权利要求1所述的双向多步deBruijn图的压缩存储和构造方法，其特征在于，所述步骤S14中的对称互补表为{A->T,C->G， G->C,T->A}。

5.根据权利要求1所述的双向多步deBruijn图的压缩存储和构造方法，其特征在于，所述步骤S14具体为，将片段t的字符串进行反转，使用对称互补表将反转的字符串中每个字符变为其互补字符，得到互补字符的字符串v，并再次使用步骤S13中的核酸编码表将字符串v进行编码，并表示为一个64位的整数b；

6.根据权利要求1所述的双向多步deBruijn图的压缩存储和构造方法，其特征在于，所述步骤S22中，若片段t没有之前或之后的片段，则对pre或者lat值赋为空或NULL。

7.根据权利要求1所述的双向多步deBruijn图的压缩存储和构造方法，其特征在于，步骤S24中正向位置映射表为{A:0，C:1，G:2，T:3}，位置查询字符为pre的最后一个字符。

8.根据权利要求1所述的双向多步deBruijn图的压缩存储和构造方法，其特征在于，步骤S25中反向位置映射表为{A:4，C:5，G:6，T:7}，位置查询字符为lat的第一个字符的互补字符。