2010-03-02 16:54:25
BLASTをスタンドアローンで動かすことにする。数千のBLAST解析をWWWサーバを使って解析するのはあまり現実的な方法ではないことがよく判ったからだ。前にMacOSXでは動かしていたことがあるのだが、Ubuntuでは初めてである。Ubuntuの場合、なんとncbiblastがリポジトリに入っているのだ。Synaptic（あるいはapt-get install）で、一発でインストールできるのだ。blast2をSynapticで選んでインストールすれば、libncbi6, ncbi-dataとともにインストールされる。/usr/binに、bl2seq、blast2、blastall、blastall_old、blastcl3、blastclust、blastpgp、formatdb、impala、megablast、rpsblast、seedtopなどが入っていることが判る。/etc/ncbiができているが、中身は空っぽである。そこをどう使ったらいいのかよく判らない。
ホームディレクトリに.ncbiというフォルダを作り、そこにdataとblastというフォルダを作る。さらにblastの中にdbを作る。もちろん、こうしなければならないという訳ではなく、たまたま私がそうしただけだ。.ncbiフォルダにncbircというファイルを作ってみた。
[NCBI] Data=/home/ynakano/.ncbi/data [BLAST] BLASTDB=/home/ynakano/.ncbi/blast/db今回は16S rRNA配列だけに限ったBLAST検索をしたいので、専用のデータベースを作るために、http://rdp.cme.msu.edu/misc/resources.jspから、release10_18_unaligned.fa.gzをダウンロード。105MB。展開してrdp10_18.fastaという名前で保存。1.27GB。これから、formatdbコマンドでBLAST用データベースを作るのである。
formatdb -i rdp10_18.fasta -p F -o Tこれで、以下のようなファイルが生成する。
formatdb.log rdp10_18.fasta.nhr rdp10_18.fasta.nsd rdp10_18.fasta.nsq rdp10_18.fasta.nin rdp10_18.fasta.nsi試しに以下のような16S rRNA配列を含むファイルを作り、test.fastaという名前で保存する。
> test1.seq AGAGTTTGATTATGGCTCAGATTG..... > test2.seq AGAGTTTGATCATGGCTCAGGACG..... > test3.seq AGAGTTTGACCTGGCTCAGGACGA.....次に、以下のようなコマンドで解析を実行する。
blastall -p blastn -d rdp10_18.fasta -i test.fasta -o test.out結果がtest.outに保存される。こんな内容である。
BLASTN 2.2.20 [Feb-08-2009]
Reference: Altschul, Stephen F., Thomas L. Madden, Alejandro A. Schaffer,
Jinghui Zhang, Zheng Zhang, Webb Miller, and David J. Lipman (1997),
"Gapped BLAST and PSI-BLAST: a new generation of protein database search
programs", Nucleic Acids Res. 25:3389-3402.
Query= test1.seq length=251 xy=3909_1517 region=4
run=R_2010_02_08_17_57_12_
(251 letters)
Database: rdp10_18.fasta
1,358,427 sequences; 1,177,075,220 total letters
Searching..................................................done
Score E
Sequences producing significant alignments: (bits) Value
S000944605 uncultured Neisseria sp.; 501C06(oral); AM420167 396 e-109
S000921397 uncultured Neisseria sp.; EHFS1_S14d; EU071525 394 e-108
S000613477 Neisseria meningitidis; CMCC(B)29019; DQ201323 387 e-106
.....
>S000944605 uncultured Neisseria sp.; 501C06(oral); AM420167
Length = 1504
Score = 396 bits (200), Expect = e-109
Identities = 242/251 (96%), Gaps = 5/251 (1%)
Strand = Plus / Plus
Query: 1 agagtttgattatggctcagattgaacgctggcggcatgctttacacatgcaagtcggac 60
|||||||||| |||||||||||||||||||||||||||||||||||||||||||||||||
Sbjct: 1 agagtttgatcatggctcagattgaacgctggcggcatgctttacacatgcaagtcggac 60
.....
これで目的のファイルを目的のデータベースで解析できることが判ったが、ちょっと使いづらい解析結果である。WWW解析のときのようにPythonで操作したいのだが、NCBIWWW.qblastのときとどうも違うのである。一日こればかりやっている訳にもいかないので、今日はこれで。もう夕方になってしまったけれど。