Skip to content

Commit d09b58a

Browse files
Enable complex contrast strings in DESeq2 (#9473)
* Enable complex contrast strings * Update docker image * Add test case with limma contrast string * Format changes and add test with shrinkage
1 parent c984ada commit d09b58a

6 files changed

Lines changed: 357 additions & 112 deletions

File tree

modules/nf-core/deseq2/differential/environment.yml

Lines changed: 1 addition & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -5,3 +5,4 @@ channels:
55
- bioconda
66
dependencies:
77
- bioconda::bioconductor-deseq2=1.34.0
8+
- bioconda::bioconductor-limma=3.50.0

modules/nf-core/deseq2/differential/main.nf

Lines changed: 2 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -4,8 +4,8 @@ process DESEQ2_DIFFERENTIAL {
44

55
conda "${moduleDir}/environment.yml"
66
container "${ workflow.containerEngine == 'singularity' && !task.ext.singularity_pull_docker_container ?
7-
'https://depot.galaxyproject.org/singularity/bioconductor-deseq2:1.34.0--r41hc247a5b_3' :
8-
'biocontainers/bioconductor-deseq2:1.34.0--r41hc247a5b_3' }"
7+
'https://community-cr-prod.seqera.io/docker/registry/v2/blobs/sha256/a1/a15f5d61792b60b6179afd885db27d3fe60eb4c42e805c8887ed0416d88cb484/data' :
8+
'community.wave.seqera.io/library/bioconductor-deseq2_bioconductor-limma:b56a0c9ddc3e87e1' }"
99

1010
input:
1111
tuple val(meta), val(contrast_variable), val(reference), val(target), val(formula), val(comparison)

modules/nf-core/deseq2/differential/templates/deseq2_differential.R

Lines changed: 20 additions & 11 deletions
Original file line numberDiff line numberDiff line change
@@ -231,6 +231,7 @@ for (file_input in c('count_file', 'sample_file')){
231231

232232
library(DESeq2)
233233
library(BiocParallel)
234+
library(limma)
234235

235236
################################################
236237
################################################
@@ -425,17 +426,25 @@ dds <- DESeq(
425426

426427
if (!is.null(opt\$contrast_string)) {
427428
coef_names <- resultsNames(dds)
428-
if (!(opt\$contrast_string %in% coef_names)) {
429-
stop(sprintf(
430-
"Contrast '%s' not in design. Available coefficients: %s",
431-
opt\$contrast_string,
432-
paste(coef_names, collapse = ", ")
433-
))
434-
}
429+
if (opt\$contrast_string %in% coef_names) {
430+
# Direct coefficient name
431+
comp.results <- run_results(name = opt\$contrast_string)
432+
if (opt\$shrink_lfc) {
433+
comp.results <- run_shrink(coef = opt\$contrast_string)
434+
}
435+
} else {
436+
# Parse as limma-style contrast expression
437+
design_mat <- model.matrix(as.formula(model), data = as.data.frame(colData(dds)))
438+
colnames(design_mat) <- make.names(colnames(design_mat))
439+
numeric_contrast <- as.numeric(
440+
limma::makeContrasts(contrasts = opt\$contrast_string, levels = colnames(design_mat))
441+
)
435442

436-
comp.results <- run_results(name = opt\$contrast_string)
437-
if (opt\$shrink_lfc) {
438-
comp.results <- run_shrink(coef = opt\$contrast_string)
443+
# Run DESeq2 results with numeric contrast
444+
comp.results <- run_results(contrast = numeric_contrast)
445+
if (opt\$shrink_lfc) {
446+
comp.results <- run_shrink(contrast = numeric_contrast)
447+
}
439448
}
440449
} else {
441450
contrast_var_tg_ref <- c(contrast_variable,
@@ -452,7 +461,7 @@ if (!is.null(opt\$contrast_string)) {
452461

453462
if (!is.null(opt\$transcript_lengths_file)){
454463
size_factors = estimateSizeFactorsForMatrix(counts(dds) / assays(dds)[["avgTxLength"]])
455-
}else {
464+
} else {
456465
size_factors = sizeFactors(dds)
457466
}
458467

Lines changed: 5 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,5 @@
1+
process {
2+
withName: 'DESEQ2_DIFFERENTIAL' {
3+
ext.args = { "--round_digits 2 --sample_id_col sample --vs_method rlog --shrink_lfc TRUE --seed 1" }
4+
}
5+
}

modules/nf-core/deseq2/differential/tests/main.nf.test

Lines changed: 96 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -57,6 +57,102 @@ nextflow_process {
5757
}
5858
}
5959

60+
test("RNAseq - Feature Counts - formula + limma contrast string - interaction") {
61+
config './contrasts_interaction.config'
62+
63+
when {
64+
process {
65+
"""
66+
input[0] = Channel.of([
67+
'id': 'genotype_WT_KO',
68+
'formula': '~ 0 + genotype',
69+
'comparison': 'genotypeWT - genotypeKO'
70+
])
71+
.map {
72+
tuple(it, it.variable, it.reference, it.target, it.formula, it.comparison)
73+
}
74+
75+
input[1] = Channel.of([
76+
[ id: 'test' ],
77+
file("https://github.com/nf-core/test-datasets/raw/differentialabundance/modules_testdata/variancepartition_dream/metadata.tsv", checkIfExists: true),
78+
file("https://github.com/nf-core/test-datasets/raw/differentialabundance/modules_testdata/variancepartition_dream/counts.tsv", checkIfExists: true)
79+
])
80+
81+
ch_spikes = [[],[]]
82+
ch_lengths = [[],[]]
83+
input[2] = ch_spikes
84+
input[3] = ch_lengths
85+
"""
86+
}
87+
}
88+
89+
then {
90+
assertAll(
91+
{ assert process.success },
92+
{ assert snapshot(
93+
process.out.results,
94+
process.out.size_factors,
95+
process.out.normalised_counts,
96+
process.out.rlog_counts,
97+
process.out.vst_counts,
98+
process.out.model,
99+
process.out.versions,
100+
file(process.out.dispersion_plot[0][1]).name,
101+
file(process.out.rdata[0][1]).name
102+
).match() },
103+
{ assert path(process.out.session_info[0][1]).text.contains("DESeq2_1.34.0") }
104+
)
105+
}
106+
}
107+
108+
test("RNAseq - Feature Counts - formula + limma contrast string - shrunken") {
109+
config './contrasts_shrink.config'
110+
111+
when {
112+
process {
113+
"""
114+
input[0] = Channel.of([
115+
'id': 'genotype_WT_KO',
116+
'formula': '~ 0 + genotype',
117+
'comparison': 'genotypeWT - genotypeKO'
118+
])
119+
.map {
120+
tuple(it, it.variable, it.reference, it.target, it.formula, it.comparison)
121+
}
122+
123+
input[1] = Channel.of([
124+
[ id: 'test' ],
125+
file("https://github.com/nf-core/test-datasets/raw/differentialabundance/modules_testdata/variancepartition_dream/metadata.tsv", checkIfExists: true),
126+
file("https://github.com/nf-core/test-datasets/raw/differentialabundance/modules_testdata/variancepartition_dream/counts.tsv", checkIfExists: true)
127+
])
128+
129+
ch_spikes = [[],[]]
130+
ch_lengths = [[],[]]
131+
input[2] = ch_spikes
132+
input[3] = ch_lengths
133+
"""
134+
}
135+
}
136+
137+
then {
138+
assertAll(
139+
{ assert process.success },
140+
{ assert snapshot(
141+
process.out.results,
142+
process.out.size_factors,
143+
process.out.normalised_counts,
144+
process.out.rlog_counts,
145+
process.out.vst_counts,
146+
process.out.model,
147+
process.out.versions,
148+
file(process.out.dispersion_plot[0][1]).name,
149+
file(process.out.rdata[0][1]).name
150+
).match() },
151+
{ assert path(process.out.session_info[0][1]).text.contains("DESeq2_1.34.0") }
152+
)
153+
}
154+
}
155+
60156
test("mouse - contrasts - matrix") {
61157

62158
config './contrasts_matrix.config'

0 commit comments

Comments
 (0)