關(guān)于GO 注釋的心得體會(huì)
目前對于GO功能注釋的思路有 以下常見的四種:
1斥滤、BLAST+InterProScan => OmicsBox
可以先通過使用 blastp 進(jìn)行nr注釋缤灵,那么問題就來了臀叙,nr 數(shù)據(jù)庫截至2018年12月已經(jīng)有120G左右的數(shù)據(jù)量了角塑,根據(jù)我的個(gè)人使用經(jīng)驗(yàn)簿透,對25000條蛋白進(jìn)行32線程的并行化比對砸讳,差不多要半個(gè)月左右(心中一萬匹cnm飄過)幕与。很顯然速兔,在計(jì)算資源有限的條件下這是不現(xiàn)實(shí)的,那么大概就會(huì)出現(xiàn)兩種思路:
(1)縮小數(shù)據(jù)庫:一種是使用Swiss-Prot等小的數(shù)據(jù)庫進(jìn)行注釋淤齐;另一種就是根據(jù)基于taxid構(gòu)建 nr子數(shù)據(jù)庫了股囊,具體的構(gòu)建可見 基于taxid構(gòu)建Blast database_bioinfomatics2medicine_新浪博客 (2)能不能使用一個(gè)更高效的比對工具呢,答案是肯定的 那就是 GitHub - bbuchfink/diamond: Accelerated BLAST compatible local sequence aligner.
其主頁介紹如下:
DIAMOND is a sequence aligner for protein and translated DNA searches, designed for high performance analysis of big sequence data. The key features are:
Pairwise alignment of proteins and translated DNA at 500x-20,000x speed of BLAST.
Frameshift alignments for long read analysis.
Low resource requirements and suitable for running on standard desktops or laptops.
Various output formats, including BLAST pairwise, tabular and XML, as well as taxonomic classification.
對的床玻,它就是那么快毁涉,相對于 blastp 的半個(gè)月,在 ‘--more-sensitive’ 模式下也僅需要4-5個(gè)小時(shí)锈死。
在nr注釋完之后,就到使用 OmicsBox mapping GO id 了∑堆撸現(xiàn)在所有的服務(wù)都要收費(fèi)了,https://www.biobam.com/omicsbox/#plans. 但價(jià)格比以前(Blast2GO)要公道很多, 如果不缺少資金可以購買待牵,但如果不想花錢其屏,那么不如就來個(gè)本地化吧,對于本地化網(wǎng)上也有很多的教程了:
blast2go本地化2017教程:blast2go本地化2017教程 - 生信技能樹
非root權(quán)限的blast2go的安裝和使用(二)· blast2go的數(shù)據(jù)和軟件準(zhǔn)備及使用:https://mp.weixin.qq.com/s?__biz=MzI1MjU5MjMzNA==&mid=2247486358&idx=1&sn=4b095e9a200a419079d947d930c38abd&chksm=e9e02237de97ab21f4227ee76eca0aa94433fe6cc374e030093eb3acc2094a546110763bb411&mpshare=1&scene=23&srcid=0101kMbzD2el0RPASBiC9pFi#rd
陳連福 centos6.9 代碼:
“# installing Blast2go Databases
mkdir /opt/biosoft/blast2go
cd /opt/biosoft/blast2go
wget http://archive.geneontology.org/latest-full/go_monthly-assocdb-data.gz
wget ftp://ftp.ncbi.nlm.nih.gov/gene/DATA/gene_info.gz ### ascp -k 1 -T -l 200M -i ~/.aspera/connect/etc/asperaweb_id_dsa.openssh --host=ftp-private.ncbi.nlm.nih.gov --user=anonftp --mode=recv /gene/DATA/gene_info.gz 缨该。/
wget ftp://ftp.ncbi.nlm.nih.gov/gene/DATA/gene2accession.gz? #### ascp -k 1 -T -l 200M -i ~/.aspera/connect/etc/asperaweb_id_dsa.openssh --host=ftp-private.ncbi.nlm.nih.gov --user=anonftp --mode=recv /gene/DATA/gene2accession.gz ./
wget ftp://ftp.pir.georgetown.edu/databases/idmapping/idmapping.tb.gz
## 對于不能使用aspera加速下載的文件可以試試 lftp -e 'pget -n NUM -c url; exit'實(shí)行多線程下載
gzip -dv go_monthly-assocdb-data.gz
gzip -dv gene_info.gz
gzip -dv gene2accession.gz
gzip -dv idmapping.tb.gz
tar zxf ~/software/local_b2g_db.tar.gz
mv local_b2g_db/* ./? && rm -rf local_b2g_db/
perl -p -i -e 's/go_201512-assocdb-data/go_2017-assocdb-data/' install_blast2goDB.sh
./install_blast2goDB.sh
#一個(gè)自己的例子 #
##下載NR 數(shù)據(jù)庫到自己的Blast+ db 中
wget ftp://ftp.ncbi.nlm.nih.gov/blast/executables/blast+/2.2.28/ncbi-blast-2.2.28+-x64-linux.tar.gz
tar -zxf ncbi-blast-2.2.28+-x64-linux.tar.gz -C /opt/biosoft/
/opt/biosoft/ncbi-blast-2.2.28+/bin/blastdbcmd -db nr -entry all -out nr.faa
### 值得注意的是這里我用的blastdbcmd是ncbi-blast-2.2.28+的偎行,輸出的 nr.faa 文件中序列名是 “>gi|66816243|ref|XP_642131.1| hypothetical protein DDB_G0277827” ,但如果用的是 ncbi-blast-2.6.0+的 blastdbcmd贰拿,輸出沒有GI號蛤袒,原因可能是:As of September 2016, the integer sequence identifiers known as "GIs" will no longer be included in the GenBank, GenPept, and FASTA formats supported by NCBI for the display of sequence records.In addition, the FASTA format will no longer include the database source abbreviation. Please refer to the NCBI News Announcement posting for more detail.
wget https://github.com/bbuchfink/diamond/releases/download/v0.9.24/diamond-linux64.tar.gz
### 同樣注意版本,因?yàn)閐iamond輸出的XML與 b2g4pipe_v2.5 不兼容而出現(xiàn) “Annotation of 0 seqs with 0 annots finished. Now searching for orfan IPRs...” 這個(gè)issues的討論膨更,你可以在到下面的鏈接去看妙真。
tar -zxf diamond-linux64.tar.gz
mv diamond ~/bin/
cp ../interpro/proteins.fasta ./
ascp -k 1 -T -l 200M -i ~/.aspera/connect/etc/asperaweb_id_dsa.openssh --host=ftp-private.ncbi.nlm.nih.gov --user=anonftp --mode=recv /pub/taxonomy/accession2taxid/prot.accession2taxid.gz ./prot.accession2taxid.gz
ascp -k 1 -T -l 200M -i ~/.aspera/connect/etc/asperaweb_id_dsa.openssh --host=ftp-private.ncbi.nlm.nih.gov --user=anonftp --mode=recv /pub/taxonomy/taxdmp.zip ./
unzip taxdmp.zip
/opt/biosoft/ncbi-blast-2.2.28+/bin/blastdbcmd -db nr -entry all -out nr.faa
diamond makedb --in nr.faa -d nr --taxonmap prot.accession2taxid.gz --taxonnodes nodes.dmp #### DIAMOND v0.9.24
diamond blastp --query proteins.fasta --more-sensitive --db nr --evalue 1e-5 --salltitles --threads 64 --outfmt 5 --out nr.xml
perl -p -e 's/diamond 0.9.24/BLASTP 2.2.26/' nr.xml > nr_new.xml
perl -p -i -e 's/^Dbacces.dbname=.*/Dbacces.dbname=b2gdb/' b2gPipe.properties
perl -p -i -e 's/^Dbacces.dbhost=.*/Dbacces.dbhost=127.0.0.1/' b2gPipe.properties
java -cp *:ext/*: es.blast2go.prog.B2GAnnotPipe -in ../nr_new.xml -out go -annot -dat -annex
# 至于整合InterProScan的結(jié)果可以使用下面的腳本,當(dāng)然也可以在java -cp *:ext/*: es.blast2go.prog.B2GAnnotPipe -ips 來實(shí)現(xiàn)荚守,或者將b2g4pipe的輸出文件加載到 Blast2GO Basic 中進(jìn)行后續(xù)可視化分析珍德。
merge_interpro_to_go.pl b2g4pipe/go.annot ../interpro/interpro.tsv > go.annot
# Anyone use diamond output xml file? for b2gpipe? · Issue #79 · bbuchfink/diamond · GitHub
# diamond to blast2go erro? · Issue #165 · bbuchfink/diamond · GitHub
2、InterProScan
對于 InterProScan 可以本地化矗漾,也可以用 interProScan5.pl 將序列發(fā)送到官方服務(wù)器進(jìn)行注釋锈候,資源允許,還是推薦本地化敞贡,網(wǎng)絡(luò)版感覺還是有點(diǎn)慢~一個(gè)本地化教程:
3泵琳、eggNOG-Mapper (包括網(wǎng)頁版和本地化)
(1)關(guān)于這個(gè)注釋方法網(wǎng)頁版的推薦自己去摸索吧(就是和NCBI提交類似)=> http://eggnogdb.embl.de/#/app/emapper
(2)本地化推薦兩個(gè)教程吧:
序列功能注釋神器:eggNOG-mapper,KEGG/COG/KOG/GO/BiGG 一網(wǎng)打盡 : 序列功能注釋神器:eggNOG-mapper嫡锌,KEGG/COG/KOG/GO/BiGG 一網(wǎng)打盡 ? Biostack.org
應(yīng)該是最好的eggnog-mapper功能注釋教程:應(yīng)該是最好的eggnog-mapper功能注釋教程 - 簡書(生信媛)
4虑稼、PANNZER2-終極注釋網(wǎng)站
此部分為2020年12月24日后加,該工具也是一個(gè)網(wǎng)頁版工具势木。相較于eggNOG-mapper蛛倦,PANNZER2注釋上的基因條目與Blast2GO差不多,而且啦桌,2.5W個(gè)基因3-4h就可以拿到結(jié)果溯壶。如果再整合上InterProScan的結(jié)果及皂,那么,這個(gè)方法便將高效且改、準(zhǔn)確集于一身验烧。
前方高能:
對于準(zhǔn)確性的個(gè)人思考,基于比對又跛,也就是基因相似性的注釋方法(blastp,diamond)可能會(huì)取決于數(shù)據(jù)庫的大小碍拆,因?yàn)樵赽last2go中有涉及到相似度cut這一步,如果數(shù)據(jù)庫不夠大慨蓝,那么相似度有的可能只為50(在cut值之上)感混,所以只使用Swiss-Prot等小數(shù)據(jù)庫注釋時(shí)可能需要提高cut值;當(dāng)然也會(huì)聽到直系同源和旁系同源的說法,所以在新版Blast2go中增加了filter GO with Taxid (可能是看到eggnog的文章慌了)這一新模塊礼烈,減少旁系同源的注釋結(jié)果弧满。那么如果要在老版的Blast2go中運(yùn)用這個(gè)模塊,或許可以換個(gè)思路此熬,在比對時(shí)庭呜,通過Taxid來過濾比對結(jié)果。
eggnog的注釋犀忱,存在兩種模式diamond和hummer募谎,前者官方推薦為存在接近物種,也就是在直系同源注釋中有更好的表現(xiàn)阴汇,hummer則反之. 同樣的近哟,只要是比對,就會(huì)有cut相似度這一說鲫寄,數(shù)據(jù)庫的大小及深度(涉及的物種),就會(huì)影響結(jié)果疯淫。其一個(gè)特點(diǎn)在于直系同源注釋地来。這點(diǎn)Blast2go中出了filter with Taxid 來補(bǔ)充。
小提示:
參考文獻(xiàn)
http://blog.sina.com.cn/s/blog_16152d7d70102xnc3.html
https://github.com/bbuchfink/diamond
https://mp.weixin.qq.com/s?__biz=MzI1MjU5MjMzNA==&mid=2247486358&idx=1&sn=4b095e9a200a419079d947d930c38abd&chksm=e9e02237de97ab21f4227ee76eca0aa94433fe6cc374e030093eb3acc2094a546110763bb411&mpshare=1&scene=23&srcid=0101kMbzD2el0RPASBiC9pFi#rd
http://www.jinciwei.cn/a342744.html
https://github.com/bbuchfink/diamond/issues/79
https://github.com/bbuchfink/diamond/issues/159
https://github.com/bbuchfink/diamond/issues/165
http://www.biostack.org/?p=698
http://www.reibang.com/p/e646c0fa6443
NGS 生物信息學(xué)分析 V6.0 陳連福 鄭越
http://www.reibang.com/p/6296385adf21