diff --git a/DIMS/GenerateViolinPlots.R b/DIMS/GenerateViolinPlots.R
index 9d1f3e7..12621a3 100644
--- a/DIMS/GenerateViolinPlots.R
+++ b/DIMS/GenerateViolinPlots.R
@@ -17,6 +17,7 @@ path_metabolite_groups <- cmd_args[3]
file_ratios_metabolites <- cmd_args[4]
file_expected_biomarkers_iem <- cmd_args[5]
file_explanation <- cmd_args[6]
+file_previous_runs <- cmd_args[7]
# load functions
source(paste0(export_scripts_dir, "generate_violin_plots_functions.R"))
@@ -32,6 +33,7 @@ expected_biomarkers_df <- expected_biomarkers_df %>%
HMDB_name = Metabolite
)
explanation_violin_plot <- readLines(file_explanation)
+data_previous_runs <- read.delim(file_previous_runs)
# Set global variables
iem_variables <- list(
@@ -79,6 +81,7 @@ make_and_save_violin_plot_pdfs(
run_name,
protocol_name,
explanation_violin_plot,
+ data_previous_runs,
number_of_metabolites
)
diff --git a/DIMS/GenerateViolinPlots.nf b/DIMS/GenerateViolinPlots.nf
index ec65a2e..2f843e7 100755
--- a/DIMS/GenerateViolinPlots.nf
+++ b/DIMS/GenerateViolinPlots.nf
@@ -22,6 +22,7 @@ process GenerateViolinPlots {
$params.path_metabolite_groups \
$params.file_ratios_metabolites \
$params.file_expected_biomarkers_IEM \
- $params.file_explanation
+ $params.file_explanation \
+ $params.file_previous_runs
"""
}
diff --git a/DIMS/export/generate_violin_plots_functions.R b/DIMS/export/generate_violin_plots_functions.R
index c314ce7..04ff040 100644
--- a/DIMS/export/generate_violin_plots_functions.R
+++ b/DIMS/export/generate_violin_plots_functions.R
@@ -131,6 +131,7 @@ calculate_zscore_ratios <- function(metabolites_ratios_df, intensities_zscores_d
#' @param run_name: string containing the run name
#' @param protocol_name: string containing the protocol name
#' @param explanation_violin_plot: vector of strings containing the explanation of the violin plots
+#' @param data_previous_runs: data from previous DIMS runs (matrix)
#' @param number_of_metabolites: list containing the number of metabolites for the top and lowest table
make_and_save_violin_plot_pdfs <- function(
zscore_patients_df,
@@ -141,6 +142,7 @@ make_and_save_violin_plot_pdfs <- function(
run_name,
protocol_name,
explanation_violin_plot,
+ data_previous_runs,
number_of_metabolites) {
# Get all patient IDs
patient_col_names <- remove_suffix_from_items(get_colnames_by_prefix(zscore_patients_df, "P"), "_Zscore")
@@ -191,7 +193,14 @@ make_and_save_violin_plot_pdfs <- function(
)
}
# generate normal violin plots
- create_pdf_violin_plots(pdf_dir, patient_id, metab_perpage, top_metabs_patient, explanation_violin_plot)
+ create_pdf_violin_plots(
+ pdf_dir,
+ patient_id,
+ metab_perpage,
+ top_metabs_patient,
+ explanation_violin_plot,
+ data_previous_runs
+ )
}
}
}
@@ -597,6 +606,28 @@ prepare_toplist <- function(patient_id, zscore_patients, num_of_highest_metaboli
return(top_metab_pt)
}
+#' Add data from previous runs (min, max, 5%, 95%) to patient_zscore_df
+#'
+#' @param patient_zscore_df: dataframe with metabolite Z-scores per patient (dataframe)
+#' @param data_previous_runs: data from previous DIMS runs (matrix)
+#'
+#' @return patient_zscore_df: dataframe with metabolite Z-scores and info from previous runs (dataframe)
+add_previous_runs <- function(patient_zscore_df, data_previous_runs) {
+ data_previous_runs <- as.data.frame(data_previous_runs)
+ for (row_nr in 1:nrow(patient_zscore_df)) {
+ # match data from previous run to data from current run based on common name
+ common_name <- trimws(patient_zscore_df$HMDB_name[row_nr])
+ find_rownr <- which(data_previous_runs$common_name == common_name)
+ # add mean, min, max and 5-95% interval of data from previous runs
+ if (length(find_rownr) == 1) {
+ patient_zscore_df[row_nr, c("mean", "min", "max", "min95", "max95")] <-
+ data_previous_runs[find_rownr, c("average_patients", "min", "max", "min95", "max95")]
+ }
+ }
+
+ return(patient_zscore_df)
+}
+
#' Create a pdf with table with metabolites and violin plots
#'
#' @param pdf_dir: location where to save the pdf file (string)
@@ -604,7 +635,8 @@ prepare_toplist <- function(patient_id, zscore_patients, num_of_highest_metaboli
#' @param metab_perpage: list of dataframes, each dataframe contains data for a page in de pdf (list)
#' @param top_metab_pt: dataframe with increased and decreased metabolites for this patient (dataframe)
#' @param explanation: text that explains the violin plots and the pipeline version (string)
-create_pdf_violin_plots <- function(pdf_dir, patient_id, metab_perpage, top_metab_pt, explanation) {
+#' @param data_previous_runs: data from previous DIMS runs, not used for dIEM plots (matrix)
+create_pdf_violin_plots <- function(pdf_dir, patient_id, metab_perpage, top_metab_pt, explanation, data_previous_runs = NULL) {
# set parameters for plots
plot_height <- 9.6
plot_width <- 6
@@ -680,6 +712,10 @@ create_pdf_violin_plots <- function(pdf_dir, patient_id, metab_perpage, top_meta
# extract original data for patient of interest (pt_name)
patient_zscore_df <- metab_zscores_df %>%
filter(Sample == patient_id)
+ # add information for 5-95% interval of data from previous runs
+ if (!is.null(data_previous_runs)) {
+ patient_zscore_df <- add_previous_runs(patient_zscore_df, data_previous_runs)
+ }
# Remove patient of interest and retain only other patient data
metab_zscores_df <- metab_zscores_df %>%
@@ -734,7 +770,34 @@ create_violin_plot <- function(metab_zscores_df, patient_zscore_df, sub_perpage,
mutate(HMDB_name = factor(HMDB_name, levels = y_order)) %>%
arrange(HMDB_name)
+ # find y-axis values for each metabolite
+ if (any(grepl("min95", colnames(patient_zscore_df)))) {
+ rectange_df <- patient_zscore_df |>
+ dplyr::distinct(HMDB_name, min95, max95) |>
+ dplyr::mutate(
+ y = match(HMDB_name, y_order),
+ ymin = y - 0.4,
+ ymax = y + 0.4
+ )
+ } else {
+ rectange_df <- patient_zscore_df
+ rectange_df[c("min95", "max95", "y", "ymin", "ymax")] <- 0
+ }
+
ggplot_object <- ggplot(metab_zscores_df, aes(x = Z_score, y = HMDB_name)) +
+ geom_rect(
+ data = rectange_df,
+ aes(
+ xmin = min95,
+ xmax = max95,
+ ymin = ymin,
+ ymax = ymax
+ ),
+ inherit.aes = FALSE,
+ fill = "blue",
+ alpha = 0.15,
+ colour = NA
+ ) +
# Make violin plots
geom_violin(scale = "width", na.rm = TRUE) +
# Add Z-score for the selected patient, shape=22 gives square for patient of interest
diff --git a/DIMS/tests/testthat/_snaps/generate_violin_plots/violin-plot-p2025m1.svg b/DIMS/tests/testthat/_snaps/generate_violin_plots/violin-plot-p2025m1.svg
index 89edec1..b0e3045 100644
--- a/DIMS/tests/testthat/_snaps/generate_violin_plots/violin-plot-p2025m1.svg
+++ b/DIMS/tests/testthat/_snaps/generate_violin_plots/violin-plot-p2025m1.svg
@@ -32,28 +32,30 @@
-
-
+
+
-
-
-
-
-Z=2.34
-Z=0.31
+
+
+
+
+
+
+Z=2.34
+Z=0.31
-metab3
-metab1
-
-
+metab3
+metab1
+
+
diff --git a/DIMS/tests/testthat/test_generate_violin_plots.R b/DIMS/tests/testthat/test_generate_violin_plots.R
index d5fe3a8..8017cca 100644
--- a/DIMS/tests/testthat/test_generate_violin_plots.R
+++ b/DIMS/tests/testthat/test_generate_violin_plots.R
@@ -508,13 +508,16 @@ testthat::test_that("create_pdf_violin_plots: Create a pdf with a table of top m
Metabolite = c("Increased", "metab1", "Decreased", "metab11"),
`Z-score` = c("", "2.45", "", "-1.51")
)
+
+ test_data_previous_runs <- NULL
expect_silent(create_pdf_violin_plots(
test_pdf_dir,
test_patient_id,
test_metab_perpage,
test_top_metab_pt,
- test_explanation
+ test_explanation,
+ test_data_previous_runs
))
out_pdf_violinplots <- file.path(test_pdf_dir, "R_P2025M1.pdf")
@@ -788,6 +791,7 @@ testthat::test_that("make_and_save_violin_plot_pdfs: Make and save violin plots
highest = 2,
lowest = 1
)
+ test_data_previous_runs <- NULL
expect_silent(make_and_save_violin_plot_pdfs(
test_zscore_patients_df,
@@ -798,6 +802,7 @@ testthat::test_that("make_and_save_violin_plot_pdfs: Make and save violin plots
test_run_name,
test_protocol_name,
test_explanation_violin_plot,
+ test_data_previous_runs,
test_number_of_metabolites
))
@@ -960,3 +965,16 @@ testthat::test_that("save_patient_no_iem: Save a list of patient IDs to a text f
expect_snapshot_file("missing_probability_scores.txt")
file.remove("missing_probability_scores.txt")
})
+
+testthat::test_that("add_previous_runs: Information for metabolites from previous runs is correctly added", {
+ test_acyl_carnitines_df <- read.delim(test_path("fixtures/", "test_acyl_carnitines_df.txt"))
+ test_patient_id <- "P2025M1"
+ test_patient_zscore_df <- test_acyl_carnitines_df %>% filter(Sample == test_patient_id)
+ test_data_previous_runs <- data.frame(matrix(1:12, ncol = 6, nrow = 2))
+ colnames(test_data_previous_runs) <- c("common_name", "average_patients", "min95", "max95", "min", "max")
+ test_data_previous_runs[, "common_name"] <- c("metab1", "metab3")
+
+ expect_equal(add_previous_runs(test_patient_zscore_df, test_data_previous_runs)$min95[1], 5)
+ expect_equal(add_previous_runs(test_patient_zscore_df, test_data_previous_runs)$max95[2], 8)
+})
+