安裝本地blast序列比對軟件，我們可以搜索一個查詢序列定制數(shù)據(jù)庫，例如想研究一個新測序的基因組，或者感興趣的一組蛋白質(zhì)序列。有時我們希望把程序插入到一個流程中，例如搜索一個大量的查詢序列，例如你的測序數(shù)據(jù)含有大量的污染片段，你想知道這些片段比對到了什么物種。
blast軟件安裝
從NCBI下載安裝包https://ftp.ncbi.nlm.nih.gov/blast/executables/blast+/LATEST/。下載2.11.0 linux版本，下載并解壓，解壓之后BLAST就安裝好了。用戶需要設(shè)置環(huán)境變量，目的是為了告訴系統(tǒng)在哪里可以找到安裝好的BLAST軟件。
wget https://ftp.ncbi.nlm.nih.gov/blast/executables/blast+/LATEST/ncbi-blast-2.11.0+-x64-linux.tar.gz
wget https://ftp.ncbi.nlm.nih.gov/blast/executables/blast+/LATEST/ncbi-blast-2.11.0+-x64-linux.tar.gz.md5
md5sum -c ncbi-blast-2.11.0+-x64-linux.tar.gz.md5
tar -zxvf ncbi-blast-2.11.0+-x64-linux.tar.gz
#加入到環(huán)境變量
export PATH=$PATH:$PWD/ncbi-blast-2.11.0+/bin
nt/nr fasta下載
從NCBI（ftp://ftp.ncbi.nlm.nih.gov/blast/db/FASTA/
）下載nt、nr fasta文件。
wget ftp://ftp.ncbi.nlm.nih.gov/blast/db/FASTA/nt.gz
wget ftp://ftp.ncbi.nlm.nih.gov/blast/db/FASTA/nr.gz
wget ftp://ftp.ncbi.nlm.nih.gov/blast/db/FASTA/nt.gz.md5
wget ftp://ftp.ncbi.nlm.nih.gov/blast/db/FASTA/nr.gz.md5
md5sum -c nt.gz.md5
md5sum -c nr.gz.md5
tar -xzvf nr.gz
tar -xzvf nt.gz
mkdir nr_db
mkdir nt_db
makeblastdb -in nr -dbtype prot -title make_nr -parse_seqids -out ./nr_db/nr -logfile make_nr.log
makeblastdb -in nt -dbtype nucl -title make_nt -parse_seqids -out ./nt_db/nt -logfile make_nt.log
或者使用NCBI處理好的db。
wget -c ftp://ftp.ncbi.nlm.nih.gov/blast/db/nt*
wget -c ftp://ftp.ncbi.nlm.nih.gov/blast/db/nr*
序列比對
- blastp：蛋白序列與蛋白庫作比對，直接比對蛋白序列的同源性。
- blastx：核酸序列與蛋白庫作比對，將核酸序列先翻譯成蛋白序列，再將其與蛋白庫作比對。
-blastn：核酸序列與核酸庫的比對，直接比對核酸序列的同源性。 - tblastn：蛋白序列對核算庫的比對，現(xiàn)將核酸庫翻譯成蛋白庫，再將蛋白序列與翻譯后的蛋白庫進行比對。
- tblastx：核酸與核酸數(shù)據(jù)庫在蛋白質(zhì)水平比較
如果是fastq先轉(zhuǎn)為fasta序列格式在進行比對，比對結(jié)果如下圖1所示。將相同序列復(fù)制到NCBI網(wǎng)站進行比對，如下圖2，兩種結(jié)果最前面的比對條目基本一致。
awk '{if(NR%4 == 1){print ">" substr($0, 2)}}{if(NR%4 == 2){print}}' test.fastq > test.fasta
blastn -query test.fasta -out test.result -db ./nt_db/nt
圖1 本地blast比對結(jié)果
圖2 NCBI blast比對結(jié)果
使用NCBI網(wǎng)站下載的db進行比對。
wget -c ftp://ftp.ncbi.nlm.nih.gov/blast/db/nr*
for i in `ls nt_db_NCBI/*gz`do;tar -zxvf $i;done
## 輸出一條最優(yōu)比對結(jié)果
blastn -query test1.fa -out test1.align -db ./nt_db_NCBI/nt -outfmt 6 -subject_besthit -num_threads 4
自建索引文件
建立索引需要輸入為序列的fasta文件，可以是蛋白或者核酸，輸出索引帶前綴名稱，最終生成比對索引數(shù)據(jù)庫。
- in 輸入文件（test.fa)
- dbtype 選擇類型，蛋白或者核酸（prot、nucl)
- out 輸出索引文件
makeblastdb -in test.fa -dbtype nucl -out ./test
