Question

perl编程处理DNA序列文件fasta格式，输出大于某数值bp的文件

如题，比如有个序列文件“dna.fasta”里面有5条序列：

>ATGCATGCCC
>ATGC
>ATGCATG
>ATG
>ATGCAT
要从中筛选出大于5bp的序列，然后输出到另外一个序列文件“dna_over5bp.fasta”中。

5这个值是可以自己设置的，以后用1000，2000也适用，这个程序脚本应该怎么编写？还有，编出来了怎么适用。。。

谢谢啦！

Answer 1

不需要你自己重新来写，单枪匹马，刀耕火种，不是一个好的办法。
正好我也在做bioinformatics这块儿，学了点bioperl，我帮你写了个小脚本，用这个工具能很快捷的搞定你的序列！！
我大致给你说下，具体自己上bioperl网站看How To。
现在假定你的序列都为fasta格式，先用
cat *.fasta > single_all_fasta.fasta
将所有fasta序列整合到一个fasta格式中。
然后编写bioperl脚本（你要事先安装bioperl，在debian和ubuntu下很简单，自带源里就有，直接sudo apt-get install bioperl即可，很方便），你参考下：
脚本使用方法为：
perl get_fasta_length.pl single_all_fasta.fasta length.result
脚本我没有测试，一些小错误你可以修改下：

#!/usr/bin/perl -w
use strict;
use Bio::SeqIO;
my $in = shift; #读取文件名称
my $out = shift; #读取输出文件名称
my $flag = 0; #计数
open WH,">> $out";
$seqio_obj= Bio::SeqIO->new (-file =>'$in',
-f =>'fasta');
while($seq_obj = $seqio_obj->next_seq() ) {
print WH "$seq_obj->desc \t $seq_obj->length\n" ;
print "正在处理第$flag条序列...\n";
$flag++;
}
print "完毕，共处理$flag条序列！结果保存至$out文件\n";
close WH;

建议：
学好bioinfo，不会perl不行，不会bioperl更不行，它不仅可以轻松批量处理序列格式，解析blast结果，甚至在脚本里做blast，而且还可以操作biosql数据库。呵呵 建议深入学习下。