写在前面
之前介绍过使用
EggNOG-mapper
对非模式物种快速完成KEGG注释的推文https://mp.weixin.qq.com/s/kIf6C2u3FID3ZeLtsB4eZQ
,不过毕竟EggNOG-mapper
是第三方软件,而且注释结果仍然存在一个Gene ID对应多个KO Number的情况,虽然对富集结果影响似乎不会太大,但是影响多少肯定是会有。这里介绍KEGG官方提供的注释工具KofamKOALA,同样能够对任意物种进行KEGG注释,注释结果中Gene ID对应唯一一个KO Number,避免了大多数第三方注释工具一个ID对应多个KO Number的情况,使富集结果更加可靠!
KofamKOALA官网
https://www.genome.jp/tools/kofamkoala/
一、网页版KofamKOALA
点击提交任务之后需要在邮件中二次确认,提交任务
注释完成之后,点击邮件中的链接
- 在结果页面中,可以直接下载注释结果,然后接着用TBtools进行富集分析
-
同时也提供了许多对注释结果的探索功能,比如查看注释上的KO Number详细信息,比对情况等。
二、本地化KofamKOALA
- 如果需要注释的序列太多,文件太大,使用网页版工具确实可能会太慢。这个时候当然是选择将其本地化!在本地根据服务器情况提高线程加速注释。
下载 KOfam 和 KofamScan
下载并解压KOfam
wget https://www.genome.jp/ftp/db/kofam/ko_list.gz
wget https://www.genome.jp/ftp/db/kofam/profiles.tar.gz
#下载完成之后解压
gunzip ko_list.gz
tar -xzvf profiles.tar.gz
下载并解压kofam_scan
wget https://www.genome.jp/ftp/tools/kofam_scan/kofam_scan-1.3.0.tar.gz
tar -xzvf kofam_scan-1.3.0.tar.gz
将kofam_scan加入环境变量
echo export PATH=/home/jiawei_li/tools/kofamscan/kofam_scan-1.3.0:\$PATH >> ~/.bashrc
source ~/.bashrc
安装依赖,KofamScan需要Ruby,HMMER和GNU Parallel
#版本需求
- Ruby >= 2.4
- HMMER >= 3.1
- GNU Parallel
conda安装Ruby,自己编译也行 但是我比较懒折腾。。。
conda install -c conda-forge ruby
- 如果喜欢手动装Ruby,也可以试试,但也许会报错,需要进行一系列操作去安装
wget https://cache.ruby-lang.org/pub/ruby/3.1/ruby-3.1.0.tar.gz
tar -xzvf ruby-3.1.0.tar.gz
cd cd ruby-3.1.0
./configure
make
make install
- 安装hmmer
conda install -c bioconda hmmer
- 安装GNU Parallel
conda install -c conda-forge parallel
修改配置文件,指定依赖软件以及KOfam的路径
首先获取相关软件的路径(注意:配置文件中不需要指定Ruby的路径,但是需要确保Ruby加入了环境变量,能够直接调用)
which parallel hmmscan
#/tools/parallel
#~/miniconda3/bin/hmmscan
使用官方模板,修改配置文件
cd kofam_scan-1.3.0
cp config-template.yml config.yml
- 注意,由于hmmscan和parallel都已经加入环境变量,可以直接调用,因此在配置文件中无需配置这两个软件的路径
# Path to your KO-HMM database
# A database can be a .hmm file, a .hal file or a directory in which
# .hmm files are. Omit the extension if it is .hal or .hmm file
profile: /home/jiawei_li/tools/kofamscan/profiles
# Path to the KO list file
ko_list: /home/jiawei_li/tools/kofamscan/ko_list
# Path to an executable file of hmmsearch
# You do not have to set this if it is in your $PATH
#hmmsearch: /home/jiawei_li/miniconda3/bin/hmmscan
# Path to an executable file of GNU parallel
# You do not have to set this if it is in your $PATH
#parallel: /tools/parallel
# Number of hmmsearch processes to be run parallelly
cpu: 8
对蛋白序列进行注释(注意:注释序列必须为蛋白序列)
exec_annotation -o test.querry2KO --cpu 8 --format mapper -E 1e-5 test.pep.fa
-
注意:
--format
参数在help中给出了四种,但主要有两种格式。
#help文档中的解释
-f, --format <format> Format of the output [detail]
detail: Detail for each hits (including hits below threshold)
detail-tsv: Tab separeted values for detail format
mapper: KEGG Mapper compatible format
mapper-one-line: Similar to mapper, but all hit KOs are listed in one line
-
--format mapper
对每一个Gene ID只保留最佳的KO Number,结果文件只含有Gene ID和KO Number的映射信息。
-
--format detail
保留所有信息,包括Gene ID对应上的每一个KO Number,比对分数,E-value以及KO Number的详细信息等。
最后,拿到query2KO文件之后,就可以用其当背景文件,使用TBtools对基因集做富集分析啦
写在后面
学期末了,祝各位老铁这学期的奋斗都有所收获!
大家都能回家安心过个好年,见见许久未见的老友哈哈