issue-21 adding validations for new functions

R/differential_expression.R

@@ -279,15 +279,17 @@ diffex_fdr <- function(term_data) {
 #' @returns a grob
 #'
 #' @examples
+#' library(dplyr)
 #' regression_significance <- diffex_mzroll(
-#'   nplug_mzroll_normalized,
-#'   "normalized_log2_abundance",
-#'   "limitation + limitation:DR + 0",
-#'   FDR_cutoff = 0.01,
-#'   feature_labels = c("UDP", "Lactate", "Serine")
-#' )
-#'
-#' plot_volcano(regression_significance, 10, 0.1)
+#' nplug_mzroll_normalized,
+#' "normalized_log2_abundance",
+#' "limitation + limitation:DR + 0") %>%
+#' dplyr::left_join(
+#' nplug_mzroll_normalized$features %>% select(groupId, compoundName),
+#' by = "groupId")
+#' 
+#' plot_volcano(regression_significance, 10, 0.1, c("Ribose-P", "acetyl-CoA", "ATP"))
+#' 
 #' @export
 plot_volcano <- function(
   regression_significance,
@@ -316,11 +318,12 @@ plot_volcano <- function(
       is_discovery = qvalue < FDR_cutoff
     ) %>%
     ggplot(aes_string(x = effect_var)) +
-    {if ("compoundName" %in% colnames(regression_significance))
-    {geom_point(aes(y = p.value.trans, color = is_discovery, name = compoundName))} 
-      else {geom_point(aes(y = p.value.trans, color = is_discovery))} 
-    } +
-    geom_text(aes(label = ifelse(compoundName %in% feature_labels, compoundName, ""), y = p.value.trans, vjust = -0.75)) +
+    {if ("compoundName" %in% colnames(regression_significance)) {
+      geom_point(aes(y = p.value.trans, color = is_discovery, name = compoundName)) +
+        geom_text(aes(label = ifelse(compoundName %in% feature_labels, compoundName, ""), y = p.value.trans, vjust = -0.75))
+      } 
+      else {geom_point(aes(y = p.value.trans, color = is_discovery))}
+      } +
     facet_wrap(~term, scales = "free_x") +
     scale_x_continuous("Effect size") +
     scale_y_continuous(expression(-log[10] ~ "pvalue")) +

R/mutate_mzroll_list.R

@@ -66,6 +66,15 @@ floor_peaks <- function(mzroll_list,
 #'@return triple omic data with imputed missing peaks
 #'
 #' @examples
+#' library(dplyr)
+#' 
+#' lod_values <- nplug_mzroll_augmented[["measurements"]] %>%
+#' dplyr::select(groupId, log2_abundance) %>%
+#' dplyr::distinct(groupId, .keep_all = TRUE)
+#' 
+#' mzroll_list <- nplug_mzroll_augmented
+#' mzroll_list$measurements <- mzroll_list$measurements %>% filter(groupId != 2 & sampleId != 1)
+#' 
 #' mzroll_list_imputed <- impute_missing_peaks(mzroll_list, lod_values, "log2_abundance", 0.15)
 #'
 #'@export
@@ -76,6 +85,18 @@ impute_missing_peaks <- function(mzroll_list,
 
   test_mzroll_list(mzroll_list)
 
+  stopifnot(colnames(lod_values) %in% c("groupId", rlang::sym(quant_var)))
+
+  if (nrow(lod_values) > nrow(lod_values %>% dplyr::distinct(groupId, .keep_all = TRUE))) {
+    stop("only one value per feature must be provided to impute missing peaks")
+  }
+
+  features <- mzroll_list$features %>% dplyr::select(groupId)
+
+  if(!all(features$groupId %in% lod_values$groupId) | nrow(features) != nrow(lod_values)) {
+    stop("groupId values of lod_value table and feature table of triple omic data must match")
+  }
+
   valid_quant_var <- setdiff(
     mzroll_list$design$measurements$variable,
     c(mzroll_list$design$feature_pk, mzroll_list$design$sample_pk)
@@ -94,12 +115,12 @@ impute_missing_peaks <- function(mzroll_list,
     )
 
   ## impute missing peaks
-  missing_peaks_imputed <- left_join(
+  missing_peaks_imputed <- dplyr::left_join(
     missing_peaks, 
     lod_values, 
     by = c("groupId")) %>%
-    rowwise() %>%
-    mutate(!!rlang::sym(quant_var) := rnorm(1, mean = !!rlang::sym(quant_var)+1, sd = imputation_sd))
+    dplyr::rowwise() %>%
+    dplyr::mutate(!!rlang::sym(quant_var) := stats::rnorm(1, mean = !!rlang::sym(quant_var)+1, sd = imputation_sd))
 
   ## merge measured peaks with imputed peaks
   completed_peaks <- dplyr::bind_rows(
@@ -125,6 +146,15 @@ impute_missing_peaks <- function(mzroll_list,
 #'@return triple omic data with imputed missing peaks
 #'
 #' @examples
+#' library(dplyr)
+#' 
+#' lod_values <- nplug_mzroll_augmented[["measurements"]] %>%
+#' dplyr::select(groupId, log2_abundance) %>%
+#' dplyr::distinct(groupId, .keep_all = TRUE)
+#' 
+#' mzroll_list <- nplug_mzroll_augmented
+#' mzroll_list$measurements <- mzroll_list$measurements %>% filter(groupId != 2 & sampleId != 1)
+#' 
 #' mzroll_list_imputed <- fill_in_missing_peaks(mzroll_list, lod_values, "log2_abundance", 0.15)
 #' mzroll_list_imputed <- fill_in_missing_peaks(mzroll_list, 12, "log2_abundance")
 #'
@@ -134,8 +164,7 @@ fill_in_missing_peaks <- function(mzroll_list,
                                   quant_var = "log2_abundance",
                                   imputation_sd = 0.15) {
   if (is.data.frame(fill_values)) {
-    stopifnot(colnames(fill_values) %in% c("groupId", rlang::sym(quant_var)))
-    output <- impute_missing_peaks(mzroll_list, lod_values, quant_var, imputation_sd)
+    output <- impute_missing_peaks(mzroll_list, fill_values, quant_var, imputation_sd)
   } else if (is.numeric(fill_values)) {
     output <- floor_peaks(mzroll_list, fill_values, quant_var)
   } else {

vignettes/NPLUG.Rmd

@@ -287,7 +287,7 @@ Two visualizations that are particularly useful for EDA are heatmaps (I'm sure y
 
 ```{r static_eda, fig.height = 8, fig.width = 8}
 samples_with_pcs <- final_processed_data %>%
-  romic::add_pca_loadings(value_var = "normalized_log2_abundance", npcs = 5)
+  romic::add_pcs(value_var = "normalized_log2_abundance", npcs = 5)
 
 romic::plot_bivariate(
   samples_with_pcs$samples,