R包clusterProfiler比较不同dataset富集结果
2018-10-30 本文已影响53人
214b3ff96d82
生物信息分析中,接触最多的莫过于基因富集分析,故在此基础上目前已经开发了很多种富集分析工具,如网页版的DAVID、KOBAS、GOEAST,WebGestalt,基迪奥平台等,本地版工具如TBtools,此外常用的R包有pathfindR、topGO、clusterProfiler等。在最后的报告中,我们通常会以各种图表的形式来展示富集的结果,常用的富集分析结果可视化的软件有REVIGO,GOView,WEGO,cytoscape插件BinGO,基迪奥富集分析动态展示等。
这些工具各有千秋,可是依然具有一定的局限性,就是完成分析后需要转换数据才能进行可视化。Y叔开发的clusterProfiler既可以轻松完成各种富集分析又可以傻瓜式出图,所以一直受到生信工作者的青睐。本文主要就是介绍这个R包,通过compareCluster函数完成不同数据集的pathway富集结果的比较。
library(org.Hs.eg.db)
library(clusterProfiler)
library(ggplot2)
setwd("C:/Users/lenovo/Desktop")
a=read.table("testgeneid.txt",header = FALSE) # 读取输入文件,主要是基因名 gene symbol
gene=as.character(a[,1]) # 转换为字符
ego <- enrichGO(gene=gene,OrgDb='org.Hs.eg.db',keyType='SYMBOL',ont= "CC",pAdjustMethod="BH",pvalueCutoff = 0.01,qvalueCutoff = 0.05) #必须指定keytype,GO富集
ego <- enrichGO(gene=gene,OrgDb='org.Hs.eg.db',keyType='SYMBOL',ont= "ALL",pAdjustMethod="BH",pvalueCutoff = 0.01,qvalueCutoff = 0.05)
x <- enrichDO(gene = gene,
ont = "DO",
pvalueCutoff = 0.01,
pAdjustMethod = "BH",
universe = names(geneList),
minGSSize = 5,
maxGSSize = 500,
qvalueCutoff = 0.05,
readable = FALSE) # DO富集
write.csv(summary(ekk),"KEGG-enrich.csv",row.names =F)
write.csv(summary(ego),"GO-enrich.csv",row.names =F)
barplot(ego, showCategory=15,title="EnrichmentGO") #条形图
dotplot(ego,title="EnrichmentGO_dot") #气泡图
enrichMap(ego, vertex.label.cex=1.2, layout=igraph::layout.kamada.kawai) #网络图
plotGOgraph(ego) #go图
gene= bitr(gene, fromType="SYMBOL", toType="ENTREZID", OrgDb="org.Hs.eg.db")
kk <- enrichKEGG(gene$ENTREZID, organism="hsa",keyType = "kegg",pvalueCutoff=0.05, pAdjustMethod="BH",qvalueCutoff=0.1)
接下来就是使用compareCluster比较两组基因富集结果,输入文件需要转换为gene id形式。
cp = list(a.gene=gene1$ENTREZID, b.gene=gene2$ENTREZID) # 合并两个数据集,并转换为列表
xx <- compareCluster(cp, fun="enrichKEGG", organism="hsa", pvalueCutoff=0.01,pAdjustMethod = "BH",qvalueCutoff = 0.05)
dotplot(xx,showCategory=10,includeAll=TRUE)
还可以改变形状
p1 + aes(shape = GeneRatio > 0.3)
还可以改变颜色
xx <- compareCluster(cp,
fun="enrichGO",
OrgDb="org.Hs.eg.db",
ont= "BP",
pvalueCutoff=0.01,
pAdjustMethod = "BH",
qvalueCutoff = 0.05)
p2 <- p1 + scale_color_continuous(low="purple",high = "green")