我使用的是Ubuntu 10.04和Perl 5.10.1。BioPerl包有一些很好的脚本,比如将文件从genbank格式转换为GFF3格式的bp_genbank2gff3.pl。
问题:在使用bp_genbank2gff3.pl时,我得到了意想不到的结果:基因特征在最后一个GFF3列中显示为"Name=“而不是"locus_tag=”。
一位亲爱的BioPerl邮件列表成员告诉我,他使用来自BioPerl存储库的最新BioPerl版本,并获得正确的结果("locus_tag=")。我得到了一份新的副本,但它对我不起作用。太奇怪了!
重新创建情况的步骤:
$ cd ~/src
$ git clone http://github.com/bioperl/bioperl-live.git
$ export PERL5LIB="$HOME/src/bioperl-live:$PERL5LIB"
$ cd /tmp
$ wget ftp://ftp.ncbi.nih.gov/genomes/Bacteria/Escherichia_coli_E24377A/NC_009789.gbk
$ ~/src/bioperl-live/scripts/Bio-DB-GFF/genbank2gff3.PLS NC_009789.gbk下面是我的结果GFF3中的第8行:
NC_009789 GenBank gene 665 781 . - 1 ID=EcE24377A_B0001;Dbxref=GeneID:5585816;Name=EcE24377A_B0001虽然这与我同事的结果是相同的:
NC_009789 GenBank gene 665 781 . - 1 ID=EcE24377A_B0001;Dbxref=GeneID:5585816;**locus_tag**=EcE24377A_B0001注意:在我的版本中(在行尾),"Name=“标签被我同事的"locus_tag=”所取代。我不知道这里发生了什么……相同的输入,可能是相同的脚本,但输出不同(我的同事得到的输出是理想的)。我们甚至对脚本(genbank2gff3.PLS)进行了diff编辑,它们是完全相同的。
有什么想法吗?有没有人能看到他得到的结果是否和我或我的同事一样?
发布于 2010-09-19 00:00:17
查看script source
#?? should gene_name from /locus_tag,/gene,/product,/transposon=xxx
# be converted to or added as Name=xxx (if not ID= or as well)
## problematic: convert_to_name ($feature); # drops /locus_tag,/gene, tags
convert_to_name($feature); 在convert_to_name中
elsif ($g->has_tag('locus_tag')) {
($gene_id) = $g->get_tag_values('locus_tag');
$g->remove_tag('locus_tag');
$g->add_tag_value('Name', $gene_id);
}所以看起来脚本正在做它应该做的事情?
https://stackoverflow.com/questions/3742296
复制相似问题