feat: qc rules

workflow/Snakefile

@@ -24,7 +24,9 @@ rule all:
         dispersion_graph = "de_analysis/dispersion_graph.svg",
         ma_graph         = "de_analysis/ma_graph.svg",
         de_heatmap       = "de_analysis/heatmap.svg",
-        lfc_analysis     = "de_analysis/lfc_analysis.csv"
+        lfc_analysis     = "de_analysis/lfc_analysis.csv",
+        QC="NanoPlot.tar.gz"
+
 
 rule build_minimap_index: ## build minimap2 index
     input:
@@ -56,6 +58,83 @@ rule map_reads:
     minimap2 -t {resources.cpus_per_task} -ax map-ont -p {params.psec} -N {params.msec} {params.opts} {input.index} {input.fastq} > {output} 2> {log}
     """
 
+#QC and metadata with NanoPlot
+if config["summary"] == "None":
+    rule plot_samples:
+        input:
+            fastq = lambda wildcards: all_samples[wildcards.sample]
+        output:
+            directory("NanoPlot")
+        log:
+            "logs/NanoPlot/{sample}.log"
+        resources:
+            ## max of 39 for our SLURM partition
+            cpus_per_task= max(8,39) #problem with max(len(input.fastq),39)
+        conda:
+            "envs/env.yml"
+        shell: 
+            "mkdir {output}/{wildcards.sample}; "
+            "NanoPlot -t {resources.cpus_per_task} --tsv_stats "
+            "--fastq {input.fastq} -o {output}/{wildcards.sample} 2> {log}"
+else:
+    rule plot_samples:
+        input:
+            summary = config["summary"]
+        output:
+            directory("NanoPlot")
+        log:
+            "logs/NanoPlot/NanoPlot.log"
+        conda:
+            "envs/env.yml"
+        shell: 
+            "mkdir {output}/{wildcards.sample}; "
+            "NanoPlot -t {resources.cpus_per_task} --barcoded --tsv_stats "
+            "--summary {input.summary} -o {output} 2> {log}"
+
+rule compress_nplot:
+    input:
+        samples=rules.plot_samples.output
+    output:
+        "NanoPlot.tar.gz"
+    log:
+        "logs/NanoPlot/compress.log"
+    shell:
+        "tar zcvf {input}.tar.gz {input} 2> {log}"
+
+rule mapping_qc:
+    input: 
+       bam=rules.map_reads.output
+    output:
+       "QC/qualimap/{sample}"
+    log:
+       "logs/qualimap/{sample}.log"
+    conda:
+       "envs/env.yml"
+    envmodules:
+       "bio/qualimap/2.3"
+    shell:"""
+        qualimap -t {resources.cpus_per_task}  {input.bam} > {output} 2> {log}
+        """
+
+rule compress_mapqc:
+    input:
+        mapping=rules.mapping_qc.output
+    output:
+        "Qualimap.tar.gz"
+    log:
+        "logs/Qualimap/compress.log"
+    shell:
+        "tar zcvf {input}.tar.gz {input} 2> {log}"
+
+#rule qc_report:
+#    input:
+#        samples=rules.compress_nplot.output#,
+#        #mapping=rules.mapping_qc.output
+#    output:
+#        "QC/QC-report.gz"
+#    shell:
+#        "mv {input.samples} > {output}" #{input.summary}
+
 rule sam_sort:
     input:
         sam = rules.map_reads.output

workflow/config/config_template.yml

@@ -19,6 +19,9 @@ transcriptome: "/lustre/miifs01/project/m2_zdvhpc/transcriptome_data/GCA_9176273
 annotation: "/lustre/miifs01/project/m2_zdvhpc/transcriptome_data/GCA_917627325.4_PGI_CHIRRI_v4_genomic.gff"
 # these samples ought to contain all samples comprising of the
 
+#NanoPlot QC with sequencing summary. leave as None for QC with .fastq files
+summary: None
+
 # Minimap2 indexing options
 minimap_index_opts: ""
 

workflow/envs/env.yml

@@ -27,4 +27,5 @@ dependencies:
     - pydeseq2=0.4.4
     - anndata>=0.8.0
     - ensureconda
+    #- qualimap=2.3
 

workflow/profile/config.yaml

@@ -3,6 +3,11 @@ default-resources:
     slurm_partition: "smp"
 
 set-resources:
+    plot_samples:
+        cpus_per_task: 8
+        mem_mb_per_cpu: 1800
+        runtime: "1h"
+
     build_minimap_index:
         cpus_per_task: 4
         mem_mb_per_cpu: 3600
@@ -14,6 +19,11 @@ set-resources:
         runtime: "3h"
         slurm_partition: "parallel"
 
+    mapping_qc:
+        cpus_per_task: 8
+        mem_mb_per_cpu: 1800
+        runtime: "1h"
+
     sam_sort:
         cpus_per_task: 8
         mem_mb_per_cpu: 1800