CN104361264A - 基因组核酸片段数量的快速计数算法 - Google Patents

Abstract

本发明公开了一种快速统计基因组核酸片段数量的算法。该算法针对基因组核酸片段高倍数的特点，将不同的核酸片段合理地分配到内存和硬盘。该算法总是优先将倍数最高的一些核酸片段放置于内存中，而将倍数较低的核酸片段存储在硬盘上。通过本发明提出的算法，任何基因组核酸片段数量的计数都可在有限内存的限制下完成。

Description

基因组核酸片段数量的快速计数算法

技术领域

本研究主要涉及到生物信息数据处理，生物大数据挖掘，基因组分析与全基因组拼接。

背景技术

k-mer，即为kilo mer，指的是将一条基因序列，连续切割，挨个碱基划动得到的一序列长度为K的核苷酸片段。在生物信息学项目的研究中，对k-核酸片段的计量是解读基因组信息的重要环节，我们可以通过k-核酸片段了解到基因组的大小，深度，频率，同时，k-核酸片段计量广泛应用于基因组拼装拼接，评估基因组测序覆盖度，测序数据的纠错，多序列比对，重复序列检测等研究中。

近年来，随着生物信息学的发展，国外研究出了一些此类算法，如Jellyfish，BFCounter等。Jellyfish，一种已有的k-核酸片段计量算法，采用了并行无锁哈希表(lock-free hash table)数据结构，在速率上比较高效，但需要非常大的内存和硬盘空间(在处理160G的人类基因组数据需要70G内存和211G硬盘)，这严重阻碍了在没有小型机等大型计算设备的普通实验室环境下的数据处理。BFCounter，采用概率型数据结构Bloom_filter，使其不需要硬盘操作，完全使用内存即可实现计量。但其消耗内存过大(处理上述数据需要56G内存)，而且速率过慢(处理上述数据时间超过40小时)所以，在有限的内存和硬盘空间限制下，尚未有可合理权衡上述要求的高效算法。同时，在数据处理中，现有方法并没有深刻发掘基因组本身信息，如频率，深度，而这些信息可以帮助我们去合理处理重复基因和错误基因，从而极大的提升算法效率。

发明内容

本发明是为了解决现有技术中的问题，而提供了一种快速、低内存需求、低硬盘需求的k-核酸片段s计量算法。

本发明为解决公知技术中存在的技术问题所采用的技术方案包括：合理分配内存和硬盘空间，采用合理的数据结构来表示k-核酸片段s，优化k-核酸片段s的索引结构和排序。算法具体步骤如下：

1)读取原始数据文件，设定参数，包括要统计k-核酸片段的长度k，限定硬盘空间D(bits)，内存空间M(bits)，数据深度(depth)。

2)遍历数据文件一次，将文件利用位运算压缩成二进制，并且建立长度为8-核酸片段的索引index。同时，获得k-核酸片段的总个数

3)利用参数划分迭代次数和k-核酸片段划分段。

4)遍历8-核酸片段索引一次，等量分配每一个区段内的k-核酸片段个数。

5)按照不同的数据结构初始化两种种存储存储空间，分配在内存，硬盘中，用于存储和计算k-核酸片段。

6)根据迭代顺序开始计量工作。

7)按照索引区index[i](i∈(划分段)读取压缩文件至内存R，直到内存R使用量达到设定内存空间的40％为止。

8)快速排序已存于R的k-核酸片段，并合并相同的k-核酸片段，进行计量，得出的结果存储于Lx(x∈1，2，3..)中。

9)迭代7，8操作，直到已使用内存超过60％。合并L1...n中的k-核酸片段，并且建立次数索引

10)遍历已统计的k-核酸片段，根据基因组的深度划分正确基因和错误基因，将重复次数少的写入硬盘，次数多的保存在内存Lx中。每次写入硬盘的文件必须按照划分段来标记，同时，通过遍历索引找到相应的k-核酸片段写入的最低次数，保证写入k-核酸片段至硬盘必须释放60％的内存。

11)迭代7-10操作，直到所有文件读完。

12)此时所有的k-核酸片段都已经存储于内存(计量次数多的k-核酸片段)和硬盘(次数少的k-核酸片段)中。根据划分段，将内存和硬盘中处于同一区段的k-核酸片段进行合并排序，合并排序的顺序遵循哈弗曼编码方法

13)完成所有的迭代次数，完成计量工作。

附图说明

图1：程序运行图1；

图2：内存CPU使用图1；

图3：程序运行图2；

图4：内存CPU使用图2；

图5：程序运行图3；

图6：内存CPU使用图3；

具体实施方式

下面结合实施例进一步描述本发明。本发明的范围不受这些实施例的限制，本发明的范围在权利要求书中提出。

本算法中具体的实施步骤如下：

1.在NCBI上下载的基因组数据，利用SRA Toolkit工具包将其转成fastq和fasta格式的数据文件。

2.如图1所示：本算法程序的运行格式为./hhz目标数据库名称k值[-m限定内存容量][-d限定硬盘容量]图中输入格式为./hhz ERR022075_2.fasta 26-m 400-d1024。

3.三行max中的数值分别为30％，60％，100％的限定内存容量，以字节为单位。在压缩数据和建立索引过程中，不需要内存参与，其CPU和内存占用率如图2所示。(289mb为***本身程序占用)

4.快速排序，写入硬盘阶段，如图3。

5.在硬盘中，根据之前计算出的划分值，分为了相应的文件夹数量，并且在快速排序写入硬盘的过程中，根据索引写入相应的文件夹中。在每一个文件夹中，每写入一次硬盘，即成为一个独立的文件。这样，每一个划分段文件夹中，都存储着各个区段的k-核酸片段，而每一个文件中的k-核酸片段，都是按照顺序排列好的，这样就为之后的合并排序提供了最好的分配策略。本阶段占用内存主要用于创建数组A和数组B，进行k-核酸片段的读取和快速排序。

6.合并排序，写入结果文件阶段。图5中的标号即为哈弗曼编码的标号，程序始终会找长度最短的两个文件进行合并排序，每完成一次合并排序后再将其结果加入哈弗曼树中，再进行下一次合并排序。

7.根据内存的走势图我们可以发现，在合并排序的时候，内存的使用率呈现逐步增加的趋势。根据满足哈弗曼编码的合并排序过程，我们可以得知，整个合并排序的过程是将硬盘中已排完序的k-核酸片段进行合并，逐步添加至内存的一个过程(如图6所示)。该阶段内存利用率最高达到80％。

由于其他k-核酸片段算法，如Jellyfish和BFCounter，都是基于大内存空间和硬盘容量作为条件，同时Jellyfish是多线程算法，所以不可能在相同的条件下进行时间的对比，于是我们将我们的算法与dsk做主要的对比。

我们选取了深度较高的几组数据，分别为大肠杆菌，酵母，线虫的基因组测试了本发明与dsk算法。运行结果显示，本发明在高深度的情况下明显优于dsk的运行速率，并且深度越高，优势越明显。大肠杆菌ERR022075序列我们的算法(341秒)相对于DSK(689秒)节省了大约一半的时间。

表1：本发明与dsk算法的比较

Claims (1)

1.一种快速统计基因组核酸片段数量的算法，该方法可在有限内存的限制条件下完成核酸片段数量的统计。本发明总是优先将高倍数的核酸片段放置于内存，而将低倍数的核酸不断移除到硬盘中。通过不断排序与迭代上述过程，尽可能地提高在内存中计算的比例。其期特征在于，

(1)读取原始数据文件，设定参数，包括要统计k-核酸片段(表示长度为k的核酸片段)的长度k，限定硬盘空间D(bits)，内存空间M(bits)，数据深度(depth)。

(2)遍历数据文件一次，将文件利用位运算压缩成二进制，并且建立长度为8-核酸片段的索引index。同时，获得k-核酸片段的总个数N。

(3)划分迭代次数和k-核酸片段区段。

(4)遍历8-核酸片段索引一次，等量分配每一个区段内的k-核酸片段个数。

(5)按照不同的数据结构初始化两种存储存储空间，分配在内存，硬盘中，用于存储和计算k-核酸片段。

(6)根据迭代顺序开始计量工作。

(7)按照索引区index[i](i∈(划分段)读取压缩文件至内存R，直到内存R使用量达到设定内存空间的40％为止。

(8)快速排序已存于R的k-核酸片段，并合并相同的k-核酸片段，进行计量，得出的结果存储于Lx(x∈1，2，3..)中。

(9)迭代7，8操作，直到已使用内存超过60％。合并L1...n中的k-核酸片段，并且建立次数索引

(10)遍历已统计的k-核酸片段，根据基因组的深度划分正确基因和错误基因，将重复次数少的写入硬盘，次数多的保存在内存Lx中。每次写入硬盘的文件必须按照划分段来标记，同时，通过遍历索引找到相应的k-核酸片段写入的最低次数，保证写入k-核酸片段至硬盘必须释放60％的内存。

(11)迭代7-10操作，直到所有文件读完。

(12)此时所有的k-核酸片段都已经存储于内存(计量次数多的k-核酸片段)和硬盘(次数少的k-核酸片段)中。根据划分段，将内存和硬盘中处于同一区段的k-核酸片段进行合并排序，合并排序的顺序遵循哈弗曼编码方法。

(13)完成所有的迭代次数，完成计量工作。