根据基因列表,从总的fasta文件中提取相关的序列,是常常遇到的问题,这个脚本很好地帮助我实现这个动做,很速度。less
用法:perl $0 gene.list *fa > outspa
#! /usr/bin/perl -w use strict; die "perl $0 <lst><fa>\n" unless @ARGV==2; my ($lst,$fa)=@ARGV; open IN,$lst||die; my %ha; map{chomp;$ha{(split)[0]}=1}<IN>; close IN; $fa=~/gz$/?(open IN,"gzip -cd $fa|"||die):(open IN,$fa||die); $/=">";<IN>;$/="\n"; my %out; while(<IN>){ my $info=$1 if(/^(\S+)/); $/=">"; my $seq=<IN>; $/="\n"; $seq=~s/>|\r|\*//g; print ">$info\n$seq" if(exists $ha{$info} && ! exists $out{$info}); $out{$info}=1; } close IN;