From 2f2bc35daf8e8ced29d4f4208fdf659093cc6abd Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?Johannes=20K=C3=B6ster?= <johannes.koester@tu-dortmund.de>
Date: Fri, 13 Oct 2017 14:23:40 +0200
Subject: [PATCH] parallelize deseq2 via config yaml.

---
 config.yaml           |  3 +++
 rules/diffexp.smk     |  2 ++
 scripts/deseq2-init.R | 12 ++++++++----
 scripts/deseq2.R      | 12 ++++++++----
 4 files changed, 21 insertions(+), 8 deletions(-)

diff --git a/config.yaml b/config.yaml
index 71d8bd5..af74bb1 100644
--- a/config.yaml
+++ b/config.yaml
@@ -23,3 +23,6 @@ params:
   star: ""
   cutadapt-se: ""
   cutadapt-pe: ""
+  deseq2:
+    # if you have many samples or conditions, you can use more than 1 thread here.
+    threads: 1
diff --git a/rules/diffexp.smk b/rules/diffexp.smk
index 271c697..1889f54 100644
--- a/rules/diffexp.smk
+++ b/rules/diffexp.smk
@@ -19,6 +19,7 @@ rule deseq2_init:
         "../envs/deseq2.yaml"
     log:
         "logs/deseq2/init.log"
+    threads: config["params"]["deseq2"]["threads"]
     script:
         "../scripts/deseq2-init.R"
 
@@ -54,5 +55,6 @@ rule deseq2:
         "../envs/deseq2.yaml"
     log:
         "logs/deseq2/{contrast}.diffexp.log"
+    threads: config["params"]["deseq2"]["threads"]
     script:
         "../scripts/deseq2.R"
diff --git a/scripts/deseq2-init.R b/scripts/deseq2-init.R
index cff6723..6e9f3c2 100644
--- a/scripts/deseq2-init.R
+++ b/scripts/deseq2-init.R
@@ -3,10 +3,14 @@ sink(log)
 sink(log, type="message")
 
 library("DESeq2")
-library("BiocParallel")
 
-# setup parallelization
-register(MulticoreParam(snakemake@threads))
+parallel <- FALSE
+if (snakemake@threads > 1) {
+    library("BiocParallel")
+    # setup parallelization
+    register(MulticoreParam(snakemake@threads))
+    parallel <- TRUE
+}
 
 # colData and countData must have the same sample order, but this is ensured
 # by the way we create the count matrix
@@ -20,6 +24,6 @@ dds <- DESeqDataSetFromMatrix(countData=cts,
 # remove uninformative columns
 dds <- dds[ rowSums(counts(dds)) > 1, ]
 # TODO optionally allow to collapse technical replicates
-dds <- DESeq(dds, parallel=TRUE)
+dds <- DESeq(dds, parallel=parallel)
 
 saveRDS(dds, file=snakemake@output[[1]])
diff --git a/scripts/deseq2.R b/scripts/deseq2.R
index cc08efb..3f7bdc4 100644
--- a/scripts/deseq2.R
+++ b/scripts/deseq2.R
@@ -3,15 +3,19 @@ sink(log)
 sink(log, type="message")
 
 library("DESeq2")
-library("BiocParallel")
 
-# setup parallelization
-register(MulticoreParam(snakemake@threads))
+parallel <- FALSE
+if (snakemake@threads > 1) {
+    library("BiocParallel")
+    # setup parallelization
+    register(MulticoreParam(snakemake@threads))
+    parallel <- TRUE
+}
 
 dds <- readRDS(snakemake@input[[1]])
 
 contrast <- c("condition", snakemake@params[["contrast"]])
-res <- results(dds, contrast=contrast, parallel=TRUE)
+res <- results(dds, contrast=contrast, parallel=parallel)
 # shrink fold changes for lowly expressed genes
 res <- lfcShrink(dds, contrast=contrast, res=res)
 # sort by p-value