前回の講義と同じように、JupyterNotebookのページを開き、Terminalを起動してください。
ここをクリックしてファイルをダウンロードしてください。 ダウンロードしたファイルを前回行ったようにJupyterにアップロードしてください。
アップロードしたあと、Terminalでlsコマンドを入力し、blast.tar.gzというファイルがあることを確認してください。
$ ls blast.tar.gz reference.fa
blast.tar.gzは複数のファイルをまとめ(書庫とかアーカイブという)、圧縮したものである。 このままでは使用できないので、解凍し、個別のファイルに展開する。
$ tar -xvzf blast.tar.gz # -xvzf の部分をオプションという blast/ # 以下、展開されてできたファイルの一覧 blast/AtSeq/ blast/AtSeq/AT1G01010.1 blast/AtSeq/AT1G01020.1 blast/AtSeq/AT1G01030.1 blast/AtSeq/AT1G01040.1 blast/AtSeq/AT1G01050.1 blast/AtSeq/AT1G01060.1 blast/AtSeq/AT1G01070.1 blast/AtSeq/AT1G01080.1 blast/AtSeq/AT1G01090.1 blast/AtSeq/AT1G01100.1 blast/db/ blast/db/AT1Gpro.fa blast/db/ch1 blast/doc/ blast/doc/README.txt blast/TfSeq/ blast/TfSeq/TfA000003.fa blast/TfSeq/TfA000005.fa blast/TfSeq/TfA000022.fa blast/TfSeq/TfA000026.fa blast/TfSeq/TfA000045.fa blast/TfSeq/TfA000050.fa blast/TfSeq/TfA000053.fa blast/TfSeq/TfA000075.fa $ ls blast blast.tar.gz reference.fa # blast というフォルダ(青色の文字で示される)が追加されていることを確認する
tarはLinuxで標準的なファイルアーカイバーで、圧縮された書庫ファイルを解凍、展開できる。 -xvzf はオプションで、一文字ごとに意味がある。xは展開する、vは作業途中の様子を画面に表示する、zは解凍する、fは目的の書庫ファイルを示す。 ここでは、tar -xvzf file_name という使い方を覚えておけば十分である。
以後、blastフォルダで作業をする。
$ cd blast # ここでbash: cd: blast: No such file or directoryのようなエラーがでた場合は、上記のblast.tar.gzの解凍ができていない $ pwd /home/jovyan/blast
blastフォルダにはAtSeq、db、doc、TfSeqの4つのフォルダがある。
$ ls AtSeq db doc TfSeq $ ls AtSeq # AtSeqフォルダの中身を表示する AT1G01010.1 AT1G01030.1 AT1G01050.1 AT1G01070.1 AT1G01090.1 AT1G01020.1 AT1G01040.1 AT1G01060.1 AT1G01080.1 AT1G01100.1 $ ls db # dbフォルダの中身を表示する AT1Gpro.fa ch1 $ ls db/ch1 # dbフォルダの中のch1ファイルを表示する db/ch1 $ ls -l db/ch1 # dbフォルダの中のch1ファイルの詳細を表示する(-lはlsコマンドの動作を少し変更するオプション) -rw-r--r--. 1 jovyan users 30427677 Dec 2 2014 db/ch1 # 30427677はファイルの大きさ(単位はバイト(Byte, B))
BLASTは高速に似た塩基配列を検索することができる。 そのためには検索対象の塩基配列からデータベースをつくっておく必要がある。
dbフォルダの中にあるch1はシロイヌナズナの1番染色体の塩基配列(約3000万塩基対)である(less db/ch1 として、少し見てみるとよい。間違ってもcat db/ch1としないこと)。 これをBLAST検索のデータベースに変換するには以下のようにする(実行ディレクトリに注意)。
$ makeblastdb -in db/ch1 -dbtype nucl Building a new DB, current time: 12/15/2022 12:40:26 New DB name: /home/jovyan/blast/db/ch1 New DB title: db/ch1 Sequence type: Nucleotide Keep MBits: T Maximum file size: 1000000000B Adding sequences from FASTA; added 1 sequences in 0.43989 seconds. $ ls db AT1Gpro.fa ch1 ch1.nhr ch1.nin ch1.nsq # ch1.nhr ch1.nin ch1.nsqの3つのファイルが増えている。これがBLASTのデータベース。
BLAST検索にはクエリー、データベースを指定する必要がある。
AtSeqフォルダの中にある10個のファイルはシロイヌナズナの遺伝子の塩基配列である。 このうちAT1G01010.1(AtSeq/AT1G01010.1)をクエリーにして検索を行う。
今回作成したデータベース(db/ch1)は塩基配列のデータベースである。
塩基配列をクエリーに、塩基配列のデータベースを検索するには、blastnコマンドを使用する。
$ cat AtSeq/AT1G01010.1 # クエリーの中身を見ておく
>AT1G01010.1
AAATTATTAGATATACCAAACCAGAGAAAACAAATACATAATCGGAGAAATACAGATTACAGAGAGCGAGAGAGATCGAC
GGCGAAGCTCTTTACCCGGAAACCATTGAAATCGGACGGTTTAGTGAAAATGGAGGATCAAGTTGGGTTTGGGTTCCGTC
(途中省略)
AATTGGAAAAAAGTATATGTAATAATAATTAGTGCATCGTTTTGTGGTGTAGTTTATATAAATAAAGTGATATATAGTCT
TGTATAAG
$ blastn -query AtSeq/AT1G01010.1 -db db/ch1 # -query に続くファイルがクエリーとなる。-dbに続いてデータベースを指定する。
BLASTN 2.9.0+
Reference: Zheng Zhang, Scott Schwartz, Lukas Wagner, and Webb
Miller (2000), "A greedy algorithm for aligning DNA sequences", J
Comput Biol 2000; 7(1-2):203-14.
Database: db/ch1
1 sequences; 30,427,671 total letters
Query= AT1G01010.1
Length=1688
Score E
Sequences producing significant alignments: (Bits) Value
ch1 856 0.0
> ch1
Length=30427671
Score = 856 bits (463), Expect = 0.0
Identities = 463/463 (100%), Gaps = 0/463 (0%)
Strand=Plus/Plus
Query 1226 AGGTGATAAGTTCGCAGAAAAGCGAATGCGAGTGGAAAATGGCTGAAGACTCGATCAAGA 1285
||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||
Sbjct 5437 AGGTGATAAGTTCGCAGAAAAGCGAATGCGAGTGGAAAATGGCTGAAGACTCGATCAAGA 5496
(以下省略)
結果をファイルに保存するにはさらに-outオプションをつける。
$ blastn -query AtSeq/AT1G01010.1 -db db/ch1 -out AT1G01010.1.txt $ less AT1G01010.1.txt # スペースキーで先に、bキーで前に移動できる。終わるときはqキーを押す。
タンパク質をクエリーに、タンパク質データベースを検索するときには blastp コマンドを使う。 データベースの作成は makeblastdb コマンドで同じように行うことができるが、-dbtype オプションが異なるので注意する。
$ makeblastdb -in db/AT1Gpro.fa -dbtype prot タンパク質のデータベースを作成する Building a new DB, current time: 12/15/2022 13:33:16 New DB name: /home/jovyan/blast/db/AT1Gpro.fa データベースの名前 New DB title: db/AT1Gpro.fa Sequence type: Protein Keep MBits: T Maximum file size: 1000000000B Adding sequences from FASTA; added 7065 sequences in 0.287537 seconds. $ blastp -query TfSeq/TfA000003.fa -db db/AT1Gpro.fa -out TfA000003.txt
TfA000003.txtを開くと以下の内容を確認できる。
BLASTP 2.9.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.
Reference for composition-based statistics: Alejandro A. Schaffer,
L. Aravind, Thomas L. Madden, Sergei Shavirin, John L. Spouge, Yuri
I. Wolf, Eugene V. Koonin, and Stephen F. Altschul (2001),
"Improving the accuracy of PSI-BLAST protein database searches with
composition-based statistics and other refinements", Nucleic Acids
Res. 29:2994-3005.
Database: db\AT1Gpro.fa
7,065 sequences; 2,923,092 total letters
Query= TfA000003
Length=252
Score E
Sequences producing significant alignments: (Bits) Value
AT1G21720.1 397 3e-141
AT1G77440.1 393 8e-140
AT1G13060.1 63.5 2e-012
AT1G53850.1 34.3 0.010
(以下省略)
blastn などには -outfmt オプションがある。 これを指定すると出力の形式を変更することができる。 最も有用なものは -outfmt 6 で、検索結果(スコアや一致度)を表形式で出力することができる。
$ blastp -query TfSeq/TfA000003.fa -db db/AT1Gpro.fa -out TfA000003.txt -outfmt 6 TfA000003 AT1G21720.1 91.176 204 18 0 49 252 1 204 3.99e-143 397 TfA000003 AT1G77440.1 91.176 204 18 0 49 252 1 204 1.30e-141 393 TfA000003 AT1G13060.1 25.926 189 138 2 56 244 57 243 1.64e-12 63.5 TfA000003 AT1G53850.1 25.140 179 113 8 56 222 34 203 0.010 34.3 TfA000003 AT1G30410.1 23.239 142 82 4 114 234 61 196 0.13 31.6 TfA000003 AT1G70530.1 30.769 78 46 2 115 191 433 503 1.7 27.7 TfA000003 AT1G72470.1 33.898 59 39 0 93 151 488 546 1.8 27.7 TfA000003 AT1G02260.1 22.628 137 91 5 63 199 154 275 4.6 26.6 TfA000003 AT1G47250.1 22.807 57 41 1 55 111 31 84 5.8 25.8 TfA000003 AT1G54090.1 40.000 35 21 0 93 127 478 512 9.5 25.4
outfmtオプションをつけなかったときの出力と比べて、それぞれの数値が何を意味するのかを考えよ。