李玉双，刘倩，张昱

（1.燕山大学 理学院，河北 秦皇岛066004；

2.石家庄邮电职业技术学院 计算机系，河北 石家庄050021）

随着生物科学技术的迅猛发展，生物信息学越来越受到人们的重视，各种研究方法相继产生［1］.近年来，数学模型被引入到该领域，对生物信息学本身而言，这是一次从量变到质变的飞跃.众所周知，数学模型在生物序列和结构的比较中起到了很好的研究效果，在理论方面给出了很好的解释，如几何表示模型［2］、字统计模型［3］和马尔科夫模型［4］等 .隐马尔科夫模型在生物信息学的一系列问题都得到成功应用，如多序列比对［5］、基因识别［6］和蛋白质二级结构预测［7］等.伴随生物研究中数学模型和算法的不断完善，产生了许多强有力的生物信息分析工具，如进化分析、聚类分析等，部分有效的分析工具极大地依赖于生物序列和结构的比较.序列和结构的比较是最重要和最常用的原始操作，是许多其他复杂操作的基础.序列的相似性分析是生物序列和结构比较中的一个重要问题.从序列分析角度，判定两条序列同源与否的一个主要依据是探寻它们之间的相似性.文献［8］提出了转移矩阵，将DNA序列看成是离散的马尔科夫链，分别以碱基A，T，C和G在序列中出现的次数作为基准来构造转移矩阵，进而刻画11个物种的β-globin基因第一个外显子编码序列的差别.本文以序列的长度作为基准，基于碱基组合在DNA序列中出现的频率，构造了DNA序列的频率矩阵.

1 碱基组合的频率矩阵

1.1 频率矩阵的定义

给定长为n的生物序列l＝l1l2l3…ln，li∈S，S＝｛A，T，C，G｝为碱基集合.记 AA在序列中出现的次数为nAA，则定义PAA＝nAA／n.同理，可分别定义PAT，PAC，PAG，PTA，PTT，PTC，PTG，PCA，PCT，PCC，PCG，PGA，PGT，PGC，PGG.这里称 AA，AT，AC，AG，TA，TT，TC，TG，CA，CT，CC，CG，GA，GT，GC，GG为碱基组合.定义该序列对应的频率矩阵P为

由此可知，对应于文献［8］刻画的11个物种的β-globin基因第一个外显子编码序列，可以分别定义相应的频率矩阵，其碱基如表1所示 .表1中：1～11个物种分别是人类（human），家山羊（goat），负鼠目（opossum），原鸡（gallus），狐猴（lemur），小鼠（mouse），兔子（rabbit），老鼠（rat），大猩猩（gorilla），牛科动物（bovine），黑猩猩（chimpanzee）.

表1 11个物种的频率矩阵Tab.1 Frequency matrix of eleven species

从表1可以看到：11个物种中TG出现的频率都是最高，其次是GG，而TA和CG频率较低.这说明在β-globin基因的编码序列中TG和GG相对来说出现频繁，而TA和CG相对出现次数较少，有些物种甚至没有出现 .从单个物种来说，opossum和gallus又有些特殊的地方，例如TG中频率较其他物种偏低，CA中频率较高.这说明了在11个物种的β-globin基因的编码序列中opossum和gallus有着特殊性.上述结果与代琦等［8］的结论基本一致.

1.2 频率矩阵的性质

根据频率矩阵的性质1），可以计算出11个物种碱基含量的向量，即

对于序列的最后一个碱基，虽然它的含量不能通过上述向量中的对应值精确体现（由于计算的是碱基组合），但由于其他3个碱基的含量恰好就是向量中的对应值，所以能够很容易得到最后一个碱基的含量.如在human中，碱基A的含量是0.184 8，碱基T的含量是0.217 3，碱基C的含量是0.206 5，则碱基G的含量是0.391 4.图1为11个物种的碱基含量分布柱状图，可以更直观的展现碱基A，T，C，G在11个物种中的分布情况.

图1 碱基在11个物种中的分布图Fig.1 Distribution of nucleotide of eleven species

观察11个碱基含量向量及图1可以看出：11个物种序列中碱基G的含量都较高，碱基A的含量分布较为均匀；相比其他物种，gallus碱基G的含量明显偏低，lemur碱基C的含量偏低，opossum碱基G的含量偏低；human和gorilla的碱基含量几乎相等.众所周知，研究DNA序列的特殊区域能为基因组的组织结构和生物作用提供更加丰富的信息.这里借助碱基含量向量及图1可以很容易的得出特殊碱基组合的含量，如GC含量.GC含量为基因组提供了数量以及性质上的重要信息，GC含量高的DNA序列要比GC含量低的DNA序列更加稳定［9］.

根据频率矩阵的性质2），可以计算出11个物种的碱基转移向量，即

通过比较碱基含量向量和碱基转移向量不难发现，每个物种的两个向量总有两个分量是相等的.因为前者忽略了序列的最后一个碱基，后者忽略了序列的第一个碱基.如果一个序列首尾碱基相同，则这个序列对应的两个向量一定相等.从这个意义上来说，碱基转移向量也能够反映出各个碱基在序列中的含量分布.此外，除首尾碱基相同的序列（注：这11个物种首尾碱基都不同），不用计算通过比较两个向量就能确定每个物种中各个碱基的含量，如human的碱基转移向量的最后一个分量即为碱基G的含量0.391 2，这与前面计算的结果一致（微小误差是由于计算时舍位引起的）.

2 序列相似性分析

由于生物序列有其进化上的生物学意义，因此比较两条生物的相似性时，不能完全使用计算机科学中的模式匹配，常会借助“距离”来反映，如向量的欧氏距离、协方差距离、夹角距离等.文中引入矩阵的2-范数对11个物种进行相似性比较.

设P1和P2为两个物种的频率矩阵，令Q＝｜P1-P2｜，则Q的2-范数计算公式为

利用2-范数的计算公式来求两个物种的相似性大小，即求得的范数越小，代表两个物种所刻画的DNA序列越相似，两个物种越接近；反之，它们刻画的DNA序列差别越大.利用2-范数的计算公式和常用的欧式距离公式计算得到的11个物种的相似性矩阵，如表2，3所示.

表2 由2-范数算得的11个物种的相似性矩阵Tab.2 Similarity matrix of eleven species based on the 2-norm

表3 由欧氏距离算得的11个物种的相似性矩阵Tab.3 Similarity matrix of eleven species based on the Euclidian distance

比较表2，3可知：2-范数法要比常用的欧氏距离法好，但从整体上看两个方法求得的结果基本一致.即human和gorilla相似性非常高，human和chimpanzee，gorilla和chimpanzee相似性也很高，goat和bovine相似性较高；相比之下，opossum和其他物种相似性较低，这与opossum是与其他哺乳动物亲缘较远的哺乳动物相符合；Gallus和其他物种相似性也较低，这与Gallus是唯一的非哺乳动物相符合.这些结论都与相关的文献结果一致［2，8］.

3 结论

介绍一种利用DNA序列碱基组合的频率矩阵来刻画物种相似性的方法 .该矩阵的每一个分量都能够反映出对应碱基组合在序列中的含量分布情况，其行和能反映每个碱基在序列中的含量分布情况，列和能反映碱基突变的情况，而所有元素值之和为定值.相较文献［8］中的转移矩阵，频率矩阵能够更好地从整体上反映出DNA序列中碱基以及碱基组合的含量分布，显示出序列碱基突变的情况.

文中引入矩阵的2-范数对11个物种进行相似性比较，结果显示该方法要优于上述常用的距离分析方法.频率矩阵的应用在物种的相似性比较方面得到了很好的体现，借助矩阵2-范数和柱状图所得到的结果对物种的进化分析有一定的参考价值.

［1］ 王勇献，王正华.生物信息学导论：面向高性能计算的算法与应用［M］.北京：清华大学出版社，2011：28-72.

［2］ XIE Guo-sen，MO Zhong-xi.Three 3Dgraphical representations of DNA primary sequences based on the classifications of DNA bases and their applications［J］.J Theor Biol，2011，269（1）：123-130.

［3］ VINGA S，GOUVEIA-OLIVEIRA R，ALMEIDA J S.Comparative evaluation of word composition distances for the recognition of SCOP relationships［J］.Bioinformatics，2004，20（2）：206-215.

［4］ PHAM T D，ZUEGG J.A probabilistic measure for alignment-free sequence comparison［J］.Bioinformatics，2004，20（18）：3455-3461.

［5］ 罗泽举，宋丽红.隐马尔可夫模型的多序列比对的研究［J］.计算机工程与应用，2010，46（7）：171-174.

［6］ 丰月姣，贺兴时.二阶隐马尔科夫模型在基因识别中的应用［J］.佳木斯大学学报，2009，27（6）：940-942.

［7］ 石峰，莫忠息，张楚瑜.隐马尔可夫模型-改进的预测蛋白质二级结构方法［J］.生物数学学报，2004，19（2）：233-237.

［8］ 代琦.生物序列、结构比较中若干数学模型研究及应用［D］.大连：大连理工大学，2009：17-71.

［9］ GAO F，ZHANG C T.GC-Profile：A web-based tool for visualizing and analyzing the variation of GC content in genomic sequences［J］.Nucleic Acids Res，2006，34：686-691.