其他
从基因名到GO注释一步到位
大部分的生物学高通量数据处理后都是得到基因集，不管是上调下调表达基因集，还是共表达的模块基因集，都是需要注释到生物学功能数据库来看基因集的意义，最常见的是GO/KEGG数据库啦，还有很多其它在MsigDB的，比如reactome和biocarta数据库等等。
library(stringr)
library(clusterProfiler)
# 我这里演示的是brown_down_gene，是WGCNA的一个模块，基因集
# 因为表达矩阵是symbol，所以需要转为ENTREZID，才能走clusterProfiler的函数。
gene.df <- bitr(brown_down_gene$symbol, fromType="SYMBOL",
toType="ENTREZID",
OrgDb = "org.Hs.eg.db")
go <- enrichGO(gene = gene.df$ENTREZID, OrgDb = "org.Hs.eg.db", ont="all")
barplot(go, split="ONTOLOGY")+ facet_grid(ONTOLOGY~., scale="free")
如果你想分开计算上下调基因的GO数据库注释
library(stringr)
library(clusterProfiler)
# 通过前面的差异分析，我们拿到了 gene_up 和 gene_down 这两个基因集
# 后面的分析，只需要 gene_up 和 gene_down 这两个变量即可
go_up <- enrichGO(gene_up,
OrgDb = "org.Hs.eg.db",
ont="all",
pvalueCutoff = 0.9,
qvalueCutoff =0.9)
go_up=DOSE::setReadable(go_up, OrgDb='org.Hs.eg.db',keyType='ENTREZID')
write.csv(go_up@result,paste0(pro,'_go_down.up.csv'))
barplot(go_up, split="ONTOLOGY",font.size =10)+
facet_grid(ONTOLOGY~., scale="free") +
scale_x_discrete(labels=function(x) str_wrap(x, width=50))+
ggsave(paste0(pro,'gene_up_GO_all_barplot.png'))
go_down <- enrichGO(gene_down,
OrgDb = "org.Hs.eg.db",
ont="all",
pvalueCutoff = 0.9,
qvalueCutoff =0.9)
go_down=DOSE::setReadable(go_down, OrgDb='org.Hs.eg.db',keyType='ENTREZID')
write.csv(go_down@result,paste0(pro,'_go_down.up.csv'))
barplot(go_down, split="ONTOLOGY",font.size =10)+
facet_grid(ONTOLOGY~., scale="free") +
scale_x_discrete(labels=function(x) str_wrap(x, width=50))+
ggsave(paste0(pro,'gene_down_GO_all_barplot.png'))
多组基因集的KEGG数据库富集
xx.formula <- compareCluster(ENTREZID~new, data=DEG, fun='enrichKEGG')
dotplot(xx.formula, x=~GeneRatio) + facet_grid(~new)
如果是多组基因集走GO数据库富集
de_gene_clusters$cluster)
# Run full GO enrichment test
formula_res <- compareCluster(
ENTREZID~cluster,
data=de_gene_clusters,
fun="enrichGO",
OrgDb="org.Mm.eg.db",
ont = "BP",
pAdjustMethod = "BH",
pvalueCutoff = 0.01,
qvalueCutoff = 0.05
)
# Run GO enrichment test and merge terms
# that are close to each other to remove result redundancy
lineage1_ego <- simplify(
formula_res,
cutoff=0.5,
by="p.adjust",
select_fun=min
)
感兴趣的可以把这个结果跟3个出名的网页工具进行比较
https://amp.pharm.mssm.edu/Enrichr/ http://www.webgestalt.org/ https://biit.cs.ut.ee/gprofiler
另外，强推Y叔clusterProfiler的一些可视化方法
barplot cnetplot dotplot emapplot gseaplot goplot upsetplot
文末友情宣传
生信爆款入门-全球听（买一得五）（第4期），你的生物信息学入门课 数据挖掘第2期（两天变三周，实力加量），医学生/临床医师首选技能提高课 生信技能树的2019年终总结 ，你的生物信息学成长宝藏 2020学习主旋律，B站74小时免费教学视频为你领路，还等什么，看啊！！！