Adding new rule to get transcript counts

config/cluster.json

@@ -31,5 +31,9 @@
     "nanostat": {
         "threads": "4",
         "mem": "16g"
+    },
+    "nanocount": {
+        "threads": "2",
+        "mem": "16g"
     }
 }

config/config.json

@@ -1,4 +1,5 @@
 {
     "options": {
+        "nanocount_em_iter": "100"
     }
 }

workflow/Snakefile

@@ -120,6 +120,12 @@ rule all:
             join(workpath, "{name}", "bams", "{name}.sorted.genome.metrics"),
             name=samples
         ),
+        # NanoCount to quantify isoform counts
+        # @imported from `rule nanocount` in rules/quant.smk
+        expand(
+            join(workpath, "{name}", "counts", "{name}.nanocount.transcripts.tsv"),
+            name=samples
+        ),
         # MultiQC to aggregate results into a report
         # @imported from `rule multiqc` in rules/qc.smk
         join(workpath, "reports", "multiqc_report.html"),
@@ -130,4 +136,5 @@ include: join("rules", "common.smk")
 include: join("rules", "download.smk")
 include: join("rules", "trim.smk")
 include: join("rules", "qc.smk")
-include: join("rules", "map.smk")
+include: join("rules", "map.smk")
+include: join("rules", "quant.smk")

workflow/rules/quant.smk

@@ -0,0 +1,28 @@
+# Quantification rules
+rule nanocount:
+    """
+    Quality-control step to gather various statistics from a BAM file.
+    This tool is supported by MultiQC.
+    Github: https://github.com/a-slide/NanoCount
+    @Input:
+        Sorted Transcriptomic BAM file (scatter)
+    @Output:
+        NanoCount metrics file
+    """
+    input:
+        bam = join(workpath, "{name}", "bams", "{name}.sorted.transcriptome.bam"),
+    output:
+        counts = join(workpath, "{name}", "counts", "{name}.nanocount.transcripts.tsv"),
+    params:
+        rname   = "nanocount",
+        em_iter = config['options']['nanocount_em_iter'],
+    conda: depending(join(workpath, config['conda']['modr']), use_conda)
+    container: depending(config['images']['modr'], use_singularity)
+    threads: int(allocated("threads", "nanocount", cluster))
+    shell: """
+    NanoCount \\
+        -i {input.bam} \\
+        --extra_tx_info \\
+        -e {params.em_iter} \\
+        -o {output.counts} 
+    """