应该是最好的eggnog-mapper功能注释教程

eggnog-mapper实现功能注释

eggNOG-Mapper介绍

通常功能注释的思路都是基于序列相似性找直系同源基因,常见的方法就是BLAST+BLAST2GO, 或者是InterProScan。eggNOG-mapper的作者认为这种方法不够可靠,毕竟你有可能找到的的是旁系同源基因。近期对多个工具的整体评估发现eggNOG(evolutionary genealogy of genes: Non-supervised Orthologous Groups)在区分旁系同源基因和直系同源基因上表现不错,因此基于eggNOG数据库开发了eggNOG-mapper工具,用于对新序列进行功能注释。

eggNOG-mapper的算法实现如下:

第一步:序列比对。首先,每条蛋白序列用HMMER3在整理的eggNOG数据库中搜索。由于每个HMM匹配都和一个功能注释的eggNOG OG对应,这一步就提供了初步的注释信息。之后,每条蛋白序列用phmmer在最佳匹配的HMM对应的一组eggNOG蛋白中进一步搜索。最后,每条序列的最佳匹配结果以 seed ortholog 形式存放,用于获取其他直系同源基因。目前eggNOG HMM数据库中拥有1,911,745个OG,覆盖了1,678种细菌,115种古细菌,238种真核物种以及352种病毒。除了HMMER3外,还而可用DIAMOND直接对所有的eggNOG蛋白序列进行搜索,它的速度更快,适合类似于宏基因组这类大数据集,或者是已有物种和eggNOG所收集的物种比较近。当然服务器性能强大的话,还是有限选择HMMER3.

应该是最好的eggnog-mapper功能注释教程

第二步:推测直系同源基因。每个用于检索的蛋白序列的最佳匹配序列会对应eggNOG的一个蛋白, 这些蛋白基于预分析的eggNOG进化树数据库会提取一组更加精细的直系同源基因。这一步还会根据bit-screo或E-value对结果进行一次过来,剔除同源性不高的结果

应该是最好的eggnog-mapper功能注释教程

第三步:功能注释。用于搜索的蛋白序列对应的直系同源基因的功能描述就是最终的注释结果。比如说GO, KEGG, COG等。

应该是最好的eggnog-mapper功能注释教程

安装

安装本体eggnog-mapper之前,需要先保证服务器上已经安装Python3.7(BioPython模块), wget, HMMER3 和/或 DIAMON,此外你还得保证70G的空间用于存放注释数据库和FASTA文件,对于真核生物至少保证90G的服务器内存。
之后从GitHub上将软件下载到本地

git clone https://github.com/jhcepas/eggnog-mapper.git

之后需要下载所需要的数据库。eggNOG提供了107个不同物种的HMM数据库(xxxNOG)以及三个优化数据库, euk对应真核生物,bact对应细菌, arch对应古细菌(Archeabacteria), 以及一个病毒数据库(viruses). 这三个优化数据库包含了属于该分类内的所有物种的HMM模型。

cd eggnog-mapper
./download_eggnog_data.py euk

此处下载的是真核生物 。下载过程中,它会反复询问你是否要下载某一类数据,我一律选择是。
如果你的服务器有多个盘,安装软件的分区不够大,可以在数据盘中创建文件夹进行软链接

cd eggnog-mapper
mkdir /data/database/EGGNOG-DB
rmdir data && ln -s /data/database/EGGNOG-DB data

具体用法

eggnog-mapper用起来非常的简单,你需要提供蛋白序列作为输入

#假如我们现在都仍在软件安装的路径下
python emapper.py -i test/p53.fa --output p53_maNOG -d euk
python emapper.py -i test/p53.fa --output p53_maNOG -d maNOG
python emapper.py -i test/p53.fa --output p53_maNOG -d maNOG --usemem --cpu 10

eggnog-mapper默认是以HMMER进行序列搜索,尽管这可以通过-m diamond更改成DIAMOND,但结果中就会缺少一些信息列。HMM的搜索参数,有--hmm_maxhits, --hmm_evalue, --hmm_score, --hmm_qcov--Z, 一般用默认就行了。

-d/--database表示搜索数据库的类型,既可以是"euk, bact, arch" 三大类的其中一种,也可以是"eggnog-mapper/data/OG_fasta"中的小类,例如"maNOG"表示的就是哺乳动物的NOG。此外也可以是自定义HMM数据库, -d /path/to/pfam.hmm

如果服务器内存比较大,线程比较多,可以用--usemem--cpu 线程数提高运行速度。

--output表示输出文件的前缀,默认输出在当前文件夹下,--output_dir可以更改为其他文件路径。--resume表示任务重启后可以跳过之前已经完成的部分, 而--override则表示覆盖原先的输出结果。

如果你计划使用同一个数据库同时对多个物种进行注释(大部分人都没有这种需求),那么你需要以服务器模式启动eggnog-mapper,然后就可以同时启动多个任务,

# 终端1
python emapper.py -d arch --cpu 10 --servermode
# 新建一个终端
python emapper.py -d arch:localhost:51600 -i test/polb.fa -o polb_arch

结果解读

eggnog-mapper会生成三个文件,

  • [project_name].emapper.hmm_hits: 记录每个用于搜索序列对应的所有的显著性的eggNOG Orthologous Groups(OG). 所有标记为"-"则表明该序列未找到可能的OG
  • [project_name].emapper.seed_orthologs: 记录每个用于搜索序列对的的最佳的OG,也就是[project_name].emapper.hmm_hits里选择得分最高的结果。之后会从eggNOG中提取更精细的直系同源关系(orthology relationships)
  • [project_name].emapper.annotations: 该文件提供了最终的注释结果。大部分需要的内容都可以通过写脚本从从提取,一共有13列。

[project_name].emapper.annotations每一列对应的记录如下:

  1. query_name: 检索的基因名或者其他ID
  2. sedd_eggNOG_ortholog: eggNOG中最佳的蛋白匹配
  3. seed_orholog_evalue: 最佳匹配的e-value
  4. seed_ortolog_evalu: 最佳匹配的bit-score
  5. predicted_gene_name: 预测的基因名,特别指的是类似AP2有一定含义的基因名,而不是AT2G17950这类编号
  6. GO_term: 推测的GO的词条, 未必最新
  7. KEGG_KO: 推测的KEGG KO词条, 未必最新
  8. BiGG_Reactions: BiGG代谢反应的预测结果
  9. Annotation_tax_scope: 对该序列在分类范围的注释
  10. Matching_OGs: 匹配的eggNOG Orthologous Groups
  11. best_OG|evalue|score: 最佳匹配的OG(HMM模式才有)
  12. COG functional categories: 从最佳匹配的OG中推测出的COG功能分类
  13. eggNOG_HMM_model_annotation: 从最佳匹配的OG中推测出eggNOG功能描述

如果打算做富集分析,用命令行的cut/awk提取对应的列,过滤掉其中未注释的部分就行了。

推荐阅读

上一篇:阿里大规模计算平台的自动化、精细化运维之路


下一篇:如何在shell脚本中控制任务投递