Merge pull request #767 from pablovgd/issue705

Add internal function to calculate beta metrics

R/AllGenerics.R

@@ -567,7 +567,7 @@ setGeneric("chromPeakSpectra", function(object, ...)
 #' columns, their names and content depend on the used parameter object. See
 #' the respective documentation above for more details.
 #'
-#' @author Pable Vangeenderhuysen, Johannes Rainer
+#' @author Pablo Vangeenderhuysen, Johannes Rainer
 #'
 #'
 setGeneric("chromPeakSummary", function(object, param, ...)

R/XcmsExperiment-functions.R

@@ -533,7 +533,8 @@
                 )
                 if ("beta_cor" %in% cn) {
                     res[i, c("beta_cor", "beta_snr")] <- .get_beta_values(
-                        vapply(xsub[nr > 0], sum, NA_real_),
+                        vapply(xsub[nr > 0], function(z) sum(z[, "intensity"]), 
+                               NA_real_),
                         rt[keep][nr > 0])
                 }
             }
@@ -542,6 +543,42 @@
     res[!is.na(res[, "maxo"]), , drop = FALSE]
 }
 
+
+#' Calculates quality metrics for a chromatographic peak. 
+#'
+#' @param x `list` of peak matrices (from a single MS level and from a single
+#'     file/sample).
+#'
+#' @param rt retention time for each peak matrix.
+#'
+#' @param peakArea `matrix` defining the chrom peak area.
+#'
+#' @author Pablo Vangeenderhuysen
+#'
+#' @noRd
+.chrom_peak_beta_metrics <- function(x, rt, peakArea, ...) {
+  res <- matrix(NA_real_, ncol = 2L, nrow = nrow(peakArea))
+  rownames(res) <- rownames(peakArea)
+  colnames(res) <- c("beta_cor","beta_snr")
+  for (i in seq_len(nrow(res))) {
+    rtr <- peakArea[i, c("rtmin", "rtmax")]
+    keep <- which(between(rt, rtr))
+    if (length(keep)) {
+      xsub <- lapply(x[keep], .pmat_filter_mz,
+                     mzr = peakArea[i, c("mzmin", "mzmax")])
+      nr <- vapply(xsub, nrow, NA_integer_)
+      res[i, c("beta_cor", "beta_snr")] <- .get_beta_values(
+          vapply(xsub[nr > 0], function(z) sum(z[, "intensity"]), NA_real_),
+          rt[keep][nr > 0])
+      }
+    }
+  res
+}
+
+
+
+
+
 #' Difference to the original code is that the weighted mean is also calculated
 #' if some of the peak intensities in the profile matrix are 0
 #'

tests/testthat/test_XcmsExperiment.R

@@ -825,6 +825,19 @@ test_that(".chrom_peak_intensity_centWave works", {
     ## pks[11, ].
 })
 
+
+## That's from XcmsExperiment-functions.R
+test_that(".chrom_peak_beta_metrics works", {
+  x <- Spectra::peaksData(spectra(xmse[2L]))
+  rt <- rtime(spectra(xmse[2L]))
+  pks <- chromPeaks(xmse)[chromPeaks(xmse)[, "sample"] == 2L, ]
+
+  res <- .chrom_peak_beta_metrics(x, rt, pks, sampleIndex = 2L,
+                                        cn = colnames(pks))
+  expect_equal(nrow(res), nrow(pks))
+
+})
+
 ## That's from XcmsExperiment-functions.R
 test_that(".chrom_peak_intensity_matchedFilter works", {
     x <- Spectra::peaksData(spectra(xmse[2L]))