5

我可以手动下载一个 FASTA 文件，如下所示：

>lcl|CR543861.1_gene_1...
ATGCTTTGGACA...
>lcl|CR543861.1_gene_2...
GTGCGACTAAAA...

通过单击“发送至”并选择“基因特征”，FASTA Nucleotide 是此页面上唯一的选项（这很好，因为这就是我想要的）。

使用这样的脚本：

#!/usr/bin/env perl
use strict;
use warnings;
use Bio::DB::EUtilities;

my $factory = Bio::DB::EUtilities->new(-eutil   => 'efetch',
                                       -db      => 'nucleotide',
                                       -id      => 'CR543861',
                                       -rettype => 'fasta');
my $file = 'CR543861.fasta';
$factory->get_Response(-file => $file);

我得到一个看起来像这样的文件：

>gi|49529273|emb|CR543861.1| Acinetobacter sp. ADP1 complete genome
GATATTTTATCCACA...

将整个基因组序列集中在一起。如何获取第一个（手动下载的）文件中的信息？

我看了其他几个帖子：

以及EUtilities Cookbook 中的这一部分。

我尝试获取并保存一个 GenBank 文件（因为我得到的 .gb 文件中的每个基因似乎都有单独的序列），但是当我使用 Bio::SeqIO 处理它时，我只会得到一个大序列。

4

1 回答 1

4

使用该登录号和返回类型，您将获得完整的基因组序列。如果您想获取单个基因序列，请指定您需要完整的 genbank 文件，然后解析出基因。这是一个例子：

#!/usr/bin/env perl

use 5.010;
use strict;
use warnings;
use Bio::SeqIO;
use Bio::DB::EUtilities;


my $factory = Bio::DB::EUtilities->new(-eutil   => 'efetch',
                                       -email   => 'foo@bar.com',
                                       -db      => 'nucleotide',
                                       -id      => 'CR543861',
                                       -rettype => 'gb');
my $file = 'CR543861.gb';
$factory->get_Response(-file => $file);

my @gene_features = grep { $_->primary_tag eq 'gene' } 
                    Bio::SeqIO->new(-file => $file)->next_seq->get_SeqFeatures;

for my $feat_object (@gene_features) {
    for my $tag ($feat_object->get_all_tags) {
        # open a filehandle here for writing each to a separate file
        say ">",$feat_object->get_tag_values($tag);
        say $feat_object->spliced_seq->seq;
        # close it!
    } 
}

这会将每个基因写入同一个文件（如果你重定向它，现在它只是写入 STDOUT），但我指出你可以在哪里做一些小改动，将它们写入单独的文件。解析 genbank 有时会有点棘手，因此阅读文档总是有帮助的，尤其是优秀的Feature Annotation HOWTO。

于 2014-02-28T16:58:54.067 回答