我有大量的Fasta格式的多序列,假如要做以下处理:
问题1 去掉特定长度的短序列
若要去除其中小于五十碱基的序列,请问如何操作;
问题2 取motif上下游一定长度的序列(含motif)
假如我要取每个序列中motif为“AAAA”及其上下游10个碱基的序列片段,并输出位置信息,请问如何编程处理
第一个问题,除去 < 50 bp 的序列,可以直接用 bioperl
use Bio::SeqIO;
my $o_seqi = Bio::SeqIO->new( -file => $infile, -format => 'fasta', );
my $o_seqo = Bio::SeqIO->new( -file => ">$outfile", -format => 'fasta', );
while (my $o_seq = $o_seqi->next_seq) { next if ($o_seq->length < 50);
$o_seqo->write_seq($o_seq); } |
|
第二个,同样用 bioperl
use Bio::SeqIO;
my $o_seqi = Bio::SeqIO->new( -file => $infile, -format => 'fasta', );
my $pattern = 'AAAA';
while (my $o_seq = $o_seqi->next_seq) { my $seq = $o_seq->seq; # 提取序列,成为一个字符串
if ($seq =~ /(.{10}$pattern.{10})/) { print $1; } } |
|