perl解决批量除去引物并生成FAST格式文件

问题：在一个文件夹下有大量的拼接好的序列，其所用引物都一致。需求除去上下游引物，得到目的基因，并附上名字，全部放入一个文件中，方便上传到网上BLAST库中。

模型：字符串比对问题，可以转换为对一个长的字符串匹配一段短的字符串得到相关位置，最后取出两个位置中间的片段即可。

算法：由于序列呈不规则并且片段也较长所以采用BM算法。

BM算法：

输入：具有n个字符的串T（文本）和具有m个字符的串P（模式）

输出：T中与P匹配的第一个子串的起始下标，或者表明P不是T的一个子串的指示

计算函数LAST

i<-m-1

j<-m-1

repeat?

if P[j]=T[i] then

if j=0 then

return i{ 匹配一次成功}

else

i<-i-1

j<-j-1

else?

i<-I+m-min(j,1+last(T[i]) {跳跃步骤}

j<-m-1

until i>n-1

return "匹配不成功"

其中算法用perl很容易实现，只是其中的函数LAST有几点需要注意 ，否则容易陷入死循环

1、注意T的字段在P中匹配成功后，从T中0位置到T字符位置的长度是否大于P中0位置到P字符位置的长度。

2、注意Last返回后生成的比对位置，是否方向错误向着反方向跳跃，而导致无限循环中。

针对这两点 我在我的LAST函数中加入了

1、判别T中当前比对的位置是否大于P中比对得到的位置。

2、判别当前P中比对得到的位置是否是大于之前P中比对位置。

整体perl程序思路是：

1、打开目标文件夹中文件

2、用BM算法比对并写入输出文件夹中result.fasta文件中。

3、关闭文件继续执行第一个步骤直到所有文件都打开。

4、用bioperl中fasta转换genebank，和genebank转换fasta程序处理下，来得到正式的fasta文件。

由于基因名字命名规则不同，所以暂时用文件名字代替，可后续按照需求更改。

#!/usr/bin/perl
use strict;
use Bio::SeqIO;
use Getopt::Long;
my $dirname="D:/sutest";#目标文件夹
my $targetFile="E:/sutest/";#输出文件夹
my $result='';
my $begin="ATTGCCCTT";#引物上游序列
my $begin2="TAACGGGA";#引物上游互补序列
my $end="AGGGCAA";#引物下游序列
my $end2="TTCCCGTT";#<span style="font-family: Arial,Helvetica,sans-serif;">引物下游互补序列</span>

my @ArrayStr;
my @ArrayBegin;
my @ArrayEnd;
my $beginFlag=0;
my $endFlag=0;
my @ArraySave;
opendir ( DIR,$dirname ) || die "Error in opening dir $dirnamen";
open(OUTFILE,">>".$targetFile."result.fasta");
while( (my $filename = readdir(DIR))){

if(!($filename=~/^./)&&$filename=~/^[a-zA-Z]/){

            print OUTFILE (">$filename"."n");
open(IN,$dirname."/".$filename)||die "Can not open!";
while(my $seq= <IN>){

if($seq=~/[CGTA]$/)
{
    $result=$result.$seq;

}
}
$beginFlag=&textCompare($result,$begin)+length($begin);

if($beginFlag eq 'none')
{
 $beginFlag=&textCompare($result,$begin2)+length($begin2);

}

$endFlag=&textCompare($result,$end)-1;

if($endFlag eq 'none')
{

 $endFlag=&textCompare($result,$end2)-1;

}


for(my $i = 0; $i < length($result); $i++){
 $ArraySave[$i] = substr($result,$i,1);
}
my $resultLast="";

for(my $i = $beginFlag; $i <= $endFlag; $i++){
      $resultLast.= $ArraySave[$i];

}

print OUTFILE ($resultLast."n");

close(IN);
}
}
close(OUTFILE);

my $informat = 'fasta';
my $outformat= 'genbank';
my $in = Bio::SeqIO->new(-format => $informat,# input handle
                         -file   => $targetFile."result.fasta");
my $out = Bio::SeqIO->new(-format => $outformat,# output handle
                          -file   => '>1_chr1NEW.gbk');
while( my $seq = $in->next_seq ) {
 $out->write_seq($seq);
}

my $infor = 'genbank';
my $outfor= 'fasta';
my $input = Bio::SeqIO->new(-format => $infor,# input handle
                         -file   => '1_chr1NEW.gbk');
my $output = Bio::SeqIO->new(-format => $outfor,# output handle
                          -file   => ">".$targetFile."result.fasta");
while( my $seq = $input->next_seq ) {
 $output->write_seq($seq);
}

unlink (('1_chr1NEW.gbk'));
closedir(DIR);

sub textCompare{

my($str,$begin)=@_;

 my @ArrayStr;
 my @ArrayBegin;
for(my $i = 0; $i < length($str); $i++){
 $ArrayStr[$i] = substr($str,1);
}

for(my $i = 0; $i < length($begin); $i++){
 $ArrayBegin[$i] = substr($begin,1);
}

my $strFlag=length($begin)-1;
my $keyFlag=length($begin)-1;

while($strFlag<=length($str)-1){
    if($ArrayStr[$strFlag]eq$ArrayBegin[$keyFlag])
    {

     if($keyFlag==0)
     {
      return $strFlag;
     last;
     }
     else{

          $strFlag--;
           $keyFlag--;
     }
    }
    else{

      $strFlag=$strFlag+length($begin)-&min(length($begin),(1+&lastFig((length($begin)-1),$ArrayStr[$strFlag],$strFlag,$keyFlag,@ArrayBegin)));

      $keyFlag=length($begin)-1;

    }

}
if($strFlag>length($str)-1)
{

return 'none';
}
  sub min{
  my($str1,$str2)=@_;

        if($str1>$str2){
        return $str2;

        }
        else{
        return $str1;
        }
  }

sub lastFig{
   my($beginLength,$arrayPoint,$totalNum,$currentPoint,@beginPoint)=@_;

  for(my $i = $beginLength;$i>=0;$i--){

  if($arrayPoint eq $beginPoint[$i])
      {

       if($totalNum<$i||$i>$currentPoint){

              next;

          }
          else{
       return $i;
       last;
      }

  }
    }
  return -1;

  }
}