Add ampliconstats process for all the files (Samtools development branch)

README.md

@@ -19,6 +19,7 @@ This Nextflow pipeline automates the ARTIC network [nCoV-2019 novel coronavirus
 
 You can also use cram file input by passing the --cram flag.
 You can also specify cram file output by passing the --outCram flag.
+You can also run ampliconstats by passing the --ampliconstats flag.
 
 For production use at large scale, where you will run the workflow many times, you can avoid cloning the scheme repository, creating an ivar bed file and indexing the reference every time by supplying both --ivarBed /path/to/ivar-compatible.bed and --alignerRefPrefix /path/to/bwa-indexed/ref.fa.
 

conf/illumina.config

@@ -48,6 +48,8 @@ params {
     // iVar minimum mapQ to call variant (ivar variants: -q)
     ivarMinVariantQuality = 20
 
+    // Run samtools ampliconstats 
+    ampliconstats = false
 }
 
 def makeFastqSearchPath ( illuminaSuffixes, fastq_exts ) {

modules/illumina.nf

@@ -102,6 +102,28 @@ process trimPrimerSequences {
         """
 }
 
+process makeAmpliconstats {
+
+    publishDir "${params.outdir}/${task.process.replaceAll(":","_")}", pattern: "*.stats", mode: 'copy'
+    publishDir "${params.outdir}/${task.process.replaceAll(":","_")}", pattern: "*.png", mode: 'copy'
+    publishDir "${params.outdir}/${task.process.replaceAll(":","_")}", pattern: "*.gp", mode: 'copy'
+
+    input:
+    file (ptrim_bams)
+    path(bedfile)
+
+    output:
+    path "nCoV-2019.amp.stats", emit: ampstats
+    path "*.png", emit: amppng
+    path "*.gp", emit: ampgp
+
+    script:
+        """
+        samtools ampliconstats -@8 -d 1,20,100 ${bedfile} *.mapped.primertrimmed.sorted.bam > nCoV-2019.amp.stats
+        plot-ampliconstats -size 1200,900 nCoV-2019-ampliconstats nCoV-2019.amp.stats
+        """
+}
+
 process callVariants {
 
     tag { sampleName }

workflows/illuminaNcov.nf

@@ -9,6 +9,7 @@ include {readTrimming} from '../modules/illumina.nf'
 include {indexReference} from '../modules/illumina.nf'
 include {readMapping} from '../modules/illumina.nf' 
 include {trimPrimerSequences} from '../modules/illumina.nf' 
+include {makeAmpliconstats} from '../modules/illumina.nf'
 include {callVariants} from '../modules/illumina.nf'
 include {makeConsensus} from '../modules/illumina.nf' 
 include {cramToFastq} from '../modules/illumina.nf'
@@ -121,6 +122,15 @@ workflow sequenceAnalysis {
 
       }
 
+      // Generate ampliconstats only if --ampliconstats parameter is passed
+      if(params.ampliconstats) {
+        // Create a channel to collect only the mapped.primertrimmed.sorted.bam output files from trimPrimerSequences
+
+        trimPrimerSequences.out.ptrim.map{sample, bam -> bam}.collect().set{ ch_ptrim_bam }
+
+        makeAmpliconstats(ch_ptrim_bam, ch_bedFile)
+      }
+
     emit:
       qc_pass = collateSamples.out
 }