diff --git a/R/AllGenerics.R b/R/AllGenerics.R index ff86fd99..a53f91bc 100644 --- a/R/AllGenerics.R +++ b/R/AllGenerics.R @@ -417,10 +417,10 @@ setGeneric("chromPeakData<-", function(object, value) #' spectra for one chromatographic peak (or a `Spectra` of length 0 if no #' spectrum was found for the respective chromatographic peak). #' -#' Parameters `chromPeakColumns` allow the user to add specific metadata +#' Parameter `chromPeakColumns` allows the user to add specific metadata #' columns from the chromatographic peaks (`chromPeaks`) to the returned -#' spectra object. This can be useful to retain information such as retention -#' time (`rt`), m/z (`mz`). The columns will be named as they is written in the +#' spectra object. This can be useful to keep information such as retention +#' time (`rt`), m/z (`mz`). The columns will be named as they are written in the #' `chromPeaks` object with the prefix `"chrom_peak_"`. The *peak ID* #' (i.e., the row name of the peak in the `chromPeaks` matrix) is always added #' to the spectra object as a metadata column named `"chrom_peak_id"`. @@ -811,7 +811,7 @@ setGeneric("featureDefinitions<-", function(object, value) #' #' The information from `featureDefinitions` for each feature can be included #' in the returned [Spectra()] object using the `featureColumns` parameter. -#' This is useful for retaining details such as the median retention time (`rtmed`) +#' This is useful for keeping details such as the median retention time (`rtmed`) #' or median m/z (`mzmed`). The columns will retain their names as specified #' in the `featureDefinitions` object, prefixed by `"feature_"` #' (e.g., `"feature_mzmed"`). Additionally, the *feature ID* (i.e., the row diff --git a/R/XcmsExperiment-functions.R b/R/XcmsExperiment-functions.R index 1b4ea2d1..cb6be689 100644 --- a/R/XcmsExperiment-functions.R +++ b/R/XcmsExperiment-functions.R @@ -797,8 +797,12 @@ chromPeakColumns = c("rt", "mz"), BPPARAM = bpparam()) { method <- match.arg(method) - pks <- .chromPeaks(x)[, c("mz", "mzmin", "mzmax", "rt", - "rtmin", "rtmax", "maxo", "sample")] + if (!chromPeakColumns %in% colnames(.chromPeaks(x))) + stop("One or more of the columns in 'chromPeakColumns' are not ", + "available in the 'chromPeaks' data.") + pks <- .chromPeaks(x)[, union(c("mz", "mzmin", "mzmax", "rt", + "rtmin", "rtmax", "maxo", "sample"), + chromPeakColumns)] if (ppm != 0) expandMz <- expandMz + pks[, "mz"] * ppm / 1e6 if (expandMz[1L] != 0) { diff --git a/R/XcmsExperiment.R b/R/XcmsExperiment.R index e5952975..69588eb6 100644 --- a/R/XcmsExperiment.R +++ b/R/XcmsExperiment.R @@ -1782,6 +1782,9 @@ setMethod( if (!hasFeatures(object)) stop("No feature definitions present. Please run ", "'groupChromPeaks' first.") + if (!featureColumns %in% colnames(featureDefinitions(object))) + stop("One or more of the requested 'featureColumns' are not ", + "present in the feature definitions.") if (hasAdjustedRtime(object)) object <- applyAdjustedRtime(object) features_all <- rownames(featureDefinitions(object)) diff --git a/man/chromPeakSpectra.Rd b/man/chromPeakSpectra.Rd index 753f2cfd..8899c15a 100644 --- a/man/chromPeakSpectra.Rd +++ b/man/chromPeakSpectra.Rd @@ -136,10 +136,10 @@ of \code{chromPeaks}. Each element of the list contains thus a \code{Spectra} wi spectra for one chromatographic peak (or a \code{Spectra} of length 0 if no spectrum was found for the respective chromatographic peak). -Parameters \code{chromPeakColumns} allow the user to add specific metadata +Parameter \code{chromPeakColumns} allows the user to add specific metadata columns from the chromatographic peaks (\code{chromPeaks}) to the returned -spectra object. This can be useful to retain information such as retention -time (\code{rt}), m/z (\code{mz}). The columns will be named as they is written in the +spectra object. This can be useful to keep information such as retention +time (\code{rt}), m/z (\code{mz}). The columns will be named as they are written in the \code{chromPeaks} object with the prefix \code{"chrom_peak_"}. The \emph{peak ID} (i.e., the row name of the peak in the \code{chromPeaks} matrix) is always added to the spectra object as a metadata column named \code{"chrom_peak_id"}. diff --git a/man/featureSpectra.Rd b/man/featureSpectra.Rd index 59753b5e..c669649a 100644 --- a/man/featureSpectra.Rd +++ b/man/featureSpectra.Rd @@ -105,7 +105,7 @@ spectra per feature). The information from \code{featureDefinitions} for each feature can be included in the returned \code{\link[=Spectra]{Spectra()}} object using the \code{featureColumns} parameter. -This is useful for retaining details such as the median retention time (\code{rtmed}) +This is useful for keeping details such as the median retention time (\code{rtmed}) or median m/z (\code{mzmed}). The columns will retain their names as specified in the \code{featureDefinitions} object, prefixed by \code{"feature_"} (e.g., \code{"feature_mzmed"}). Additionally, the \emph{feature ID} (i.e., the row diff --git a/tests/testthat/test_XcmsExperiment.R b/tests/testthat/test_XcmsExperiment.R index 57703518..10d0699e 100644 --- a/tests/testthat/test_XcmsExperiment.R +++ b/tests/testthat/test_XcmsExperiment.R @@ -939,6 +939,9 @@ test_that("chromPeakSpectra works", { expect_error(chromPeakSpectra(xmse, peaks = "other"), "out of bounds") pks <- c("CP242", "CP007", "CP123") + + expect_error(chromPeakSpectra(xmse, peaks = pks, + chromPeakColumns = "other"), "not available") res <- chromPeakSpectra(xmse, peaks = pks) expect_s4_class(res, "Spectra") expect_equal(length(res), 0) @@ -1081,6 +1084,8 @@ test_that("filterFeatureDefinitions works", { test_that("featureSpectra works", { expect_error(featureSpectra(xmse), "No feature definitions") + expect_error(featureSpectra(xmseg, + featureColumns = "other"), "not present") res_all <- featureSpectra(xmseg, msLevel = 1L) expect_s4_class(res_all, "Spectra") expect_true(all(rownames(featureDefinitions(xmseg)) %in%