CN102332064B - 基于基因条形码的生物物种识别方法 - Google Patents

Abstract

本发明公开了一种基于基因条形码的生物物种识别方法，步骤如下：1.基因条形码图像与基因条形码图像数据库的制作：从http://www.ncbi.nlm.nih.gov/网站下载617条原核生物体的DNA核苷酸序列建立数据库；根据现有技术制作要识别物种的基因条形码图像。2.对基因条形码图像预处理：采用灰度拉伸将灰度为[0，L]的基因条码图像拉伸到[0，255]；再应用灰度增强，加强基因条码图像的对比度。3.对基因条形码图像中外来基因片段的检索：基因条形码图像的纵向分割；横向外来基因片段的搜索。4.物种识别：确定两个物种之间的相似性度量即确定它们之间的空间距离；根据相似性度量进行物种识别及输出结果。

Description

基于基因条形码的生物物种识别方法

技术领域

本发明涉及一种生物物种的识别方法，更具体地说，本发明涉及一种基于基因条形码的自动物种分类识别方法。

背景技术

物种鉴定和分类隶属于传统生物分类学范畴，对物种进行鉴定分类是对其进行科学研究的基础和前提。虽然科学技术发展在维持生物多样性以确保生物种群安全、保护物种及预防传染性疾病流行等方面解决了许多重要的生物问题，但若想成功鉴定和描述物种，仍然是一个待解决的科研难题。2002年4月，德国科学协会提出了建立一个以DNA为基础的包含所有生物群体在内的分类系统的构想。并首次提出了“生物条形码”的概念，即利用一段短的DNA序列作为物种快速鉴定的标记，并希望以此建立起物种名称(条形码)和生物实体之间一一对应的关系。2003年美国科学院院士卡琳和他的同事以k-mer频率分布的分析为依据研究了各种基因组的性能，并指出，通过对整个基因组高达50K碱基对片段二核苷酸相对丰度的测量，标准化的k-mer频率对于单聚体的频率是总体稳定的。他们甚至提出，这种标准化的k-mer的频率分布可以作为基因组标签。2008年，周峰峰等利用条形码为所有的基因的DI-Mer频率进行灰度排序的思路，发现基因自身具有条形码的特征，并将其命名为基因条形码，该基因条形码的主要优点是为基因组提供了一种直观化、全局化的信息，使各种基因组的特征变得一目了然。他们通过研究发现基因条形码有很多有意义的属性：

1.所有的染色体基因组都有显著稳定的k-mer频率分布，从而产生一致的跨每个条码垂直条带灰度值；

2.片段的一小部份有明显不同和异常，与其他基因组的条形码(条形码中横条纹)的对比说明它有2-3类特殊的基因；

3.同一生物体的多个染色体一般有高度相似的条形码，但它们都有其异常片段的独特模式；

4.条形码的相似程度往往与‘系统发育的近似性’成正比。

发明内容

本发明所要解决的技术问题是首次提出用图形处理结合基因条形码进行物种种属鉴定识别，提供了一种基于基因条形码的生物物种识别方法。

为解决上述技术问题，本发明是采用如下技术方案实现的：所述的基于基因条形码图像的生物物种识别方法包括如下步骤：

1.基因条形码图像与基因条形码图像数据库的制作：

1）从http://www.ncbi.nlm.nih.gov/网站下载617条原核生物体的DNA核苷酸序列建立数据库。

2）基因序列是由4种核苷酸A、T、G、C构成的线性序列，首先将物种的全基因组序列切割成长度为M个碱基点的不重叠DNA片段，然后在每一段DNA片段中计算该k-mer频率与它们的反向互补k-mer频率组合而成的综合频率，如当k参数取4时即是4-mer频率，其频率为所有4个碱基点与其反向互补序列的频率的商，将计算所得的综合频率赋予一个灰度值，则基因构成了一幅灰度图，即基因条形码图像，每个基因组条形码图像最终是一个矩阵，矩阵中的每一个元素代表相应序列片段内的相应频率，本发明中取k-mer频率的k参数为4，基因片段的长度M=3000来制作基因条形码图像。

2.对基因条形码图像的预处理：

1）采用灰度拉伸将基因条码图像拉伸到[0,255]。

2）再应用灰度增强，加强基因条码图像的对比度。

3.对基因条形码图像中外来基因片段的检索：

1）基因条形码图像的纵向分割，所述的基因条形码图像的纵向分割包括如下步骤：

a.假设基因条形码图像为I，我们创建二维特征描述向量作为输入向量：

Z=(I(x,y),|▽I(x,y)|)

其中：I(x,y)和|▽I(x,y)|分别为基因条形码图像上第x行y列的灰度值和梯度值；

b.对输入向量中的灰度值I(x,y)进行处理

用灰度统计法确定阈值T，把灰度值I(x,y)按灰度阈值T分成两类即c₀与c₁，其中，c₀为灰度偏黑的点的集合，c₁为灰度偏白的点的集合，L-1为基因条形码图像的最大灰度；

$I_1(x,y)=\left\{\begin{array}{cc}{c}_0& I(x,y)\∈[0,T]\\ c_1& I(x,y)\∈[T+1,L-1]\end{array}\right.$

c.对输入向量中的梯度值|▽I(x,y)|进行处理

利用像素灰度的一阶导数算子在灰度迅速变化处得到极值来进行突变点的检测，寻找其中梯度变化较大的点即方差突变点的集合；用基因条形码图像的一阶差分直接代替基因条形码图像函数的导数，经过处理后的离散基因条形码图像函数I₁(x,y)在X方向的差分

和Y方向的差分

可通过以下公式计算：

$I_1^x{=I}_1(x,y)*(-\mathrm{1,1})$

$I_1^y=I_1(x,y)*\left(\begin{array}{c}-1\\ 1\end{array}\right)$

其中：“*”是卷积运算，由于基因条形码图像整体呈垂直的条纹状，所以只需要对图像I₁(x,y)计算X方向的差分

即可，所得到基因条形码图像像素值为： $I_1^x=|\▿I(x,y)|=|I(x,y)-I(x-1,y)|$

令输出向量Z满足

则当Z＝1时为突变点，取这些点所在基因条形码图像的坐标位置将预处理后的基因条形码图像纵向分割为H段灰度条，即：I＝{I₁,I₂,...,I_H},每份作为一个子类。

2）横向外来基因片段的搜索,步骤如下：

采用基因条形码图像的纵向分割方法把基因条形码图像切割成H份I＝{I₁,I₂,...,I_H}，每份作为一个子类，根据以下步骤进行横向搜索外来DNA片段；

a.应用灰度统计法计算每个子类I_k的阈值其中k＝1,2,…,H；

b.定义聚类标示矩阵为{r_i,j}_m*n来标识观测样本中每个成员的归属问题，这样定义{r_i,j}_m*n的元素值：对每个子类I_k若

则认为图像偏白，我们需要搜索的目标点为偏黑的点 $r_{i,j}=\left\{\begin{array}{cc}0& r_{i,j}\≥\stackrel{\‾}{T_k}\\ 1& \mathrm{else}\end{array}\right.;$ 若 $\stackrel{\‾}{T_k}\≤\mathrm{mean}\left(I_k(x,y)\right),$ 则认为图像偏黑，我们需要搜索的目标点为偏白的点 $r_{i,j}=\left\{\begin{array}{cc}1& r_{i,j}\≥\stackrel{\‾}{T_k}\\ 0& \mathrm{else}\end{array}\right.;$

c.取矩阵{r_i,j}_m*n中r_i,j＝1的点为目标点进行聚类，切割每个类中目标点所在区域的最大域，即是需要过滤的外来插入DNA片段。

4.物种识别：

（1）采用的是4阶Minkovsky距离来测算物种基因间的相似性，设物种集为{X_i,i＝1,...,N}，X_i＝(X_i1,...,X_id)^T，d为特征向量的维数，N为特征向量的个数，Minkovsky距离由下式计算：

$S(X_i,X_j)={\left[\underset{l=1}{\overset{d}{\mathrm{\Σ}}}{|X_{\mathrm{il}}-X_{\mathrm{jl}}|}^q\right]}^{\frac{1}{q}}$ 其中q=1、2、3、4

（2）根据相似性度量进行物种识别及输出结果。

与现有技术相比本发明的有益效果是：

1.本发明所述的基于基因条形码的生物物种识别方法主要是通过生物物种的核苷酸序列来制作能独特而唯一的代表该物种的条码图，并利用该物种的条码图特征来进行物种识别，这样做的主要优势在于为各种物种的种类识别提供了一种直观，快捷的方法，使各种物种的基因组特征识别变得一目了然。

2.本发明所述的基于基因条形码的生物物种识别方法通过切除该生物物种中的外来物种插入基因片段，得到更为准确的表达该物种基因特征的条形码，从而提高了物种识别的精度，速度和稳定性。

3.本发明所述的基于基因条形码的生物物种识别方法尚属国内首次提出用图形处理结合基因条形码进行物种种属鉴定识别这一方法。在大规模物种基因条形码数据库上的测试结果证实了该算法的优良性能。在生物物种种属鉴定识别中，经过该算法过滤后的基因条形码比不经过过滤的基因条形码精度更高。由于基因条形码图像具有来源于不同生理位置其条形码主框架相同，但所插入的外来物种基因片段不同这一属性。我们还选了一些原核生物体细胞核，质粒体的核苷酸序列，按本发明所提算法进行物种识别，发现识别精度也很高，所以，该算法能够有效地区分不同的生物物种，是一种高效、稳定的物种鉴定识别方法。

附图说明

下面结合附图对本发明作进一步的说明：

图1是本发明所述的基于基因条形码的生物物种识别方法的流程框图；

图2-a是本发明所述的基于基因条形码的生物物种识别方法实施例中采用大肠杆菌k12的原始基因条形码图像；

图2-b是本发明所述的基于基因条形码的生物物种识别方法实施例中将图2-a中的原始基因条形码图像经过灰度拉伸后的图像；

图2-c是本发明所述的基于基因条形码的生物物种识别方法实施例中将图2-b中所得到的基因条形码图像进行增强后的条形码图像；

图3-a是本发明所述的基于基因条形码的生物物种识别方法实施例中将图2-c的基因条形码图像进行纵向划分偏白、偏黑区域的标志图；

图3-b是本发明所述的基于基因条形码的生物物种识别方法实施例中图3-a的基因条形码图像的纵向突变点的标示图；

图3-c是本发明所述的基于基因条形码的生物物种识别方法实施例中对图3-b的基因条形码图进行纵向分割后的条形码图像；

图3-d是本发明所述的基于基因条形码的生物物种识别方法实施例中的基因条形码图像进行横向切除外来插入基因片段后的图像。

具体实施方式

下面结合附图对本发明作详细的描述：

本发明所述的基于基因条形码的生物物种识别方法是在前人的研究的基础上，按照他们所提出的k-mer频率灰度映射制作基因条形码图像的方法制作条形码，然后通过图像处理方法，将基因条形码图像中的外来物种插入基因进行过滤，最后应用概率统计算法测算每个条形码之间的差别距离，从而实现对生物体所属物种的鉴别、分类。

所述的基于基因条形码的生物物种识别方法包括如下步骤：

1.基因条形码图像与基因条形码图像数据库的制作

1）从http://www.ncbi.nlm.nih.gov/网站下载617条原核生物体的DNA核苷酸序列建立数据库。

2）参阅图2-a，我们使用文献[Zhou,F.et al.(2008)Barcodes forgenomes and applications.BMC Bioinformatics,9,546.]中所提出的方法来制作所要识别的生物物种的基因条形码图像：基因序列是由4种核苷酸A、T、G、C构成的线性序列，首先，把每个基因序列划分为一系列非重叠、等尺寸（尺寸为M个碱基点这么长）的片段；然后计算该k-mer频率，如当k-mer=4时，其频率为所有4个碱基点与其反向互补序列的频率的商:GGTA/TACC,CGAA/TTGC,GGTC/GACC……；最后将所计算的每个K-mer频率赋予一个灰度值，则基因构成了一幅灰度图，即基因条形码图像。每个基因组条形码图像最终是一个矩阵，矩阵中的每一个元素代表相应序列片段内的相应的k-mer频率。本发明中取参数k-mer=4，M=3000来制作基因条形码图像。

本发明所述的基于基因条形码的生物物种识别方法中取大肠杆菌k12为例。

2.对基因条形码图像的预处理

为了让基因条形码图像信息灰度级丰富而动态范围大，要对基因条码图像进行预处理。

1）参阅图2-b，首先采用灰度拉伸将基因条码图像拉伸到[0,255]。

2）参阅图2-c，再应用灰度增强，加强基因条码图像的对比度。

这样预处理使基因条形码图像各个纵向纹理间的对比度增强，为后续精确搜索外来插入基因区域奠定基础。

3.对基因条形码图像中外来基因片段的检索

我们用以下方法确定插入基因片段的位置：

1）基因条形码图像的纵向分割

a.假设基因条形码图像为I，我们创建以下二维特征描述向量作为输入向量：

Z=(I(x,y),|▽I(x,y)|)

其中：I(x,y)和|▽I(x,y)|分别为基因条形码图像上第x行y列的灰度值和梯度值。

b.对输入向量中的灰度值I(x,y)进行处理

用灰度统计法确定阈值T，假设基因条形码图像为I(x,y)，把像素按灰度阈值T分成两类即c₀与c₁。其中c₀为灰度偏黑的点的集合，c₁为灰度偏白的点的集合。如图3-a所示，红色所标记的点为纵向区域偏白的点,L-1为基因条形码图像的最大灰度。

$I_1(x,y)=\left\{\begin{array}{cc}{c}_0& I(x,y)\∈[0,T]\\ c_1& I(x,y)\∈[T+1,L-1]\end{array}\right.$

c.对输入向量中的梯度值|▽I(x,y)|进行处理

这一步主要是利用像素灰度的一阶导数算子在灰度迅速变化处得到极值来进行突变点的检测。寻找其中梯度变化较大的点（即方差突变点）的集合。我们用基因条形码图像的一阶差分直接代替基因条形码图像函数的导数。经过a步骤处理后的离散基因条形码图像函数I₁(x,y)在X方向的差分

和Y方向的差分

可通过以下公式计算：

$I_1^x{=I}_1(x,y)*(-\mathrm{1,1})$

$I_1^y=I_1(x,y)*\left(\begin{array}{c}-1\\ 1\end{array}\right)$

其中“*”是卷积运算，由于基因条形码图像整体呈垂直的条纹状，所以我们只需要对图像I₁(x,y)计算X方向的差分

即可。

所得到基因条形码图像像素值 $I_1^x=|\▿I(x,y)|=|I(x,y)-I(x-1,y)|$

令输出向量Z满足

则当Z＝1时为突变点，如图3-b所示，红色标记点为纵向区域灰度发生突变的点的位置。取这些点所在基因条形码图像的坐标位置进行纵向基因条形码图像分割，分割效果如图3-c所示。

2）横向外来基因片段的搜索

根据以下步骤进行横向搜索外来基因片段：

a.用纵向分割方法把基因条形码图像切割成H份即：I＝{I₁,I₂,·,I_H}，每份作为一个子类。

b.根据灰度统计法计算每个子类I_k的样本阈值其中k＝1,2,…,H。

c.定义聚类标示矩阵{r_i,j}_m*n来标识观测样本中每个成员的归属问题，这样定义{r_i，j}_m*n的元素值：对每段灰度条I_k若

则认为图像偏白，我们需要搜索的目标点为偏黑的点 $r_{i,j}=\left\{\begin{array}{cc}0& r_{i,j}\≥\stackrel{\‾}{T_k}\\ 1& \mathrm{else}\end{array}\right.;$ 若 $\stackrel{\‾}{T_k}\≤\mathrm{mean}\left(I_k(x,y)\right),$ 则认为图像偏黑，我们需要搜索的目标点为偏白的点 $r_{i,j}=\left\{\begin{array}{cc}1& r_{i,j}\≥\stackrel{\‾}{T_k}\\ 0& \mathrm{else}\end{array}\right.;$

d.取矩阵{r_i,j}_m*n中r_i,j＝1的点为目标点进行聚类，切割每个类中目标点所在区域的最大域，即是我们需要过滤的外来插入DNA片段。如图3-d所示，为我们横向切除外来插入基因片段后的基因条形码图像。

4.物种识别

1）确定两个物种之间的相似性度量即确定它们之间的空间距离：

将基因条形码图像的外来插入基因片段过滤掉之后，进一步就要考虑如何利用基因条形码图像识别生物物种。这一过程的核心是确定两个物种之间的相似性度量，作为两个物种之间的相似性，最直观的度量就是确定它们之间的空间距离，本发明采用的是4阶Minkovsky距离来测算物种基因间的相似性。设物种集为{X_i,i＝1,...,N}，X_i＝(X_i1,...,X_id)^T，d为特征向量的维数，N为特征向量的个数，Minkovsky距离由下式计算：

$S(X_i,X_j)={\left[\underset{l=1}{\overset{d}{\mathrm{\Σ}}}{|X_{\mathrm{il}}-X_{\mathrm{jl}}|}^q\right]}^{\frac{1}{q}}$ 其中q=1、2、3、4。

2）根据相似性度量进行物种识别及输出结果。

将要识别的物种与数据库中存在的物种分别计算它们之间的Minkovsky距离，其中距离最小者就是物种辨识的结果。

结果及分析

1.生物体基因条形码图像数据库

我们从http://www.ncbi.nlm.nih.gov/网站下载617条原核生物体的DNA核苷酸序列，利用文献[Zhou,F.et al.(2008)Barcodes for genomes andapplications.BMC Bioinformatics,9,546.]所提的k-mer综合频率制作基因条码图像的方法制作基因条形码图像。其中取参数k-mer=4，M=3000来制作基因条形码图像，用于建立数据库。

2.生物种属识别精度

首先预处理录入的DNA序列基因条形码图像，然后对该基因条形码图像的外来插入基因片段过滤切除，得到基因条形码图像特征。将该特征与原核生物数据库中已处理后的各种生物物种条形码图像特征进行匹配，计算其Minkovsky距离，其中距离最小者就是物种辨识的结果。从所建立的原核生物物种数据库中随机取出一副基因条形码图像作为模板，其余条形码图像作为测试样本，本发明所提出的算法的精确度可达99.15%，而不过滤切除外来插入基因片段的识别精度为96.63%，显然，在生物物种识别中，本发明所提方法提高了识别精度，具有一定优越性。

3.处理速度

对联机生物识别系统来说，处理速度至关重要，所有的实验是在IntelPentium2处理器（500MHz）的个人电脑上用matlab编程实现，表一列出了本文所提算法的预处理，特征提取以及匹配所需的平均时间。从该表可知，进行一次物种种属辨识需要大约（268+157+1.2*N）/1000秒（N是模板的数目），若N取1000，则进行一次物种辨识所需时间小于2秒，这完全能够满足一个实时生物识别系统的要求。

表一：算法处理速度

识别阶段	预处理	过滤外来插入基因片段	物种识别
				所需时间（毫秒）	268	157	1.2

实施例1：空肠弯曲杆菌亚种NCTC11168杆菌的鉴定

我们对空肠弯曲杆菌家族的亚种NCTC11168进行鉴定，根据本发明所提出的方法对空肠弯曲杆菌亚种NCTC11168杆菌的基因条形码图像进行制作，预处理，纵向分割，横向检索、切割外来插入物种基因片段，最后与原核生物数据库里的物种计算相识性，得到空肠弯曲杆菌亚种NCTC11168杆菌与同族细菌的相似性度量如表二所示，空肠弯曲杆菌亚种NCTC11168杆菌与其它异族细菌的相似性度量如表三所示：

表二：NCTC11168杆菌与同族细菌的相似性度量

空肠弯曲杆菌亚种NCTC11168杆菌	0
		空肠弯曲杆菌亚种81-176杆菌	30
空肠弯曲杆菌亚种81116杆菌	23
		空肠弯曲杆菌RM1221杆菌	33.4
空肠弯曲杆菌亚种269.97杆菌	24.5

表三：NCTC11168杆菌与异族细菌的相似性度量

隐藏嗜酸菌JF-5	515.78
		氧化亚铁硫杆菌ATCC23270	1033.2
鲍氏不动杆菌AB307-0294	738.18
		蜡样芽胞杆菌AH187	677.74
脆弱拟杆菌YCH46	970.2
		伯氏疏螺旋体菌B31	924.55
砂眼披衣菌D/UW-3/CX	794.48
		番茄溃疡病菌NCPPB382	1394.6

由表二可见，空肠弯曲杆菌亚种NCTC11168杆菌与其自身的相似性度量为零，即最相似，与其同族的其他杆菌距离均小于50。由表三可见与其不同族的细菌相似性度量均大于500。由此可见本发明所提出的方法能明显的将物种自身及其所属家族识别出来。

实施例2：砂眼披衣菌D/UW-3/CX的鉴定

我们选定砂眼披衣菌D/UW-3/CX做第二组物种鉴定识别实验，根据本发明所提出的方法对砂眼披衣菌D/UW-3/CX的基因条形码图像进行制作、处理、最后与原核生物数据库里的物种计算相识性，得到砂眼披衣菌D/UW-3/CX与同族细菌的相似性度量如表四所示，砂眼披衣菌D/UW-3/CX与其它异族细菌的相似性度量如表五所示：

表四：砂眼披衣菌D/UW-3/CX与同族细菌的相似性度量

砂眼披衣菌D/UW-3/CX	0
		砂眼披衣菌434/Bu	28.94
砂眼披衣菌A/HAR-13	31.462
		砂眼披衣菌L2b/UCH-1/proctitis	31.292

表五：砂眼披衣菌D/UW-3/CX与异族细菌的相似性度量

隐藏嗜酸菌JF-5	712.6
		氧化亚铁硫杆菌ATCC23270	658.12
鲍氏不动杆菌AB307-0294	1093.9
		蜡样芽胞杆菌AH187	1210.3
脆弱拟杆菌YCH46	526.98
		伯氏疏螺旋体菌B31	838.13
空肠弯曲杆菌亚种NCTC11168	792.59
		番茄溃疡病菌NCPPB382	937.88

由表四可见，砂眼披衣菌D/UW-3/CX与其自身的相似性度量为零，即最相似，与其同族的其他杆菌距离均小于50。由表五可见与其不同族的细菌相似性度量均大于500。由此可见本发明所提出的方法能明显的将物种自身及其所属家族识别出来。

实施例3：流感嗜血杆菌Rd KW20

我们选定流感嗜血杆菌Rd KW20做第三组物种鉴定识别实验，根据本发明所提出的方法对流感嗜血杆菌Rd KW20的基因条形码图像进行制作、处理、最后与原核生物数据库里的物种计算相似性，得到流感嗜血杆菌Rd KW20与同族细菌的相似性度量如表四所示，与其它异族细菌的相似性度量如表五所示：

表六：流感嗜血杆菌Rd KW20与同族细菌的相似性度量

流感嗜血杆菌Rd KW20	0
		流感嗜血杆菌86-028NP	23.79
流感嗜血杆菌PittEE	24.51
		流感嗜血杆菌PittGG	24.86

表七：流感嗜血杆菌Rd KW20与异族细菌的相似性度量

海栖热袍菌MSB	964.24
		生殖器支原体G37	572.82
蓝藻PCC6803	820.9
		肺炎支原体M129	521.64

大肠杆菌K12	807.72
		幽门螺杆菌26695	593.18
超嗜热菌VF5	892.52
		梅毒螺旋体亚种	963.6

由表六可见，流感嗜血杆菌Rd KW20与其自身的相似性度量为零，即最相似，与其同族的其他细菌距离均小于50。由表七可见与其不同族的细菌相似性度量均大于500。由此可见本发明所提出的方法能明显的将物种自身及其所属家族识别出来。

Claims (1)

1.一种基于基因条形码的生物物种识别方法，其特征在于，所述的基于基因条形码的生物物种识别方法包括如下步骤：

1）基因条形码图像与基因条形码图像数据库的制作：

（1）从http://www.ncbi.nlm.nih.gov/网站下载617条原核生物体的DNA核苷酸序列建立数据库；

（2）基因序列是由4种核苷酸A、T、G、C构成的线性序列，首先将物种的全基因组序列切割成长度为M个碱基点的不重叠DNA片段，然后在每一段DNA片段中计算k-mer频率与它们的反向互补k-mer频率组合而成的综合频率，如当k参数取4时即是4-mer频率，其频率为所有4个碱基点与其反向互补序列的频率的商，将计算所得的综合频率赋予一个灰度值，则基因构成了一幅灰度图，即基因条形码图像，每个基因组条形码图像最终是一个矩阵，矩阵中的每一个元素代表相应序列片段内的相应频率，取k-mer频率的k参数为4，基因片段的长度M=3000来制作基因条形码图像；

2）对基因条形码图像的预处理：

（1）采用灰度拉伸将基因条码图像的灰度级拉伸到[0，255]区间；

（2）再应用灰度增强，加强基因条码图像的对比度；

3）对基因条形码图像中外来基因片段的检索：

（1）基因条形码图像的纵向分割，所述的基因条形码图像的纵向分割包括如下步骤：

a.假设基因条形码图像为I，创建二维特征描述向量作为输入向量：

$Z=(I(x,y),|\▿I(x,y)\left|\right)$

其中：I(x,y)和

分别为基因条形码图像上第x行y列的灰度值和梯度值；

b.对输入向量中的灰度值I(x,y)进行处理

用灰度统计法确定阈值T，把灰度值I(x,y)按灰度阈值T分成两类即c₀与c₁，其中，c₀为灰度偏黑的点的集合，c₁为灰度偏白的点的集合，L-1为基因条形码图像的最大灰度；

$I_1(x,y)=\left\{\begin{array}{cc}{c}_0& I(x,y)\∈[0,T]\\ c_1& I(x,y)\∈[T+1,L-1]\end{array}\right.$

c.对输入向量中的梯度值

进行处理

利用像素灰度的一阶导数算子在灰度迅速变化处得到极值来进行突变点的检测，寻找其中梯度变化较大的点即方差突变点的集合；用基因条形码图像的一阶差分直接代替基因条形码图像函数的导数，经过处理后的离散基因条形码图像函数I₁(x,y)在X方向的差分

和Y方向的差分

可通过以下公式计算：

$I_1^x=I_1(x,y)*(-\mathrm{1,1})$

$I_1^y=I_1(x,y)*\left(\begin{array}{c}-1\\ 1\end{array}\right)$

其中：“*”是卷积运算，由于基因条形码图像整体呈垂直的条纹状，所以只需要对图像I₁(x,y)计算X方向的差分

即可，所得到基因条形码图像像素值为： $I_1^x=|\▿I(x,y)|=|I(x,y)-I(x-1,y)|$

令输出向量Z满足

$Z=\left\{\begin{array}{cc}1& (I(x,y)!=0)\&\&\left(\right|\▿I(x,y)|!=0)\\ 0& \mathrm{else}\end{array}\right.$

则当Z=1时为突变点，取这些点所在基因条形码图像的坐标位置将预处理后的基因条形码图像纵向分割为H段灰度条，即：I={I₁,I₂,...,I_H},每份作为一个子类；

（2）横向外来基因片段的搜索,步骤如下：

采用基因条形码图像的纵向分割方法把基因条形码图像切割成H份I={I₁,I₂,...,I_H}，每份作为一个子类，根据以下步骤横向搜索外来DNA片段；

a.应用灰度统计法计算每个子类I_k的阈值

，其中k=1,2,...,H；

b.定义聚类标示矩阵为{r_i,j}_m*n来标识观测样本中每个成员的归属问题，这样定义{r_i,j}_m*n的元素值：对每个子类I_k若

则认为图像偏白，搜索的目标点为偏黑的点 $r_{i,j}=\left\{\begin{array}{cc}0& r_{i,j}\≥\stackrel{\‾}{T_k}\\ 1& \mathrm{else}\end{array}\right.;$ 若 $\stackrel{\‾}{T_k}\≤\mathrm{mean}\left(I_k(x,y)\right),$ 则认为图像偏黑，搜索的目标点为偏白的点 $r_{i,j}=\left\{\begin{array}{cc}1& r_{i,j}\≥\stackrel{\‾}{T_k}\\ 0& \mathrm{else}\end{array}\right.;$

c.取矩阵{r_i,j}_m*n中r_i,j=1的点为目标点进行聚类，切割每个类中目标点所在区域的最大域，即是需要过滤的外来插入DNA片段;

4）物种识别：

（1）采用的是4阶Minkovsky距离来测算物种基因间的相似性，设物种集为{X_i,i=1,...,N}，X_i=(X_i1,...,X_id)^T，d为特征向量的维数，N为特征向量的个数，Minkovsky距离由下式计算：

$S(X_i,X_j)={\left[\underset{l=1}{\overset{d}{\mathrm{\Σ}}}{|X_{\mathrm{il}}-X_{\mathrm{jl}}|}^q\right]}^{\frac{1}{q}}$     其中q=1、2、3、4

（2）根据相似性度量进行物种识别及输出结果。

Images