nf-core · njspix · Dec 8, 2022 · Dec 12, 2022 · Dec 12, 2022 · Dec 12, 2022
@@ -49,6 +49,24 @@
 - 🐛 fix `ignore_3prime_r2` param #299
 - 🐛 removed unused directory #297
 
+## Working list...
+
+### Pipeline Updates
+
+- Add `Biscuit` aligner as a separate, 3rd workflow
+  - Add `biscuit_index` parameter
+- Add a `--save_merged` option to save merged (`cat`'d) fastq files for the same sample
+
+### Refactoring
+
+- 🧹 Refactoring of parameters, especially those that generalize across aligners:
+  - `--non_directional` is now under 'Alignment options' rather than 'Bismark options', as both the Bismark and Biscuit aligners support directional and non-directional alignment.
+  - New 'Methylation calling options' section
+    - `--merge_cg` replaces the former `--cytosine_report` option. Biscuit and bwa-meth produce stranded methylation calls by default, while bismark does not. This flag abstracts away all the necessary flags for each tool, and defaults to producing merged (unstranded) methylationc calls.
+    - `--meth_cutoff`, `--comprehensive`, and the `--ignore` parameters apply to all 3 aligners
+    - `--no_overlap` applies to biscuit and bismark but NOT bwa-meth
+  - `--nomeseq` has been moved to 'Special Library Types' as both `bismark` and `biscuit` support NOMe-seq.
+
 ## [v2.3.0](https://github.com/nf-core/methylseq/releases/tag/2.3.0) - 2022-12-16
 
 ### Pipeline Updates

@@ -22,23 +22,23 @@ On release, automated continuous integration tests run the pipeline on a full-si
 
 ## Pipeline Summary
 
-The pipeline allows you to choose between running either [Bismark](https://github.com/FelixKrueger/Bismark) or [bwa-meth](https://github.com/brentp/bwa-meth) / [MethylDackel](https://github.com/dpryan79/methyldackel).
-Choose between workflows by using `--aligner bismark` (default, uses bowtie2 for alignment), `--aligner bismark_hisat` or `--aligner bwameth`.
-
-| Step                                         | Bismark workflow | bwa-meth workflow     |
-| -------------------------------------------- | ---------------- | --------------------- |
-| Generate Reference Genome Index _(optional)_ | Bismark          | bwa-meth              |
-| Merge re-sequenced FastQ files               | cat              | cat                   |
-| Raw data QC                                  | FastQC           | FastQC                |
-| Adapter sequence trimming                    | Trim Galore!     | Trim Galore!          |
-| Align Reads                                  | Bismark          | bwa-meth              |
-| Deduplicate Alignments                       | Bismark          | Picard MarkDuplicates |
-| Extract methylation calls                    | Bismark          | MethylDackel          |
-| Sample report                                | Bismark          | -                     |
-| Summary Report                               | Bismark          | -                     |
-| Alignment QC                                 | Qualimap         | Qualimap              |
-| Sample complexity                            | Preseq           | Preseq                |
-| Project Report                               | MultiQC          | MultiQC               |
+The pipeline allows you to choose between running either [Bismark](https://github.com/FelixKrueger/Bismark), [bwa-meth](https://github.com/brentp/bwa-meth) / [MethylDackel](https://github.com/dpryan79/methyldackel), or [biscuit](https://huishenlab.github.io/biscuit/).
+Choose between workflows by using `--aligner bismark` (default, uses bowtie2 for alignment), `--aligner bismark_hisat`, `--aligner bwameth`, or `--aligner biscuit`.
+
+| Step                                         | Bismark workflow | bwa-meth workflow     | Biscuit workflow |
+| -------------------------------------------- | ---------------- | --------------------- | ---------------- |
+| Generate Reference Genome Index _(optional)_ | Bismark          | bwa-meth              | Biscuit `index`  |
+| Merge re-sequenced FastQ files               | cat              | cat                   | cat              |
+| Raw data QC                                  | FastQC           | FastQC                | FastQC           |
+| Adapter sequence trimming                    | Trim Galore!     | Trim Galore!          | Cutadapt         |
+| Align Reads                                  | Bismark          | bwa-meth              | biscuit `align`  |
+| Deduplicate Alignments                       | Bismark          | Picard MarkDuplicates | samblaster       |
+| Extract methylation calls                    | Bismark          | MethylDackel          | Biscuit `pileup` |
+| Sample report                                | Bismark          | -                     | Biscuit `qc`     |
+| Summary Report                               | Bismark          | -                     | Biscuit `qc`     |
+| Alignment QC                                 | Qualimap         | Qualimap              | Qualimap         |
+| Sample complexity                            | Preseq           | Preseq                | Preseq           |
+| Project Report                               | MultiQC          | MultiQC               | MultiQC          |
 
 ## Usage
 

@@ -34,6 +34,15 @@ process {
         ]
     }
 
+    withName: CAT_FASTQ {
+        publishDir = [
+            path: { "${params.outdir}/merged_fastq" },
+            mode: params.publish_dir_mode,
+            pattern: '*.gz',
+            enabled: params.save_merged
+        ]
+    }
+
     withName: FASTQC {
         ext.args = '--quiet'
         publishDir = [
@@ -209,11 +218,11 @@ process {
             params.comprehensive   ? ' --comprehensive --merge_non_CpG' : '',
             params.meth_cutoff     ? " --cutoff ${params.meth_cutoff}" : '',
             params.nomeseq         ? '--CX' : '',
-            params.ignore_r1 > 0   ? "--ignore ${params.ignore_r1}" : '',
-            params.ignore_3prime_r1 > 0   ? "--ignore_3prime ${params.ignore_3prime_r1}" : '',
-            meta.single_end ? '' : (params.no_overlap           ? ' --no_overlap'                         : '--include_overlap'),
-            meta.single_end ? '' : (params.ignore_r2        > 0 ? "--ignore_r2 ${params.ignore_r2}"       : ""),
-            meta.single_end ? '' : (params.ignore_3prime_r2 > 0 ? "--ignore_3prime_r2 ${params.ignore_3prime_r2}": "")
+            meta.single_end ? '' : (params.no_overlap           ? ' --no_overlap'          : '--include_overlap'),
+            meta.single_end ? '' : (params.ignore           > 0 ? "--ignore ${params.ignore}"               : ""),
+            meta.single_end ? '' : (params.ignore_3prime    > 0 ? "--ignore_3prime ${params.ignore_3prime}" : ""),
+            meta.single_end ? '' : (params.ignore_r2        > 0 ? "--ignore_r2 ${params.ignore_r2}"         : ""),
+            meta.single_end ? '' : (params.ignore_3prime_r2 > 0 ? "--ignore_r2 ${params.ignore_3prime_r2}"  : "")
         ].join(' ').trim() }
         publishDir = [
             [
@@ -304,6 +313,173 @@ process {
         ]
     }
 
+    withName: BISCUIT_INDEX {
+        ext.args = ''
+        publishDir = [
+            path: { "${params.outdir}/${params.aligner}/reference_genome" },
+            saveAs: { it =~ /.*\.yml/ ? null : it },
+            mode: params.publish_dir_mode,
+            enabled: params.save_reference
+        ]
+    }
+
+    withName: BISCUIT_ALIGN {
+        ext.args = { [
+            // Use directional alignment (align read 1 to parent, read 2 to daughter).
+            // Note  that `pbat` libaries are directional, but 'backwards' (read 1 aligns to daughter);
+            // this 'inversion' is handled in the `biscuit` subworkflow by reversing read 1 and read 2 if `pbat` is set.
+            ( params.pbat || params.single_cell || params.zymo || (!params.non_directional) ) ? '-b 1' : '',
+        ].join(' ').trim() }
+        ext.args2 = ''
+        publishDir = [
+            [
+                path: { "${params.outdir}/${params.aligner}/alignments/raw" },
+                mode: params.publish_dir_mode,
+                pattern: "*.bam.*",
+                enabled: (params.save_align_intermeds || params.skip_deduplication || params.rrbs)
+            ]
+        ]
+    }
+
+    withName: BISCUIT_BLASTER {
+        ext.args = { [
+            // Use directional alignment (align read 1 to parent, read 2 to daughter).
+            // Note  that `pbat` libaries are directional, but 'backwards' (read 1 aligns to daughter);
+            // this 'inversion' is handled in the `biscuit` subworkflow by reversing read 1 and read 2 if `pbat` is set.
+            ( params.pbat || params.single_cell || params.zymo || (!params.non_directional) ) ? '-b 1' : '',
+        ].join(' ').trim() }
+        ext.args2 = { [
+            meta.single_end ? "--ignoreUnmated" : ""
+        ].join(' ').trim() }
+        ext.args3 = ''
+        publishDir = [
+            [
+                path: { "${params.outdir}/${params.aligner}/alignments/dup_marked" },
+                mode: params.publish_dir_mode,
+                pattern: "*.bam*",
+                enabled: (params.save_align_intermeds || params.skip_deduplication || params.rrbs)
+            ]
+        ]
+    }
+
+    withName: BISCUIT_BSCONV {
+        ext.prefix = { "${meta.id}_bsconv" }
+        ext.args = { "-f ${params.bs_conv_filter}" }
+        publishDir = [
+                path: { "${params.outdir}/${params.aligner}/alignments/bisulfite_conversion_filtered" },
+                mode: params.publish_dir_mode,
+                pattern: "*.bam",
+                enabled: params.save_align_intermeds
+            ]
+    }
+
+    withName: BISCUIT_PILEUP {
+        errorStrategy = 'retry'
+        maxRetries = 3
+        ext.args = { [
+            params.no_overlap ? '' : '-d',
+            params.nomeseq ? '-N' : '',
+            params.ignore        > 0 ? "-5 ${params.ignore}"        : '',
+            params.ignore_3prime > 0 ? "-3 ${params.ignore_3prime}" : '',
+        ].join(' ').trim() }
+        publishDir = [
+            path: { "${params.outdir}/${params.aligner}/snp_data/" },
+            mode: params.publish_dir_mode,
+            enabled: params.save_align_intermeds,
+            pattern: "*.vcf.gz"
+        ]
+    }
+
+    withName: BISCUIT_VCF2BED {
+        ext.args = { [
+            (params.meth_cutoff) ? "-k ${params.meth_cutoff}" : '-k 1',
+            (params.comprehensive) ? '-t c' : '',
+            (params.nomeseq) ? '-t hcg' : ''
+        ].join(' ').trim() }
+        publishDir = [
+            path: { "${params.outdir}/${params.aligner}/methylation_calls/" },
+            mode: params.publish_dir_mode,
+            enabled: !params.merge_cg,
+            pattern: "*.bed.gz"
+        ]
+    }
+
+    withName: BISCUIT_VCF2BED_NOME {
+        ext.args = { [
+            (params.meth_cutoff) ? "-k ${params.meth_cutoff}" : '-k 1',
+            (params.nomeseq) ? '-t gch' : ''
+        ].join(' ').trim() }
+        publishDir = [
+            path: { "${params.outdir}/${params.aligner}/accessibility_data/" },
+            mode: params.publish_dir_mode,
+            enabled: !params.merge_cg,
+            pattern: "*.bed.gz"
+        ]
+    }
+
+    withName: BISCUIT_MERGECG {
+        ext.prefix = { "${meta.id}_mergecg" }
+        ext.args = ''
+        publishDir = [
+            path: { "${params.outdir}/${params.aligner}/methylation_calls/" },
+            mode: params.publish_dir_mode,
+            enabled: params.merge_cg,
+            pattern: "*.bed.gz"
+        ]
+    }
+
+    withName: BISCUIT_MERGECG_NOME {
+        ext.prefix = { "${meta.id}_mergecg" }
+        ext.args = '-N'
+        publishDir = [
+            path: { "${params.outdir}/${params.aligner}/accessibility_data/" },
+            mode: params.publish_dir_mode,
+            enabled: params.merge_cg,
+            pattern: "*.bed.gz"
+        ]
+    }
+
+    withName: BISCUIT_QC {
+        ext.args = ''
+        publishDir = [
+            [
+                path: { "${params.outdir}/${params.aligner}/qc/mbias/cpg" },
+                mode: params.publish_dir_mode,
+                pattern: "*_CpGRetentionByReadPos.txt"
+            ],
+            [
+                path: { "${params.outdir}/${params.aligner}/qc/mbias/cph" },
+                mode: params.publish_dir_mode,
+                pattern: "*_CpHRetentionByReadPos.txt"
+            ],
+            [
+                path: { "${params.outdir}/${params.aligner}/qc/dup_report/" },
+                mode: params.publish_dir_mode,
+                pattern: "*_dup_report.txt"
+            ],
+            [
+                path: { "${params.outdir}/${params.aligner}/qc/insert_size/" },
+                mode: params.publish_dir_mode,
+                pattern: "*_isize_table.txt"
+            ],
+            [
+                path: { "${params.outdir}/${params.aligner}/qc/mapping_qual/" },
+                mode: params.publish_dir_mode,
+                pattern: "*_mapq_table.txt"
+            ],
+            [
+                path: { "${params.outdir}/${params.aligner}/qc/mapping_strand/" },
+                mode: params.publish_dir_mode,
+                pattern: "*_strand_table.txt"
+            ],
+            [
+                path: { "${params.outdir}/${params.aligner}/qc/conversion/" },
+                mode: params.publish_dir_mode,
+                pattern: "*_totalReadConversionRate.txt"
+            ],
+        ]
+    }
+
     withName: PICARD_MARKDUPLICATES {
         ext.args = "--ASSUME_SORTED true --REMOVE_DUPLICATES false --VALIDATION_STRINGENCY LENIENT --PROGRAM_RECORD_ID 'null' --TMP_DIR tmp"
         ext.prefix = { "${meta.id}.markdup.sorted" }
@@ -381,6 +557,18 @@ process {
         ]
     }
 
+    withName: SAMTOOLS_INDEX_BSCONV {
+        ext.args = ""
+        publishDir = [
+            [
+                path: { "${params.outdir}/${params.aligner}/alignments/bisulfite_conversion_filtered" },
+                mode: params.publish_dir_mode,
+                pattern: "*.bai",
+                enabled: params.save_align_intermeds
+            ]
+        ]
+    }
+
     withName: SAMTOOLS_INDEX_DEDUPLICATED {
         ext.args = ""
         publishDir = [
@@ -407,12 +595,17 @@ process {
     }
 
     withName: METHYLDACKEL_EXTRACT {
-        ext.args = [
-            params.comprehensive ? ' --CHG --CHH' : '',
-            params.ignore_flags ? " --ignoreFlags" : '',
-            params.methyl_kit ? " --methylKit" : '',
-            params.min_depth > 0 ? " --minDepth ${params.min_depth}" : ''
-        ].join(" ").trim()
+        ext.args = {
+            def simple_args = [
+                params.comprehensive ? ' --CHG --CHH' : '',
+                params.ignore_flags ? " --ignoreFlags" : '',
+                params.methyl_kit ? " --methylKit" : '',
+                params.meth_cutoff ? " --minDepth ${params.meth_cutoff}" : '',
+            ].join(" ")
+            def ignore_bases = [params.ignore, params.ignore_3prime, params.ignore_r2, params.ignore_3prime_r2].join(",")
+            def ignore_args = ["--nOT", "--nOB", "--nCTOT", "--nCTOB"].collect { "$it $ignore_bases" }.join(" ")
+            return [simple_args, ignore_args].join(" ").trim()
+        }
         publishDir = [
             [
                 path: { "${params.outdir}/methyldackel" },

@@ -20,6 +20,7 @@ nextflow.enable.dsl = 2
 params.fasta = WorkflowMain.getGenomeAttribute(params, 'fasta')
 params.fasta_index = WorkflowMain.getGenomeAttribute(params, 'fasta_index')
 params.bismark_index = WorkflowMain.getGenomeAttribute(params, 'bismark')
+params.biscuit_index = WorkflowMain.getGenomeAttribute(params, 'biscuit')
 params.bwa_meth_index = WorkflowMain.getGenomeAttribute(params, 'bwa_meth')
 
 /*