remove (AssaySupport) legacy functions and add more tests

NAMESPACE

@@ -38,7 +38,6 @@ export(loadSpectra)
 export(loadSpectraMzML)
 export(normalize)
 export(normalizeByFactor)
-export(peaks2df)
 export(plotCurves)
 export(plotPeak)
 export(shiftMassAxis)
@@ -66,7 +65,6 @@ importFrom(MALDIquant,referencePeaks)
 importFrom(MALDIquant,removeBaseline)
 importFrom(MALDIquant,snr)
 importFrom(MALDIquant,totalIonCurrent)
-importFrom(MALDIquant,trim)
 importFrom(MALDIquant,warpMassPeaks)
 importFrom(MALDIquant,warpMassSpectra)
 importFrom(MALDIquantForeign,importBrukerFlex)
@@ -123,6 +121,7 @@ importFrom(nplr,getYcurve)
 importFrom(nplr,nplr)
 importFrom(purrr,map)
 importFrom(purrr,map_dbl)
+importFrom(purrr,map_vec)
 importFrom(scales,comma)
 importFrom(scales,scientific)
 importFrom(stats,median)

R/MALDIassay_Class.R

@@ -22,6 +22,79 @@ setClass("MALDIassay",
   )
 )
 
+.validMALDIassay <- function(object) {
+  if (length(object@avgSpectra) < 1) {
+    return("Length of avgSpectra can't be 0.")
+  }
+
+  if (!isMassSpectrumList(object@avgSpectra)) {
+    return("avgSpectra must be a list of class MALDIquant::MassSpectrum objects.")
+  }
+
+  if (length(object@avgPeaks) < 1) {
+    return("Length of avgPeaks can't be 0.")
+  }
+
+  if (!isMassPeaksList(object@avgPeaks)) {
+    return("avgPeaks must be a list of class MALDIquant::MassPeaks objects.")
+  }
+
+  if (length(object@singlePeaks) < 1) {
+    return("Length of singlePeaks can't be 0.")
+  }
+
+  if (!isMassPeaksList(object@singlePeaks)) {
+    return("singlePeaks must be a list of class MALDIquant::MassPeaks objects.")
+  }
+
+  if (length(object@avgSpectra) != length(object@avgPeaks)) {
+    return(paste0("avgSpectra (", length(object@avgSpectra),
+                  ") and avgPeaks (", length(object@avgPeaks),
+                  ") must have the same length."))
+  }
+
+  if (length(object@singlePeaks) != length(object@included_specIdx)) {
+    return(paste0("singlePeaks (", length(object@singlePeaks),
+                  ") and included_specIdx (", length(object@included_specIdx),
+                  ") must have the same length."))
+  }
+
+  if (!(length(object@singlePeaks) == length(object@mzShifts) || length(object@mzShifts) == 1)) {
+    return(paste0("singlePeaks (", length(object@singlePeaks),
+                  ") and mzShifts (", length(object@mzShifts),
+                  ") must have the same length or mzShifts must have length 1."))
+  }
+
+  if (!(length(object@singlePeaks) == length(object@normFactors) || length(object@normFactors) == 1)) {
+    return(paste0("singlePeaks (", length(object@singlePeaks),
+                  ") and normFactors (", length(object@normFactors),
+                  ") must have the same length or normFactors must have length 1."))
+  }
+
+  if (!(length(object@singlePeaks) == length(object@singleSpecSpots) || length(object@singleSpecSpots) == 1)) {
+    return(paste0("singlePeaks (", length(object@singlePeaks),
+                  ") and singleSpecSpots (", length(object@singleSpecSpots),
+                  ") must have the same length or singleSpecSpots must have length 1."))
+  }
+
+  if (length(object@fits) != length(unique(object@stats[["mz"]]))) {
+    return(paste0("fits (", length(object@fits),
+                  ") and the unique m/z values in stats (", length(unique(object@stats[["mz"]])),
+                  ") must have the same length."))
+  }
+
+  TRUE
+}
+
+setValidity("MALDIassay", method=.validMALDIassay)
+
+#' Print summary of MALDIassay object
+#'
+#' @param object object of class MALDIassay
+#'
+#' @return nothing, prints to console
+#'
+#' @noRd
 show_MALDIassay <- function(object) {
   mz <- round(getNormMz(object), digits = 2)
   varFilterMethod <- getVarFilterMethod(object)

R/fitCurve.R

@@ -13,7 +13,7 @@
 #' @param SNR                 Numeric, signal to noise ratio for peak detection.
 #' @param halfWindowSize      2ction. See `MALDIquant::detectPeaks()`.
 #' @param allowNoMatches      Logical, if normMz can not be found in a spectrum, proceed and exclude spectrum or stop
-#' @param normMeth            Character, normalization method. Can either be "TIC", "PQM", "median" or "mz". If "mz" then the normMz is used. If none no normalization is done.
+#' @param normMeth            Character, normalization method. Can either be "TIC", "median" or "mz". If "mz" then the normMz is used. If none no normalization is done.
 #' @param SinglePointRecal    Logical, perform single point recalibration to normMz
 #' @param verbose             Logical, print logs to console.
 #'
@@ -53,7 +53,7 @@ fitCurve <- function(spec,
                      SNR = 3,
                      halfWindowSize = 3,
                      allowNoMatches = TRUE,
-                     normMeth = c("mz", "TIC", "PQN", "median", "none"),
+                     normMeth = c("mz", "TIC", "median", "none"),
                      SinglePointRecal = TRUE,
                      verbose = TRUE) {
 

R/getMzShift.R

@@ -1,70 +1,19 @@
-#' Get mass shift for target mz - internal
-#'
-#' @param peaksdf       data.frame with peaks information as generated by peaks2df()
-#' @param tol           Numeric, tolerance around targetMz
-#' @param targetMz      Numeric, target mass
-#' @param tolppm        Logical, is the tolerance provided in ppm (TRUE) or Dalton (FALSE)
-#' @param allowNoMatch  Logical, stop if targetMz is not found in single spectrum?
-#'                      If TRUE spectra without targetMz match will be excluded.
-#'
-#' @return             List with two entries:
-#'                                       MzShift The mass shift for each spectrum
-#'                                       specIdx The index of the spectra with a match for targetMz
-#'
-#' @importFrom dplyr %>% mutate filter group_by arrange pull .data
-#' @noRd
-
-.getMzShift <- function(peaksdf, tol, targetMz, tolppm = TRUE, allowNoMatch = TRUE) {
-  plot_Idx <- sort(unique(peaksdf$plotIdx))
-
-  if(tolppm) {
-    tol <- tol * targetMz / 1e6
-  }
-
-
-  f_resdf <- peaksdf %>%
-    mutate(match = .data$mz > targetMz - tol & .data$mz < targetMz + tol) %>%
-    filter(match) %>%
-    mutate(mz.diff = round(targetMz - .data$mz, 4)) %>%
-    group_by(.data$plotIdx) %>%
-    filter(abs(.data$mz.diff) == min(abs(.data$mz.diff))) %>%
-    arrange(.data$plotIdx)
-
-
-
-  if (!all(plot_Idx %in% (f_resdf %>% pull(.data$plotIdx)))) {
-    if (!allowNoMatch) {
-      stop("Could not find ", targetMz, " for all spectra! Consider adjusting tol.\n")
-    }
-    warning("Could not find ", targetMz, " in spectrum ", paste(which(!(plot_Idx %in% (f_resdf %>% pull(.data$plotIdx)))), collapse = ", "), ".\n")
-    specIdx <- sort(which(plot_Idx %in% (f_resdf %>% pull(.data$plotIdx))))
-  } else {
-    specIdx <- plot_Idx
-  }
-  if (length(specIdx) < 1) {
-    stop("Could not find targetMz in any spectrum! Consider adjusting tol.\n")
-  }
-
-  return(list(
-    mzshift = pull(f_resdf, .data$mz.diff),
-    specIdx = specIdx
-  ))
-}
-
-
 #' Get mass shift for target mz
 #'
-#' @param peaks    List of MALDIquant::MassPeak
-#' @param targetMz Numeric, target mass
-#' @param tol      Numeric, tolerance around targetMz
-#' @param tolppm   Logical, tolerance supplied in ppm
-#' @param verbose  Logical, print logs to the console.
+#' @param peaks         List of MALDIquant::MassPeak
+#' @param targetMz      Numeric, target mass
+#' @param tol           Numeric, tolerance around targetMz
+#' @param tolppm        Logical, tolerance supplied in ppm
+#' @param verbose       Logical, print logs to the console.
+#' @param allowNoMatch  Logical, allow no matches.
 #'
 #' @return
 #' List with two entries:
 #' `MzShift` The mass shift for each spectrum
 #' `specIdx` The index of the spectra with a match for targetMz
 #' @export
+#' @importFrom MALDIquant match.closest
+#' @importFrom purrr map_vec
 #' 
 #' @examples
 #' data(Blank2022peaks)
@@ -73,38 +22,62 @@ getMzShift <- function(peaks,
                        targetMz,
                        tol,
                        tolppm = FALSE,
-                       verbose = TRUE) {
+                       verbose = TRUE,
+                       allowNoMatch = TRUE) {
   stopifnot(isMassPeaksList(peaks))
   nm <- names(peaks)
   stopifnot(!is.null(nm))
   stopifnot(all(!is.na((as.numeric(nm)))))
-
-  # perform single point mass re-calibration
-  mzShift <- .getMzShift(
-    peaksdf = peaks2df(peaks),
-    tol = tol,
-    targetMz = targetMz,
-    tolppm = tolppm,
-    allowNoMatch = TRUE
-  )
+
+  allIdx <- 1:length(peaks)
+
+  if(tolppm) {
+    tol <- tol * targetMz / 1e6
+  }
+
+  mzShift <- map_vec(peaks, 
+          function(x) {
+            mz <- mass(x)
+
+            idx <- match.closest(targetMz, 
+                          table = mz, 
+                          tolerance = tol, 
+                          nomatch = NA_integer_)
+            mzdiff <- targetMz-mz[idx]
+            return(mzdiff)
+          })
+
+  specIdx <- which(!is.na(mzShift))
+  mzShift <- na.omit(mzShift)
+
+  if (!all(allIdx %in% specIdx)) {
+    warning("Could not find ", targetMz, " in spectrum ", paste(which(!(allIdx %in% specIdx)), collapse = ", "), ".\n")
+    specIdx <- sort(which(allIdx %in% specIdx))
+  } else {
+    specIdx <- allIdx
+  }
+  if (length(specIdx) < 1) {
+    stop("Could not find targetMz in any spectrum! Consider adjusting tol.\n")
+  }
+
   if(verbose) {
-    cat("found mz", targetMz, "in", length(mzShift$specIdx), "/",
+    cat("found mz", targetMz, "in", length(specIdx), "/",
         length(peaks), "spectra\n")
-    cat(timeNow(), "mzshift was", mean(mzShift$mzshift),
-        "in mean and", max(abs(mzShift$mzshift)), " abs. max.\n")  
+    cat(timeNow(), "mzshift was", mean(mzShift),
+        "in mean and", max(abs(mzShift)), " abs. max.\n")  
   }
-
 
-  if(length(unique(nm)) != length(unique(nm[mzShift$specIdx]))) {
+  if(length(unique(nm)) != length(unique(nm[specIdx]))) {
     # stop if a single condition got filtered completely
     u_nm <- unique(nm)
-    u_fil <- unique(nm[mzShift$specIdx])
+    u_fil <- unique(nm[specIdx])
     label_removed <- u_nm[which(!(u_nm %in% u_fil))]
 
     stop("Could not find ", targetMz, " in any spectrum with label ",
          paste0(label_removed, collapse = ", "),
          ".\n Consider increasing tol.\n")
   }
-  return(mzShift)
+  return(list(mzshift = mzShift,
+              specIdx = specIdx))
 }
 

R/getNormFactors.R

@@ -1,6 +1,6 @@
 #' Get normalization factors from peak data.frame
 #'
-#' @param peaksdf       data.frame with peaks information as generated by peaks2df()
+#' @param peaks         List of Object of class MALDIquant::MassPeaks 
 #' @param targetMz      Numeric, target mass
 #' @param tol           Numeric, tolerance around targetMz
 #' @param tolppm        Logical, is the tolerance provided in ppm (TRUE) or Daltion (FALSE)
@@ -16,41 +16,42 @@
 #' 
 #' @examples
 #' data(Blank2022peaks)
-#' getNormFactors(peaks2df(Blank2022peaks), targetMz = 760.585, tol = 0.1, tolppm = FALSE)
-getNormFactors <- function(peaksdf, targetMz, tol, tolppm = TRUE, allowNoMatch = TRUE) {
-  plot_Idx <- sort(unique(peaksdf$plotIdx))
-
+#' getNormFactors(Blank2022peaks, targetMz = 760.585, tol = 0.1, tolppm = FALSE)
+getNormFactors <- function(peaks, targetMz, tol, tolppm = TRUE, allowNoMatch = TRUE) {
+  allIdx <- 1:length(peaks)
+  
   if (tolppm) {
-    tol <- (tol / 1e6)
-    resdf <- peaksdf %>%
-      mutate(match = .data$mz > targetMz - .data$mz * tol & .data$mz < targetMz + .data$mz * tol)
-  } else {
-    resdf <- peaksdf %>%
-      mutate(match = .data$mz > targetMz - tol & .data$mz < targetMz + tol)
-  }
-
-  f_resdf <- resdf %>%
-    filter(match) %>%
-    mutate(mz.diff = round(targetMz - .data$mz, 4)) %>%
-    group_by(.data$plotIdx) %>%
-    filter(abs(.data$mz.diff) == min(abs(.data$mz.diff))) %>%
-    arrange(.data$plotIdx)
-
-
-  if (!all(plot_Idx %in% (f_resdf %>% pull(.data$plotIdx)))) {
+    tol <- tol * targetMz / 1e6
+  } 
+
+  norm_factor <- map_vec(peaks, 
+                         function(x) {
+                           mz <- mass(x)
+
+                           idx <- match.closest(targetMz, 
+                                                table = mz, 
+                                                tolerance = tol, 
+                                                nomatch = NA_integer_)
+
+                           return(intensity(x)[idx])
+                         })
+  specIdx <- which(!is.na(norm_factor))
+  norm_factor <- na.omit(norm_factor)
+
+  if (!all(allIdx %in% specIdx)) {
     if (!allowNoMatch) {
       stop("Could not find ", targetMz, " for all spectra! Consider adjusting tol.\n")
     }
-    warning("Could not find ", targetMz, " in spectrum ", paste(which(!(plot_Idx %in% (f_resdf %>% pull(.data$plotIdx)))), collapse = ", "), ".\n")
-    specIdx <- which(plot_Idx %in% (f_resdf %>% pull(.data$plotIdx)))
+    warning("Could not find ", targetMz, " in spectrum ", paste(which(!(allIdx %in% specIdx)), collapse = ", "), ".\n")
+    specIdx <- which(allIdx %in% specIdx)
   } else {
-    specIdx <- plot_Idx
+    specIdx <- allIdx
   }
   if (length(specIdx) < 1) {
     stop("Could not find targetMz in any spectrum! Consider adjusting tol.\n")
   }
   return(list(
-    norm_factor = pull(f_resdf, .data$int),
+    norm_factor = norm_factor,
     specIdx = specIdx
   ))
 }

R/getters.R

@@ -323,13 +323,12 @@ getPeakStatistics <- function(object, summarise = FALSE) {
 #' data(Blank2022res)
 #' getRecalibrationError(Blank2022res)
 getRecalibrationError <- function(object) {
-  peaks <- peaks2df(getAvgPeaks(object))
-  mzdev <- .getMzShift(
-    peaksdf = peaks,
+  mzdev <- getMzShift(
+    peaks = getAvgPeaks(object),
     tol = getNormMzTol(object),
     targetMz = getNormMz(object),
-    tolppm = FALSE,
-    allowNoMatch = TRUE
+    tolppm = FALSE, 
+    verbose = FALSE
   )
 
   res_df <- tibble(

R/helpers.R

@@ -24,7 +24,7 @@ transformConc2Log <- function(conc) {
 }
 
 
-#' Check if object is of class MALDIassay
+#' Check if object if of class MALDIassay
 #'
 #' @param object object to text
 #'