CN110310709A - 一种基于参考序列的基因压缩方法 - Google Patents

Abstract

本发明公开了一种基于参考序列的基因压缩方法，首先任意选取一个基因序列作为参考序列。其次，获取参考序列的小写字符和ACGT，并以二元组表示小写字符。然后，读取参考文件，获得参考文件的头部、换行信息、小写字符、N字符、碱基信息和其他字符，并将换行长度、小写字符、N字符和其他字符表示成二元组。接着，匹配参考序列和待压缩序列的小写字符二元组。最后匹配Hash值。解压缩过种采用压缩过程相反的步骤。采用本压缩方法的压缩比高，压缩速度快，而且二元组编码与基因次序无关，有利于分布式存储和分析基因序列。

Description

一种基于参考序列的基因压缩方法

技术领域

本发明涉及一种基因序列压缩的方案，主要用于解决基因序列数据过大，存储和传输成本高等问题，属于压缩算法领域。

背景技术

基因是DNA上有遗传效应的片断，人类的生老病死等都与基因有关。基因数据研究可以获得对生命运行机制和疾病机理等的深入研究，在生物医药学和相关生物技术产业，如制药、农、林、牧、渔、环保等的发展发挥越来越重要的作用，对于推动精准医疗，助力解决三大民生问题之一的医疗问题，具有重要的作用。因此，基因数据因其重要的社会价值和科研价值受到国际社会的广泛重视。自1990年正式启动的国际人类基因组计划以来，随着基因测序技术的不断进步，基因测序成本的不断降低，测序速度不断提高，众多国家和组织纷纷启动基因工程计划。2017年12月28日，我国启动“中国十万人基因组计划”，这是我国在人类基因组研究领域实施的首个重大国家计划，也是目前世界最大规模的人类基因组计划。随着各种测序项目的展开，产生的序列数据量呈指数规模增长，而且未来增长速度会更快。基因数据增长的速度大大超过了存储和传输带宽增长的速度，给存储和传输带来了很大的压力。如何以更高的效率存储基因数据，减轻存储和传输压力，在基因研究和应用中着十分重要的作用。

DNA序列数据具有与其他数据截然不同的特性，DNA序列是仅由A、G、C、T四个符号构成的超长序列，构成种类简单但是序列长度巨大。很大一部分DNA序列至今无法确定其用途，如果数据压缩过程中出现丢失，可能造成不可估量的损失，所以DNA序列必须保证无损压缩。另外，DNA序列中碱基对的排列并不是随机的，具有特定的概率分布和规律性。而且，DNA序列具有高度的相似性。首先，不同物种间的DNA序列相似度很高，同一物种间的DNA序列相似性更为明显。其次，同一个体内的不同片断的DNA序列也存在着许多精确重复。利用DNA这些信息特点，工业界和学术界提出了众多利用DNA序列特征的DNA序列压缩方法。经过对现在技术的文献检索发现，2000年T Matsumoto和K Sadakane在Genome Informatics上的“Biological sequence compression algorithms”提出了CTW+LZ方法，将上下文树加权(Context Tree Weighting,CTW)方法和LZ压缩方法相结合，使用多个编码模型对DNA序列的不同片断进行压缩。2002年，X Chen和M Li在Bioinformatics上的“DNACompress:fastand effective DNA sequence compression”提出了DNACompress压缩方法，使用了Pattern Hunter工具搜索DNA序列的重复与近似重复片断，提高了方法的整体速度。2005年，G Korodi和I Tabus在ACM Transactions on Information Systems上的“AnEffective Normalized Maximum Likelihood Algorithm for DNA SequenceCompression”提出了GeNML方法，对具有不同数据特点的DNA片断使用不同的编码策略和概率模型进行压缩。2013年，Sebastian Wandelt and Uif Leser在IEEE/ACM Transactionson Computational Biology and Bioinformatics的“FRESCO:Referential Compressionof Highly Similar Sequences”提出了一种叫FRESCO的快速基因压缩方法，它采用了一种用参考基因来表示被压缩基因的方法。2015年，Xiaojing Xie,Shuigeng Zhou和JihongGuan在IEEE/ACM Transactions on Computational Biology and Bioinformatics的“CoGI：Towards Compressing Genomes as an Image”上提出了一种用图模型来表示基因数据，从而可以利用图压缩技术来压缩基因模型的方法。总结这些DNA序列压缩方法可以分为两大类：基于非参考序列的DNA序列压缩方法和基于参考序列的DNA序列压缩方法，这些方法都有效的提高了压缩比和压缩效率。但总体而言，DNA序列的压缩技术仍处于起步阶段，组成基因片断的生物信息特征及片断内部的细节重复特性并没有被充分发挥利用。基因序列之间的特征也还没有被充分挖掘。

发明内容

发明目的：为了克服现有技术中存在的不足，本发明提供一种面向多条序列的以二元组形式表示基因序列的基于参考序列的基因压缩方法，该方法能显著提高压缩效率，实现高效存储。首先选取参考序列，然后根据参考序列将待压缩序列的碱基信息、小写字母、N字符和其他字符的数据以二元组数据结构形式表示。通过对比参考基因文件和待压缩基因文件，根据基因具有高度相似性特点，在压缩文件中以二元组形式写入两个基因文件中相同的部分，并将不同的部分单独写入压缩文件，比较的元素包括碱基信息、小写字符、N字符、换行信息、以及其他字符。既充分利用了基因文件的高度相似性特点，又规避了单独文件压缩的局限性，大大的提高了基因压缩的效率，具有较高的研究意义。

技术方案：为实现上述目的，本发明采用的技术方案为：

一种基于参考序列的基因压缩方法，首先遍历参考基因文件，提取小写字符和碱基信息；然后遍历待压缩基因文件，提取碱基信息、小写字符、换行信息、N字符和其他字符；最后通过比较，将生成的二元组写入文件，实现基因压缩，具体包括以下步骤：

步骤1)、遍历参考基因文件；将其中小写字符写入内存，此小写字符信息包括碱基、N字符和其他字符在内的所有小写字符信息；遍历过程中遇到大写的碱基字符时，将其写入内存；遇到小写的碱基字符时，将其转成大写字符，然后写入内存，并且使用开始位置与长度二元组记录参考基因文件小写字符的位置信息，得到参考基因文件小写字符二元组；

步骤2)、参考基因文件碱基信息Hash值计算；使用k-mer方法将基因分段，对碱基ACGT分别编码，并计算Hash值，使得每一个Hash值代表一个基因片段；

步骤3)、遍历待压缩基因文件；将其中小写字符写入内存，此小写字符信息包括碱基、N字符和其他字符在内的所有小写字符信息；遍历过程中遇到大写的碱基字符、N字符和其他字符时，将其写入内存；遇到小写的碱基字符、N字符和其他字符时，将其转成大写字符，然后写入内存；同时将换行信息也写入内存；并且使用开始位置与长度二元组记录压缩基因文件小写字符的位置信息，得到压缩基因文件小写字符二元组；

步骤4)、小写字符写入压缩文件；对比已经生成的小写字符二元组，将待压缩基因文件的小写字符二元组对比参考基因文件的小写字符二元组，进一步将相同的部分用二元组表示，并写入压缩文件；不同的部分直接写入压缩文件；

步骤5)、换行信息、N字符和其他字符写入压缩文件；将带压缩文件中换行信息、N字符和其他字符以二元组形式表示并将二元组写入压缩文件；

步骤6)、碱基信息写入压缩文件；使用k-mer方法对带压缩基因文件每一段求取Hash值，并根据Hash值匹配得到与参考基因文件相同的k-mer段，记录该k-mer段的开始位置与连续重复最大长度为一个新二元组，直至遍历待压缩基因文件小写字符二元组完成，生成新二元组写入压缩文件。

优选的：步骤2)中使用k-mer方法将基因分段，对碱基ACGT分别编码，并计算Hash值，使得每一个Hash值代表一个基因片段的方法：

首先将A、C、G、T信息编码，分别使用0、1、2、3表示；其次从首位开始，选取连续的11位字符，将第一位的数字与4的10次方相乘，将第二位的数字与4的9次方相乘，以此类推，直至最后一位数字与4的0次方相乘，然后将各个结果相加得到最终的结果即为Hash值，并存入Hash表中。

优选的：步骤2)中Hash值计算方法：

对于一段字符(x₁₀x₉x₈x₇x₆x₅x₄x₃x₂x₁x₀)，其中x_i对应的编码为code(x_i)，则计算Hash值的方式为

HashValue表示Hash值。

本发明相比现有技术，具有以下有益效果：

本发明提出的基于二元组模型的基因压缩方法，提高了基因压缩的效率。基因序列转换成二元组形式后，一个对顺序有严格要求的基因序列变成与顺序无关，有利于利用分布式存储和计算提升基因压缩和分析的效率。

附图说明

图1本发明方法的压缩/解压缩流程图

图2 ref信息提取流程图

图3 Hash表构建流程图

图4 tar信息提取流程图

具体实施方式

下面结合附图和具体实施例，进一步阐明本发明，应理解这些实例仅用于说明本发明而不用于限制本发明的范围，在阅读了本发明之后，本领域技术人员对本发明的各种等价形式的修改均落于本申请所附权利要求所限定的范围。

一种基于参考序列的基因压缩方法，，主要用于解决基因序列数据过大，存储和传输成本高等问题。对基因文件单独压缩具有一个压缩率瓶颈，而基因文件之间存在着高度的相似性，利用参考文件对基因文件进行压缩能够大大地提升基因压缩率。首先从众多待压缩基因序列中选取参考序列，然后将待压缩基因与压缩基因进行运算，达到存储中去除冗余数据的目的，最后对参考序列等进行单独压缩，如图1所示，具体包括以下步骤：

步骤1)、遍历参考基因文件；如图2所示，将其中小写字符写入内存，此小写字符信息包括碱基、N字符和其他字符在内的所有小写字符信息；遍历过程中遇到大写的碱基字符时，将其写入内存；遇到小写的碱基字符时，将其转成大写字符，然后写入内存，并且使用开始位置与长度二元组记录参考基因文件小写字符的位置信息，得到参考基因文件小写字符二元组；

参考序列只需要获取其全部的碱基信息即可，需要将文件中小写字符a、c、g、t表示的碱基信息转换成大写字符A、C、G、T表示，并且使用开始位置与长度二元组记录小写字符的位置信息。

步骤2)、参考基因文件碱基信息Hash值计算；如图3所示，使用k-mer方法将基因分段，对碱基ACGT分别编码，并计算Hash值，使得每一个Hash值代表一个基因片段；

首先将A、C、G、T信息编码，分别使用0、1、2、3表示。其次从首位开始，选取连续的11位字符，将第一位的数字与4的10次方相乘，将第二位的数字与4的9次方相乘，以此类推，直至最后一位数字与4的0次方相乘，然后将各个结果相加得到最终的结果即为Hash值，并存入Hash表中。

步骤3)、遍历待压缩基因文件；如图4所示，将其中小写字符写入内存，此小写字符信息包括碱基、N字符和其他字符在内的所有小写字符信息；遍历过程中遇到大写的碱基字符、N字符和其他字符时，将其写入内存；遇到小写的碱基字符、N字符和其他字符时，将其转成大写字符，然后写入内存；同时将换行信息也写入内存；并且使用开始位置与长度二元组记录压缩基因文件小写字符的位置信息，得到压缩基因文件小写字符二元组；

从待压缩序列中得到小写字符信息、N字符信息、未测明字符信息以及A、C、G、T信息。其中小写字符信息和A、C、G、T信息的获取和参考序列一样取得二元组，N字符信息同样记录开始位置与长度二元组信息，而未测明字符信息由于数量微少，记录其存在的位置以及使用0、1、2、3表示的碱基二元组信息。

步骤4)、小写字符写入压缩文件；对比已经生成的小写字符二元组，将待压缩基因文件的小写字符二元组对比参考基因文件的小写字符二元组，进一步将相同的部分用二元组表示，并写入压缩文件；不同的部分直接写入压缩文件；

参考序列和待压缩序列的碱基信息的小写字符连续相同匹配。小写字符序列在基因序列中也存在着大量的重复，匹配连续重复的小写字符序列可以降低存储。使用顺序和倒序两次匹配的方法，寻求最大连续匹配长度。若不能匹配，则写入(0,1)，如匹配长度为1，则同样写入(0,1)。匹配成功和匹配不成功的第一个序列均放入到diff_pos_vec序列中。

步骤5)、换行信息、N字符和其他字符写入压缩文件；将带压缩文件中换行信息、N字符和其他字符以二元组形式表示并将二元组写入压缩文件；

步骤6)、碱基信息写入压缩文件；使用k-mer方法对带压缩基因文件每一段求取Hash值，并根据Hash值匹配得到与参考基因文件相同的k-mer段，记录该k-mer段的开始位置与连续重复最大长度为一个新二元组，直至遍历待压缩基因文件小写字符二元组完成，生成新二元组写入压缩文件。

从待压缩序列第一组开始计算Hash值并与参考序列的Hash表匹配，若有多组可以匹配，则选取长度最大的一组。并将不能匹配的少许组直接写入到文件中。

1)小写字符二元组生成，参考基因文件和待压缩基因文件都需要生成此二元组，如图2所示。首先设置length＝0并开始遍历基因文件，通过函数判断该字符是否为小写。若不是小写字符，则执行length++，并进行下一个字符判断；若是小写字符，则将length赋值给二元组的第一个位置，同时执行length＝0并迭代，若字符仍为小写，则length++，直至出现下一个大写字符，将length赋值给二元组的第二个位置，同时执行length＝0。

2)碱基信息提取，参考基因文件和带压缩基因文件都需要将碱基信息提取到内存中。遍历基因文件，若是大写字符的碱基信息，则直接写入到ref_seq_code或tar_seq_code中；若是小写字符的碱基信息，则将其转化为对应的大写字符，并写入到ref_seq_code或tar_seq_code中，直至文件读取完成。

3)Hash值计算，如图3所示，将碱基信息A、C、G、T分别编码为0、1、2、3，取k-mer每段长度为11个字符。对于一段字符(x₁₀x₉x₈x₇x₆x₅x₄x₃x₂x₁x₀)，其中x_i对应的编码为code(x_i)，则计算Hash值的方式为

4)小写字符二元组进一步压缩。经过研究表明，基因文件之间小写字符二元组的相似性也非常地高。因此，将待压缩基因文件的小写字符二元组根据参考基因文件的二元组使用新二元组表示，可以进一步降低压缩文件大小，增大压缩率。比较待压缩基因文件第一个二元组，在参考基因文件二元组中遍历查找符合条件，并通过迭代比较出连续重复长度最大的二元组，记录该二元组的开始位置与连续重复最大长度为新二元组。直至遍历待压缩基因文件小写字符二元组完成，将新二元组写入文件。

5)碱基信息匹配。利用计算的Hash值实现碱基信息的匹配，比较待压缩基因文件第一个k-mer段Hash值，在参考基因文件Hash表中遍历查找符合条件的匹配Hash值，并通过迭代比较出连续重复长度最大的k-mer段，记录该k-mer段的开始位置与连续重复最大长度为一个二元组。直至遍历待压缩基因文件小写字符二元组完成，将新二元组写入文件。

压缩的过程具体步骤如下所示：

1.任意选取一个基因序列作为参考序列。

2.读取参考序列，并提取信息，如图2所示。

3.计算Hash值，构建Hash表，如图3所示。

4.待压缩基因序列读取，提取相关信息，如图4所示。

5.匹配小写字符二元组，将大量重复的小写字符二元组表示为连续重复最长的开始位置与长度二元组。

6.匹配Hash值，将大量重复的字符中连续重复最长的开始位置与长度写入文件，并把不能匹配的字符的连续长度以及相关字符直接写入文件。

本发明首先任意选取一个基因序列作为参考序列。其次，获取参考序列的小写字符和ACGT，并以二元组表示小写字符。然后，读取参考文件，获得参考文件的头部、换行信息、小写字符、N字符、碱基信息和其他字符，并将换行长度、小写字符、N字符和其他字符表示成二元组。接着，匹配参考序列和待压缩序列的小写字符二元组，利用他们大量重复的特点再次压缩，将待压缩序列二元组表示成开始重复的位置和长度二元组。最后匹配Hash值，由于ACGT大量重复，只在文件中写入最长重复开始位置和长度，并将未能匹配成功的字符直接编码写入。解压缩过种采用压缩过程相反的步骤。采用本压缩方法的压缩比高，压缩速度快，而且二元组编码与基因次序无关，有利于分布式存储和分析基因序列。

以上所述仅是本发明的优选实施方式，应当指出：对于本技术领域的普通技术人员来说，在不脱离本发明原理的前提下，还可以做出若干改进和润饰，这些改进和润饰也应视为本发明的保护范围。

Claims (3)

1.一种基于参考序列的基因压缩方法，其特征在于，包括以下步骤：

步骤1)、遍历参考基因文件；将其中小写字符写入内存，此小写字符信息包括碱基、N字符和其他字符在内的所有小写字符信息；遍历过程中遇到大写的碱基字符时，将其写入内存；遇到小写的碱基字符时，将其转成大写字符，然后写入内存，并且使用开始位置与长度二元组记录参考基因文件小写字符的位置信息，得到参考基因文件小写字符二元组；

步骤2)、参考基因文件碱基信息Hash值计算；使用k-mer方法将基因分段，对碱基ACGT分别编码，并计算Hash值，使得每一个Hash值代表一个基因片段；

步骤3)、遍历待压缩基因文件；将其中小写字符写入内存，此小写字符信息包括碱基、N字符和其他字符在内的所有小写字符信息；遍历过程中遇到大写的碱基字符、N字符和其他字符时，将其写入内存；遇到小写的碱基字符、N字符和其他字符时，将其转成大写字符，然后写入内存；同时将换行信息也写入内存；并且使用开始位置与长度二元组记录压缩基因文件小写字符的位置信息，得到压缩基因文件小写字符二元组；

步骤4)、小写字符写入压缩文件；对比已经生成的小写字符二元组，将待压缩基因文件的小写字符二元组对比参考基因文件的小写字符二元组，进一步将相同的部分用二元组表示，并写入压缩文件；不同的部分直接写入压缩文件；

步骤5)、换行信息、N字符和其他字符写入压缩文件；将带压缩文件中换行信息、N字符和其他字符以二元组形式表示并将二元组写入压缩文件；

步骤6)、碱基信息写入压缩文件；使用k-mer方法对带压缩基因文件每一段求取Hash值，并根据Hash值匹配得到与参考基因文件相同的k-mer段，记录该k-mer段的开始位置与连续重复最大长度为一个新二元组，直至遍历待压缩基因文件小写字符二元组完成，生成新二元组写入压缩文件。

2.根据权利要求1所述基于参考序列的基因压缩方法，其特征在于：步骤2)中使用k-mer方法将基因分段，对碱基ACGT分别编码，并计算Hash值，使得每一个Hash值代表一个基因片段的方法：

首先将A、C、G、T信息编码，分别使用0、1、2、3表示；其次从首位开始，选取连续的11位字符，将第一位的数字与4的10次方相乘，将第二位的数字与4的9次方相乘，以此类推，直至最后一位数字与4的0次方相乘，然后将各个结果相加得到最终的结果即为Hash值，并存入Hash表中。

3.根据权利要求2所述基于参考序列的基因压缩方法，其特征在于：步骤2)中Hash值计算方法：

对于一段字符(x₁₀x₉x₈x₇x₆x₅x₄x₃x₂x₁x₀)，其中x_i对应的编码为code(x_i)，则计算Hash值的方式为

HashValue表示Hash值。