序列比对正文

一前言

DNA蛋白质序列比对是生物信息学中的基本手段，深入了解比对算法的核心内容是学习掌握生物信息学的必需内容。本专题通过 Perl 程序分别实现动态规划算法与模糊匹配的方式（即字符串相似度算法），处理比较多条序列，从而比较两种算法的区别和优劣。进而通过实际操作（程序运行）了解动态规划算法与其它算法之间的区别，掌握动态规划算法的特点和优势。

二本论

动态规划算法

动态规划解决序列比对问题的基本思想:使用迭代法计算出两个序列的相似分值, 并存入一个得分矩阵中, 根据得分矩阵回溯寻找最优的比对序列.

序列比对是动态规划的一个重要应用。序列比对问题通常是使用编辑操作（替换、插入、删除一个要素等）进行序列转换。每次操作对应不同成本，目标是找到编辑序列的最低成本。可以很自然地想到使用递归解决这个问题，序列A到B的最优编辑通过以下措施之一实现：插入B的第一个字符，对A和B的剩余序列进行最优比对；删去A的第一个字符，对A和B进行最优比对；用B的第一个字符替换A的第一个字符，对A的剩余序列和B进行最优比对。局部比对可在矩阵中列表表示，单元(i,j)表示A[1..i]到b[1..j]最优比对的成本[1]。单元(i,j)的成本计算可通过累加相邻单元的操作成本并选择最优解实现。

两条序列分别作为横坐标和纵坐标放置，组成一个路径矩阵，即得分矩阵, 矩阵元素(i,j)值为比对的得分值。在得分矩阵中到达位置为(i,j)的某一个元素有三种可能的路径：通过位置i- 1,j- 1 的对角方向, 没有空位罚分：通过列 j 的垂直方向, 通过行i 的水平方向, 空位罚分的值取决于插入空格的个数。序列比对序列同源(homology)指的是序列来自相同的祖先, 意味着这些序列具有相同的进化历史, 而序列的相似性(similarity)指的是两序列在某参数条件下的相像, 它可以用相同残基的百分比或是其他的方法来表示[2]。列之间的相似度是可以量化的参数，列是否同源需要有进化事实的验证, 显著的相似性通常意味着同源是运用某种特定的数学模型或算法, 找出两个或多个序列之间的最大匹配碱基或残基数, 比对算法的结果在很大程度上反映了序列之间的相似性程度以及它们的生物学特征。

动态规划序列比对算法

动态规划算法比较序列相似度程序如下：

open(INFILE, 'C:\Users\Administrator\Desktop\生信比对算法应用专题比对序列.txt) or die "Can't open\生信比对算法应用专题比对序列.txt $!";

o pen(OUTFILE, '>C:\Users\Administrator\Desktop\') or die "Can't open : $!";

while()

{ if ( $_ =~ / SEQUENCE1:/ )

{ $_=~s/ SEQUENCE1:和.0" 字符串单元存放的数组是从1 开始计数的单元中，开始扫描! 和" 字符串" 如果遇到字符串单元相同的时候权值相乘最后得到一个/92"!./#.0"/92 即为匹配度" 根据/92 的值进行冒泡排序" 权值由大到小排列[4]。

模糊匹配序列比对算法程序

下面是采用模糊匹配方法对 2 个文本中的序列进行的相似度处理程序：

#!/usr/bin/perl

use strict;

use Tie::File;

my $count=0;

tie my @array ,'Tie::File',"数据1",memory=>0 ;

#将数据绑定成数组的形式，给内存为250MB，控制程序的内存利用，维持稳定和性能open OUT1,"数据2" or die "$!";

open IN,">./字符串相似的数据.txt"; #将相似的数据导入到这个表中

my @ip=@array;

my @data=;

for(my $i=0;$i<=$#ip;$i++){

for (@data){ my $s2=$_; chomp($s2); chomp($ip[$i]);

if (leng($ip[$i],$s2)) {

if ((levenshtein($ip[$i],$s2)/leng11($ip[$i],$s2))<{ #我查询的80%的相似度

#print levenshtein($ip[$i],$s2)."\n";

print IN $ip[$i].'和这个很相似'.$s2."\n";

} else {

print $i."\n";

}

print IN1 $i."\n";

} } close OUT1; close IN;

sub levenshtein($$){

my @A=split @B);

my ($i, $j, $cur, $next);

for $i (0..$#A){ $cur=$i+1;

for $j (0..$#B){ $next=min( $W[$j+1]+1, $cur+1,($A[$i] ne $B[$j])+$W[$j]);

$W[$j]=$cur;$cur=$next;} $W[@B]=$next;}

return $next;}

sub min($$$){

if ($_[0] < $_[2]){ pop @_; } else { shift @_; }

return $_[0] < $_[1]? $_[0]:$_[1];

}

sub leng #当2个数据进行模糊匹配的数据，用最长字符的数据底数相除

{ my ($s1,$s2) = @_;

my $len1 = length $s1;

my $len2 = length $s2;

my $ss1=$len1>$len2?$len1:$len2;

my $ss2=$len1<$len2?$len1:$len2;

if ($ss2/$ss1> #只有最大字符和最小字符长度的比较大与80%的时候才调用算法，这样可以提高速度

{ return ($s1,$s2) }

}

sub leng11

{ my ($s1,$s2) = @_;

my $len1 = length $s1;

my $len2 = length $s2;