Fix rule thread assignment at parse time

Author: lczech

workflow/Snakefile

@@ -126,14 +126,11 @@ def make_resource_fn(rule_name, kind):
     return _fn
 
 
-def get_cpus(rule_name):
-    return int(resources_config.get(rule_name, {}).get("cpus", resources_config["default"]["cpus"]))
-
-
 # Set the resources for all rules automatically,
 # without having to specify this for all of them individually.
+# We also re-set the threads here, just in case.
 # Cannot name the iteration variable `rule` here, as that conflicts...
 for wf_rule in workflow.rules:
     wf_rule.resources["mem_mb"]  = make_resource_fn(wf_rule.name, "mem")
     wf_rule.resources["runtime"] = make_resource_fn(wf_rule.name, "time")
-    wf_rule.threads = get_cpus(wf_rule.name)
+    wf_rule.threads = get_rule_threads(wf_rule.name)

workflow/rules/annotation.smk

@@ -90,6 +90,7 @@ rule snpeff:
     params:
         # For finding the chromosome names used by snpeff, add `-v` here
         extra=config["params"]["snpeff"]["extra"],
+    threads: get_rule_threads("snpeff")
     conda:
         "../envs/snpeff.yaml"
     wrapper:
@@ -224,6 +225,7 @@ rule vep:
         # Plugin args can be added as well, e.g. via an entry "MyPlugin,1,FOO", see docs.
         plugins=config["params"]["vep"]["plugins"],
         extra=config["params"]["vep"]["extra"],
+    threads: get_rule_threads("vep")
     log:
         "logs/annotation/vep-annotate.log",
     conda:

workflow/rules/calling-bcftools-combined.smk

@@ -41,7 +41,7 @@ rule call_variants:
         mpileup=get_mpileup_params,
         # Optional parameters for bcftools call (except -v, -o, -m).
         call=config["params"]["bcftools"]["call"],
-    threads: 1 # Dummy, but will be overwritten by our automatic resources
+    threads: get_rule_threads("call_variants")
     log:
         "logs/calling/bcftools-call/{contig}.log",
     benchmark:

workflow/rules/calling-bcftools-individual.smk

@@ -40,7 +40,7 @@ rule call_variants:
         mpileup=get_mpileup_params,
         # Optional parameters for bcftools call (except -v, -o, -m).
         call=config["params"]["bcftools"]["call"],
-    threads: 1 # Dummy, but will be overwritten by our automatic resources
+    threads: get_rule_threads("call_variants")
     log:
         "logs/calling/bcftools-call/{sample}-{contig}.log",
     benchmark:

workflow/rules/calling-freebayes.smk

@@ -50,15 +50,16 @@ rule call_variants:
     output:
         # touch("calling/called/{contig}.vcf.done"),
         pipe("calling/called/{contig}.vcf"),
-    log:
-        "logs/calling/freebayes/{contig}.log",
-    benchmark:
-        "benchmarks/calling/freebayes/{contig}.log"
     params:
         # Optional extra parameters.
         extra=config["params"]["freebayes"]["extra"] + know_variants_extra(),
         # Reference genome chunk size for parallelization (default: 100000)
         chunksize=config["params"]["freebayes"]["chunksize"],
+    threads: get_rule_threads("call_variants")
+    log:
+        "logs/calling/freebayes/{contig}.log",
+    benchmark:
+        "benchmarks/calling/freebayes/{contig}.log"
     group:
         "call_variants"
     # wrapper:
@@ -83,7 +84,7 @@ rule compress_vcf:
         ),
         # protected("calling/called/{contig}.vcf.gz")
         touch("calling/called/{contig}.vcf.gz.done"),
-    threads: 1 # Dummy, but will be overwritten by our automatic resources
+    threads: 1  # Dummy, but will be overwritten by our automatic resources
     log:
         "logs/calling/compress-vcf/{contig}.log",
     group:

workflow/rules/calling-haplotypecaller.smk

@@ -63,16 +63,17 @@ rule call_variants:
             else temp("calling/called/{sample}-{contig}.g.vcf.gz.tbi")
         ),
         done=touch("calling/called/{sample}-{contig}.g.vcf.gz.done"),
-    log:
-        "logs/calling/gatk-haplotypecaller/{sample}-{contig}.log",
-    benchmark:
-        "benchmarks/calling/gatk-haplotypecaller/{sample}-{contig}.log"
     params:
         # The intervals param here is where the contig variable is propagated to haplotypecaller.
         # Contigs are used as long as no restrict-regions are given in the config file.
         intervals=get_gatk_intervals,
         extra=config["params"]["gatk"].get("HaplotypeCaller-extra", ""),
         java_opts=config["params"]["gatk"].get("HaplotypeCaller-java-opts", ""),
+    threads: get_rule_threads("call_variants")
+    log:
+        "logs/calling/gatk-haplotypecaller/{sample}-{contig}.log",
+    benchmark:
+        "benchmarks/calling/gatk-haplotypecaller/{sample}-{contig}.log"
     group:
         "call_variants"
     conda:
@@ -142,10 +143,6 @@ rule genomics_db_import:
     output:
         db=directory("calling/genomics_db/{contig}"),
         done=touch("calling/genomics_db/{contig}.done"),
-    log:
-        "logs/calling/gatk-genomicsdbimport/{contig}.log",
-    benchmark:
-        "benchmarks/calling/gatk-genomicsdbimport/{contig}.log"
     params:
         # Here, we actually use the intervals to provide them to the wrapper.
         intervals=get_gatk_intervals,
@@ -157,6 +154,11 @@ rule genomics_db_import:
         + " "
         + config["params"]["gatk"].get("GenomicsDBImport-extra", ""),
         java_opts=config["params"]["gatk"].get("GenomicsDBImport-java-opts", ""),
+    threads: get_rule_threads("genomics_db_import")
+    log:
+        "logs/calling/gatk-genomicsdbimport/{contig}.log",
+    benchmark:
+        "benchmarks/calling/gatk-genomicsdbimport/{contig}.log"
     resources:
         tmpdir=config["params"]["gatk"].get("GenomicsDBImport-temp-dir", ""),
     conda:
@@ -204,6 +206,7 @@ rule combine_calls:
             else ""
         ),
         java_opts=config["params"]["gatk"]["CombineGVCFs-java-opts"],
+    threads: get_rule_threads("combine_calls")
     log:
         "logs/calling/gatk-combine-gvcfs/{contig}.log",
     benchmark:
@@ -274,6 +277,7 @@ rule genotype_variants:
         + " "
         + config["params"]["gatk"]["GenotypeGVCFs-extra"],
         java_opts=config["params"]["gatk"]["GenotypeGVCFs-java-opts"],
+    threads: get_rule_threads("genotype_variants")
     log:
         "logs/calling/gatk-genotype-gvcfs/{contig}.log",
     benchmark:
@@ -325,6 +329,7 @@ rule merge_variants:
             if platform.system() == "Darwin"
             else ""
         ),
+    threads: get_rule_threads("merge_variants")
     log:
         "logs/calling/picard-merge-genotyped.log",
     benchmark:

workflow/rules/duplicates-dedup.smk

@@ -71,7 +71,7 @@ rule sort_reads_dedup:
         extra=config["params"]["samtools"]["sort"],
         tmp_dir=config["params"]["samtools"]["temp-dir"],
     # Samtools takes additional threads through its option -@
-    threads: 1  # This value - 1 will be sent to -@. Weird flex, but okay.
+    threads: get_rule_threads("sort_reads_dedup")
     log:
         "logs/mapping/dedup-samtools-sort/{sample}.log",
     benchmark:

workflow/rules/duplicates-picard.smk

@@ -23,10 +23,6 @@ rule mark_duplicates:
         ),
         metrics="qc/dedup/{sample}.metrics.txt",
         done=touch("mapping/dedup/{sample}.bam.done"),
-    log:
-        "logs/mapping/picard-markduplicates/{sample}.log",
-    benchmark:
-        "benchmarks/mapping/picard-markduplicates/{sample}.log"
     params:
         # Take the params from the config.
         # On MacOS (we experienced it with 10.16, 11, and 12 so far), there is an issue between Java
@@ -40,6 +36,11 @@ rule mark_duplicates:
             else ""
         ),
         java_opts=config["params"]["picard"]["MarkDuplicates-java-opts"],
+    threads: get_rule_threads("mark_duplicates")
+    log:
+        "logs/mapping/picard-markduplicates/{sample}.log",
+    benchmark:
+        "benchmarks/mapping/picard-markduplicates/{sample}.log"
     group:
         "mapping_extra"
     conda:

workflow/rules/filtering-gatk-variantfiltration.smk

@@ -29,6 +29,7 @@ rule gatk_hard_filter_calls:
         filters=get_filter,
         extra=config["params"]["gatk-variantfiltration"]["extra"],
         java_opts=config["params"]["gatk-variantfiltration"]["java-opts"],
+    threads: get_rule_threads("gatk_hard_filter_calls")
     log:
         "logs/calling/gatk-variantfiltration/{vartype}.log",
     benchmark:

workflow/rules/filtering-gatk-vqsr.smk

@@ -101,6 +101,7 @@ rule gatk_variant_recalibrator:
         # Extras
         extra=get_variant_recalibrator_extra,
         java_opts=config["params"]["gatk-vqsr"]["variantrecalibrator-java-opts"],
+    threads: get_rule_threads("gatk_variant_recalibrator")
     log:
         "logs/calling/gatk-variantrecalibrator/{vartype}.log",
     benchmark:
@@ -132,15 +133,16 @@ rule gatk_apply_vqsr:
             else temp("calling/filtered/all.{vartype}.recalibrated.vcf.gz")
         ),
         done=touch("calling/filtered/all.{vartype}.recalibrated.vcf.gz.done"),
-    log:
-        "logs/calling/gatk-applyvqsr/{vartype}.log",
-    benchmark:
-        "benchmarks/calling/gatk-applyvqsr/{vartype}.log"
     params:
         # set mode, must be either SNP, INDEL or BOTH
         mode="{vartype}",
         extra=get_apply_vqsr_extra,
         java_opts=config["params"]["gatk-vqsr"]["applyvqsr-java-opts"],
+    threads: get_rule_threads("gatk_apply_vqsr")
+    log:
+        "logs/calling/gatk-applyvqsr/{vartype}.log",
+    benchmark:
+        "benchmarks/calling/gatk-applyvqsr/{vartype}.log"
     conda:
         # We overwrite the original yaml, as this wrapper here (version 0.85.0) and the one above
         # for the variantrecalibrator (also 0.85.0) use different GATK versions originally...