CN111063394B - 基于基因序列的物种快速查找及建库方法、系统和介质 - Google Patents

Abstract

本发明公开了一种基于基因序列的物种快速查找方法，包括以下步骤：1)输入包含待查找的未知物种的特征信息的集合A，所述特征信息包括未知物种的采样基因序列子序列及其反向互补序列中字典序的最小值K^* _j，i1、采样基因序列子序列在所属基因序列中的起始位置P_j，i1以及判断K^* _j，i1是否为采样基因序列子序列反向互补序列的逻辑符号D_j，i1；2)将集合A分别与预设的已知物种数据库匹配计算相似性；3)根据相似性最大的匹配结果确定待查找基因序列的物种查找结果。本发明通过设置未知物种的采样基因序列子序列及其反向互补序列中字典序的最小值K^* _j，i1，将待查询基因序列中反向互补的子序列调整为序列相同，减少了需要匹配的子序列种类，从而在进行查询时提升了效率。

Description

基于基因序列的物种快速查找及建库方法、系统和介质

技术领域

本发明涉及一种物种查找方法，尤其涉及一种基于基因序列的物种快速查找及建库方法、系统和介质。

背景技术

将基因序列快速定位到多个物种的参考序列在生物信息中有诸多应用场景：(1)对原始基因测序数据进行快速质控，确定物种来源以判定是否存在样本污染；(2)对原始基因测序数据进行有参无损压缩，通常需要快速确定待压缩基因序列的物种来源以选择合适的参考基因序列，进而获得更好的压缩性能；(3)对原始基因数据进行归档管理以及对已归档但物种信息丢失的数据进行快速的测序物种来源鉴定，这些场景均需要对原始序列进行快速的物种判别分析。

传统的基于序列比对的方法如BLAST，通过空间种子序列索引的方法找到查询序列在参考基因序列上的位置(所谓空间种子序列是指对一段连续的基因序列，将其中部分序列进行掩码屏蔽，查找时只关心未屏蔽的部分与参考序列是否完全匹配，如ACGTAGT，屏蔽其第二位和第四位的碱基，则对应种子序列为ANGNAGT,其中N可以是ACGT四种碱基中的任意一种，故与其匹配的序列共有16种)，根据参考序列的查找位置对种子序列两端进行扩展并对扩展序列进行碱基级别的局部动态规划比对，其计算复杂度高；另一方面为了提高查找的灵敏度，通常需要对不同形状的空间种子构建查找索引表，单个形状的种子序列查找表尺寸通常在100G以上，对内存的要求太高，非常不适合对上述场景下的需求做出快速应答。

另一种基于kmer(所谓kmer是指任意基因序列的一段长度为k的子序列)的基因序列物种识别的方法基于签名序列的频谱，所述签名序列通常为k比较小的kmer，通常取k＝4,其基本思想在于不同物种的签名序列的频数分布存在差异，通过对待识别的未知测序基因序列进行频谱统计并与参考物种的频谱向量计算某种相似性距离(如欧几里得距离)，通过设定一定阈值将测序基因序列划分到相应物种，这种方法的优点在于能比较高效地对单个测序序列进行物种识别，但由于第二代测序技术产生的测序序列读长通常小于300bp，因此通过少量序列统计得到的频谱具有很大误差，为此经常需要对整个文件中的测序序列进行统计，在实际应用中不可避免地会遇到IO瓶颈，也无法满足上述场景的需求。

发明内容

本发明要解决的技术问题就在于：针对现有技术存在的技术问题，本发明提供一种查询效率高、占用存储空间小的基于基因序列的物种快速查找方法。

为解决上述技术问题，本发明提出的技术方案为：

一种基于基因序列的物种快速查找方法，包括以下步骤：

1)输入包含待查找的未知物种的特征信息的集合A，所述特征信息包括未知物种的采样基因序列子序列及其反向互补序列中字典序的最小值K^* _j，i1，采样基因序列子序列在所属基因序列中的起始位置P_j，i1以及判断K^* _j，i1是否为采样基因序列子序列反向互补序列的逻辑符号D_j，i1；

2)将集合A分别与预设的已知物种数据库匹配计算相似性；

3)根据相似性最大的匹配结果确定待查找基因序列的物种查找结果。

优选的，步骤1)之前还包括获取待查找的未知物种的特征信息的集合A的步骤，具体为：

对未知物种的基因序列集多点采样，得到采样包，所述采样包包括指定数量的采样基因序列，提取采样包中所有采样基因序列的子序列并进行规范化处理得到未知物种的特征信息的集合A。

优选的，步骤2)具体包括以下步骤

2.1)遍历选取一个已知物种数据库作为当前已知物种数据库；

2.2)遍历选取一个集合A中的特征信息作为当前特征信息，根据当前特征信息匹配当前物种数据库中的内容，并根据匹配结果获得被选取的已知物种参考基因序列中的参考基因段，计算当前特征信息对应的采样基因序列与参考基因段的相似性，所述已知物种参考基因序列为同一已知物种的基因序列首尾相接，并删除基因序列中被屏蔽的碱基后的基因序列；

2.3)判断集合A中的特征信息是否已经遍历完毕，是则选取每条采样基因序列与当前物种参考基因段相似性的最大值并进入下一步，否则跳转步骤2.2)；

2.4)判断已知物种数据库是否已经遍历完毕，是则进入步骤3)，否则跳转步骤2.1)。

优选的，步骤2.2)包括以下步骤：

2.2.1)选取集合A中的一个特征信息作为当前特征信息；

2.2.2)获取当前已知物种数据库对应的位图向量，根据当前特征信息中K^* _j，i1的位图索引结果验证当前已知物种数据库中是否存在匹配当前特征信息的内容，是则进入下一步，否则跳转步骤2.3)；

2.2.3)根据当前特征信息中K^* _j，i1对应的位图向量位置和当前特征信息所占字节数计算偏移量，获取当前已知物种数据库偏移量位置的对应内容，将内容解码获得被选取的已知物种参考基因序列中的参考基因段起始位置，

2.2.4)根据起始位置和当前特征信息对应的采样基因序列的长度在已知物种参考基因序列上截取参考基因段，计算当前特征信息对应的采样基因序列与参考基因段的相似性。

优选的，步骤3)具体为：针对预设的已知物种，分别将同一物种对应的所有相似性结果进行加权计算，使得加权计算结果最大时的已知物种为待查找基因序列的物种查找结果。

本发明还提出一种基于基因序列的物种建库方法，包括以下步骤：

S1)输入包含已知物种的特征信息的集合B，所述特征信息包括已知物种参考基因序列子序列及其反向互补序列中字典序的最小值K^* _i2、已知物种参考基因序列子序列在所属基因序列中的起始位置P_i2以及判断K^* _i2是否为已知物种参考基因序列子序列反向互补序列的逻辑符号D_i2；

S2)设置位图向量，建立位图索引，对集合B中特征信息编码并计算每个特征信息对应的偏移量；

S3)建立与已知物种对应的已知物种数据库，将编码结果按照对应偏移量保存在已知物种数据库的对应位置。

优选的，步骤S1)之前还包括获取包含已知物种的特征信息的集合B的步骤，具体为：

获取一个已知物种的基因序列集，将基因序列首尾相接，并删除基因序列中被屏蔽的碱基，得到已知物种参考基因序列，提取参考基因序列唯一出现的子序列并进行规范化处理得到已知物种的特征信息的集合B。

优选的，步骤S2)包括以下步骤：

S2.1)设置位图向量并初始化位图向量为零向量；

S2.2)遍历选取集合B中的一个特征信息作为当前特征信息，根据当前特征信息中K^* _i2查找位图向量对应位置，若位图向量对应位置元素为0，将该元素置1，进入下一步，若位图向量对应位置元素为1，跳转步骤S2.4)；

S2.3)根据当前特征信息中P_i2及D_i2和当前特征信息所占字节数对当前特征信息编码，并根据当前特征信息中K^* _i2对应的位图向量位置和当前特征信息所占字节数计算编码结果对应的偏移量；

S2.4)判断集合B中的特征信息是否选取完毕，是则跳转步骤S3)，否则返回步骤S2.2)。

本发明还提出一种基于基因序列的物种快速查找系统，包括计算机设备，该计算机设备被编程或配置以执行上述基于基因序列的物种快速查找方法的步骤，或被编程或配置以执行上述基于基因序列的物种建库方法的步骤，或被编程或配置以执行上述基于基因序列的物种快速查找方法和上述基于基因序列的物种建库方法的步骤。

本发明还提出一种计算机可读存储介质，该计算机可读存储介质上存储有被编程或配置以执行上述基于基因序列的物种快速查找方法的计算机程序，或存储有被编程或配置以执行上述基于基因序列的物种建库方法的计算机程序。

与现有技术相比，本发明的优点在于：

(1)本发明通过特征信息匹配预设的已知物种数据库计算相似性，特征信息中包括未知物种的采样基因序列子序列及其反向互补序列中字典序的最小值K^* _j，i1，将待查询基因序列中反向互补的子序列调整为序列相同，减少了需要匹配的子序列种类，从而在进行查询时提升了效率；

本发明还具有以下优点：

(2)本发明通过对未知物种的基因序列集多点采样，得到采样包，采样包中的基因序列数量小于未知物种的基因序列集，根据采样包建立特征信息的集合进行查询，有效减少了当前物种数据库查询量，解决了传统方法需要通过整个基因序列集进行查询时遇到的IO瓶颈。

附图说明

图1为本实施例的基于基因序列的物种快速查找方法的步骤图。

图2为本实施例的基于基因序列的物种快速查找方法的预处理步骤流程图。

图3为本实施例的基于基因序列的物种快速查找方法步骤2)的流程图。

图4为本实施例的基于基因序列的物种建库方法的步骤图。

图5为本实施例的基于基因序列的物种建库方法的预处理步骤流程图。

图6为本实施例的基于基因序列的物种建库方法步骤S2)的流程图。

具体实施方式

以下结合说明书附图和具体优选的实施例对本发明作进一步描述，但并不因此而限制本发明的保护范围。

如图1所示，本发明的基于基因序列的物种快速查找方法，包括以下步骤：

1)输入包含待查找的未知物种的特征信息的集合A，所述特征信息包括未知物种的采样基因序列子序列及其反向互补序列中字典序的最小值K^* _j，i1、采样基因序列子序列在所属基因序列中的起始位置P_j，i1以及判断K^* _j，i1是否为采样基因序列子序列反向互补序列的逻辑符号D_j，i1；

2)将集合A分别与预设的已知物种数据库匹配计算相似性；

3)根据相似性最大的匹配结果确定待查找基因序列的物种查找结果。

如图2所示，本实施例中，步骤1)之前还包括获取待查找的未知物种的特征信息的集合A的步骤，具体为：

获取未知物种的基因序列集S，设置n个采样包以及与采样包一一对应的采样器，设定采样包中最大基因序列数量为N，设定采样频率h，打开n个文件句柄，将文件指针定位到基因序列集S偏移量分别为

位置处；

对于每个采样器，从对应的偏移量位置开始顺序读取基因序列集的基因序列，并将读取到的基因序列作为采样基因序列保存在采样包中，为提高采样质量，本实施例在采样器采样的同时产生均匀分布伪随机数μ，若μ＜h，则采样器将读取的当前序列写入采样包，否则丢弃该序列并读取下一条基因序列，当采样包中采样基因数目大于或等于预设的最大基因序列数量N时，对应的采样器终止采样，本实施例采用逆转法产生均匀分布伪随机数μ，作为产生伪随机数的方法，逆转法相比于直接法和接受拒绝法原理更简单、编程更方便、适用性更广；

采样结束后，丢弃采样基因序列数目小于N的采样包，得到满足条件的采样包的集合，记采样包的集合为DS＝{ds₁,ds₂,...,ds_n}，若所有采样包都不满足条件，则减少采样器个数，重新定位文件指针并开始采样直到得到满足条件的采样包的集合。

选取采样包的集合DS中的一个采样包ds_i作为当前采样包，其中i为当前采样包在采样包的集合DS中的序号，1≤i≤n，对于该采样包中的所有采样基因序列S_i,j，1≤j≤N，可以有以下处理方式：

分别单独对采样基因序列进行处理，即选取采样包中一条采样基因序列作为当前采样基因序列，对当前采样基因序列处理后进入步骤1开始进行物种查找，查找完毕后选取采样包中下一条采样基因序列进行相同处理，直到采样包中的采样基因序列选取完毕；

统一对采样基因序列进行处理，即分别对采样包中的采样基因序列进行处理，得到所有采样基因序列的处理结果后再进入步骤1开始进行物种查找。

本实施例中，为了提高工作效率，防止资源闲置，采用统一对采样基因序列进行处理的方式，分别对于当前采样包中的采样基因序列，以步长1提取长度为k的子序列，考虑到计算机内存和采样基因序列的大小，本实施例中k＝15，并对所有子序列进行规范化处理得到与子序列一一对应的特征信息，每个特征信息分别对应一个采样基因序列，所有特征信息组成当前采样包对应的未知物种的特征信息的集合A。

本实施例对子序列规范化处理得到特征信息的过程具体为：

分别对子序列及其在所属采样基因序列起始位置建立键值对＜K_j，i1,P_j，i1＞，其中K_j，i1为当前采样包序号为j的采样基因序列中第i1个子序列，P_j，i1为子序列在所属采样基因序列的起始位置，其中1≤i1≤n_i，j，n_i，j为当前采样包序号为j的采样基因序列键值对的数目，其值为

其中L_i，j为当前采样包序号为j的采样基因序列长度，s为子序列采样步长，k为子序列长度；

设置包含多种特征信息的元素组，本实施例中为三元组，计算键值对＜K_j，i1,P_j，i1＞中第

～～一个元素K_j，i1的反向互补序列K_j，i1，取K_j，i1和K_j，i1字典序中的较小者，标记为K^* _j，i1，将K^* _j，i1作为三元组的第一个元素，将待查询基因序列中反向互补的子序列调整为序列相同，减少了需要匹配的子序列种类，从而在进行查询时提升了效率，将键值对＜K_j，i1,P_j，i1＞中第二个元素P_j，i1作为三元组的第二个元素，将判断K^* _j，i1是否为子序列反向互补序列的逻辑符号D_j，i1作为三元组的第三个元素，K^* _j，i1为子序列时，D_j，i1＝0，K^* _j，i1为子序列反向互补序列时，D_j，i1＝1，得到与每个键值对一一对应的三元组

如图3所示，本实施例的步骤2)具体包括以下步骤

2.1)遍历选取一个已知物种数据库作为当前已知物种数据库；

2.2)遍历选取一个集合A中的三元组作为当前三元组，根据当前三元组匹配当前物种数据库中的内容，并根据匹配结果获得被选取的已知物种参考基因序列中的参考基因段，计算当前三元组对应的子序列所属采样基因序列与参考基因段的相似性，本实施例中的已知物种参考基因序列为预设的已知物种数据库对应的预处理后的物种基因序列，预处理的过程具体为：将同一已知物种的基因序列首尾相接，并删除基因序列中被屏蔽的碱基，得到长度为L_i0的基因序列；

2.3)判断集合A中的三元组是否已经遍历完毕，是则选取每条采样基因序列与当前物种参考基因段相似性的最大值并进入下一步，否则跳转步骤2.2)；

2.4)判断已知物种数据库是否已经遍历完毕，是则进入步骤3)，否则跳转步骤2.1)。

本实施例的步骤2.2)包括以下步骤：

2.2.1)选取集合A中的一个三元组作为当前三元组，本实施例以三元组

为例继续后续步骤；

2.2.2)获取当前已知物种数据库对应的位图向量，根据当前三元组中K^* _j，i1的位图索引结果验证当前已知物种数据库中是否存在匹配当前三元组的内容，是则进入下一步，否则跳转步骤2.3)，本步骤中的验证过程具体为：

2.2.2a)初始化位图偏移值idx为0，将K^* _j，i1的每个字符映射为对应的自然数，本实施例中的映射方式为将K^* _j，i1的每个字符右移一位，并分别将每个右移后的字符与自然数3做按位与运算得到对应的自然数，其目的是将字符集{A,C,G,T}恒等映射成自然数{0,1,3,2}；

2.2.2b)按照从左往右的顺序遍历选取K^* _j，i1上的字符对应的字符作为当前字符，将idx左移2位，与当前字符对应的自然数做按位或运算，并将计算所得结果作为新的idx；

2.2.2c)判断K^* _j，i1上的字符是否遍历完毕，是则根据idx的值找到位图向量的对应位置，否则返回上一步；

2.2.2d)判断位图向量的对应位置的元素是否为1，本实施例中预设的物种数据库在建库时需要建立位图索引，将相关内容存入物种数据库时，位图向量对应位置元素置1作为标记，因此位图向量的对应位置的元素若为1，说明当前物种数据库中存在当前三元组匹配结果，进入步骤2.2.3)，否则说明当前物种数据库中不存在当前三元组匹配结果，跳转步骤2.3)；

2.2.3)获取当前已知物种数据库偏移量位置的对应内容，将内容解码获得被选取的已知物种参考基因序列中的参考基因段起始位置，本实施例中根据当前三元组中K^* _j，i1对应的位图向量位置和当前三元组所占字节数计算偏移量，本步骤具体为：

2.2.3a)统计位图向量中位于偏移值idx之前为1的元素个数x_i1，获取当前三元组所占字节空间δ_j，偏移量的函数表达式如下

上式中，L_i，j为当前采样包序号为j的采样基因序列长度；

2.2.3b)获取预设的物种数据库中偏移量位置对应内容并解码为元组<Pos,Dir>，其中Pos为起始位置参考值，Dir为判断当前三元组对应的采样基因序列是否为原序列反向互补序列的逻辑符号，计算已知物种参考基因序列中的参考基因段起始位置，函数表达式如下：

上式中，||S_i,j||为当前采样包序号为j的采样基因序列长度，k为子序列的长度，t为当前三元组对应的子序列在序号为j的采样基因序列的起始位置；

2.2.4)根据起始位置和当前三元组对应的采样基因序列的长度在已知物种参考基因序列上截取参考基因段，计算当前三元组对应的采样基因序列与参考基因段的相似性系数，函数表达式如下：

上式中，||S_i,j||为当前采样包序号为j的采样基因序列长度，hd为当前三元组对应的采样基因序列原序列与参考基因段的海明距离,hd’为当前三元组对应的采样基因序列原序列的反向互补序列与参考基因段的海明距离。通过海明距离计算，规避了耗时的动态规划比对，有助于物种识别的速度提升。

由于每个三元组分别对应一个采样基因序列，因此每个相似性系数也分别对应一个采样基因序列，本实施例步骤2.3)中，选取每条采样基因序列与当前物种参考基因段相似性的最大值具体为：根据每个相似性系数对应的采样基因序列对所有相似性系数分组，每组中所有相似性系数均对应同一采样基因序列，选取每一组中的最大相似性系数，每一组的最大相似性系数即为该组对应的采样基因序列与当前物种参考基因段相似性的最大值。

本实施例步骤3)具体为：针对预设的已知物种，分别将同一物种对应的所有相似性结果进行加权计算，使得加权计算结果最大时的已知物种为待查找基因序列的物种查找结果。本实施例的加权计算结果为加权最大相似性指数，加权最大相似性指数的函数表达式如下：

上式中，c_τ为当前已知物种的权重，τ为当前物种编号，m为物种编号最大值，c_k'为编号为k’的已知物种的权重，θ_j,k'为当前采样包中编号为j的采样基因序列与编号为k’的物种的参考基因段的相似性最大值。

为了得到更加精确的结果，本实施例可以在步骤3)之前继续选取采样包的集合DS中的其他采样包作为当前采样包，提取当前采样包所有采样基因序列的子序列并对子序列规范化处理后得到未知物种特征信息的集合，然后重复步骤1)至步骤2)的步骤，直到采样包的集合DS中的采样包选取完毕，最后计算加权最大相似性指数，函数表达式如下：

上式中，c_τ为当前已知物种的权重，τ为当前物种编号，m为物种编号最大值，c_k”为编号为k”的已知物种的权重，θ_i,j,k”为第i个采样包中编号为j的采样基因序列与编号为k”的物种的参考基因段的相似性最大值。

如图4所示，本发明还提出一种基于基因序列的物种建库方法，包括以下步骤：

S1)输入包含已知物种的特征信息的集合B，所述特征信息包括已知物种参考基因序列子序列及其反向互补序列中字典序的最小值K^* _i2、已知物种参考基因序列子序列在所属基因序列中的起始位置P_i2以及判断K^* _i2是否为已知物种参考基因序列子序列反向互补序列的逻辑符号D_i2；

S2)设置位图向量，建立位图索引，对集合B中特征信息编码并计算每个特征信息对应的偏移量；

S3)建立与已知物种对应的已知物种数据库，将编码结果按照对应偏移量保存在已知物种数据库的对应位置。

如图5所示，本实施例的步骤S1)之前还包括获取包含已知物种的特征信息的集合B的步骤，具体为：

获取一个已知物种的基因序列集，将基因序列首尾相接，并删除基因序列中被屏蔽的碱基，得到已知物种参考基因序列，记已知物种参考基因序列长度为L_i0，对已知物种参考基因序列从起始位置开始以步长s提取长度为k的子序列，考虑到计算机内存和常见物种参考基因序列大小，具体实施取k＝15,s＝4；

对于每个子序列，删除与其字符顺序相同的重复子序列，使得处理后的子序列各不相同，成为唯一出现的子序列。

对唯一出现的子序列进行规范化处理，分别对子序列及其在已知物种参考基因序列起始位置建立键值对＜K_i2,P_i2＞，其中K_i2为已知物种参考基因序列中第i2个子序列，P_i2为子序列在所属采样基因序列的起始位置，其中1≤i2≤n_i0，n_i0为当前采样包序号为j的采样基因序列键值对的数目，其值为

其中L_i0为已知物种参考基因序列长度，s为子序列采样步长，k为子序列长度；

设置包含多种特征信息的元素组，本实施例中为三元组，计算键值对＜K_i2,P_i2＞中第一个元素K_i2的反向互补序列

取K_i2和/>

字典序中的较小者，标记为K^* _i2，将K^* _i2作为三元组的第一个元素，将键值对＜K_i2,P_i2＞中第二个元素P_i2作为三元组的第二个元素，将判断K^* _i2是否为子序列反向互补序列的逻辑符号D_i2作为三元组的第三个元素，K^* _i2为子序列时，D_i2＝0，K^* _i2为子序列反向互补序列时，D_i2＝1，得到与每个键值对一一对应的三元组

每一个三元组均对应当前已知物种；

将所有唯一出现的子序列对应的三元组组成当前已知物种的特征信息的集合B。

如图6所示，本实施例的步骤S2)包括以下步骤：

S2.1)设置位图向量的长度为2^k并初始化位图向量为零向量；

S2.2)遍历选取集合B中的一个三元组作为当前三元组，本实施例以三元组

为例继续后续步骤，根据当前特征信息中/>

查找位图向量对应位置，若位图向量对应位置元素为0，将该元素置1，进入下一步，若位图向量对应位置元素为1，跳转步骤S2.4)，为避免重复选取，本实施例中还删除集合B中的当前三元组，本步骤中根据当前三元组中/>

查找位图向量对应位置的具体步骤为：

2.2.1)初始化位图偏移值idx’为0，将

的每个字符映射为对应的自然数，本实施例中的映射方式为将/>

的每个字符右移一位，并分别将每个右移后的字符与自然数3做按位与运算得到对应的自然数，其目的是将字符集{A,C,G,T}恒等映射成自然数{0,1,3,2}；

2.2.2)按照从左往右的顺序遍历选取

上的字符对应的字符作为当前字符，将idx’左移2位，与当前字符对应的自然数做按位或运算，并将计算所得结果作为新的idx’；

2.2.3)判断

上的字符是否遍历完毕，是则根据idx’的值找到位图向量的对应位置，否则返回上一步。

S2.3)根据当前三元组中P_i2及D_i2和当前三元组所占字节数对当前三元组编码，并根据当前特征信息中

对应的位图向量位置和当前三元组所占字节数计算编码结果对应的偏移量，具体步骤如下：

S2.3.1)将当前三元组中P_i2左移一位并与D_i2做按位或运算得到x；

S2.3.2)获取当前三元组所占字节空间δ_i，截取x低位的8*δ_i位作为编码后的值δ_ki，δ_i的函数表达式如下：

上式中，L_i0为已知物种参考基因序列长度；

由于

与位图向量建立了恒等映射，故无需对/>

进行编码存储，进一步节省存储空间，减少了查找时的内存消耗，在查询过程中，本实施例所建立的数据库中的编码结果δ_ki解码为元组<Pos,Dir>的方法为：将δ_ki转化为二进制后左移一位得到Pos，将δ_ki转化为二进制后与1做按位与运算得到Dir；

S2.3.3)统计位图向量中位于偏移值idx’之前为1的元素个数x_i2，根据当前三元组所占字节空间δ_i，编码结果δ_ki对应偏移量的函数表达式如下：

上式中，L_i0为已知物种参考基因序列长度；

S2.4)判断集合B中的三元组是否选取完毕，是则跳转步骤S3)，否则返回步骤S2.2)。

本实施例的步骤S3)中，当前已知物种对应的物种数据库d_l中包括对应的位图查询索引和集合B中所有三元组的编码结果，l为当前已知物种的编号，1≤l≤m，m为已知物种编号最大值，对于所有已知物种，可以依次选取一个已知物种作为当前已知物种，将当前已知物种的基因序列预处理为参考基因序列后，进入步骤1开始建库，直到所有已知物种均建立对应的物种数据库。

本实施例的步骤S3)之后还包括计算已知物种权重的步骤，具体为：

根据当前已知物种参考基因序列中子序列起始位置P_i2以及子序列长度k计算所有子序列覆盖参考基因序列的长度l_i2，l_i2的函数表达式如下：

上式中，k为子序列的长度，l_i2的初始值为k，l_i2-1为截止上一子序列的所有子序列覆盖参考基因序列的长度，P_i2-1为上一子序列在当前已知物种参考基因序列的起始位置；

根据所有子序列覆盖参考基因序列的长度l_i2计算当前已知物种的权重，函数表达式如下：

c_l＝l_i2/L_i0 (9)

上式中，l为当前已知物种的编号，L_i0为当前已知物种参考基因序列长度。

本发明还提出一种基于基因序列的物种快速查找系统，包括计算机设备，该计算机设备被编程或配置以执行上述基于基因序列的物种快速查找方法的步骤，或被编程或配置以执行上述基于基因序列的物种建库方法的步骤，或被编程或配置以执行上述基于基因序列的物种快速查找方法和上述基于基因序列的物种建库方法的步骤。

本发明还提出一种计算机可读存储介质，该计算机可读存储介质上存储有被编程或配置以执行上述基于基因序列的物种快速查找方法的计算机程序，或存储有被编程或配置以执行上述基于基因序列的物种建库方法的计算机程序。

上述只是本发明的较佳实施例，并非对本发明作任何形式上的限制。虽然本发明已以较佳实施例揭露如上，然而并非用以限定本发明。因此，凡是未脱离本发明技术方案的内容，依据本发明技术实质对以上实施例所做的任何简单修改、等同变化及修饰，均应落在本发明技术方案保护的范围内。

Claims (9)

1.一种基于基因序列的物种快速查找方法，其特征在于，包括以下步骤：

1)输入包含待查找的未知物种的特征信息的集合A，所述特征信息包括未知物种的采样基因序列子序列及其反向互补序列中字典序的最小值K^* _j，i1，采样基因序列子序列在所属基因序列中的起始位置P_j，i1以及判断K^* _j，i1是否为采样基因序列子序列反向互补序列的逻辑符号D_j，i1；

2)将集合A分别与预设的已知物种数据库匹配计算相似性，包括：

2.2)遍历选取一个集合A中的特征信息作为当前特征信息，根据当前特征信息匹配当前物种数据库中的内容，并根据匹配结果获得被选取的已知物种参考基因序列中的参考基因段，计算当前特征信息对应的采样基因序列与参考基因段的相似性，所述已知物种参考基因序列为同一已知物种的基因序列首尾相接，并删除基因序列中被屏蔽的碱基后的基因序列，包括：

2.2.1)选取集合A中的一个特征信息作为当前特征信息；

2.2.2)获取当前已知物种数据库对应的位图向量，根据当前特征信息中K^* _j，i1的位图索引结果验证当前已知物种数据库中是否存在匹配当前特征信息的内容，是则进入下一步，否则跳转步骤2.3)，包括：

2.2.2a)初始化位图偏移值idx为0，将K^* _j，i1的每个字符映射为对应的自然数；

2.2.2b)按照从左往右的顺序遍历选取K^* _j，i1上的字符作为当前字符，将idx左移2位，与当前字符对应的自然数做按位或运算，并将计算所得结果作为新的idx；

2.2.2c)判断K^* _j，i1上的字符是否遍历完毕，是则根据idx的值找到位图向量的对应位置，否则返回上一步；

2.2.2d)判断位图向量的对应位置的元素是否为1，若为1，进入步骤2.2.3)，否则跳转步骤2.3)；

2.2.3)根据当前特征信息中K^* _j，i1对应的位图向量位置和当前特征信息所占字节数计算偏移量，获取当前已知物种数据库偏移量位置的对应内容，将内容解码获得被选取的已知物种参考基因序列中的参考基因段起始位置，

2.2.4)根据起始位置和当前特征信息对应的采样基因序列的长度在已知物种参考基因序列上截取参考基因段，计算当前特征信息对应的采样基因序列与参考基因段的相似性；

3)根据相似性最大的匹配结果确定待查找基因序列的物种查找结果。

2.根据权利要求1所述的基于基因序列的物种快速查找方法，其特征在于，步骤1)之前还包括获取待查找的未知物种的特征信息的集合A的步骤，具体为：

对未知物种的基因序列集多点采样，得到采样包，所述采样包包括指定数量的采样基因序列，提取采样包中采样基因序列的子序列并进行规范化处理得到未知物种的特征信息的集合A。

3.根据权利要求1所述的基于基因序列的物种快速查找方法，其特征在于，步骤2)具体包括以下步骤：

步骤2.2)之前包括：

2.1)遍历选取一个已知物种数据库作为当前已知物种数据库；

步骤2.2)之后包括：

2.3)判断集合A中的特征信息是否已经遍历完毕，是则选取每条采样基因序列与当前物种参考基因段相似性的最大值并进入下一步，否则跳转步骤2.2)；

2.4)判断已知物种数据库是否已经遍历完毕，是则进入步骤3)，否则跳转步骤2.1)。

4.根据权利要求1所述的基于基因序列的物种快速查找方法，其特征在于，步骤3)具体为：针对预设的已知物种，分别将同一物种对应的所有相似性结果进行加权计算，使得加权计算结果最大时的已知物种为待查找基因序列的物种查找结果。

5.一种基于基因序列的物种建库方法，应用于权利要求1所述的基于基因序列的物种快速查找方法的已知物种数据库，其特征在于，包括以下步骤：

S1)输入包含已知物种的特征信息的集合B，所述特征信息包括已知物种参考基因序列子序列及其反向互补序列中字典序的最小值K^* _i2、已知物种参考基因序列子序列在所属基因序列中的起始位置P_i2以及判断K^* _i2是否为已知物种参考基因序列子序列反向互补序列的逻辑符号D_i2；

S2)设置位图向量，建立位图索引，对集合B中特征信息编码并计算每个特征信息对应的偏移量；

S3)建立与已知物种对应的已知物种数据库，将编码结果按照对应偏移量保存在已知物种数据库的对应位置。

6.根据权利要求5所述的基于基因序列的物种建库方法，其特征在于，步骤S1)之前还包括获取包含已知物种的特征信息的集合B的步骤，具体为：

获取一个已知物种的基因序列集，将基因序列首尾相接，并删除基因序列中被屏蔽的碱基，得到已知物种参考基因序列，提取参考基因序列中唯一出现的子序列并进行规范化处理得到已知物种的特征信息的集合B。

7.根据权利要求5所述的基于基因序列的物种建库方法，其特征在于，步骤S2)包括以下步骤：

S2.1)设置位图向量并初始化位图向量为零向量；

S2.2)遍历选取集合B中的一个特征信息作为当前特征信息，根据当前特征信息中K'^*查找位图向量对应位置，若位图向量对应位置元素为0，将该元素置1，进入下一步，若位图向量对应位置元素为1，跳转步骤S2.4)；

S2.3)根据当前特征信息中P_i2及D_i2和当前特征信息所占字节数对当前特征信息编码，并根据当前特征信息中K^* _i2对应的位图向量位置和当前特征信息所占字节数计算编码结果对应的偏移量；

S2.4)判断集合B中的特征信息是否选取完毕，是则跳转步骤S3)，否则返回步骤S2.2)。

8.一种基于基因序列的物种快速查找系统，包括计算机设备，其特征在于，该计算机设备被编程或配置以执行权利要求1～4中任意一项所述基于基因序列的物种快速查找方法的步骤，或被编程或配置以执行权利要求5～7中任意一项所述基于基因序列的物种建库方法的步骤，或被编程或配置以执行权利要求1～4中任意一项所述基于基因序列的物种快速查找方法和权利要求5～7中任意一项所述基于基因序列的物种建库方法的步骤。

9.一种计算机可读存储介质，其特征在于，该计算机可读存储介质上存储有被编程或配置以执行权利要求1～4中任意一项所述基于基因序列的物种快速查找方法的计算机程序，或存储有被编程或配置以执行权利要求5～7中任意一项所述基于基因序列的物种建库方法的计算机程序。

Images