add multiqc testing

.github/workflows/pipeline.yaml

@@ -130,6 +130,13 @@ jobs:
           conda activate default;
           snakemake --profile profiles/gh-actions --use-conda --conda-frontend mamba test_iqtree;
           conda deactivate
+      # Test IQTREE
+      - name: test multiqc
+        shell: bash -l {0}
+        run:
+          conda activate default;
+          snakemake --profile profiles/gh-actions --use-conda --conda-frontend mamba test_multiqc;
+          conda deactivate
       # Test Report Generation
       - name: test report
         shell: bash -l {0}

workflow/Snakefile

@@ -29,9 +29,10 @@ envs_dir = os.path.join(project_dir, "workflow", "envs")
 logs_dir = os.path.join(project_dir, "workflow", "logs")
 report_dir = os.path.join(project_dir, "workflow", "report")
 rules_dir = os.path.join(project_dir, "workflow", "rules")
+config_dir = os.path.join(project_dir, "config")
 
 # Include the config file
-configfile: os.path.join(project_dir, "config", "snakemake.yaml")
+configfile: os.path.join(config_dir, "snakemake.yaml")
 
 # Resources
 
@@ -52,9 +53,10 @@ report: report_dir + "/workflow.rst"
 
 rule all:
     """
-    The default pipeline target.
+    The default pipeline targets.
     """
     input:
+        results_dir + "/multiqc/multiqc_report.html",
         results_dir + "/iqtree/iqtree.core-filter" + str(config['snippy_missing_data']) + ".treefile",
 
 # -----------------------------------------------------------------------------#

workflow/envs/qc.yaml

@@ -7,3 +7,4 @@ dependencies:
   - samtools=1.9
   - qualimap=2.2.2d
   - multiqc=1.9
+  - pandoc=2.10.1

workflow/report/multiqc.rst

@@ -0,0 +1 @@
+MultiQC Report.

workflow/rules/alignment.smk

@@ -43,7 +43,7 @@ rule eager:
     os.path.join(envs_dir,"eager.yaml")
   log:
     html = os.path.join(logs_dir, "eager_{reads_origin}","{biosample}.html"),
-    #txt = os.path.join(logs_dir, "eager_{reads_origin}","{biosample}.log"),
+    txt = os.path.join(logs_dir, "eager_{reads_origin}","{biosample}.log"),
   shell:
     "cd {results_dir}/eager_{wildcards.reads_origin}; "
     "nextflow run nf-core/eager -r {config[eager_rev]} \
@@ -63,7 +63,7 @@ rule eager:
         --bam_unmapped_type discard \
         --max_cpus {threads} \
         --max_memory {resources.mem_mb}.MB \
-        -resume; "
+        -resume 1> {log.txt}; "
     "{scripts_dir}/eager_cleanup.sh {results_dir} {wildcards.reads_origin} {wildcards.biosample}; "
 
 # -----------------------------------------------------------------------------#
@@ -111,7 +111,7 @@ rule snippy_pairwise_bam:
         snp_txt = results_dir + "/snippy_pairwise_{reads_origin}/{sample}/{sample}_snippy.txt",
         snippy_aln = results_dir + "/snippy_pairwise_{reads_origin}/{sample}/{sample}_snippy.aligned.fa",
     threads:
-        workflow.cores,
+        (workflow.cores / 2) if (workflow.cores > 1) else workflow.cores
     log:
         os.path.join(logs_dir, "snippy_pairwise_{reads_origin}","{sample}.log")
     conda:

workflow/rules/phylogeny.smk

@@ -33,6 +33,6 @@ rule iqtree:
             -o {config[iqtree_outgroup]} \
             -seed {params.seed} \
             --runs {config[iqtree_runs]} \
-            -fconst `$(snp-sites -C {input.full_aln})` \
+            -fconst `snp-sites -C {input.full_aln}` \
             {config[iqtree_other]} \
             -pre {results_dir}/iqtree/iqtree.core-filter{config[snippy_missing_data]} 1>{log}"

workflow/rules/qc.smk

@@ -7,29 +7,38 @@ rule qualimap:
     Run qualimap metrics on the output of snippy pairwise.
     """
     input:
-        snippy_dir = results_dir + "/snippy_pairwise_assembly/{sample}",
+        snippy_dir = results_dir + "/snippy_pairwise_{reads_origin}/{sample}",
     output:
-        bamq = results_dir + "/qualimap/{sample}/{sample}.bam",
-        html = results_dir + "/qualimap/{sample}/qualimapReport.html"
+        dir = directory(results_dir + "/qualimap_{reads_origin}/{sample}"),
+        bamq = results_dir + "/qualimap_{reads_origin}/{sample}/{sample}.bam",
+        html = results_dir + "/qualimap_{reads_origin}/{sample}/qualimapReport.html",
     conda:
         os.path.join(envs_dir,"qc.yaml")
     log:
-        os.path.join(logs_dir, "qualimap/{sample}.log")
+        os.path.join(logs_dir, "qualimap_{reads_origin}/{sample}.log")
     threads:
-        workflow.cores
+        1
     shell:
         "samtools view -b -q {config[snippy_map_qual]} {input.snippy_dir}/{wildcards.sample}_snippy.bam > {output.bamq}; "
-        "qualimap bamqc -bam {output.bamq} --skip-duplicated -c -outformat 'HTML' -outdir {results_dir}/qualimap/{wildcards.sample} -nt {threads} 1> {log}; "
+        "qualimap bamqc -bam {output.bamq} --skip-duplicated -c -outformat 'HTML' -outdir {output.dir} -nt {threads} 1> {log}; "
 
 
 rule multiqc:
     """
     Run multiqc on miscellaneous data files.
     """
     input:
-        qualimap_asm_dir = expand(results_dir + "/qualimap/{sample}", sample=identify_assembly_sample()),
+        qualimap_asm_dir = expand(results_dir + "/qualimap_assembly/{sample}", sample=identify_assembly_sample()),
+        qualimap_local_dir = expand(results_dir + "/qualimap_local/{sample}", sample=identify_local_sample()),
+    output:
+        report(results_dir + "/multiqc/multiqc_report.html",
+                caption=os.path.join(report_dir,"multiqc.rst"),
+                category="Quality Control",
+                subcategory="MultiQC"),
+        dir = directory(results_dir + "/multiqc/"),
     conda:
         os.path.join(envs_dir,"qc.yaml")
+    log:
+        os.path.join(logs_dir, "multiqc/multiqc.log")
     shell:
-        "echo {input}; "
-        "multiqc --force {input.qualimap_asm_dir}"
+        "multiqc --pdf -c {config_dir}/multiqc.yaml --outdir {output.dir} --force {input.qualimap_asm_dir} {input.qualimap_local_dir} 2> {log}"