好久不见啦,办公室的工作有点繁忙鸭~今天遇到了一个挺简单的问题,但是用awk/sed/grep来解决好像有点难,就是提取fasta文件中第一条序列。因为很多fasta文件的序列部分都是自动换行的,所以我决定写一个简单的perl脚本来处理!
用法是:perl split_fasta.pl <dirname> <no>
其中dirname
是fasta文件所在目录路径,便于批量处理;no
是想要截取的第几条序列。
split_fasta.pl
脚本内容如下:
use strict;
use warnings;
my $argnum;
my $dirname;
my $no;
my @files;
$argnum = @ARGV;
if($argnum != 2 ){
die "usage: perl split_fasta.pl <dirname> <no>";
};
$dirname = $ARGV[0];
$no = $ARGV[1];
@files = glob("$ARGV[0]/*.fasta");
@ARGV = @files;
foreach my $i(@ARGV){
my $a = 0;
my $fasta = $i;
open IN, "$fasta" or die "fail to open $fasta:$!\n";
open OUT,">$fasta.seq$no.out" or die "fail to create $fasta.seq$no.out:$!\n";
while (<IN>){
chomp;
if ($_ =~ /^>/){
$a ++;
}
if ($a == $no){
print OUT $_, "\n";
}
}
close IN;
close OUT;
}
原理很简单啦,算是给自己一个小小的perl练习!
网友评论