Contents

1 Overview

Gene Set Enrichment Analysis (GSEA) is a very powerful and interesting computational method that allows an easy correlation between differential expressed genes and biological processes. Unfortunately, although it was designed to help researchers to interpret gene expression data it can generate huge amounts of results whose biological meaning can be difficult to interpret.

Many available tools rely on the hierarchically structured Gene Ontology (GO) classification to reduce reundandcy in the results. However, due to the popularity of GSEA many more gene set collections, such as those in the Molecular Signatures Database (MSigDB), are emerging. Since these collections are not organized as those in GO, their usage for GSEA do not always give a straightforward answer or, in other words, getting all the meaninful information can be challenging with the currently available tools. For these reasons, GSEAmining was born to be an easy tool to create reproducible reports to help researchers make biological sense of GSEA outputs.

Given the results of GSEA, GSEAmining clusters the different gene sets collections based on the presence of the same genes in the leadind edge (core) subset. Leading edge subsets are those genes that contribute most to the enrichment score of each collection of genes or gene sets. For this reason, gene sets that participate in similar biological processes should share genes in common and in turn cluster together. After that, GSEAmining is able to identify and represent for each cluster:

In each case, positive and negative enrichments are shown in different colors so it is easy to distinguish biological processes or genes that may be of interest in that particular study.

2 Installation

You can install GSEAmining from Bioconductor:

if (!requireNamespace("BiocManager", quietly = TRUE))
    install.packages("BiocManager")

BiocManager::install("GSEAmining")

Or directly from GitHub:

install.packages("devtools") # If you have not installed "devtools" package
library(devtools)
devtools::install_github("oriolarques/GSEAmining")

3 Input data format: gm_filter

By default GSEAmining is designed to accept the resuls of the GSEA function from the clusterProfiler package. For more information about this function click here. However, it is not mandatory to use this package.

In this example the data corresponds to GSEA analysis of differential expressed genes from treated versus control samples in HGPalmer-PDX-P30 experiment. Differential gene expression (tableTop) was obtained using the oligo and limma R packages.

# A geneList contains three features:
# 1.numeric vector: fold change or other type of numerical variable
# 2.named vector: every number has a name, the corresponding gene ID
# 3.sorted vector: number should be sorted in decreasing order
tableTop_p30 <- as.data.frame(tableTop_p30)
geneList = tableTop_p30[,2]
names(geneList) = as.character(tableTop_p30[,1])
library(clusterProfiler)
# Read the .gmt file from MSigDB
gmtC2<- read.gmt("gmt files/c2.all.v7.1.symbols.gmt")
gmtC5<- read.gmt('gmt files/c5.all.v7.1.symbols.gmt')
gmtHALL <- read.gmt('gmt files/h.all.v7.1.symbols.gmt')

# Merge all the gene sets
gmt_all <- rbind(gmtC2, gmtC5, gmtHALL)
GSEA_p30<-GSEA(geneList, TERM2GENE = gmt_all, nPerm = 1000, pvalueCutoff = 0.5)

# Selection of gene sets with a specific thershold in terms of NES and p.adjust
genesets_sel <- GSEA_p30@result

Data should be in a data.frame with, at least three columns labelled as follows:

# Structure of the data included in the package
data('genesets_sel', package = 'GSEAmining')
tibble::glimpse(genesets_sel)
#> Rows: 52
#> Columns: 4
#> $ ID              <chr> "WATANABE_COLON_CANCER_MSI_VS_MSS_DN", "GO_RNA_BINDIN…
#> $ NES             <dbl> 2.511796, 2.234696, 2.200913, 2.167304, 2.127964, 2.1…
#> $ p.adjust        <dbl> 0.02463671, 0.02716912, 0.03017176, 0.02612857, 0.022…
#> $ core_enrichment <chr> "DACH1/DUOX2/GRM8/CDHR1/SPINK1/CEL/CYP2B6/C10orf99/SL…

gm_filter: Filter the input data

GSEA outputs normally presents tens or hundreds of genesets but many of them may not meet the thresholds for considering them significantly enriched. For that it is better to filter the data for a better visualization.

The gm_filter() function of GSEAmining allows this process to be very easy, especially if the format of your data meets the aforementioned criteria.

library(GSEAmining)
data("genesets_sel", package = 'GSEAmining')
gs.filt <- gm_filter(genesets_sel, 
                     p.adj = 0.05, 
                     neg_NES = 2.6, 
                     pos_NES = 2)

4 Clustering

4.1 gm_clust: Creation of a gm_clust object

Using the gm_clust() function, we can create an object that will contain the hierarchical clustering of the gene sets according to their genes in the core_enrichment column. This function accepts the data frame created with gm_filter. In the process, a distance matrix is calculated using the binary method (from the dist() function in stats) and then a cluster with the complete method (from hclust() function in stats) is created.

# Create an object that will contain the cluster of gene sets.
gs.cl <- gm_clust(gs.filt)

4.2 gm_dendplot: Plot the cluster

The gm_dendplot() function uses (i) the filtered data frame and (ii) the gm_clust object to plot the cluster dendrogram. It shows which gene sets are positively or negatively enriched coloring them in red or blue respectively. Additionally by default the function will highlight every other cluster with a rectangle to facilitate the differentiation of all the clusters.

gm_dendplot(gs.filt, 
            gs.cl)

It is possible to tune the height of the dendrogram, the colors for the enrichment and the use of rectangles (and their sizes) for highlighting clusters.

gm_dendplot(gs.filt, 
            gs.cl, 
            col_pos = 'orange', 
            col_neg = 'black', 
            rect = TRUE,
            dend_len = 20, 
            rect_len = 2)