富集分析

clusterProfiler
https://yulab-smu.top/biomedical-knowledge-mining-book/index.html
https://www.jianshu.com/p/d04f3bb09e16
https://www.jianshu.com/p/bb7442bb8cad

GO enrichment analysis

rm(list = ls())
library("clusterProfiler")
library("org.Hs.eg.db")
library("enrichplot")
library("ggplot2")
rt=read.csv("nrDEG_edgeR_signif.csv")
## ID转换
gene <- bitr(rt$ID, 
             fromType="SYMBOL", toType="ENTREZID", 
             OrgDb="org.Hs.eg.db")
## GO 
gene <- gene$ENTREZID
ego <- enrichGO(gene = gene,
                OrgDb = org.Hs.eg.db,
                pAdjustMethod = "BH",
                pvalueCutoff =0.05, 
                qvalueCutoff = 0.05,
                ont="BP",
                readable =T)
head(ego)
bp <- data.frame(ego)
write.csv(ego, file="GO-bp.csv", quote=F, row.names = F)

KEGG enrichment analysis

R.utils::setOption( "clusterProfiler.download.method",'auto' )
kk <- enrichKEGG(gene = gene, 
                 organism = "hsa", 
                 pvalueCutoff =0.05, 
                 qvalueCutoff =0.05)
head(kk)
kegg <- data.frame(kk)
write.csv(kk, file="KEGG.csv", row.names = F)

Gene Set Enrichment Analysis

GSEA()函数

hallmarks <- read.gmt("KEGG gene sets as Gene Symbols.gmt")
gsea <- GSEA(geneList, 
             minGSSize = 10,
             maxGSSize = 500,
             eps = 1e-10,
             pvalueCutoff = 0.05,
             pAdjustMethod = "BH",
             TERM2GENE =hallmarks,
             by = "fgsea")
res <- data.frame(gsea)
gseaplot2(res,"KEGG_PROGESTERONE_MEDIATED_OOCYTE_MATURATION",
          color = "#0072b5", pvalue_table = T)

GO GSEA

rm(list = ls())
library("clusterProfiler")
library("org.Hs.eg.db")
library("enrichplot")
library("ggplot2")

rt=read.csv("nrDEG_edgeR_signif.csv")

## ID转换
gene <- bitr(rt$ID, 
             fromType="SYMBOL", toType="ENTREZID", 
             OrgDb="org.Hs.eg.db")

## 去重合并
gene <- dplyr::distinct(gene, SYMBOL, .keep_all=TRUE)
gene_df <- data.frame(logFC=rt$logFC, SYMBOL = rt$ID)
gene_df <- merge(gene_df, gene, by="SYMBOL")

## 排序很重要
geneList <- gene_df$logFC
names(geneList) = gene_df$ENTREZID
geneList = sort(geneList, decreasing = TRUE)

## GO GSEA
ego3 <- gseGO(geneList     = geneList,
              OrgDb        = org.Hs.eg.db,
              ont          = "BP",
              minGSSize    = 100,
              maxGSSize    = 500,
              pvalueCutoff = 0.05,
              verbose      = FALSE)

head(ego3)
write.csv(ego3, file="gseGObp.csv", quote=F, row.names = F)

KEGG GSEA

rm(list = ls())
library("clusterProfiler")
library("org.Hs.eg.db")
library("enrichplot")
library("ggplot2")

rt=read.csv("nrDEG_edgeR_signif.csv")

## ID转换
gene <- bitr(rt$ID, 
             fromType="SYMBOL", toType="ENTREZID", 
             OrgDb="org.Hs.eg.db")

## 去重合并
gene <- dplyr::distinct(gene, SYMBOL, .keep_all=TRUE)
gene_df <- data.frame(logFC=rt$logFC, SYMBOL = rt$ID)
gene_df <- merge(gene_df, gene, by="SYMBOL")

## 排序很重要
geneList <- gene_df$logFC
names(geneList) = gene_df$ENTREZID
geneList = sort(geneList, decreasing = TRUE)

## KEGG GSEA
kk2 <- gseKEGG(geneList     = geneList,
               organism     = 'hsa',
               minGSSize    = 120,
               pvalueCutoff = 0.05,
               verbose      = FALSE)
head(kk2)
write.csv(kk2, file="gseaKEGG.csv", row.names = F)