差异分析后可以得到基因表达的logFC值,p值。设定差异基因的筛选标准后,选取基因,可以进一步进行GO分析,KEGG 分析,网络分析。主要使用clusterprofiler package进行分析。
加载相关的包
library(DOSE)
library(GO.db)
library(org.Hs.eg.db)
library(GSEABase)
library(clusterProfiler)
其实加载package, 主要的是clusterprofiler
基因名称的转换
基因名称的转换可以使用多种方法,在线工具可以使用DAVID工具,此处仍然使用 Y叔的clusterprofiler package,进行基因名称转换。clusterprofiler有相关的示例。
gene = bitr(gene, fromType="SYMBOL", toType="ENTREZID", OrgDb="org.Hs.eg.db")
head(gene)
GO 分析
GO 分析有三种结果,细胞组分CC,生物过程BP, 分子功能 MF
先展示分开分析三个过程的结果
ego_CC <- enrichGO(gene = gene$ENTREZID,
OrgDb= org.Hs.eg.db,
ont = "CC",
pAdjustMethod = "BH",
minGSSize = 1,
pvalueCutoff = 0.01,
qvalueCutoff = 0.01,
readable = TRUE)
ego_BP <- enrichGO(gene = gene$ENTREZID,
OrgDb= org.Hs.eg.db,
ont = "BP",
pAdjustMethod = "BH",
minGSSize = 1,
pvalueCutoff = 0.01,
qvalueCutoff = 0.01,
readable = TRUE)
ego_MF <- enrichGO(gene = gene$ENTREZID,
OrgDb= org.Hs.eg.db,
ont = "MF",
pAdjustMethod = "BH",
minGSSize = 1,
pvalueCutoff = 0.01,
qvalueCutoff = 0.01,
readable = TRUE)
barplot(ego_CC, showCategory=20,title="EnrichmentGO_CC")#条状图,按p从小到大排的
dotplot(ego_BP,title="EnrichmentGO_BP_dot")#点图,按富集的数从大到小的。
将以上的三个结果统一绘制于一幅图中,并保存结果。
go <- enrichGO(gene = de, OrgDb = "org.Hs.eg.db", ont="all")
library(ggplot2)
p <- dotplot(go, split="ONTOLOGY") +facet_grid(ONTOLOGY~., scale="free")
p
x <- go
y <- setReadable(x, OrgDb = org.Hs.eg.db, keyType="ENTREZID")
write.csv(y@result,file = 'commonDiff-enrichment-readabla.csv')
KEGG分析
kk <- enrichKEGG(gene = gene$ENTREZID,
organism ="human",
pvalueCutoff = 0.01,
qvalueCutoff = 0.01,
minGSSize = 1,
#readable = TRUE ,
use_internal_data =FALSE)
barplot(kk)
dotplot(kk)
保存GO分析结果和KEGG 分析结果
EGG <- enrichKEGG(gene= gene$ENTREZID,
organism = 'hsa',
pvalueCutoff = 0.05)
x <- EGG
yEGG <- setReadable(x, OrgDb = org.Hs.eg.db, keyType="ENTREZID")
write.csv(yEGG@result,file = 'commonDiff-KEGG-readable.csv')
dotplot(EGG)
test <- data.frame(EGG)
browseKEGG(EGG, 'hsa04110')
save(go,EGG, file = 'commonDiff-enrichment-analysis.Rdata')
write.csv(go@result, file = "commonDiff-go.csv")
write.csv(EGG@result, file = 'commonDiff-KEGG.csv')
其实按照简单的生物信息分析思路,这里已经有了雏形,有了差异分析,有了富集分析。差异分析可以是按照设定的cutoff值进行筛选,也可以进行GSEA分析。富集分析中得到基因的分布。其实对于一个临床医生来说。这还远远不够,有了这些富集分析,那怎么进行下一步的分析,得到更有价值的信息或者分析,这应该是临床医生做生信分析因该考虑的。