使用clusterProfiler对非模式植物进行注释

❝
非模式植物这个经常有人问，这里提到了clusterProfiler + GOMAP.
然后小伙伴给出了clusterProfiler + eggNOG的解决方法。使用eggNOG的话，就不限于植物了，「其它各种真核生物都能搞」，所以啊，还在等什么，clusterProfiler做富集分析，搞起来！
BTW:非常感谢小伙伴们的支持，继续征集clusterProfiler在各种非模式生物上的应用场景和实例。
PS: 如果你在学习生信的路上，在考研考博和读书的路上，有什么故事和体会，也欢迎来投稿，因为关注群体学生比较多，经历可以给后来者参考，也能鼓励后来者。我自己各种考研考博没人要，以及遇到变态老板的故事，估计也鼓励了不少人。
❞

流程

使用eggNOG对基因组蛋白序列进行注释
解析eggNOG结果文件
使用clusterProfiler进行富集分析

1.使用eggNOG对基因组进行注释

要想进行富集分析，首先要有背景数据集的GO注释和KEGG注释，这里选用eggNOG进行注释。原因主要有两个：

从我个人的使用经验来看，这个软件的注释结果优于InterProScan等软件
是在线服务器，点点鼠标上传就能注释，无需复杂配置。今天看见公众号截图中有人建议GOMAP，我没用过GOMAP，但是从该软件说明来看还是需要部署到本地的。
eggNOG虽然是web server，但一次最多可以注释10万条序列，应该是完全可以满足需求的，且一般几小时就能注释完。
我一般使用的参数如下

2.解析eggNOG文件

首先需要去GO下载GO的obo文件,这里我使用go-basic.obo然后我写了个脚本可以把obo文件解析为如下格式：

python parse_go_obofile.py -i go-basic.obo \
                           -o go.tb

处理后的文件：

GO	Description	level
GO:0000001	mitochondrion inheritance	biological_process
GO:0000007	low-affinity zinc ion transmembrane transporter activity	molecular_function
...	...	...

处理eggNOG结果需要用到这一文件

等eggNOG注释完后，使用另一个parse_eggNOG.py脚本处理eggNOG的结果。「参数说明」
「-i」 eggNOG的注释结果
「-g」 上一步根据obo解析出来的文件
「-O」 参考物种(只用于KEGG注释，使用KEGG三字母物种缩写表示).设置这个参数的原因是我做KEGG富集的时候发现有的基因会出现在非常荒唐的通路上，比如某个植物基因富集到了癌症的相关通路，后来发现原因是有的比较基础的KO可能与癌症通路有关，如果不使用参考物种，直接用KO去寻找map的话就会出现上述的情况。这里使用参考物种可以把没有出现在参考物种中的通路给过滤掉。植物我选择拟南芥和水稻作为参考，同样的如果做非模式动物的话，可以考虑设置一些动物物种来排除富集到植物的通路上
「-o」 输出结果文件夹。会在该文件夹生成GOannotation.tsv和KOannotation.tsv两个文件

python parse_eggNOG.py -i panax_ginseng.annotations \
                       -g go.tb \
                       -O ath,osa \
                       -o D:\RProject\MedicalPlantDB

处理的结果文件有两个：「GOannotation.tsv」和「KOannotation.tsv」 分别对应GO注释和KO注释,使用这两个文件就可以进行富集分析了

处理后的文件：GOannotation.tsv

Gene	GO	level
Pg_S3686.2	GO:0000165	biological_process
Pg_S3686.2	GO:0003674	molecular_function
...	...	...

KOannotation.tsv

Gene	KO	pathway	decription
Pg_S6540.1	K05907	map00920	Sulfur metabolism
Pg_S6540.1	K05907	map01100	Metabolic pathways
...	...	...

「ps」：上述脚本在投稿前已确认可以正常使用

3.富集分析

总的来说就是利用clusterProfiler中的enricher这个通用函数进行富集

library(clusterProfiler)
KOannotation <- read.delim('D:/RProject/MedicalPlantDB/KOannotation.tsv', stringsAsFactors=FALSE)
GOannotation <- read.delim('D:/RProject/MedicalPlantDB/GOannotation.tsv', stringsAsFactors=FALSE)
GOinfo <- read.delim('D:/RProject/MedicalPlantDB/go.tb', stringsAsFactors=FALSE)
# 前面获取gene list的过程略
gene_list<- # 你的gene list
# GO富集
## 拆分成BP，MF，CC三个数据框
GOannotation = split(GOannotation, with(GOannotation, level))
## 以MF为例
enricher（gene_list,
          TERM2GENE=GOannotation[['molecular_function']][c(2,1)],
          TERM2NAME=GOinfo[1:2]）
# KEGG富集
enricher（gene_list,
          TERM2GENE=KOannotation[c(3,1)],
          TERM2NAME=KOannotation[c(3,4)])

优缺点

优点

理论上针对所有有基因组蛋白序列的物种都可以注释，甚至没有基因组但是有参考转录本也可以注释
相对来说没有用到特别多的依赖工具，两个python脚本也只使用了最基本的包。

缺点

中间需要解析eggNOG的结果文件，我个人用的python写的脚本，一是与R衔接不连贯，二是速度较慢，主要是我考虑做一个背景数据集能用很久，所以偷懒没有对脚本的性能进行优化。
「PS」：本人使用i5-7300U，16G内存笔记本，解析go-basic.obo文件需要5~6分钟；解析eggNOG结果文件（约2w基因左右）需要几十秒。
无法一张表完成基因、转录本等多层次的富集信息。主要是因为我做的植物很多时候基因组质量不高，根本没有多个转录本的注释，所以就没考虑这一问题。

其他

其实很多人知道非模式植物可以用AnnotationHub包里面带的一些数据集，但是一是徐洲更说过这个数据集有时并不可靠。二是实际上我自己使用AnnotationForge根据eggNOG注释结果试图打包物种数据集的时候发现，物种数据集里的GO注释是依赖GO.db的，而GO.db的更新往往是滞后的。此外还有一个非常有趣的事是eggNOG原始的注释文件使用AnnotationHub打包成注释包后体积非常大，水平有限没找出原因，所以综上我选择了现在的这个方案。

本站仅提供存储服务，所有内容均由用户发布，如发现有害或侵权内容，请点击举报。