Home > バイオインフォマティクスツール・ソフトウェア

バイオインフォマティクスツール・ソフトウェア Archive

カテゴリ： バイオインフォマティクスツール・ソフトウェア の記事一覧

| | 古い記事一覧へ≫

Ads by Google

上記の広告は１ヶ月以上更新のないブログに表示されています。
新しい記事を書く事で広告が消せます。
  • Comments (Close): -
  • TrackBack (Close): -

ClustalWでアライメント

Bio::AlignIO->newFh()を使うことにより、アウトプットのフォーマットを指定することができます。

use Bio::Tools::Run::Alignment::Clustalw;
use Bio::AlignIO;
use strict;
my $infile = 'test2.fas';
my $outfile = $infile.'.out';
my @params = ('ktuple' => 2, 'type' => 'protein', 'matrix' => 'BLOSUM',
'outfile' => $outfile);
my $factory = Bio::Tools::Run::Alignment::Clustalw->new(@params);
my $aln = $factory->align($infile);
my $string = Bio::AlignIO->newFh('-format' => 'clustalw');
print $string $aln;

01:use Bio::Tools::Run::Alignment::Clustalw;
02:use Bio::AlignIO;
03:use strict;
04:my $infile = 'test2.fas';
05:my $outfile = $infile.'.out';
06:my @params = ('ktuple' => 2, 'type' => 'protein', 'matrix' => 'BLOSUM',
07:'outfile' => $outfile);
08:my $factory = Bio::Tools::Run::Alignment::Clustalw->new(@params);
09:my $aln = $factory->align($infile);
10:my $string = Bio::AlignIO->newFh('-format' => 'clustalw');
11:print $string $aln;

Mac OS Xユーザーのためのバイオインフォマティクスソフト

eBiotoolsとBioXを入れることで、GUI上で

o EBIの検索
o Multiple ClustalW
o 2つの配列間のGlobal alignment
o Dotplot
o Searching features in DNA
o Restriction map
o Staden package

などを使うことが出来ます。

http://www.ebioinformatics.org/からeBiotoolsとBioXをそれぞれダウンロードします。

eBiotoolsはインストーラー形式になっています。BioXは展開してアプリケーションフォルダーあたりに入れておけばよいでしょう。

blast、EMBOSS、PHYLIP、clustalw、t-coffee、ViennaRNA、wwwblastなどが入ったみたいです。

ClustalW,EMBOSSでの系統樹の書き方

ClustalW編

treeオプションで以下のように指定する。

% clustalw -infile="test2.aln" -tree -outputtree=phylip

この結果、test2.phというファイルができます。

EMBOSS編

EMBOSSで系統樹を書く方法はいくつかあります。

例えば、fdnadistとffitchを使ってみます。

% fdnadist -sequence=test2.aln -outfile=test2.fdnadist -method=f

% ffitch -datafile=test2.dst -intreefile=test2.ph -outfile=test2.ffitch

F84は距離計算の方法です。他にもkimuraの方法などがあります。

f (F84 distance model)

k (Kimura 2-parameter distance)

j (Jukes-Cantor distance)

l (LogDet distance)

s (Similarity table)

ffitchでのtreeファイルの読み込みはオプションです。

読み込ませない場合は-nointreefileと指定します。

インターラクティブモードでも実行できます。

結果としてtest2.treefileというファイルができます。

BioPerlでFASTAファイルから指定したIDの配列を抜き出す

BioPerlでFASTAファイルから指定したIDの配列を抜き出すことをしたいと思います。
use Bio::DB::Fastaを使って、FASTAファイル中から指定したIDの配列をFASTA形式で出力させます。

use Bio::SeqIO;
use Bio::DB::Fasta;
$seq_out = Bio::SeqIO->new(-fh => \*STDOUT, '-format' => 'fasta');
$file = "sample.fasta";
$id = $ARGV[0];
$db = Bio::DB::Fasta->new($file);
$seqobj = $db->get_Seq_by_id($id);
$seq_out->write_seq($seqobj);

解説
01:use Bio::SeqIO;
02:use Bio::DB::Fasta;
03:$seq_out = Bio::SeqIO->new(-fh => \*STDOUT, '-format' => 'fasta');
04:$file = "sample.fasta";
05:$id = $ARGV[0];
06:$db = Bio::DB::Fasta->new($file);
07:$seqobj = $db->get_Seq_by_id($id);
08:$seq_out->write_seq($seqobj);

01:〜02:おまじない
03:
04:
05:
06:
07:
08:

ClustalWで大量の配列をアライメントする

ClustalWで大量の配列のアライメントをしようと思うと、かなりの時間がかかりますが、FASTA(FastAlignment)を実行すると、早くアライメントすることが出来る。

以下そのコード。
% clustalw -output=PIR -infile=test.fasta -align -quicktree -ktuple=2 -topdiags=4 -window=4 -pairgap=5 -score= -type=DNA -transweight=0.5 -gapopen=15 -gapext=6.66 -maxdiv=30 -dnamatrix=clustalw >test.clwlog&

BioPerlでデータベースID検索

いままではデータベースからIDのわかっているタンパク質の配列などをとってこさせていましたが、
実際に検索する際には、IDなどはわかっていないので、検索することになります。

今回はキーワード検索した結果を一覧表示するプログラムです。

use Bio::DB::GenBank;
use Bio::DB::Query::GenBank;
$query = "Homo sapiens[ORGN] AND CYP2C[TITL] AND 100:3000[SLEN]";
$query_obj = Bio::DB::Query::GenBank->new(-db => 'nucleotide', -query => $query );
$gb_obj = Bio::DB::GenBank->new;
$stream_obj = $gb_obj->get_Stream_by_query($query_obj);
while ($seq_obj = $stream_obj->next_seq) {
$n++;
print $n, ". ",$seq_obj->display_id, " (",$seq_obj->length, ") ",
substr($seq_obj->desc,0,50), "...\n", ;
}

解説
01:use Bio::DB::GenBank;
02:use Bio::DB::Query::GenBank;
03:$query = "Homo sapiens[ORGN] AND CYP2C[TITL] AND 100:3000[SLEN]";
04:$query_obj = Bio::DB::Query::GenBank->new(-db => 'nucleotide', -query => $query );
05:$gb_obj = Bio::DB::GenBank->new;
06:$stream_obj = $gb_obj->get_Stream_by_query($query_obj);
07:while ($seq_obj = $stream_obj->next_seq) {
08:$n++;
09:print $n, ". ",$seq_obj->display_id, " (",$seq_obj->length, ") ",
10substr($seq_obj->desc,0,50), "...\n", ;
11:}

01:おまじない
02:おまじない
03:query="□□□□"で検索ワードを設定する。
　　[ORGN]は生物種
[TITL]はタイトル
を指定しています。
04:
05:
06:
07:
08:
09:
10:
11:


では次に検索の対象をタンパク質に絞って、その一覧をfasta形式のファイルに出力してみます。

use Bio::SeqIO;
use Bio::DB::GenBank;
use Bio::DB::Query::GenBank;
$query = "Homo sapiens[ORGN] AND TNFa[TITL]";
$query_obj = Bio::DB::Query::GenBank->new(-db => 'protein', -query => $query );
$gb_obj = Bio::DB::GenBank->new;
$stream_obj = $gb_obj->get_Stream_by_query($query_obj);
$seqio_obj = Bio::SeqIO->new(-file => '>sequence.fas', -format => 'fasta' );
while ($seq_obj = $stream_obj->next_seq) {
$seqio_obj->write_seq($seq_obj);
}

解説
01:use Bio::SeqIO;
02:use Bio::DB::GenBank;
03:use Bio::DB::Query::GenBank;
04:$query = "Homo sapiens[ORGN] AND TNFa[TITL]";
05:$query_obj = Bio::DB::Query::GenBank->new(-db => 'protein', -query => $query );
06:$gb_obj = Bio::DB::GenBank->new;
07:$stream_obj = $gb_obj->get_Stream_by_query($query_obj);
08:$seqio_obj = Bio::SeqIO->new(-file => '>sequence.fas', -format => 'fasta' );
09:while ($seq_obj = $stream_obj->next_seq) {
10:$seqio_obj->write_seq($seq_obj);
11:}

01:〜03:おまじない
04:
05:
06:
07:
08:
09
10:
11:09の{を閉じた。

BioPerlで塩基配列をアミノ酸配列に翻訳してみよう

この実習では塩基配列をアミノ酸配列に翻訳してみます。
これも今までどおり翻訳しろという命令があるので難しくありません^^

翻訳を行う命令は、
translate()
です。　
これもやはりSeq Objectが返ってくるので配列だけを表示するためには、
seq()
を追加する必要があります。

use Bio::Perl;
$test_dna = get_sequence('genbank', "X12671");
print $test_dna->translate()->seq(),"\n";


実行結果をみると、配列が塩基配列ではA,T,G,Cの４文字だったのがいろいろとアルファベットが出ているので翻訳されているようですね。

ここでアルゴリズムの話になるのですが、アミノ酸配列に翻訳する時に考えなければならないこととしてDNA配列の1文字目から読み始めるか、2文字目、あるいは3文字から読み始めるかを考える必要があります。

これを指定する方法が　translace() に準備されています。

2文字目から開始されるフレームで翻訳を行う場合には、translate(undef,undef,2)　と、
3文字目からなら、translate(undef,undef,3) になります。

このオプションをつかったコードが以下です。実行してみましょう。

use Bio::Perl;
$test_dna = get_sequence('genbank', "X12671");
$short_dna = $test_dna->trunc(10,20);
print "\n", $short_dna->seq(),"\n Translates to";

$translate = $short_dna->translate;
print "\nFrame1: " ,$translate->seq();
$translate = $short_dna->translate(undef,undef,1);
print "\nFrame2: " ,$translate->seq();
$translate = $short_dna->translate(undef,undef,2);
print "\nFrame3: " ,$translate->seq(), "\n";


ちゃんと、翻訳されてるか、コドン表をみて確認するといいでしょう。

BioPerlで相補鎖をつくろう！

相補配列を作る命令は
revcom()
です。ただし相補鎖の文字列ではなく、Seq Objectを返してくるから配列だけを返すように指定しなればいけません。

それでまずは相補鎖のSeq Objectを返すコードを実行しましょう。
use Bio::Perl;
$test_dna = get_sequence('genbank', "X12671");
$rev_dna = $test_dna->revcom();
print "\n", $test_dna->seq(), "\n";
print "reverse\n", $rev_dna->seq(),"\n";


コード解説
01:use Bio::Perl;
02:$test_dna = get_sequence('genbank', "X12671");
03:$rev_dna = $test_dna->revcom();
04:print "\n", $test_dna->seq(), "\n";
05:print "reverse\n", $rev_dna->seq(),"\n";

03:このrevcom() が、相補鎖を計算する命令にになっています。

そして、相補鎖を計算した結果である　$rev_dna は、Seq Object、すなわち、データベースからとってきた配列情報以外も含んでいます。
そこで５行目が大事です。
05:Seq Objectから相補鎖の配列だけを表すため、相補鎖が入っている変数　$rev_dna に　配列だけを出力させる命令　seq() を組み合わせています。

このままでは、長すぎて見えないので、出力を　more につなげて見てみましょう。

相補鎖になっているのは確認できましたか？

more をつかっても、やっぱり長すぎて、困るって人のために、
trunc() という関数を紹介します。

trunc(10,20) という命令は、　subseq(10,20) と同じように、部分配列を取り出す命令です。
しかし、文字列を返すsuqbseq(10,20) と異なり、　seq オブジェクトを返します。

それでは使用法を見てみましょう。
use Bio::Perl;
$test_dna = get_sequence('genbank', "X12671");
$truncated_dna = $test_dna->trunc(10,20);
print "\n", $truncated_dna->seq(), "\n";
print "reverse...\n", $truncated_dna->revcom()->seq(),"\n";


01:use Bio::Perl;
02:$test_dna = get_sequence('genbank', "X12671");
03:$truncated_dna = $test_dna->trunc(10,20);
04:print "\n", $truncated_dna->seq(), "\n";
05:print "reverse...\n", $truncated_dna->revcom()->seq(),"\n";

03:10番目から20番目の文字列だけを取り出してきています。

04:部分文字列を、　5行目で部分文字列の相補鎖を表示しています。

以上でBioPerlを使った相補鎖の作り方を終了します。

BioPerlで配列を表示してみよう

いままでの実習ではとってきた配列をみる方法として拡張子を.txtに変更してテキストエディタでみるという流れでした。ただいちいち拡張子を変えていたら面倒なのでプログラミングをして表示させましょう。

まず今までの実習が完璧な方は

use Bio::Perl;
$test_dna = get_sequence('genbank', "X12671");
print $test_dna;
というコードを使えばいいじゃん！　っておもうかもしれません。
実際私も思っていました。
ってことはだめなのか！やってられねーPerlといわずにまずtest5.plなどとして実行してみてください。

またそのファイルの場所は
XPの場合 C:\Documents and Settings\ユーザー名に置きましょう。
アクセス手順はマイコンピューター→Cドライブ→Documents and Settings→ユーザー名とクリックしていけばたどり着きます。
Vistaの場合　C:\Users\ユーザー名に置きましょう。
アクセス手順はコンピューター→ローカルディスク(C)→ユーザー→ユーザー名でたどり着きます。

それで実行した結果は


あれ？配列がでないぞ　となりますが、これはプログラムのエラーではありません。
ただ配列がでてないのでエラーと思っていただいてかまいません。

このようにわけわからないことが起きた理由を考えてみましょう。
このプログラムの命令を日本語で言うならば、
BioPerlを使ってGenbankからIDがX12671のDNA配列をとってきてそれを表示させろ
です。
やりたいこととコードがあってるじゃないか　と思いますよね。
ただ普通にインターネットブラウザーでGenbankにアクセスしX12671を検索した際の結果をみればすぐわかっていただけると思います。

X12671での検索結果

配列情報だけでなく、ほかの関連情報が表示されています。
$test_dna には、こういった配列以外も格納したので、print文で$test_dna を表示しろといっても変になっちゃうのです。

で、そこで配列だけを取り出す命令を紹介します。
下のコードをtest6.plなどと保存して実行してみてください。
use Bio::Perl;
$test_dna = get_sequence('genbank', "X12671");
$seq = $test_dna->seq();
print $seq;

いきなり画面いっぱいに、塩基配列が表示されましたか。

表示が確認できましたらコードの説明に入りたいと思います。

01:use Bio::Perl;
02:$test_dna = get_sequence('genbank', "X12671");
03:$seq = $test_dna->seq();
04:print $seq;

01:02:今まで通り。

03:ここがきもです。先ほど確認したように$test_dna というのは、Genbankからとってきた配列情報以外も含む情報です。この情報をSeq Object(シーク　オブジェクト)といいます。
そして、Seq Objectから配列だけを取ってこいという命令が、seq() です。
つなげ方は-> をつかいます。前回の実習で分子量を計算した際の　-> と同じ意味です。
この場合、　$test_dna　（に格納されているデータの）の　seq() (配列)　という意味で　-> が使われています。

この、$test_dnaのあとにつけるシリーズを紹介します。
・seq() : 配列を文字列としてとってこいという命令

・subseq(A,B) : (A,Bは整数)配列のA文字目からB文字目をとってこいという命令。

・display_id : とってきた配列の名前の文字列を返します。

・accession_number() : 　アクセッション　番号の文字列を返します。

・alphabet() : DNA,RNA,Protein といった、何の配列であるかの情報を示す文字列を返します。

他にも、　annotation() とか　speces() とかもあるのですが、そのへんはまたあとでとします。

それでは、このシリーズをつかった命令を使ってみましょう。
use Bio::Perl;
$test_dna = get_sequence('genbank', "X12671");
$subseq = $test_dna->subseq(20,30);
print $subseq,"\n";

$id = $test_dna->display_id;
print $id, "\n";

$acnum = $test_dna->accession_number();
print $acnum,"\n";

$alpha = $test_dna->alphabet();
print $alpha,"\n";

1行目に配列が２０文字目から30文字目表示されて
2行目にIDが表示されて
3行目にアクセッション　番号の文字列が表示されて
4行目に何の配列かが表示
されれば成功です。
ちなみに"\n"というのは改行しろという命令です。
改行で情報が読み取りやすくなるので便利な命令です＾＾

BioPerlで塩基配列を扱ってみよう

この実習では塩基配列をアミノ酸配列への翻訳を行います。
いままでずっと、アミノ酸配列を扱ってきたので、気分転換に塩基配列を扱ってみたいと思います。
基本的には変わりません。というか流れは全く変わりません。変わったのは扱う配列の種類です。

まずは塩基配列のデータをgenbank からとってきましょう。
コードは以下。テキストエディタなどにコピペしてファイル名はtest4.plとしましょう。

use Bio::Perl;
use Bio::Tools::SeqStats;
$test_dna = get_sequence('genbank', "X12671");
$weight = Bio::Tools::SeqStats->get_mol_wt($test_dna);
print "\nMolecular weight of the sequence is between ",
$$weight[0], " and " , $$weight[1], "\n";


またそのファイルは
XPの場合 C:\Documents and Settings\ユーザー名に置きましょう。
アクセス手順はマイコンピューター→Cドライブ→Documents and Settings→ユーザー名とクリックしていけばたどり着きます。
Vistaの場合　C:\Users\ユーザー名に置きましょう。
アクセス手順はコンピューター→ローカルディスク(C)→ユーザー→ユーザー名でたどり着きます。

それでですが、配列をとってくるだけでは実習１と変わらないので、前回BioPerlで分子量を計算してみようで行った分子量の計算方法を復習がてらしてみましょう。

解説
01:use Bio::Perl;
02:use Bio::Tools::SeqStats;
03:$test_dna = get_sequence('genbank', "X12671");
04:$weight = Bio::Tools::SeqStats->get_mol_wt($test_dna);
05:print "\nMolecular weight of the sequence is between ",
06:$$weight[0], " and " , $$weight[1], "\n";

かなり前回のコードと似通っているので違いだけ説明します。

詳しい解説を必要とする方はBioPerlで分子量を計算してみようをごらんになってください。

前回との違い
一つ目 03が、swissprot からgenbank に変わった。（IDも）
二つ目03と04の変数が$test_seqから$test_dnaに変わった。

まぁ一つ目はただ単に配列を取りに行くデータベースがかわっただけで、二つ目も変数の名前がかわっただけです。本当に前回と同じことをやっています。

BioPerlで分子量を計算してみよう

この実習では、タンパク質の分子量の計算を行います。
分子量を計算するためには、use Bio::Perlに加えて、use Bio::Tools::SeqStats;の記述をします。

use Bio::Tools::SeqStatsの詳細が知りたい方は
perldoc Bio::Tools::SeqStats
でマニュアルを見てください。（なんとも自分勝手な筆者を許してくださいｗ）
ではプログラムコードを示します。
テキストエディタを開いて以下のコードをコピペしてください。
use Bio::Perl;
use Bio::Tools::SeqStats;
$test_seq = get_sequence('swissprot', "roa1_human");
$weight = Bio::Tools::SeqStats->get_mol_wt($test_seq);
print "\nMolecular weight of the sequence is between ",
$$weight[0], " and " , $$weight[1], "\n";

ファイル名はtest3.plとして保存します。

またそのファイルは
XPの場合 C:\Documents and Settings\ユーザー名に置きましょう。
アクセス手順はマイコンピューター→Cドライブ→Documents and Settings→ユーザー名とクリックしていけばたどり着きます。
Vistaの場合　C:\Users\ユーザー名に置きましょう。
アクセス手順はコンピューター→ローカルディスク(C)→ユーザー→ユーザー名でたどり着きます。

それでは、test3.pl の解説を行っていきます。

01:use Bio::Perl;
01:use Bio::Tools::SeqStats;
03:$test_seq = get_sequence('swissprot', "roa1_human");
04:$weight = Bio::Tools::SeqStats->get_mol_wt($test_seq);
05:print "\nMolecular weight of the sequence is between ",
06:$$weight[0], " and " , $$weight[1], "\n";

01:C言語で言う#inculde
BioPerlを使ったプログラムを書きますよ〜というおまじないだと思ってください。
BioPerlを使うときはいつでも書くと覚えておいてください。

02:BioToolをつかいますよー　というおまじないです。
BioToolを使うときは必ず必要な記述です。

03:02:Perl では、$○○　と、＄で始まる単語は変数といいます。
変数は値(データ)を入れる（記憶する）入れ物です。

すなわち、$test_seq というのは変数、すなはち入れ物を用意したわけです。

その入れ物に何かを入れろ！という命令 ＝以下に来ています。

get_sequence(□□□□)　は、名前の通り、sequence をとってこい(getしてこい）という命令です。

□□□□の部分には、sequenceを取ってくるデータベース名と、　そのデータベースのID（アクセッションナンバー）を記します。
今回の場合、　Genbankというデータベースの、　roa1_human というタンパク質のデータを取ってこいと命令しています。

つまりまとめると02:では、Genbankからroa1_human というタンパク質のデータを取ってきて、$test_seq という変数（入れ物）の中に入れろ、という命令をしています。

04:この行で、分子量の計算を行っています。
命令は、Bio::Tools::SeqStats->get_mol_wt( ) という長めの命令になっています。

この命令の最初に注目してもらいたいのですが、02:と同じ言葉すなわち、Bio::Tools::SeqStats
が出てきていますね。

つまり04はBio::Tools::SeqStatsに含まれるget_mol_wt( )っていう命令を実行しなさいと書いてあります。日本語で言う"に含まれる"をPerl では、"->"であらわします。
私はなんとなく左の方が大きい感じがするから">"と覚えました。
正式な意味はわかりません><

そうしてget_mol_wt() の引数（括弧の中身）である$test_seq が、　何の分子量を計算するかを表しています。
05.06:長い文ですが、2行で一つの表示命令になっています。
C言語でいうprintfです。Perlではprintで表します。
つまりこの部分は計算結果を画面に表示しろという命令だけです。
長い割に内容が薄かったですねｗ

行をまたがるのはperlでは命令の区切りは、セミコロン（;）なので問題ありません。
セミコロンが出てくるまでが１つの命令なので注意してください。

ちなみに出力結果をみると、　between 38715 and 38715 と同じ数字がならんでいます。

それは分子量が曖昧なアミノ酸配列があるからなのです。
よくわからないと思いますので実際にアミノ酸コード表をみてみましょう。
アミノ酸コード表
上から順にみていくと・・・・

最後の方にB(Asx)やZ (Glx）が見つかると思います。これは２つのアミノ酸のどちらか一方をあらわしています。それで分子量が□から○といったように曖昧となるのです。

BioPerlのマニュアル、チュートリアルをみよう

BioPerlのマニュアルをみて勉強を進めましょう。
コマンドプロンプトを起動したら
perldoc Bio::Perl
を入力して実行してみましょう。
このperldocという命令はその後に、マニュアルを見たいモジュール名を書くことで、そのモジュールのマニュアルを見ることができます。

実行すると、次のような画面がでてきます。　

（スペースを押すと次の画面にスクロールします。）
残念ながら英語です。
ただ易しめな英語でかかれてるから、読みやすいと思います。　

より詳しい使い方は他のマニュアルも参照のこと。
"for more info see Bio::Seq" とかって書いてあったら、　
perldoc Bio::Seq
で見れます。

このマニュアル以外に、チュートリアルも存在します。

チュートリアルは、　
コマンドプロンプトを起動して
perldoc bptutorial.pl
を実行すれば見れます。
この記事で好奇心を刺激された方はマニュアル、チュートリアルで自学を進めていくことをおすすめします＾＾

BioPerlでBLASTを使ってみよう

今回の実習では、データベースから取ってきたタンパク質にBLASTをかけて、その結果を　roa1.blast というファイルに出力するということを行いたいと思います。

まずはテキストエディタを開いて
use Bio::Perl;
$test_seq = get_sequence('swissprot', "roa1_human");
$blast_result = blast_sequence($test_seq);
write_blast(">roa1.blast",$blast_result);
をコピペしてファイル名はtest2.plで保存してください。

またそのファイルは
XPの場合 C:\Documents and Settings\ユーザー名に置きましょう。
アクセス手順はマイコンピューター→Cドライブ→Documents and Settings→ユーザー名とクリックしていけばたどり着きます。
Vistaの場合　C:\Users\ユーザー名に置きましょう。
アクセス手順はコンピューター→ローカルディスク(C)→ユーザー→ユーザー名でたどり着きます。

それでコマンドプロンプトを起動したら
perl test2.pl
と入力しエンターを押しましょう。
そうするとNCBI(BLASTを提供する団体)のサーバーが相同性検索した結果がroa1.blastという名前で出力されます。この作業にはNCBIのサーバーの混み具合によってかかる時間が変わりますので、タイムアウトした場合には再実行してください。

そうするとtest2.plをおいたフォルダにroa1.blastというファイルができます。
中身をみるには拡張子の.blastを.txtに変更してテキストエディタなどでみてください。
データの見方についてはBLASTの解説で行いたいと思います。
少々お待ちください＾＾

プログラムコードの解説
01:use Bio::Perl;
02:$test_seq = get_sequence('swissprot', "roa1_human");
03:$blast_result = blast_sequence($test_seq);
04:write_blast(">roa1.blast",$blast_result);

実習１データベースから配列をとってこようから続いてきた方は01,02のコードについては解説済みです。03からごらんになってください。
01:C言語で言う#inculde
BioPerlを使ったプログラムを書きますよ〜というおまじないだと思ってください。
BioPerlを使うときはいつでも書くと覚えておいてください。

02:Perl では、$○○　と、＄で始まる単語は変数といいます。
変数は値(データ)を入れる（記憶する）入れ物です。

すなわち、$test_seq というのは変数、すなはち入れ物を用意したわけです。

その入れ物に何かを入れろ！という命令 ＝以下に来ています。

get_sequence(□□□□)　は、名前の通り、sequence をとってこい(getしてこい）という命令です。

□□□□の部分には、sequenceを取ってくるデータベース名と、　そのデータベースのID（アクセッションナンバー）を記します。
今回の場合、　Genbankというデータベースの、　roa1_human というタンパク質のデータを取ってこいと命令しています。

つまりまとめると02:では、Genbankからroa1_human というタンパク質のデータを取ってきて、$test_seq という変数（入れ物）の中に入れろ、という命令をしています。

03:blast_sequence(□□□□) で、引数□□□□（括弧の中の文字列や引数のこと、この場合、$test_seq ）のタンパク質に対してBLASTを実行しろと命令しています。

つまり03ではデータベースから取ってきたタンパク質　$test_seq に対してBLASTを実行しろと命令しています。

04:BLASTの結果をファイルに出力するのには、write_blastという命令を使います。

write_blast の引数には、
１．出力するファイル名、
２．BLASTの実行結果
を指定します。

これで解説は終了です。

ちなみにこのBLASTをかける方法は前回ローカルでBLAST解析 windowsXP編によって自分のパソコン内だけでBLASTを実行する方法をあげました。
今回のBioPerlでBLASTはPerlでネットワークを通じたBLAST解析を行う方法でした。　

またまたちなみにblast_sequence() は、　
$blast_result = blast_sequence('MFVEGGTFASEDDDSASAEDE');

のように、直接アミノ酸配列を渡すことも可能です。
しかし、DNA配列に対しては、実行することができません。

データベースから配列をとってこよう

まずは、Genbankのデータベースからタンパク質のデータを取ってくるという実習をしましょう。
それではテキストエディタを開いて下記をコピペ（コピーして貼り付け）してください。
use Bio::Perl;
$test_seq = get_sequence('genbank', "roa1_human");
write_sequence(">roa1_human.fasta",'fasta',$test_seq);

そしたらファイル名はtest1.plとしましょう。
またそのファイルは
XPの場合 C:\Documents and Settings\ユーザー名に置きましょう。
アクセス手順はマイコンピューター→Cドライブ→Documents and Settings→ユーザー名とクリックしていけばたどり着きます。
Vistaの場合　C:\Users\ユーザー名に置きましょう。
アクセス手順はコンピューター→ローカルディスク(C)→ユーザー→ユーザー名でたどり着きます。

そうしたらコマンドプロンプトを起動して、
perl test1.pl
を入力してエンターを押してください。
少しするとtest1.plのファイルがある場所にroa_human.fastaというファイルが作成されます。
そうしたら成功です。ちなみにファイルの中身はファイルの拡張子を.fastaから.txtに変更してテキストエディタでみることができます。
きっと以下のようなタンパク質のデータが並んでいると思います。
確認できたら以下の解説を読んで理解を深めましょう。
>P09651 RecName: Full=Heterogeneous nuclear ribonucleoprotein A1; Short=hnRNP core protein A1; AltName: Full=Helix-destabilizing protein; AltName: Full=Single-strand RNA-binding protein.
MSKSESPKEPEQLRKLFIGGLSFETTDESLRSHFEQWGTLTDCVVMRDPNTKRSRGFGFV
TYATVEEVDAAMNARPHKVDGRVVEPKRAVSREDSQRPGAHLTVKKIFVGGIKEDTEEHH
LRDYFEQYGKIEVIEIMTDRGSGKKRGFAFVTFDDHDSVDKIVIQKYHTVNGHNCEVRKA
LSKQEMASASSSQRGRSGSGNFGGGRGGGFGGNDNFGRGGNFSGRGGFGGSRGGGGYGGS
GDGYNGFGNDGGYGGGGPGYSGGSRGYGSGGQGYGNQGSGYGGSGSYDSYNNGGGRGFGG
GSGSNFGGGGSYNDFGNYNNQSSNFGPMKGGNFGGRSSGPYGGGGQYFAKPRNQGGYGGS
SSSSSYGSGRRF

ではtest1.plの解説をしてきます。

01:use Bio::Perl;
02:$test_seq = get_sequence('genbank', "roa1_human");
03:write_sequence(">roa1_human.fasta",'fasta',$test_seq);

01:C言語で言う#inculde
BioPerlを使ったプログラムを書きますよ〜というおまじないだと思ってください。
BioPerlを使うときはいつでも書くと覚えておいてください。

02:Perl では、$○○　と、＄で始まる単語は変数といいます。
変数は値(データ)を入れる（記憶する）入れ物です。

すなわち、$test_seq というのは変数、すなはち入れ物を用意したわけです。

その入れ物に何かを入れろ！という命令 ＝以下に来ています。

get_sequence(□□□□)　は、名前の通り、sequence をとってこい(getしてこい）という命令です。

□□□□の部分には、sequenceを取ってくるデータベース名と、　そのデータベースのID（アクセッションナンバー）を記します。
今回の場合、　Genbankというデータベースの、　roa1_human というタンパク質のデータを取ってこいと命令しています。

つまりまとめると02:では、Genbankからroa1_human というタンパク質のデータを取ってきて、$test_seq という変数（入れ物）の中に入れろ、という命令をしています。

03:write_sequece（□□□□） という命令は、　データをファイルに出力しろという命令です。

□□□□には
どこに？どうなふうに？なにを？出力するかを書いていきます。
今回は
１．　">roa1_human.fasta" の部分で、出力するファイル名を指定しています。

２．　'fasta' で、fasta 形式で出力することを決めています。

３．　$test_seq で、出力したいデータを決めています。

つまり03:をまとめるとGenbankからとってきたroa1_humanのデータをroa1_humanという名前でfasta形式でファイルに出力しろという命令をしたことになります。

ついでに指定できるデータベースの一例と出力形式の一例を示しておきます。
参考になったら私自身うれしいです＾＾
使えるデータベースの一例

"swissprot", "genbank", "embl"
$seq_object = get_sequence('swissprot',"ROA1_HUMAN");
$seq_object = get_sequence('embl',"AI129902");
$seq_object = get_sequence('genbank',"AI129902");

使える出力形式の一例

fasta, genbank, swissprot
write_sequence(">roa1_human.fasta",'fasta',$test_seq);
write_sequence(">roa1_human.gb",'genbank',$test_seq);
write_sequence(">roa1_human.sp",'EMBL',$test_seq);
＊データベースの指定にシングルクオート（’’）をつかっていてIDにはダブルクオート（””）を使っています。その違いは、変数を展開するか、展開しないかが違います。ダブルクオートの場合が展開されるほうです。　
test1.plの例だと、変数を使っていないので、シングルクオートでもダブルクオートでもOKです。　
わかりにくかった場合は、引数に文字列そのものを渡す場合はシングルクオートを、変数を渡す時にはダブルクオート又はクオートなしにすればいいと思ってください。

windowsXPにBioPerlをインストールしよう

windowsでBioPerlをインストールするには２段階ステップを踏む必要があります。
ちなみにこの方法によるBioPerlインストールはWindows XP SP3,Vistaで成功できました。

１） ActivePerl のインストール
２） BioPerl のインストール
です。

１） ActivePerl のインストール
activeperl1.gif
Nextをクリック
activeperl2.gif
I accept the terms in the License AgreementにチェックをいれてNextをクリック
activeperl3.gif
ここもそのままNextをクリックでかまいません。
ちなみに設定内容はインストールする機能を取捨選択できます。
デフォルトではすべての機能をインストールしているのでそれでかまわないと思います。
activeperl4.gif
この設定もいじらずこのままNextでかまいません。
説明すると、1個目のオプションでは、Perlコマンドをどこからでも実行できるように環境変数に加えています。これからのコマンドプロンプトを使った作業では必須です。
2個目のオプションは.pl拡張子のファイルをActivePerlに関連づけています。
これも便利ですのでいじらずに大丈夫です。
activeperl5.gif
Installをクリックしていよいよインストール開始。
activeperl6.gif
少し時間がかかります。
activeperl7.gif
ついにAcitivePerlのインストールを終えました。
このあと一回再起動しましょう。
無事インストールできたかためすには、再起動後コマンドプロンプトを起動して
perl -v
と入力してエンターを押してみましょう。(perlスペース-v)
インストールされたPerlのバージョンが出てくれば成功です。

*コマンドプロンプトの起動手順は
「スタート」→「すべてのプログラム」→「アクセサリ」→「コマンドプロンプト」
これからPerlを使う際には毎回起動するプログラムです。

インストールが確認できたら次のステップであるBioPerlインストールに向かいましょう。

２） BioPerl のインストール
BioPerlってのは、Ｐｅｒｌっていうプログラム言語を通して使う便利ツールみたいなものだと思ってください。
このperl が使えるパソコンにBioPerlｌをインストールすることによって使えるようになります。
ではインストールに入りましょう
まずは、コマンドプロンプトを起動してください。

起動が終ったら、以下のコマンドを打ち込んでください。

ppm

エンターを押してしばらくするとPerl Package Manegerが起動します。

perl -MCPAN -e "install Bundle::BioPerl"
途中、いろいろと聞かれますが、基本的に全部　Enter　キーを押していればOKです。

このコマンドが終了するまで、１０分くらいかかります。
また、止まってるかと不安になる場面も何度かありますが大丈夫です。
辛抱強く待ちましょう。

BioPerl がインストールされていることの確認。
さて、ここでBioPerl がちゃーんとインストールできたかどうかの確認をしたいと思います。
コマンドプロンプトから次のように入力してください。
perl
use Bio::Perl
^D
^D は　Ctrlキーを押しながら ｄ を入力してください。　^D の後に　enter キーを入力するのを忘れずに。
インストールに成功していると何も問題なく、終了します。
インストールに失敗していると、　Can't locate Bio/Perl.pm in ... と表示されます。

では　次にBioPerlを使った実習をしましょう！

KNOB-Windows でバイオインフォマティクスを簡単に始める方法

Windows でバイオインフォマティクスを始めるにはまず環境を構築する必要があります。
Perl や BLAST等のソフト、プログラムのインストールのことです。

てっとり早くバイオインフォマティクスを始めてみたい人には KNOB ( Knoppix for Bio ) がお勧めです。KNOB はノブと読みます。

KNOB

この KNOB を使えばバイオインフォマティクスに必要な環境を手に入れることができます。
Knoppix for Bio は、バイオインフォマティクスの解析や学習用にカスタマイズされた 1 CD bootable Linux です。ベースには、Knoppix-jp を用いています。

KNOB のメリット
・インストールされている Windows には一切影響がない。
　　KNOB はハードディスクに一切の影響を与えません。ハードディスクの中身は閲覧できますが、書き込むことはできません。KNOBはメモリ上につくられた仮想的なディスクにデータを読み書きします。KNOB 起動時に、Windows を同時に起動することはできません。

・データベースは、online のものを利用するので、CD には持たなくてもよい。
　　EMBOSS からは、emboss.defalut を利用。BioRuby, BioPerl は、BioRegistry? (/etc/bioinformatics/seqdatabase.ini)を利用することで、GenBank? などの配列を自由に使えます。また、Windowsのパーティションにデータベースを入れておけば、それを利用することも可能です。

・一枚でみんな幸せ
　　DNA配列解析環境を購入しても、ライセンスの問題で、一台のマシンに入れることしかできないので、他の誰かが使っていたら使えない。その点、KNOPPIX for Bio は、SSH でリモートログインできるので、スペックの許す限り、何人ででも利用できます。VNC Server も同梱されています。

・コピペ系Bioinformaticsで十分な人も幸せ
SSH, VNC を利用しなくても、EMBOSS + EMBOSS::GUI を利用して、自分の windows マシンの In　　ternet Explorer から、コピペ系配列解析ができる!


KNOBにインストールされているソフト一覧
GNUPLOT 4.0.0-2.1
Ruby 1.8.4
mysql-ruby(MySQL)
libdbi-ruby(DBI)
libdb4.3(BDB)
libreadline-ruby
libdbm-ruby
libcurses-ruby
libbigdecimal-ruby
irb
librexml-ruby
libopenssl-ruby
libiconv-ruby
liberb-ruby
libdrb-ruby
libdl-ruby
libdbm-ruby
soap4r
rdoc
libzlib-ruby
libyaml-ruby
libxmlrpc-ruby
libtest-unit-ruby
libstrscan-ruby
ri
Perl v5.8.7
GD
SOAP::Lite
Bundle::XML
MedBio 0.10.1
Primer3 1.1.8-1
EMBOSS 3.0.0
BioRuby 1.0
Bioperl 1.4-1
G-language 1.0.1 BETA
R 4.0.0-2.1
BioConductor 1.7
QTLCart 1.01-9-1
FASTA2 v2.1u00
FASTA3 version 3.4t25
Graphviz 2.2.1-1sarge1
SIM4 2003.09.21-3
meme-mast 3.0.4-3
meta-meme 3.2-2
Sun JRE+JDK(一部) 1.5.0_06
Taverna 1.4

※本来ならばソフトの役割（塩基配列の相同性を検索するソフトなど）を説明すべきかもしれませんが、今回はKNOBがどれだけバイオインフォマティクス系のソフトを収録しているかを強調する形にしました。

化学構造描画ソフト

化学構造を表示出来るソフトウェア

SketchEl
SketchEl
Javaで起動できるフリーウェアで、どのOSでも使えます。手軽で分かりやすいインターフェースなので、簡単な用途には便利です。

JChemPaint (CDK project)
JChemPaintダウンロードサイト
化学構造（２次元）を描画するフリーソフトです。これもJavaアプリケーションです。

JChemPaint
BKChem
BKChem
　※Pythonライブラリが必要。

Chemul [けむり]
http://hossy.net/chemul/
日本の方が開発された、化学構造を閲覧できるプログラムで、Flashで閲覧できるところが特徴的。上記サイトで、アミノ酸などの多くの有機分子を閲覧できます。

タンパク質立体構造表示ツール

・Jmol（英語サイト）
化学物質やタンパク質の三次元構造ビューワーです。

・CueMol (日本語サイト)
日本発のタンパク質構造モデリングソフト。フリーウェア（Windows）。ユーザーインターフェースがわかりやすく、商用ソフトにひけをとりません。

・Chime（英語サイト）

・Rasmol（英語サイト）

・Visual Molecular Dynamics

プラスミド描画ツール

PlasMapper
　PlasMapper

PlasMapperは、カナダのアルバータ大学で開発されたWebベースのプログラムで、無料で使用できます。プラスミドDNA配列を入力するだけで、自動的にプラスミドマップを作成します。20kbpまでのプラスミドマップを描画できます。また、描画した図はPNG,JPG,SVG,SVGZ形式で出力することができます。

このツールを使うと、本当に簡単に作図できます。

Reference → PubMed:15415271

その他
・プラスミド描画ソフトβ（日本語サイト）
ブラウザ上でプラスミドマップを作成できる。SVGによるベクトル描画により、曲線が滑らか。

・pDRAW32 DNA analysis software（英語サイト）
プラスミドマップを描画できるフリーウェア。

RNAやDNAの構造

mfold
mfoldはRNAやDNAの二次構造を予測するプログラムです。

・Michael Zuker's Home Page
・RNA & DNA Folding Applications

SimTK RNA Folding
SimTK RNA Folding
SimTKのプロジェクトの中の、ToRNADoはPDBコードやファイルからRNA二次構造を表示・作画できるフリーソフトです。

・National NIH Center for Physics-based Simulation of Biological Structures

Home > バイオインフォマティクスツール・ソフトウェア

スポンサード リンク
最新記事
カテゴリ
月別アーカイブ
最新記事はRSSで

 RSSリーダーで購読

Google Readerへ追加

My Yahooへ追加

gooへ追加

はてなRSSへ追加

Livedoorへ追加

コメント
最新トラックバック
リンク
ブログで苗木 グリムス
blog parts

track feed

Powered by SEO Stats

ページランク

個人サイト検索 あったかNAVI

この日記のはてなブックマーク数

フィードメーター - バイオ　やばいぉ・・・

tracker あわせて読みたいブログパーツ
スポンサードリンク
火災保険
www.sure-support.com
中古トラック
www.s-auto.jp
バインダー
www.sty.co.jp
看板
www.omoshiro.co.jp
パソコン修理
pc119.unit.ne.jp

Page Top