1、加载需要的包
library(DOSE)
library(GO.db)
library(org.Hs.eg.db)
library(topGO)
library(GSEABase)
library(clusterProfiler)
2、导入数据
gene<-read.csv("./data.csv",header = F)
head(gene)
V1
1 ABCA4
2 ABCA6
3 ABCC9
4 ABR
5 ACSL6
6 ADAM6
我这里的数据其实就是一列基因名,是通过其它分析挑选出来的几百个基因,想要看看这些基因是否能够富集到某些特殊的通路
3、ID转化
做GO分析是不能直接用基因名的,必须得先转化成entre id
symbol=as.character(gene[,1])
eg = bitr(symbol, fromType="SYMBOL", toType="ENTREZID", OrgDb="org.Hs.eg.db")
id = as.character(eg[,2])
head(id)
[1] "24" "23460" "10060" "29" "23305" "8755"
也可以这样,同时转换成两种:
ids <- bitr(symbol, fromType="SYMBOL", toType=c("UNIPROT", "ENSEMBL"), OrgDb="org.Hs.eg.db")
head(id)
SYMBOL UNIPROT ENSEMBL
1 ABCA4 P78363 ENSG00000198691
2 ABCA4 Q6AI28 ENSG00000198691
3 ABCA6 Q8N139 ENSG00000154262
4 ABCC9 A0A024RAV7 ENSG00000069431
5 ABCC9 O60706 ENSG00000069431
6 ABR B7Z2X0 ENSG00000159842
或者转化KEGG的id
eg2np <- bitr_kegg(symbol, fromType='kegg', toType='ncbi-proteinid', organism='hsa')
bitr_kegg("Z5100", fromType="kegg", toType='ncbi-geneid', organism='ece')
bitr_kegg("Z5100", fromType="kegg", toType='uniprot', organism='ece')
4、GO分析
这里做了BP(生物学过程)和MF(分子功能)分析
ego <- enrichGO(gene = id,
OrgDb = org.Hs.eg.db,
ont = "MF",
pAdjustMethod = "BH",
pvalueCutoff = 0.05,
qvalueCutoff = 0.05,
readable = TRUE)
做BP则改为ont="BP"即可
分析的结果大概是这样子的:
总共9列,是我手里的这些基因能够富集到的MF,接下来就是可视化
dotplot(ego)
barplot(ego, showCategory=15)
enrichMap(ego)
plotGOgraph(ego)
5、KEGG通路分析(物种为人类)
kk <- enrichKEGG(gene = id,
organism = 'hsa', #小鼠的是mmu
pvalueCutoff = 0.05,
pAdjustMethod = "BH",
qvalueCutoff = 0.05)
结果如下:
可视化:
dotplot(kk)
另外,如果想看看富集到的通路图可以直接用:
browseKEGG(kk, 'hsa04974')
会打开浏览器调到KEGG数据库的这条通路上,并且富集的基因会以不同的颜色标出
6、DO(disease ontology)分析
也可以看看这些基因是否与某些疾病相关
x <- enrichDO(gene = id,
ont = "DO",
pvalueCutoff = 0.05,
pAdjustMethod = "BH",
minGSSize = 1,
maxGSSize = 500,
qvalueCutoff = 1,
readable = TRUE)
7、其它
如果还是觉得DAVID更可靠一点,也不用在网页上用,也可以直接使用这个包,但是要先去注册一个邮箱
david <- enrichDAVID(gene = gene,
idType = "ENTREZ_GENE_ID",
listType = "Gene",
annotation = "KEGG_PATHWAY",
david.user = "clusterProfiler@hku.hk")