clusterProfiler包做GO分析

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"即可
分析的结果大概是这样子的：

EGO.PNG

总共9列，是我手里的这些基因能够富集到的MF，接下来就是可视化

dotplot(ego)

dotplot.png

barplot(ego, showCategory=15)

barplot.png

enrichMap(ego)

map.jpg

plotGOgraph(ego)

tree.png

5、KEGG通路分析(物种为人类)

kk <- enrichKEGG(gene = id,
                 organism = 'hsa',
                 pvalueCutoff = 0.05,
                 pAdjustMethod = "BH",
                 qvalueCutoff = 0.05)

结果如下：

kk.PNG

可视化：

dotplot(kk)

另外，如果想看看富集到的通路图可以直接用：

browseKEGG(kk, 'hsa04974')

会打开浏览器调到KEGG数据库的这条通路上，并且富集的基因会以不同的颜色标出

hsa04974.png

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")