2022-06-14eggNOG注释及cluster Profi

eggNOG网站注释蛋白序列得到文件query_seqs.fa.emapper.annotations

python环境

#git clone https://github.com/Hua-CM/HuaSmallTools.git下载注释包(parse_go_obofile.py 和 parse_eggNOG.py)

#wget http://purl.obolibrary.org/obo/go/go-basic.obo下载Go数据库的obo文件

输入数据准备

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

$python parse_go_obofile.py -i go-basic.obo -o go.tbold

$less go.tbold |awk 'BEGIN{FS="\t";OFS="\t"}{print $2,$1,$3}' > tmp && mv tmp go.tb (必不可少的一步，不转化parse_eggNOG.py会报错, go.tb为parse_eggNOG.py用临时文件，富集时还是用go.tbold文件)

$python parse_eggNOG.py -i query_seqs.fa.emapper.annotations \

                      -g go.tb \

                      -O ath,osa \

                      -o ./  #生成GO和KO数据库

--------------------------参数说明----------------------------------------------------------

「-i」 eggNOG的注释结果

「-g」 上一步根据obo解析出来的文件

「-O」 参考物种(只用于KEGG注释，使用KEGG三字母物种缩写表示).设置这个参数的原因是我做KEGG富集的时候发现有的基因会出现在非常荒唐的通路上，比如某个植物基因富集到了癌症的相关通路，后来发现原因是有的比较基础的KO可能与癌症通路有关，如果不使用参考物种，直接用KO去寻找map的话就会出现上述的情况。这里使用参考物种可以把没有出现在参考物种中的通路给过滤掉。植物我选择拟南芥和水稻作为参考，同样的如果做非模式动物的话，可以考虑设置一些动物物种来排除富集到植物的通路上

「-o」 输出结果文件夹。会在该文件夹生成GOannotation.tsv和KOannotation.tsv两个文件

--------------------------------------------------------------------------------------------------

#(在界面直接跑python parse_eggNOG.py -i query_seqs.fa.emapper.annotations -g go.tb -O ath,osa -o ./) 处理的结果文件有两个：「GOannotation.tsv」和「KOannotation.tsv」 分别对应GO注释和KO注释,使用这两个文件就可以进行富集分析了

富集分析

#R环境

library(clusterProfiler)

KOannotation <- read.delim("KOannotation.tsv", stringsAsFactors=FALSE)

GOannotation <- read.delim("GOannotation.tsv", stringsAsFactors=FALSE)

GOinfo <- read.delim("go.tbold", stringsAsFactors=FALSE)

#根据具体比较组来，上面步骤不需要重做。

gene_list <- read.csv("test.txt",header=FALSE) ###提供差异基因数据集，不要表头

gene_list <- gene_list$V1 ####提取基因所在列，如gene_list <- gene_list$Row.names，根据实际情况提取；此处为类型转换，因为enricher识别的gene_list为facter，charater类型，而不识别data.drame；转化为character类型的话#gene_list <- as.character(gene_list$V

Go富集

GOannotation = split(GOannotation, with(GOannotation, level))

df_GO <- enricher(gene_list,TERM2GENE=GOannotation[['MF']][c(2,1)],TERM2NAME=GOinfo[1:2])

df_GO_df <- as.data.frame(df_GO) ###易读模式，查看GO富集结果

dotplot(df_GO)

####柱状图用barplot(df_GO, showCategory = 10)，10表示选择前10个term###

###做有向无环图plotGOgraph(df_GO)###要先安装topGo####

####输出表格write.csv(df_GO, "0_up1.csv")###

df_GO <- enricher(gene_list,TERM2GENE=GOannotation[['CC']][c(2,1)],TERM2NAME=GOinfo[1:2])

dotplot(df_GO)

df_GO <- enricher(gene_list,TERM2GENE=GOannotation[['BP']][c(2,1)],TERM2NAME=GOinfo[1:2])

dotplot(df_GO)

KEGG富集

enricher(gene_list,TERM2GENE=KOannotation[c(3,1)],TERM2NAME=KOannotation[c(3,4)])

df_KO <- enricher(gene_list,TERM2GENE=KOannotation[c(3,1)],TERM2NAME=KOannotation[c(3,4)])

dotplot(df_KO)

优缺点

优点

理论上针对所有有基因组蛋白序列的物种都可以注释，甚至没有基因组但是有参考转录本也可以注释

相对来说没有用到特别多的依赖工具，两个python脚本也只使用了最基本的包。

缺点

中间需要解析eggNOG的结果文件，我个人用的python写的脚本，一是与R衔接不连贯，二是速度较慢，主要是我考虑做一个背景数据集能用很久，所以偷懒没有对脚本的性能进行优化。

无法一张表完成基因、转录本等多层次的富集信息。主要是因为我做的植物很多时候基因组质量不高，根本没有多个转录本的注释，所以就没考虑这一问题。

引用：https://mp.weixin.qq.com/s/Mr3YLoc_-Y1WeLKJku1TzQ