バイオインフォマティクス初心者の日常

バイオインフォマティクス関連の研究・勉強などなど

Bio::Perlを利用したAccession NumberからのGenbankファイルのダウンロードとGenbankからタンパク質FASTAへの変換方法

目的の生物種全ゲノムのGenbankファイルとそれに対応する全タンパク質FASTAファイルを手動でダウンロードするのが面倒だったので、Perlで自動化してみた。
Bio::Perlを利用したAccession NumberからのGenbankファイルのダウンロード及びGenbankファイルからの全タンパク質FASTAファイルへの変換のPerl自動化スクリプトのメモです。

Accession NumberからGenbankファイルをダウンロードするプログラム(acc2gbk.pl)

#!/usr/bin/perl

# Command usage:
#   $ perl acc2gbk.pl [$acc] [$gbk_outfile]
#
# Function:
#   Convert accession number to genbank file
# 
# Example usage:
#   $ perl acc2gbk.pl NC_000911.1 NC_000911.1.gbk 

use strict;
use warnings;

use Bio::Perl;
use Bio::DB::GenBank;

my ($acc, $gbk_outfile) = @ARGV;
acc2gbk($acc, $gbk_outfile);

sub acc2gbk{
	my ($acc, $gbk_outfile) = @_;
	
	my $gbk = Bio::DB::GenBank->new;
	my $seq = $gbk->get_Stream_by_acc($acc)->next_seq;
	write_sequence(">$gbk_outfile", 'genbank', $seq);
}

Genbankファイルから全タンパク質FASTAファイルへ変換するプログラム(gbk2faa.pl)

#!/usr/bin/perl

# Command usage:
#   $ perl gbk2faa.pl [$gbk_infile] [$faa_outfile]
# 
# Function:
#   Convert genbank to protein fasta

use strict;
use warnings;

use Bio::SeqIO;

my ($gbk_infile, $faa_outfile) = @ARGV;
gbk2faa($gbk_infile, $faa_outfile);

sub gbk2faa{
	my ($gbk_infile, $faa_outfile) = @_;
	
	my $gbk_all = Bio::SeqIO->new( -file => $gbk_infile, -format => 'genbank' );
	
	my $faa_outinfo;
	while(my $gbk = $gbk_all->next_seq){
		# Read all CDS features in genbank record
		my @cds_features = grep{ $_->primary_tag eq 'CDS' } $gbk->get_SeqFeatures;
		
		foreach my $cds (@cds_features){
			next if( $cds->has_tag('protein_id') ne "1" || $cds->has_tag('translation') ne "1" );
			
			# Get Accession number & Protein sequence
			my $acc = ($cds->get_tag_values('protein_id'))[0];
			my $seq = ($cds->get_tag_values('translation'))[0];
			
			$faa_outinfo .= ">$acc\n$seq\n";
		}
	}
	# Output protein fasta
	open(NEW,">$faa_outfile");
	print NEW $faa_outinfo;
	close(NEW);
}

実行コマンド例

$ perl acc2gbk.pl NC_000911.1 NC_000911.1.gbk
$ perl gbk2faa.pl NC_000911.1.gbk NC_000911.1.faa

上記コマンドを実行することで、Accession NumberがNC_000911.1である'Synechocystis sp. PCC 6803'の全ゲノムGenbankファイル(NC_000911.1.gbk)とそれに対応する全タンパク質FASTAファイル(NC_000911.1.faa)を取得することができる。