GEO数据挖掘系列,第15篇学习笔记:学习富集分析基础概念及运行代码
为了无缝衔接上一篇学习笔记,该文中的序号将接着上一篇来标注
7.富集分析
7.1entrezID
使用的数据:差异基因的
ENTREZ
ID
常说的基因名用symbol表示,富集分析指定用ENTREZID;两者的转换并非一一对应,损失/增加部分基因是正常现象
图1
7.2富集分析背景数据库
7.2.1 KEGG数据库
做富集分析需要背景数据库,常见的有GO和KEGG
图2
KEGG里面有很多通路如信号通路、代谢通路等等,在富集分析时使用KEGG就相当于将基因进行归类。例如哪些基因负责DNA复制,哪些基因又负责细胞循环……将这些基因归类成若干个类别,每个类别都有一个编号如hsa03030,这个类别或者说这条通路里面的基因则都是负责DNA复制的
7.2.2 GO数据库
GO数据库能够用来找细胞的功能
图3
7.3富集分析结果
如图4是进行富集分析后得到的表格,以及其对应的含义
图4
总结图4,可以简单概括富集分析到底是在做什么:
图5
用图4中的
bgRatio
、
geneRatio
来判断是否富集得足够多
举例理解:
假设数据库中收录了1000个基因,其中的50个是通路1的,那么遇到通路1的基因的概率是5%
当随机取200个,200×5%=10,取到10个左右则可以看
成
这10个基因确实是偶然筛序进来的
但如果筛到了通路1的基因25个,这比例就比5%大得多了,那么就可以理解为通路1确实对自己的研究起作用
这就叫做超几何分布检验
图6
7.4富集分析代码
以下代码将实现用GO和KEGG做的富集分析,输入数据即上一篇笔记中得到的deg
rm(list = ls())
load(file = 'step4output.Rdata')
library(clusterProfiler)# Y叔的R包
library(ggthemes)
library(org.Hs.eg.db)
library(dplyr)
library(ggplot2)
library(stringr)
library(enrichplot)
#(1)输入数据
gene_diff = deg$ENTREZID[deg$change != "stable"] # 提取差异基因
#(2)富集
# KEGG
ekk organism = 'hsa')# 如果用魔法这一步容易联网失败
ekk OrgDb = org.Hs.eg.db,
keyType = "ENTREZID")# 这个函数将ENTREZID转为SYMBOL
# GO
ego OrgDb= org.Hs.eg.db,
ont = "ALL",# 细胞组分、分子功能、生物过程全要,更改参数可只要其中一个,例如把ALL改成BP
readable = TRUE)
#setReadable和readable = TRUE都是把富集结果表格里的基因名称转为symbol
class(ekk)
图7
#(3)可视化
# KEGG
barplot(ekk)
#或者是dotplot
barplot()
横坐标是count数,纵坐标是通路名字/介绍,颜色根据p.adjust值上色
图8
dotplot()
可显示geneRatio比例数
图9
看自己研究的通路在两张图中的哪个图位置比较靠上,就可以用哪张;但其实也无所谓,两张都可以用
# GO
barplot(ego, split = "ONTOLOGY") +
facet_grid(ONTOLOGY ~ ., space = "free_y",scales = "free_y") #ggplot中画分面的函数
# ONTOLOGY ~ .是横着切,切成图10
# ~ ONTOLOGY 是竖着切