Mac で local BLAST を走らせる

UB3/informatics/bioinformatics/blast_local_mac

このページの最終更新日: 2023/06/10

  1. インストール
  2. BLAST の実行方法
  3. 解析例 1: Local BLAST の結果を fasta として得る
    • seqkit を使う方法
    • Batch Entrez を使う方法
  4. エラーの記録

広告

インストール

Homebrew を使ってインストールできる。たぶん brew install blast だったはず。

2017 年 9 月は、次のように NCBI からダウンロードしてインストールしたので、一応記録しておく。

  1. NBCI のページ から、LATEST にゲストとして FTP で接続。
  2. ウィンドウが現れるので、ncbi-blast-2.6.0+.dmg をダブルクリックして開く。フォルダ内に .pkg という拡張子の installation package があるので、これを開き指示通りにインストールする。
  3. デフォルトの設定は、blast プログラムは usr/local/ncbi/blast というフォルダにインストールされる。よって path を通したり しなくてもよい。
  4. または、.tar.gz で終わるファイルをダウンロードして解凍すると、そのまま blast のフォルダになる。そのフォルダにパスを通すか、コマンドのパスを指定して実行する。

なお、Linux の場合は apt-get でインストールできる。sudo apt-get install ncbi-blast+ でよい。

blast --version でバージョンが見られるかと思ったが、このコマンドは使えない。ただし、エラーメッセージとともにバージョンが表示されたので、結果的に OK だった。

BLAST の実行方法

BLAST のバージョンによって微妙に違うかもしれない。以下は少なくとも 2020 年 4 月時点、BLAST 2.10.0 で有効なコマンドである。


BLAST database の構築

file.fasta といいうファイルから、塩基配列のデータベースを構築する。アミノ酸配列の場合、dbtype は prot である。

makeblastdb -in file.fasta -parse_seqids -dbtype nucl -out databasename

makeblastdb コマンドのオプションをまとめておく。網羅はしていないので、公式マニュアル も参照のこと。

in

データベースを作る fasta ファイルの指定。

out

データベース名の指定。何も指定しない場合は、input する fasta file の名前 test.fasta などがデータベース名になる。

dbtype

塩基 nucl またはタンパク質 prot を指定。

parse_seqids

parse は「解析する」という意味の単語。-parse_seqids は、もともとの fasta file に含まれている sequence id を解析せよというオプションのようで、よくわからないが 常に入れておく必要がある ようだ。

Ref. 1 より: Please note that this option should be provided consistently among the various applications involved in creating BLAST databases. For instance, the filtering applications as well as convert2blastmask should use this option if makeblastdb uses it also.

公式マニュアルにはこのように書かれているが、いまいち意味がわからない。The identifier should begin right after the ">" sign on the definition line and contain no spaces and the -parse_seqids flag should be used.


うまくいくと、6 つのファイルが作られる。塩基配列の場合、ファイル名は file.fasta.nhr, file.fasta.nin, file.fasta.nog, file.fasta.nsd, file.fasta.nsi, and file.fasta.nsq であり、アミノ酸データベースでは、file.fasta.nhr, file.fasta.nin, file.fasta.nsq, file.fasta.phr, file.fasta.pin, file.fasta.psq になるようである。


まず script コマンドを実行しておくのが良いだろう。これで、ターミナルの内容がログとして記録されることになる (参考: コマンドの一覧)。

query.fasta というファイルを、database.fasta という BLAST database に対して検索する以下のコマンドが基本形。e value のオプション、output format、output file の名前を指定している。

blastn -query query.fasta -db database.fasta -evalue 1e-3 -outfmt 2 -out result.out

blast 系コマンドのオプション。基本的には、blastn、blastp など blast の種類を問わず使えるはずである。

query

query を指定。必須オプション。

db

データベースを指定。これも必須。

out

アウトプットファイルを指定。これがないと、Terminal 上に結果が出力される。

evalue

evalue の閾値を指定。

outfmt

-outfmt 2 のようにしてアウトプットフォーマットを指定。以下の解析例 1 で、6 を指定して結果を fasta file として得る方法を解説している。

max_target_seqs

-max_target_seqs 5 などのようにして、ヒット数を指定 (Ref)。BLAST 結果からトップヒットのみを抽出するスクリプトも検索すると見つかるが、こっちの方が簡単だろう。数が 5 未満だと警告が出るが、実行には問題ない。


解析例 1: Local BLAST の結果を fasta として得る

Local BLAST のデフォルトのアウトプットは、次のように blast の種類の宣言と reference から始まっている。貴重な情報ではあろうが、blast の結果を次の解析にそのまま引き渡す場合、これは少々不便である。

ここでは、BLAST の結果を fasta として得る方法を 2 つ解説する。この他にも、awk を使う方法 など、いくつも方法があるようである。

TBLASTN 2.10.0+

Reference: Stephen F. Altschul, 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.


seqkit を使う方法

seqkit がない場合は、homebrew を使ってインストールする。

outfmt で 6 を指定する。6 は便利なオプションで、結果をテーブルとして出力でき、さらに何をテーブルに含めるかを指定できる。

このページ に specifier の一覧がある。今回必要なのは、subject の sequence id なので sseqid である。つまり実行するコマンドは

blastn -query query.fasta -db database.fasta -evalue 1e-3 -outfmt "6 sseqid" -out result.txt

のようになる。

result.txt には、subject sequence id が順に保存されている。これを while read line で seqkit grep に入れていけば良い。

seqkit grep -nrp ${line} original.fasta >> output.fasta

original.fasta は、blast のデータベースを作った fasta ファイルと同じである場合が多いだろう。


Batch Entrez を使う方法

かなり昔の解析例なので色々と至らない部分があるが、一応ここに載せておく。grep/findstr で開くセッション番号を探し、Batch Entrez でダウンロードする方法である。

たとえば Notch という遺伝子について、データベースにある類似の遺伝子をまとめてダウンロードしたいときの手順。

1. Query とデータベースの準備

  1. NCBI から特定のキーワードで遺伝子のリストをダウンロードする。たとえば "Notch" で検索し、数千個の Notch gene の配列を notch.fasta という FASTA 形式のファイルとしてダウンロードしたとする。
  2. 未知の Notch gene の配列を探すため、NCBI に登録されている全ての配列をダウンロードする。このページ の Databases から、ゲノム、トランスクリプトームなど様々なデータをダウンロードできる。既に BLAST database の形式になっているので、makeblastdb をする必要はない。

2. 検索結果から番号を抽出

  1. BLAST を走らせ、結果をファイルに保存する。検索条件にもよるが、結果のファイルは数 GB になり、アプリケーションでは開けない場合も。
  2. デフォルトの output format の場合、> の次に accession 番号があり、その次に alignment が続くという形式になっているはずである。
  3. Windows なら findstr で、Mac なら grep で > を検索し、これを含む行だけを他のファイルに保存する。

>NM_0168731 Mus musculus gene abc, transcript
>NM_0169871 Homo sapiens gene abc, mRNA
>NM_0489412 Mus musculus gene abc, transcript

のように、番号とその他の情報がひたすら並ぶ テキストファイル になっているはずである。

この段階で、さらに特定のキーワードを使って情報を絞ることが可能だと思う。

Accession 番号とその他の情報の間にはスペースがあるので、これを Excel でスペース区切りテキストとして開くと、番号だけが最初のカラムに来る。> を置換で消去すれば、

NM_0168731
NM_0169871
NM_0489412

のように番号のみが並ぶファイルになる。


3. Batch Entrez で一括 download

  1. このページ にテキストファイルをアップロードすれば、配列を一括でダウンロードすることができる。

エラーの記録

Error: Too many positional arguments (1), the offending value.

ハイフンとダッシュの問題。(キーボードによるかもしれないが)、ゼロの隣にある - を使っていることを確認する。多くの場合、どこかからコピペしたときに見られるエラー。手打ちすれば OK になるはず。


広告

References

  1. BLAST+ manual. Pdf file: Last access 2020/04/26.

コメント欄

サーバー移転のため、コメント欄は一時閉鎖中です。サイドバーから「管理人への質問」へどうぞ。

アップデート前、このページには以下のようなコメントを頂いていました。ありがとうございました。

2018/04/24 08:26 Linux firefoxだと、リロードするたびに書き込まれるというエラーになる

2018/04/24 08:24 良いページをみつけたのでメモ  http://d.hatena.ne.jp/aaikmyz/20110320/1300605293