在生物信息學(xué)，BLAST(基本的局部比對搜索工具)是一種算法，用于比較生物的序列信息，諸如蛋白質(zhì)的氨基酸序列、DNA或RNA核苷酸序列。BLAST搜索使研究人員能夠?qū)⒁粋€目標(biāo)蛋白或核苷酸序列(稱為查詢)與一個數(shù)據(jù)庫進(jìn)行比較，并識別某個特定閾值以上的與目標(biāo)序列相似的庫序列。在線的BLAST功能只需要將序列輸入，然后BLAST，最后會給出相似的庫序列(按照同源性高低排列)。在線BLAST明顯存在兩個缺陷：1.目標(biāo)序列BLAST后的得到的庫序列是多個物種的，如果我們只是要特定物種的話，還需要再自己找。2.如果目標(biāo)序列有很多，則需要執(zhí)行很多次，效率低。而本地化的BLAST+能很好的解決上面兩個問題。

1.BLAST+的下載與安裝

BLAST+下載地址: ftp://ftp.ncbi.nlm.nih.gov/blast/executables/blast+/LATEST/

針對不同的平臺(linux,mac,windows)選擇性的進(jìn)行安裝，如圖一。
下面以mac為例進(jìn)行介紹：

1.下載:ncbi-blast-2.9.0+-x64-macosx.tar.gz
2.下載完成后，進(jìn)行解壓縮。
3.打開終端，輸入以下代碼：
#將目標(biāo)地址設(shè)為環(huán)境變量
open ~/.bash_profile       #打開環(huán)境變量配置窗口。
在彈出的窗口中輸入如下：
export PATH=$PATH:/Users/bcl/Desktop/software_plugin/ncbi-blast-2.9.0+/bin    #將blast文件夾中的bin地址設(shè)為環(huán)境變量。
#驗證是否成功配置為環(huán)境變量
blastn                    #出現(xiàn)圖二說明OK。

圖一

圖二

2.BLAST+的使用

下面以一個實例來說明怎么使用。比如說你想把番茄的一部分序列對應(yīng)到桃基因組數(shù)據(jù)庫中，并按照同源性高低得到與目標(biāo)序列相似的庫序列(根據(jù)不同的篩選標(biāo)準(zhǔn)，得到的序列條數(shù)會不一樣)，應(yīng)該怎么做呢？

2.1 序列的下載

首先你得從相應(yīng)的網(wǎng)站下載番茄一部分序列(因為我們以blastp為例，所以只下載了它的蛋白序列)，以及下載得到桃的基因組文件(只用其中的蛋白序列)。下載得到的文件放在一個文件夾中。

2.2 庫的構(gòu)建

利用下載得到的桃的基因組文件進(jìn)行構(gòu)建蛋白庫。在終端中輸入如下代碼：

cd /Users/bcl/Desktop/blastp        #cd到存放文件的地址。
makeblastdb -in Prunus_persica_v2.0.a1.allTrs.pep.fa -dbtype prot  
#下面為注釋：
-in:代表用來構(gòu)庫的基因組文件的地址和文件名
-dbtype:代表構(gòu)建的是什么庫，prot代表蛋白庫。

構(gòu)建成功會出現(xiàn)圖三以及文件夾中會多出三個文件(圖四)。

圖三

圖四

2.3 目標(biāo)序列與庫序列進(jìn)行比對

終端下輸入如下代碼：

blastp -help     #查看幫助文檔
blastp -query JA_regulated_genes.fasta -db Prunus_persica_v2.0.a1.allTrs.pep.fa -out JA_regulated_genes_blastp.txt -evalue 1e-5 -outfmt 6        
#下面為注釋。
-query：代表的是要進(jìn)行blastp的番茄一部分序列文件
-db：為前面用于建庫桃的的基因組文件
-out：為輸出文件的文件名
-evalue：為篩選標(biāo)準(zhǔn)(evalue越低，相似性越高)
-outfmt：輸出文件格式，6代表表格形式。

輸出的結(jié)果按照同源性從高到低排列(目標(biāo)序列與庫序列進(jìn)行同源比對時根據(jù)evalue篩選出來的可能不止一個)，那么如何選出每個最匹配的呢?

圖5

3.R studio中提取輸出文件中同源性最高的序列

打開R studio，輸入以下代碼，即可實現(xiàn)。

#讀取文件
txt <- read.csv('/Users/bcl/Desktop/blastp/JA_regulated_genes_blastp.txt',header=FALSE,sep='\t',encoding='UTF-8')   #第一個參數(shù)為上一步比對后生成文件的地址和文件名，后面的參數(shù)不需要做調(diào)整。
#構(gòu)建一個布爾向量，索引。
index<- duplicated(txt$V1) #duplicated()在數(shù)據(jù)框函數(shù)中應(yīng)用廣泛
index
#篩選數(shù)據(jù)(兩種)
            #method1
txt<-txt[which(index==FALSE),]  #選擇非重復(fù)數(shù)據(jù)
txt
            #method2
txt<-txt[!index,]                              #選擇非重復(fù)數(shù)據(jù)
txt
#數(shù)據(jù)寫入到txt中，把行名和列名設(shè)為flase,并去掉引號
write.table(txt,file='/Users/bcl/Desktop/blastp/JA_regulated_genes_select.txt',row.names=FALSE,col.names = FALSE,quote=F)    #第二個參數(shù)與讀取的時候同

下圖為經(jīng)過選擇后的文件可發(fā)現(xiàn)已提取出同源性最高的序列。

圖6

本文由博客一文多發(fā)平臺 OpenWrite 發(fā)布！