Perl基础实例---多序列的最长公共子序列LCS

Posted 2021-04-17 EasyPerl

在生物应用中， 经常需要来比较两个或者多个不同生物dna， 并且求其的最长公共子序列。

最长公共子序列问题（简称LCS），是给定N个序列X1，X2，。。。Xn 求其最长的公共子序列。

求2个序列的最长公共子序列， 一般情况下最佳的方式是动态规划， 今天我们不讲动态规划。

我们有3个序列X，Y，Z， 今天通过暴力搜索的方式来求其最大公共子序列， 该方式不一定是最佳的解决方案， 但一定是可行的解决方案。 

我们定义3个序列：

my $x = 'ATCGGAGCAGT';my $y = 'TTCCAGAGA';my $z = 'TGAG';my @all = sort {length $a <=> length $b} ($x, $y, $z);my $shortest = shift @all;my $length = length $shortest;my $length_tmp = length $shortest;

我们先取出长度最短的一个序列， 以及计算其长度备用。

由于我们需要取最长公共子序列， 意味着这个结果$result必然同时是所有序列的子序列， 那么最长的情况也就是结果$result等于最短的序列$shortest。我们采用暴力破解的方法， 即一一枚举所有的可能， 这样肯定可以得到我们需要的公共子序列。由于我们只需要最大的子序列， 不需要所有的子序列， 因此我们从最长的可能（即$result = $shortest）开始枚举， 依次减短长度， 只要我们得到一个符合条件的， 那么其必须是最大公共子序列（也可能是同长度的几个情况之一）。

设置循环熔断条件， 当所有情况枚举完， 仍然没有结果， 那么就推出：

while($length > 0){ 

  依次从长开始取$shortest序列的子序列，不成功则尝试下一个：

  for my $index(0...$length_tmp - $length){ my $flag = 1;

 我们设置一个标志位， 用来记录中间匹配结果

 由于我们有3个序列， 因此需要匹配多次， 只有所有序列有成功才能算成功

    

 通过substr函数来取子序列

 my $sub = substr($shortest, $index, $length);  for my $a(@all){

 只要有一个序列匹配不成功， 那么意味着不是公共子序列

  则马上退出， 并尝试下一个子序列， 这样可以提高效率

 unless($a =~ /$sub/){ $flag = 0; last; } }
 if($flag == 1){ #一旦成功就退出，因为只需要一个结果 print "$sub\n"; exit; } } #同长度的子序列都失败， 则长度-1， 以此类推  $length--;}