From 2486e2529f70e2029a5faf1ab629d02dcb5ade9f Mon Sep 17 00:00:00 2001 From: Gregor Sturm Date: Mon, 11 May 2026 16:20:58 +0200 Subject: [PATCH 1/4] Simplify dream contrast handling --- .../variancepartition/dream/templates/dream.R | 84 ++++++++----------- 1 file changed, 36 insertions(+), 48 deletions(-) diff --git a/modules/nf-core/variancepartition/dream/templates/dream.R b/modules/nf-core/variancepartition/dream/templates/dream.R index d8dcb493922f..b14d87a2d6de 100644 --- a/modules/nf-core/variancepartition/dream/templates/dream.R +++ b/modules/nf-core/variancepartition/dream/templates/dream.R @@ -272,8 +272,39 @@ if (as.logical(opt\$apply_voom)) { vobjDream <- countMatrix } -# Fit the DREAM model with ddf and reml options -fitmm <- dream(vobjDream, form, metadata, ddf = opt\$ddf, reml = opt\$reml, BPPARAM = bp) +# If contrast_string is provided, use that for makeContrast +if (!is.null(opt\$contrast_string)) { + contrast_string <- opt\$contrast_string +} else { + # Construct the contrast_string + treatment_target <- make.names(paste0(opt\$contrast_variable, opt\$contrast_target)) + treatment_reference <- make.names(paste0(opt\$contrast_variable, opt\$contrast_reference)) + contrast_string <- paste0(treatment_target, "-", treatment_reference) +} + +# define contrasts outside DREAM +design = model.matrix(reformulas::nobars(form), metadata) + +# print diagnostic output of design matrix +head(design, 3) +cat("Raw coefficient names:\n") +print(colnames(design)) + +# Use "normalized" coefficient names to compute contrast vector +normalized_coef_names = make.names(colnames(design)) +cat("Coefficient names used for contrasts:", paste(normalized_coef_names, collapse = ", "), "\n") + +contrast_mat = makeContrasts(contrast = contrast_string, levels = normalized_coef_names) +# however, reset the coefficient names in the contrast matrix to use the "original" contrast names used by DREAM +# DREAM will check if the coefficient names match its internal design matrix, so there is no danger that +# design() and dream() will do things differently without noticing. +rownames(contrast_mat) = colnames(design) +stopifnot("exactly one contrast defined"=ncol(contrast_mat) == 1) + +# Fit the DREAM model with ddf and reml options. +# Specifying contrast matrix directly here saves time refitting the model and avoids downstream issues +# with repeated calling of eBayes +fitmm <- dream(vobjDream, form, metadata, ddf = opt\$ddf, reml = opt\$reml, BPPARAM = bp, L=contrast_mat) # Parse stdev_coef_lim and winsor_tail_p into numeric vectors stdev_coef_lim_vals <- as.numeric(unlist(strsplit(opt\$stdev_coef_lim, ","))) @@ -285,52 +316,9 @@ fitmm <- eBayes(fitmm, proportion = opt\$proportion, trend = opt\$trend, robust = opt\$robust, winsor.tail.p = winsor_tail_p_vals) -# Display design matrix -head(fitmm\$design, 3) -cat("Raw coefficient names from dream():\n") -print(colnames(fitmm\$design)) - -# Normalize coefficient names consistently before building contrasts -normalized_coef_names <- make.names(colnames(fitmm\$design)) -colnames(fitmm\$design) <- normalized_coef_names - -if (!is.null(colnames(fitmm\$coefficients))) { - colnames(fitmm\$coefficients) <- normalized_coef_names -} -if (!is.null(colnames(fitmm\$stdev.unscaled))) { - colnames(fitmm\$stdev.unscaled) <- normalized_coef_names -} -if (!is.null(fitmm\$cov.coefficients)) { - rownames(fitmm\$cov.coefficients) <- normalized_coef_names - colnames(fitmm\$cov.coefficients) <- normalized_coef_names -} - -cat("Coefficient names used for contrasts:", paste(normalized_coef_names, collapse = ", "), "\n") - -# If contrast_string is provided, use that for makeContrast -if (!is.null(opt\$contrast_string)) { - contrast_string <- opt\$contrast_string -} else { - # Construct the contrast_string - treatment_target <- make.names(paste0(opt\$contrast_variable, opt\$contrast_target)) - treatment_reference <- make.names(paste0(opt\$contrast_variable, opt\$contrast_reference)) - contrast_string <- paste0(treatment_target, "-", treatment_reference) -} - -# Use makeContrasts if contrast_string exists -if (is_valid_string(contrast_string)) { - cat("Using contrast string:", contrast_string, "\n") - - contrast_matrix <- makeContrasts(contrast = contrast_string, levels = normalized_coef_names) - fit2 <- contrasts.fit(fitmm, contrast_matrix) - fit2 <- eBayes(fit2, proportion = opt\$proportion, - stdev.coef.lim = stdev_coef_lim_vals, - trend = opt\$trend, robust = opt\$robust, - winsor.tail.p = winsor_tail_p_vals) - results <- topTable(fit2, number = Inf, - adjust.method = opt\$adjust.method, - p.value = opt\$p.value, lfc = opt\$lfc, confint = opt\$confint) -} +results <- topTable(fit2, coef=colnames(contrast_mat), number = Inf, + adjust.method = opt\$adjust.method, + p.value = opt\$p.value, lfc = opt\$lfc, confint = opt\$confint) results\$gene_id <- rownames(results) results <- results[, c("gene_id", setdiff(names(results), "gene_id"))] From d995b36e19d7a8fd4f845356dc866ab126f41aae Mon Sep 17 00:00:00 2001 From: Gregor Sturm Date: Mon, 11 May 2026 16:40:10 +0200 Subject: [PATCH 2/4] fix pre-commit --- modules/nf-core/variancepartition/dream/templates/dream.R | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/modules/nf-core/variancepartition/dream/templates/dream.R b/modules/nf-core/variancepartition/dream/templates/dream.R index b14d87a2d6de..827f5a871823 100644 --- a/modules/nf-core/variancepartition/dream/templates/dream.R +++ b/modules/nf-core/variancepartition/dream/templates/dream.R @@ -296,8 +296,8 @@ cat("Coefficient names used for contrasts:", paste(normalized_coef_names, collap contrast_mat = makeContrasts(contrast = contrast_string, levels = normalized_coef_names) # however, reset the coefficient names in the contrast matrix to use the "original" contrast names used by DREAM -# DREAM will check if the coefficient names match its internal design matrix, so there is no danger that -# design() and dream() will do things differently without noticing. +# DREAM will check if the coefficient names match its internal design matrix, so there is no danger that +# design() and dream() will do things differently without noticing. rownames(contrast_mat) = colnames(design) stopifnot("exactly one contrast defined"=ncol(contrast_mat) == 1) From f30d5733f7ba1e7dd242585122ed1521f9ba37e1 Mon Sep 17 00:00:00 2001 From: Gregor Sturm Date: Mon, 11 May 2026 17:03:13 +0200 Subject: [PATCH 3/4] Use lme4 as source for nobars --- modules/nf-core/variancepartition/dream/templates/dream.R | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/modules/nf-core/variancepartition/dream/templates/dream.R b/modules/nf-core/variancepartition/dream/templates/dream.R index 827f5a871823..f5e286d48514 100644 --- a/modules/nf-core/variancepartition/dream/templates/dream.R +++ b/modules/nf-core/variancepartition/dream/templates/dream.R @@ -283,7 +283,7 @@ if (!is.null(opt\$contrast_string)) { } # define contrasts outside DREAM -design = model.matrix(reformulas::nobars(form), metadata) +design = model.matrix(lme4::nobars(form), metadata) # print diagnostic output of design matrix head(design, 3) From 3e79fb3d24de31e4865ae8e05d65f0d184888ee8 Mon Sep 17 00:00:00 2001 From: Gregor Sturm Date: Mon, 11 May 2026 17:15:25 +0200 Subject: [PATCH 4/4] Apply suggestion from @delfiterradas Co-authored-by: Delfina Terradas <155591053+delfiterradas@users.noreply.github.com> --- modules/nf-core/variancepartition/dream/templates/dream.R | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/modules/nf-core/variancepartition/dream/templates/dream.R b/modules/nf-core/variancepartition/dream/templates/dream.R index f5e286d48514..65a88aac8f99 100644 --- a/modules/nf-core/variancepartition/dream/templates/dream.R +++ b/modules/nf-core/variancepartition/dream/templates/dream.R @@ -316,7 +316,7 @@ fitmm <- eBayes(fitmm, proportion = opt\$proportion, trend = opt\$trend, robust = opt\$robust, winsor.tail.p = winsor_tail_p_vals) -results <- topTable(fit2, coef=colnames(contrast_mat), number = Inf, +results <- topTable(fitmm, coef=colnames(contrast_mat), number = Inf, adjust.method = opt\$adjust.method, p.value = opt\$p.value, lfc = opt\$lfc, confint = opt\$confint)