## load required packages
library(Seurat)
library(cowplot)
library(dplyr)
library(ggplot2)
library(DT)
library(paletteer)
library(forcats)
  • Cluster marker genes GSEA using Gene Ontology
    • BP
    • MF
    • CC

GSEA by cluster marker genes

Each tab includes the followings:

  • GSEA (BP)
  • GSEA (MF)
  • GSEA (CC)
  • GSEA (data table)

GSEA function

library(clusterProfiler)
library(org.Hs.eg.db)

gene.ego = function(cluster_number) {
  geneset = mks %>%
    filter(cluster == !!cluster_number) %>%
    dplyr::select(gene) %>% 
    pull()

  genes_to_convert <- tryCatch({
    bitr(geneset, fromType = "SYMBOL", toType = "ENTREZID", OrgDb = "org.Hs.eg.db")
  }, error = function(e) {
    message("Error converting gene symbols: ", e)
    return(NULL)
  })

  if (is.null(genes_to_convert) || nrow(genes_to_convert) == 0) {
    message("No valid gene IDs were found.")
    return(NULL)
  }

  gene.ego <- enrichGO(gene     = genes_to_convert$ENTREZID, 
                       ont = 'ALL',
                       OrgDb = "org.Hs.eg.db",
                       pvalueCutoff = 0.01, qvalueCutoff = 0.01, readable = TRUE)
  return(gene.ego)
}
custom_barplot = function(df, ont, n){
  df = df@result 
  p= df %>%
    filter(ONTOLOGY == ont) %>%
    slice_max(order_by = Count, n = n) %>%
    mutate(Description = forcats::fct_reorder(Description, -Count, .desc = TRUE)) %>%
    ggplot(aes(x = Count, y = Description, fill = -qvalue, size=Count)) +
    geom_point(shape = 21) +
    scale_size_area(max_size = 5) +
    scale_fill_gradient(low = "ivory", high = "salmon") +
    labs(fill = "q-value(-)") +
    theme_minimal() +
    theme(plot.background = element_rect(color = "black", fill = NA, size = 0.2),
          panel.border = element_rect(color = "black", fill = NA, size = 0.2)) +
    ylab("") +
    ggtitle(ont)
  return(p)
}

Cluster 0

Run GSEA

gene.ego.out = gene.ego(cluster_number = 0)
gene.ego.out0 = gene.ego.out

GSEA plot

custom_barplot(df=gene.ego.out, ont="BP", n=20)

custom_barplot(df=gene.ego.out, ont="MF", n=20)

custom_barplot(df=gene.ego.out, ont="CC", n=20)

GSEA info table

gene.ego.out@result %>% DT::datatable(extensions = "Buttons", 
                                      width="800px",
                options = list(scrollX = TRUE,
                               dom="Bfrtip", buttons=c("csv","excel")))

Cluster 1

Run GSEA

gene.ego.out = gene.ego(cluster_number = 1)
gene.ego.out1= gene.ego.out

GSEA plot

custom_barplot(df=gene.ego.out, ont="BP", n=20)

custom_barplot(df=gene.ego.out, ont="MF", n=20)

custom_barplot(df=gene.ego.out, ont="CC", n=20)

GSEA info table

gene.ego.out@result %>% DT::datatable(extensions = "Buttons", 
                                      width="800px",
                options = list(scrollX = TRUE,
                               dom="Bfrtip", buttons=c("csv","excel")))

Cluster 2

Run GSEA

gene.ego.out = gene.ego(cluster_number = 2)
gene.ego.out1= gene.ego.out

GSEA plot

custom_barplot(df=gene.ego.out, ont="BP", n=20)

custom_barplot(df=gene.ego.out, ont="MF", n=20)

custom_barplot(df=gene.ego.out, ont="CC", n=20)

GSEA info table

gene.ego.out@result %>% DT::datatable(extensions = "Buttons", 
                                      width="800px",
                options = list(scrollX = TRUE,
                               dom="Bfrtip", buttons=c("csv","excel")))