knitr::opts_chunk$set(
  collapse = TRUE,
  comment = "#>"
)

This R-package aims to perform metabolite set enrichment analysis (MSEA) on single-cell metabolomics datasets. In contrast to bulk-metabolomics, metabolite annotation is often more ambiguous with fully resolved molecular structures. That means, annotations are vectors of isomeric (and/or isobaric) molecules, complicating downstream MSEA. This package uses a boostrapping approach by performing enrichment analyses many times with random sampling of the isomers/isobars.

options(stringsAsFactors = FALSE, warn = -1)

## install devtools if not installed
if(!("devtools" %in% rownames(installed.packages()))){
  install.packages("devtools",  repos = c(CRAN = "http://cran.rstudio.com"))
  }

## install bmetenrichr if not installed
if(!("bmetenrichr" %in% rownames(installed.packages()))){
   devtools::install_github(repo = "martijnmolenaar/bmetenrichr", build_vignettes = TRUE)
}

library(bmetenrichr)

The package contains example data from Rappez et al., 2021 (https://doi.org/10.1038/s41592-021-01198-0).

data("Rappez_et_al")

## the main input is a single-cell metabolomics matrix with molecules as rows and cells as columns
Rappez_et_al$sc_matrix[1:10,1:10]

## for the molecules, a vector with molecular formulas plus adduct is required
rownames(Rappez_et_al$sc_matrix)[1:10]

## a conditions vector is required to define to which condition a given cell belongs
Rappez_et_al$conditions[1:10]

## in this analysis, only specific annotations should be included as others are extracellular 
Rappez_et_al$cellular[1:10]

enrichment analysis with isomers

The main enrichment object can be generated as follows. By default, bmetenrichr uses LION (Molenaar et al., 2019, GigaScience, https://doi.org/10.1093/gigascience/giz061) as metabolite sets.

myTestRun <-
  initEnrichment(scmatrix = Rappez_et_al$sc_matrix,
                 annotations = rownames(Rappez_et_al$sc_matrix),
                 conditions = Rappez_et_al$conditions,
                 include = Rappez_et_al$cellular,
                 condition.x = "U",
                 condition.y = "F"                    )

First, metabolites are ranked by rankScore()

## rank metabolites, in this case by t.test statistic

myTestRun <- rankScore(myTestRun, ranking.by = 't.test')

Then, perform enrichment analysis with n = 100 bootstraps

myTestRun <- calcEnrichment(myTestRun, n = 100)

The enrichment analysis can be visualized with plotEnrichment(), here with enrichment score (ES) on x-axis

plotEnrichment(myTestRun, min.annotations = 5, q.value.cutoff = .1, by.statistic = "ES")

Plots can also be arranged with q.values on x-axis, and with LION IDs

plotEnrichment(myTestRun, min.annotations = 5, q.value.cutoff = .05, plotIDs = T, 
               by.statistic = "q.value")

The enrichment analysis can also be exported as data.frame:

enrichmentTable(myTestRun)[1:10,]

To compare other conditions, use setConditions():

## now, let's test FIT vs F

myTestRun <-  setConditions(object = myTestRun, condition.x = 'F', condition.y = 'FIT')
myTestRun

## rank metabolites, in this case by t.test statistic

myTestRun <- rankScore(myTestRun, ranking.by = 't.test')

## and perform enrichment analysis

myTestRun <- calcEnrichment(myTestRun, n = 100)

plotEnrichment(myTestRun, min.annotations = 5, q.value.cutoff = .05)

enrichment analysis with isomers and isobars

By default, only isomers are included. With isobars = TRUE, it's also possible to include isobars within a set m/z range:

## create object

myTestRun <-
  initEnrichment(scmatrix = Rappez_et_al$sc_matrix, 
                 isobars = TRUE,                      ## to include isobars (default is FALSE)
                 mass_range_ppm = 3,                  ## mass range to define isobars
                 polarization_mode = "positive",      ## mode is important to include the right adducts
                 annotations = rownames(Rappez_et_al$sc_matrix), 
                 conditions = Rappez_et_al$conditions,
                 include = Rappez_et_al$cellular,
                 condition.x = "U",
                 condition.y = "F"                    )

Downstream, the same workflow can be used:

## rank metabolites, in this case by t.test statistic

myTestRun <- rankScore(myTestRun, ranking.by = 't.test')

## perform enrichment analysis with n = 100 bootstraps

myTestRun <- calcEnrichment(myTestRun, n = 100)

## example of the annotations, that now also include isobars

myTestRun$annotations[[
  sample(which(sapply(myTestRun$isobars_list, length) > 1), size = 1)]][1:10]

## plot enrichment analysis, with q.values on x-axis, and with LION IDs

plotEnrichment(myTestRun, min.annotations = 5, q.value.cutoff = .05, 
               by.statistic = "q.value")


martijnmolenaar/bmetenrichr documentation built on Aug. 16, 2024, 12:20 a.m.