Read PASEF DDA MS2 precursor information

DESCRIPTION

@@ -1,6 +1,6 @@
 Package: MsBackendTimsTof
 Title: Mass spectrometry Data Backend for Bruker TimsTOF Files
-Version: 0.1.2
+Version: 0.2.0
 Authors@R:
     c(person(given = "Johannes", family = "Rainer",
              email = "[email protected]",
@@ -10,6 +10,10 @@ Authors@R:
              email = "[email protected]",
              role = "aut",
              comment = c(ORCID = "0000-0001-9438-6909")),
+      person(given = "Roger", family = "Gine Bertomeu", 
+             email = "[email protected]",
+             comment = c(ORCID = "0000-0003-0288-9619"),
+             role = "aut"),
       person(given = "Steffen", family = "Neumann",
              email = "[email protected]",
              role = "ctb",
@@ -33,7 +37,9 @@ Imports:
     methods,
     opentimsr (>= 1.0.11),
     MsCoreUtils,
-    S4Vectors
+    S4Vectors,
+    dplyr,
+    utils
 Suggests:
     testthat,
     knitr (>= 1.1.0),
@@ -48,4 +54,4 @@ BugReports: https://github.com/RforMassSpectrometry/MsBackendTimsTof/issues
 URL: https://github.com/RforMassSpectrometry/MsBackendTimsTof
 biocViews: Infrastructure, MassSpectrometry, Metabolomics, DataImport
 Roxygen: list(markdown=TRUE)
-RoxygenNote: 7.1.2
+RoxygenNote: 7.3.0

NAMESPACE

@@ -4,12 +4,16 @@ export(MsBackendTimsTof)
 exportClasses(MsBackendTimsTof)
 exportMethods(length)
 importFrom(BiocParallel,bplapply)
+importFrom(BiocParallel,bpmapply)
 importFrom(BiocParallel,bpparam)
 importFrom(IRanges,NumericList)
+importFrom(MsCoreUtils,between)
 importFrom(MsCoreUtils,i2index)
 importFrom(MsCoreUtils,rbindFill)
 importFrom(S4Vectors,DataFrame)
 importFrom(Spectra,coreSpectraVariables)
+importFrom(dplyr,left_join)
+importFrom(dplyr,rename)
 importFrom(methods,"slot<-")
 importFrom(methods,as)
 importFrom(methods,new)
@@ -19,11 +23,20 @@ importFrom(opentimsr,OpenTIMS)
 importFrom(opentimsr,opentims_set_threads)
 importFrom(opentimsr,query)
 importFrom(opentimsr,setup_bruker_so)
+importFrom(opentimsr,table2df)
 importFrom(stats,setNames)
 importFrom(utils,capture.output)
+importFrom(utils,packageVersion)
 importMethodsFrom(Spectra,backendInitialize)
+importMethodsFrom(Spectra,collisionEnergy)
 importMethodsFrom(Spectra,dataStorage)
+importMethodsFrom(Spectra,isolationWindowLowerMz)
+importMethodsFrom(Spectra,isolationWindowTargetMz)
+importMethodsFrom(Spectra,isolationWindowUpperMz)
 importMethodsFrom(Spectra,msLevel)
 importMethodsFrom(Spectra,peaksVariables)
+importMethodsFrom(Spectra,precursorCharge)
+importMethodsFrom(Spectra,precursorIntensity)
+importMethodsFrom(Spectra,precursorMz)
 importMethodsFrom(Spectra,spectraData)
 importMethodsFrom(Spectra,spectraVariables)

NEWS.md

@@ -1,3 +1,14 @@
+# MsBackendTimsTof 0.2
+
+## MsBackendTimsTof 0.2.0
+
+- Read PASEF MS/MS data from the TDF tables.
+
+- Add methods for `precursorMz`, `precursorCharge`, `precursorIntensity`, `collisionEnergy`,   `isolationWindowLowerMz`, `isolationWindowTargetMz` and `isolationWindowUpperMz`.
+
+- Update version info automatically in the `version` slot.
+
+
 # MsBackendTimsTof 0.1
 
 ## MsBackendTimsTof 0.1.2

R/MsBackendTimsTof-functions.R

@@ -6,7 +6,7 @@
 #'
 #' @importFrom MsCoreUtils rbindFill
 #'
-#' @importFrom opentimsr OpenTIMS CloseTIMS
+#' @importFrom opentimsr OpenTIMS CloseTIMS table2df
 #'
 #' @return initialized `MsBackendTimsTOF` object
 #'
@@ -33,6 +33,7 @@
         x@frames$polarity <- .format_polarity(x@frames$polarity)
     x@indices <- do.call(rbind, lapply(L, "[[", 2))
     x@fileNames <- setNames(seq_len(length(file)), file)
+
     x
 }
 
@@ -288,6 +289,11 @@ MsBackendTimsTof <- function() {
 # can we assume that tms@all_coulmns is the same for all the TimsTOF?
 .TIMSTOF_COLUMNS <- c("mz", "intensity", "tof", "inv_ion_mobility")
 
+.TIMSTOF_MS2_COLUMNS <- c("precursorMz", "precursorCharge",
+                          "precursorIntensity", "collisionEnergy",
+                          "isolationWindowLowerMz", "isolationWindowTargetMz",
+                          "isolationWindowUpperMz")
+
 #' @importFrom methods as
 #'
 #' @importFrom S4Vectors DataFrame
@@ -349,6 +355,22 @@ MsBackendTimsTof <- function() {
         res[["dataStorage"]] <- dataStorage(x)
         core_cols <- core_cols[core_cols != "dataStorage"]
     }
+    if (any(core_cols %in% .TIMSTOF_MS2_COLUMNS)) {
+        res_ms2_data <- .calculate_core_ms2_information(
+            x, 
+            columns = core_cols[core_cols %in% .TIMSTOF_MS2_COLUMNS]
+        )
+        for (col in core_cols[core_cols %in% .TIMSTOF_MS2_COLUMNS]) {
+            if (col == "precursorCharge"){
+                res[[col]] <- as.integer(
+                    res_ms2_data[, col]
+                )
+            } else {
+                res[[col]] <- res_ms2_data[, col]    
+            }
+            core_cols <- core_cols[core_cols != col]
+        }
+    }
     if (length(core_cols)) {
         res[core_cols] <- lapply(coreSpectraVariables()[core_cols],
                                  function(z, n) rep(as(NA, z), n), length(x))
@@ -357,3 +379,85 @@ MsBackendTimsTof <- function() {
         as(res, "DataFrame")
     else DataFrame()
 }
+
+#' @importFrom BiocParallel bpmapply bpparam
+.calculate_core_ms2_information <- function(x, columns = .TIMSTOF_MS2_COLUMNS){
+    if(!length(columns)){
+        stop("Some columns must be selected. Valid options are: ",
+             paste(.TIMSTOF_MS2_COLUMNS, collapse = " "))
+    }
+    if(!all(columns %in% .TIMSTOF_MS2_COLUMNS)){
+        stop("Invalid column/s selected: ", 
+             columns[!columns %in% .TIMSTOF_MS2_COLUMNS], "\n",
+             "Valid options are: ", paste(.TIMSTOF_MS2_COLUMNS, collapse = " "))
+    }
+    if (is.null(names(x@fileNames))){
+        out <- data.frame(matrix(NA_real_, nrow = 0, ncol = length(columns)))
+        colnames(out) <- columns
+        return(out)
+    }
+    tbls <- bpmapply(FUN = .do_calculate_core_ms2_information,
+                     x = names(x@fileNames),
+                     indices = split(as.data.frame(x@indices[,1:2]),
+                                     f = x@indices[,3]),
+                     MoreArgs = list(columns = columns),
+                     USE.NAMES = FALSE,
+                     SIMPLIFY = FALSE,
+                     BPPARAM = bpparam())
+    output <- do.call(rbind, tbls)
+    rownames(output) <- NULL
+    output[order(order(x@indices[,3])), , drop=FALSE]
+}
+
+#' @importFrom opentimsr OpenTIMS CloseTIMS table2df
+#' @importFrom MsCoreUtils between
+#' @importFrom dplyr left_join rename
+.do_calculate_core_ms2_information <- function(x, indices,
+                                               columns = .TIMSTOF_MS2_COLUMNS){
+    tms <- opentimsr::OpenTIMS(x)
+    on.exit(opentimsr::CloseTIMS(tms))
+    if (missing(indices)) {
+        indices <- unique(.query_tims(tms, frames = tms@frames$Id,
+                                      columns = c("frame", "scan")))
+    }
+    indices <- as.data.frame(indices)
+    if(!all(c("frame", "scan") %in% colnames(indices))){
+        stop("Indices must have both frame and scan columns\n",
+             "Current columns are: ", paste(colnames(indices), collapse = " "))
+    }
+
+    tbl <- opentimsr::table2df(tms, c("PASEFFrameMsMsInfo", "Precursors"))
+
+    ## Join MS/MS frame details with precursor information
+    joined <- left_join(tbl[[1]], tbl[[2]], by = c("Precursor" = "Id"))
+
+    ## Create needed MS/MS info columns
+    joined <- rename(joined, 
+                     isolationWindowTargetMz = IsolationMz,
+                     collisionEnergy = CollisionEnergy,
+                     precursorMz = MonoisotopicMz,
+                     precursorCharge = Charge,
+                     precursorIntensity = Intensity)
+
+    joined$isolationWindowLowerMz <- 
+        joined$isolationWindowTargetMz - joined$IsolationWidth / 2
+    joined$isolationWindowUpperMz <- 
+        joined$isolationWindowTargetMz + joined$IsolationWidth / 2
+
+    ## Expand to get a row for each MS/MS scan, then join with the indices
+    full <- joined[rep(seq_len(nrow(joined)),
+                       times = (joined$ScanNumEnd - joined$ScanNumBegin + 1)),]
+    full$scan <- do.call(c, 
+                         mapply(function(x,y){seq(x,y)}, 
+                                joined$ScanNumBegin, joined$ScanNumEnd,
+                                SIMPLIFY = FALSE)
+                         )
+    output <- left_join(as.data.frame(indices), full,
+                        c("frame" = "Frame", "scan"))
+    if(!all(columns %in% colnames(output))){
+        stop("Could not find the following columns in the MS/MS data tables: ", 
+             columns[!columns %in% colnames(output)])
+    }
+    rownames(output) <- NULL
+    output[, columns, drop = FALSE]
+}

R/MsBackendTimsTof.R

@@ -97,6 +97,8 @@
 #' @rdname MsBackendTimsTof
 #'
 #' @exportClass MsBackendTimsTof
+#' 
+#' @importFrom utils packageVersion
 #'
 #' @examples
 #'
@@ -137,7 +139,8 @@ setClass("MsBackendTimsTof",
                                                                   "file"))),
                                fileNames = integer(),
                                readonly = TRUE,
-                               version = "0.1"))
+                               version = as.character(
+                                   packageVersion("MsBackendTimsTof"))))
 
 #' @importFrom methods validObject
 setValidity("MsBackendTimsTof", function(object) {
@@ -214,7 +217,7 @@ setMethod("rtime", "MsBackendTimsTof", function(object) {
     .get_frame_columns(object, "rtime")
 })
 
-#' @importFrom methods "slot<-"
+#' @importFrom methods slot<-
 #'
 #' @importFrom MsCoreUtils i2index
 #'
@@ -288,6 +291,61 @@ setMethod("msLevel", "MsBackendTimsTof", function(object, ...) {
     .get_msLevel(object)
 })
 
+#' @importMethodsFrom Spectra precursorMz
+#'
+#' @rdname MsBackendTimsTof
+setMethod("precursorMz", "MsBackendTimsTof", function(object, ...) {
+    .calculate_core_ms2_information(object, columns = "precursorMz")[[1]]
+})
+
+#' @importMethodsFrom Spectra precursorCharge
+#'
+#' @rdname MsBackendTimsTof
+setMethod("precursorCharge", "MsBackendTimsTof", function(object, ...) {
+    .calculate_core_ms2_information(object,
+                                    columns = "precursorCharge")[[1]] |>
+    as.integer()
+})
+
+#' @importMethodsFrom Spectra precursorIntensity
+#'
+#' @rdname MsBackendTimsTof
+setMethod("precursorIntensity", "MsBackendTimsTof", function(object, ...) {
+    .calculate_core_ms2_information(object,
+                                    columns = "precursorIntensity")[[1]]
+})
+
+#' @importMethodsFrom Spectra collisionEnergy
+#'
+#' @rdname MsBackendTimsTof
+setMethod("collisionEnergy", "MsBackendTimsTof", function(object, ...) {
+    .calculate_core_ms2_information(object, columns = "collisionEnergy")[[1]]
+})
+
+#' @importMethodsFrom Spectra isolationWindowLowerMz
+#'
+#' @rdname MsBackendTimsTof
+setMethod("isolationWindowLowerMz", "MsBackendTimsTof", function(object, ...) {
+    .calculate_core_ms2_information(object, 
+                                    columns = "isolationWindowLowerMz")[[1]]
+})
+
+#' @importMethodsFrom Spectra isolationWindowTargetMz
+#'
+#' @rdname MsBackendTimsTof
+setMethod("isolationWindowTargetMz", "MsBackendTimsTof", function(object, ...) {
+    .calculate_core_ms2_information(object, 
+                                    columns = "isolationWindowTargetMz")[[1]]
+})
+
+#' @importMethodsFrom Spectra isolationWindowUpperMz
+#'
+#' @rdname MsBackendTimsTof
+setMethod("isolationWindowUpperMz", "MsBackendTimsTof", function(object, ...) {
+    .calculate_core_ms2_information(object,
+                                    columns = "isolationWindowUpperMz")[[1]]
+})
+
 #' @rdname MsBackendTimsTof
 setMethod("$", "MsBackendTimsTof", function(x, name) {
     if (!any(spectraVariables(x) == name))
@@ -297,3 +355,4 @@ setMethod("$", "MsBackendTimsTof", function(x, name) {
     else
         spectraData(x, name)[, 1L]
 })
+

tests/testthat/test_MsBackendTimsTof-functions.R

@@ -145,3 +145,18 @@ test_that(".inv_ion_mobility works", {
     res_2 <- .inv_ion_mobility(be_2)
     expect_equal(res_2, res[idx])
 })
+
+test_that("calculate_core_ms2_information works", {
+    tbl <- .calculate_core_ms2_information(be)
+    expect_true(all(apply(tbl, 2, is.numeric)))
+    expect_equal(ncol(tbl), length(.TIMSTOF_MS2_COLUMNS))
+    expect_equal(nrow(tbl), nrow(be@indices))
+    expect_true(!all(is.na(tbl)))
+
+    #Check that it also works with subsets of scans
+    be_subset <- be[c(120:130, 17000:17010)] #From both file 1 and 2
+    tbl <- .calculate_core_ms2_information(be_subset)
+    expect_equal(ncol(tbl), length(.TIMSTOF_MS2_COLUMNS))
+    expect_equal(nrow(tbl), nrow(be_subset@indices))
+})
+