有参转录组实战8-基因功能注释_GO_KEGG_swisspro

2024-01-23  本文已影响0人  啊辉的科研

#进行功能注释时,我们只用到蛋白文件,就是上一期提取序列的文件“Ptri.protein.fa”。

#使用命令“grep -c ">" Ptri.protein.fa”统计下“>”的个数,发现有52400个。

#新建文件夹“swissprot”

wget https://ftp.uniprot.org/pub/databases/uniprot/current_release/knowledgebase/complete/uniprot_sprot.fasta.gz

gunzip -c uniprot_sprot.fasta.gz >uniprot_sprot.fasta#解压

conda install diamond#安装

diamond makedb --in uniprot_sprot.fasta --db uniprot_sprot.fasta#建索引

nohup diamond blastp -d uniprot_sprot.fasta -q Ptri.protein.fa --max-target-seqs 1 --outfmt 6 --evalue 1e-5 > blastp.out &#注释

#查看文件blatp.out,十二列解释看表头

#新建文件夹“pfam”

wget http://ftp.ebi.ac.uk/pub/databases/Pfam/releases/Pfam34.0/Pfam-A.hmm.gz#下载

gunzip -c Pfam-A.hmm.gz > Pfam-A.hmm#解压

conda install hmmer#安装

hmmpress Pfam-A.hmm#构索引

nohup hmmscan --domtblout pfam.domtblout Pfam-A.hmm

Ptri.protein.fa &#注释

#查看文件pfam.domtblout

#Plant TFDB网站预测转录因子http://planttfdb.gao-lab.org/prediction.php

#预测了3835个转录因子。自己下载整理。

#iTAK预测转录因子、调控因子、激酶http://itak.feilab.net/cgi-bin/itak/online_itak.cgi

#接着GO注释和KEGG注释。以下重点参考https://zhuanlan.zhihu.com/p/475588763教程。

#http://eggnog-mapper.embl.de/

#打开邮箱

#开始工作

#等约半小时,打开链接下载结果

#只要这一个

#使用TBtools的这个功能

#放入注释文件

#得到几个txt文件,后面的富集会用到。

#最后自己用excel整理下,可以得到基因的各种注释信息。

#赛博朋克边缘行者

上一篇 下一篇

猜你喜欢

热点阅读