BLASTのつかいかた(Web版とローカル版)

バイオロジーにおいて最もよく使用される配列検索ツール(BLAST)の解説をします。

例題として以下の遺伝子配列が何かというのを探索します。

例題

以下の配列が何かを同定せよ。ただし、ヒトゲノム由来の配列であることがわかっているとする。

> seq1
atgagcacagcaggaaaagtaatcaaatgcaaagcagctgtgctatgggagttaaagaaa
cccttttccattgaggaggtggaggttgcacctcctaaggcccatgaagttcgtattaag
atggtggctgtaggaatctgtggcacagatgaccacgtggttagtggtaccatggtgacc
ccacttcctgtgattttaggccatgaggcagccggcatcgtggagagtgttggagaaggg
gtgactacagtcaaaccaggtgataaagtcatcccactcgctattcctcagtgtggaaaa
tgcagaatttgtaaaaacccggagagcaactactgcttgaaaaacgatgtaagcaatcct
caggggaccctgcaggatggcaccagcaggttcacctgcaggaggaagcccatccaccac
ttccttggcatcagcaccttctcacagtacacagtggtggatgaaaatgcagtagccaaa
attgatgcagcctcgcctctagagaaagtctgtctcattggctgtggattttcaactggt
tatgggtctgcagtcaatgttgccaaggtcaccccaggctctacctgtgctgtgtttggc
ctgggaggggtcggcctatctgctattatgggctgtaaagcagctggggcagccagaatc
attgcggtggacatcaacaaggacaaatttgcaaaggccaaagagttgggtgccactgaa
tgcatcaaccctcaagactacaagaaacccatccaggaggtgctaaaggaaatgactgat
ggaggtgtggatttttcatttgaagtcatcggtcggcttgacaccatgatggcttccctg
ttatgttgtcatgaggcatgtggcacaagtgtcatcgtaggggtacctcctgattcccaa
aacctctcaatgaaccctatgctgctactgactggacgtacctggaagggagctattctt
ggtggctttaaaagtaaagaatgtgtcccaaaacttgtggctgattttatggctaagaag
ttttcattggatgcattaataacccatgttttaccttttgaaaaaataaatgaaggattt
gacctgcttcactctgggaaaagtatccgtaccattctgatgttttga

BLASTのつかいかた(Web版)


1. 公式ページへ飛ぶ

公式ページ:https://blast.ncbi.nlm.nih.gov/Blast.cgi


2. クエリとリファレンスの配列によって使うBLASTを選ぶ

BLASTの使い方 - BLASTの選択


3. クエリ配列を入れた上で、オプション選んで検索

  • Database

配列がヒト/マウス由来とわかっている場合はそれを、配列が何由来かわからない場合はOthersを選べばOK。

ただし、大きいDBほど機能未知の配列なども混じっているため注意すること。

オプション説明
Human genomic + transcriptヒトゲノムなど
Mouse genomic + transcriptマウスゲノムなど
Others (nr etc.)nrは全ての配列が入っているため、生物種不明のときはnrが無難。
  • Program Selection

基本的にデフォルトのmegablastでOK。配列が見つからなかったら下のオプションを選ぶ。

オプション説明
Optimize for Highly similar sequences (megablast)類似配列があるとわかっている場合
Optimize for More dissimilar sequences (discontiguous megablast)マウスゲノムなど
Optimize for Somewhat similar sequences (blastn)nrは全ての配列が入っているため、生物種不明のときはnrが無難。



BLASTの使い方 - オプション選択


4. 結果を眺める

結果を眺めます。主に気をつけた方がいいのは、一致率(Per. Ident)とカバー率(Query Cover)です。

カバー率が低い場合、部分的な配列しか一致していません。

それでも問題ない場合もありますが、そこを考慮してください。

BLASTの使い方 - 結果の見方


BLASTのつかいかた(ローカル版)

DBも自分で作った場合や、大量の配列検索をする場合はローカル版のBLAST(+)のインストールが必要になります。


0. ヒトゲノム配列の準備

ローカル版のBLASTの場合、Refrenceとなる配列も自分で準備することがあります。

今回はWeb版と同じく、ヒトゲノム遺伝子配列を準備します。

NCBIにあるヒトゲノムリファレンス配列のGRCh38(Genome Research Consortium human build 38)を使用します。

https://www.ncbi.nlm.nih.gov/genome/guide/human/

wget ftp://ftp.ncbi.nlm.nih.gov/refseq/H_sapiens/annotation/GRCh38_latest/refseq_identifiers/GRCh38_latest_rna.fna.gz
gunzip GRCh38_latest_rna.fna.gz

1. BLASTのインストール

公式ページのいうとおりにインストールしてみます。

BLASTのバージョンなどは変更する必要がありますが、以下のようなコマンドでインストールできます。

# $HOME/src にBLASTのtar.gzファイルなどを、
# $HOME/src/bin にコマンドをインストールします。

# mkdir -p $HOME/src
# mkdir -p $HOME/src/bin
# cd $HOME/src

# wget
wget ftp://ftp.ncbi.nlm.nih.gov/blast/executables/blast+/LATEST/ncbi-blast-2.9.0+-src.tar.gz

# 解凍
tar -zxvf ncbi-blast-2.9.0+-src.tar.gz

# 移動
cd ncbi-blast-2.9.0+-src/c++

# ./configure && make && make install
./configure --prefix=$HOME/src
make
make install

2. DBのインデックス作成

BLASTで検索をかけるまえに、データベース配列の整理(インデックス作成)をする必要があります。

# アミノ酸配列の場合は -dbtype proc にする
makeblastdb -in GRCh38_latest_rna.fna -dbtype nucl -hash_index

GRCh38_latest_rna.fna のあとに3文字くらいついたファイル群ができていればOKです。

私の環境では以下のファイル群が作成されました。

$ ls -1 GRCh38_latest_rna.fna.*
GRCh38_latest_rna.fna.nhd
GRCh38_latest_rna.fna.nhi
GRCh38_latest_rna.fna.nhr
GRCh38_latest_rna.fna.nin
GRCh38_latest_rna.fna.nog
GRCh38_latest_rna.fna.nsd
GRCh38_latest_rna.fna.nsi
GRCh38_latest_rna.fna.nsq

3. 検索

上記の塩基配列をinput.faに入れ、blastしてみます。

# -outfmt 6 : タブ区切りで出力
blastn -db GRCh38_latest_rna.fna -outfmt 6 -query input.fa -out output.tsv

結果は以下のようになりました。

seq1    NM_000667.4 100.000 1128    0   0   1   1128    72  1199    0.0 2084
seq1    NM_000668.6 94.947  1128    57  0   1   1128    71  1198    0.0 1768
seq1    NM_001286650.2  94.789  1113    58  0   16  1128    208 1320    0.0 1735
seq1    NM_000669.5 93.883  1128    69  0   1   1128    72  1199    0.0 1701
seq1    NR_133005.2 94.133  784 46  0   1   784 72  855 0.0 1194
seq1    NR_133005.2 92.691  301 22  0   828 1128    855 1155    6.64e-120   435
seq1    NR_037884.1 100.000 143 0   0   119 261 3985    3843    9.22e-69    265
seq1    NM_001102470.2  79.420  379 74  2   376 752 470 846 9.22e-69    265
seq1    NM_000672.4 79.420  379 74  2   376 752 470 846 9.22e-69    265
seq1    NR_132990.2 83.069  189 31  1   564 752 300 487 2.07e-40    171

ベストヒット配列のIDで検索をかけると、Webの結果と同じくアルコールデヒドロゲナーゼであることがわかります。

$ grep "NM_000667.4" GRCh38_latest_rna.fna
>NM_000667.4 Homo sapiens alcohol dehydrogenase 1A (class I), alpha polypeptide (ADH1A), mRNA

参考

  • 旧BLASTと新BLASTのオプション対応表

http://togotv.dbcls.jp/20110225.html

  • E-valueの計算方法(DDBJ)

https://www.ddbj.nig.ac.jp/faq/ja/how-is-e-value-calculated.html

  • E-valueの計算方法(NCBI・英語)

https://www.ncbi.nlm.nih.gov/BLAST/tutorial/Altschul-1.html

https://blast.ncbi.nlm.nih.gov/Blast.cgi?CMD=Web&PAGE_TYPE=BlastDocs&DOC_TYPE=FAQ

  • E-valueの計算方法(わかりやすい) https://note.mu/tatanzaw/n/n82ed08c96c77