diff --git a/DESCRIPTION b/DESCRIPTION old mode 100644 new mode 100755 index 3812877..134f28c --- a/DESCRIPTION +++ b/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 = "Johannes.Rainer@eurac.edu", @@ -10,6 +10,10 @@ Authors@R: email = "andrea.vicini@eurac.edu", role = "aut", comment = c(ORCID = "0000-0001-9438-6909")), + person(given = "Roger", family = "Gine Bertomeu", + email = "roger.gine@iispv.cat", + comment = c(ORCID = "0000-0003-0288-9619"), + role = "aut"), person(given = "Steffen", family = "Neumann", email = "sneumann@ipb-halle.de", role = "ctb", @@ -33,7 +37,8 @@ Imports: methods, opentimsr (>= 1.0.11), MsCoreUtils, - S4Vectors + S4Vectors, + dplyr Suggests: testthat, knitr (>= 1.1.0), @@ -48,4 +53,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.1 diff --git a/NAMESPACE b/NAMESPACE old mode 100644 new mode 100755 index 4c3e182..9ce7187 --- a/NAMESPACE +++ b/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,19 @@ 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) 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) diff --git a/NEWS.md b/NEWS.md old mode 100644 new mode 100755 index 4462bc7..0fd75c5 --- a/NEWS.md +++ b/NEWS.md @@ -1,3 +1,13 @@ +# 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`. + + + # MsBackendTimsTof 0.1 ## MsBackendTimsTof 0.1.2 diff --git a/R/MsBackendTimsTof-functions.R b/R/MsBackendTimsTof-functions.R old mode 100644 new mode 100755 index c76fd3b..3f10026 --- a/R/MsBackendTimsTof-functions.R +++ b/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] +} diff --git a/R/MsBackendTimsTof.R b/R/MsBackendTimsTof.R old mode 100644 new mode 100755 index 0fae96b..2a50dd0 --- a/R/MsBackendTimsTof.R +++ b/R/MsBackendTimsTof.R @@ -97,7 +97,7 @@ #' @rdname MsBackendTimsTof #' #' @exportClass MsBackendTimsTof -#' +#' #' @examples #' #' ## Load the opentimsr package to retrieve the required shared library @@ -214,7 +214,7 @@ setMethod("rtime", "MsBackendTimsTof", function(object) { .get_frame_columns(object, "rtime") }) -#' @importFrom methods "slot<-" +#' @importFrom methods slot<- #' #' @importFrom MsCoreUtils i2index #' @@ -288,6 +288,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 +352,4 @@ setMethod("$", "MsBackendTimsTof", function(x, name) { else spectraData(x, name)[, 1L] }) + diff --git a/README.md b/README.md old mode 100644 new mode 100755 diff --git a/man/MsBackendTimsTof.Rd b/man/MsBackendTimsTof.Rd old mode 100644 new mode 100755 index 88040a1..6bc5be7 --- a/man/MsBackendTimsTof.Rd +++ b/man/MsBackendTimsTof.Rd @@ -18,6 +18,13 @@ \alias{spectraData,MsBackendTimsTof-method} \alias{show,MsBackendTimsTof-method} \alias{msLevel,MsBackendTimsTof-method} +\alias{precursorMz,MsBackendTimsTof-method} +\alias{precursorCharge,MsBackendTimsTof-method} +\alias{precursorIntensity,MsBackendTimsTof-method} +\alias{collisionEnergy,MsBackendTimsTof-method} +\alias{isolationWindowLowerMz,MsBackendTimsTof-method} +\alias{isolationWindowTargetMz,MsBackendTimsTof-method} +\alias{isolationWindowUpperMz,MsBackendTimsTof-method} \alias{$,MsBackendTimsTof-method} \title{TimsTOF data backend} \usage{ @@ -49,6 +56,20 @@ MsBackendTimsTof() \S4method{msLevel}{MsBackendTimsTof}(object, ...) +\S4method{precursorMz}{MsBackendTimsTof}(object, ...) + +\S4method{precursorCharge}{MsBackendTimsTof}(object, ...) + +\S4method{precursorIntensity}{MsBackendTimsTof}(object, ...) + +\S4method{collisionEnergy}{MsBackendTimsTof}(object, ...) + +\S4method{isolationWindowLowerMz}{MsBackendTimsTof}(object, ...) + +\S4method{isolationWindowTargetMz}{MsBackendTimsTof}(object, ...) + +\S4method{isolationWindowUpperMz}{MsBackendTimsTof}(object, ...) + \S4method{$}{MsBackendTimsTof}(x, name) } \arguments{ diff --git a/tests/testthat.R b/tests/testthat.R old mode 100644 new mode 100755 diff --git a/tests/testthat/test_MsBackendTimsTof-functions.R b/tests/testthat/test_MsBackendTimsTof-functions.R old mode 100644 new mode 100755 index c657f5c..51fe854 --- a/tests/testthat/test_MsBackendTimsTof-functions.R +++ b/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)) +}) + diff --git a/tests/testthat/test_MsBackendTimsTof.R b/tests/testthat/test_MsBackendTimsTof.R old mode 100644 new mode 100755 index 1d0fa64..a16e0fe --- a/tests/testthat/test_MsBackendTimsTof.R +++ b/tests/testthat/test_MsBackendTimsTof.R @@ -163,6 +163,7 @@ test_that("spectraData,MsBackendTimsTof works", { expect_equal(res$inv_ion_mobility, res_2$inv_ion_mobility) ## Random order + set.seed(1234) idx <- sample(seq_along(be)) be_2 <- be[idx] @@ -184,6 +185,48 @@ test_that("msLevel,MsBackendTimsTof works", { expect_identical(res, match(.get_frame_columns(be, "MsMsType"), c(0L, 8L))) }) +test_that("precursorMz,MsBackendTimsTof works", { + expect_identical(precursorMz(MsBackendTimsTof()), numeric(0)) + res <- precursorMz(be) + expect_identical(length(res), length(be)) +}) + +test_that("precursorCharge,MsBackendTimsTof works", { + expect_identical(precursorCharge(MsBackendTimsTof()), integer(0)) + res <- precursorCharge(be) + expect_identical(length(res), length(be)) +}) + +test_that("precursorIntensity,MsBackendTimsTof works", { + expect_identical(precursorIntensity(MsBackendTimsTof()), numeric(0)) + res <- precursorIntensity(be) + expect_identical(length(res), length(be)) +}) + +test_that("collisionEnergy,MsBackendTimsTof works", { + expect_identical(collisionEnergy(MsBackendTimsTof()), numeric(0)) + res <- collisionEnergy(be) + expect_identical(length(res), length(be)) +}) + +test_that("isolationWindowLowerMz,MsBackendTimsTof works", { + expect_identical(isolationWindowLowerMz(MsBackendTimsTof()), numeric(0)) + res <- isolationWindowLowerMz(be) + expect_identical(length(res), length(be)) +}) + +test_that("isolationWindowTargetMz,MsBackendTimsTof works", { + expect_identical(isolationWindowTargetMz(MsBackendTimsTof()), numeric(0)) + res <- isolationWindowTargetMz(be) + expect_identical(length(res), length(be)) +}) + +test_that("isolationWindowUpperMz,MsBackendTimsTof works", { + expect_identical(isolationWindowUpperMz(MsBackendTimsTof()), numeric(0)) + res <- isolationWindowUpperMz(be) + expect_identical(length(res), length(be)) +}) + test_that("$,MsBackendTimsTof works", { res_all <- spectraData(be)