Md cram

.github/workflows/ci.yml

@@ -134,7 +134,7 @@ jobs:
           - "mutect2"
           - "prepare_recalibration"
           - "recalibrate"
-          - "save_bam_mapped"
+          - "save_mapped"
           - "save_output_as_bam"
           - "skip_markduplicates"
           - "skip_qc"

CHANGELOG.md

@@ -27,6 +27,7 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0
 - [#731](https://github.com/nf-core/sarek/pull/731) - Run the somatic test as default on `-profile test_full`, the germline can be tested with `-profile test_full_somatic`
 - [#733](https://github.com/nf-core/sarek/pull/733) - Add description for params.cf_chrom_len
 - [#734](https://github.com/nf-core/sarek/pull/734) - nf-core modules update -a
+- [#735](https://github.com/nf-core/sarek/pull/735) - GATK Markduplicates now natively supports CRAM output
 - [#736](https://github.com/nf-core/sarek/pull/736) - More extensive CI for default test
 
 ### Fixed

conf/modules.config

@@ -386,12 +386,13 @@ process {
             // However if it's skipped, reads need to be coordinate-sorted
             // Only name sort if Spark for Markduplicates + duplicate marking is not skipped
             ext.args2        = { params.use_gatk_spark && params.use_gatk_spark.contains('markduplicates') && (!params.skip_tools || (params.skip_tools && !params.skip_tools.split(',').contains('markduplicates'))) ? '-n' : '' }
-            ext.prefix       = { params.split_fastq > 1 ? "${meta.id}".concat('.').concat(reads.get(0).name.tokenize('.')[0]) : "" }
+            ext.prefix       = { params.split_fastq > 1 ? "${meta.id}".concat('.').concat(reads.get(0).name.tokenize('.')[0]) : "${meta.id}.sorted" }
             publishDir       = [
                 mode: params.publish_dir_mode,
                 path: { "${params.outdir}/preprocessing/" },
                 pattern: "*bam",
-                saveAs: { (params.save_bam_mapped || (params.skip_tools && params.skip_tools.split(',').contains('markduplicates'))) && (meta.size * meta.numLanes == 1) ? "mapped/${meta.id}/${it}" : null }
+                // Only save if (save_mapped OR (no_markduplicates AND save_output_as_bam)) AND only a single BAM file per sample
+                saveAs: { (params.save_output_as_bam && (params.save_mapped || params.skip_tools && params.skip_tools.split(',').contains('markduplicates'))) && (meta.size * meta.numLanes == 1) ? "mapped/${meta.id}/${it}" : null }
             ]
         }
 
@@ -404,32 +405,47 @@ process {
     withName: 'MERGE_BAM|INDEX_MERGE_BAM' {
         publishDir       = [
             mode: params.publish_dir_mode,
-            path: { "${params.outdir}/preprocessing/mapped/${meta.id}/" },
-            pattern: "*{bam,bai}"
+            path: { "${params.outdir}/preprocessing/" },
+            pattern: "*{bam,bai}",
+            // Only save if (save_output_as_bam AND (no_markduplicates OR save_mapped ))
+            saveAs: { (params.save_output_as_bam && (params.save_mapped || params.skip_tools && params.skip_tools.split(',').contains('markduplicates'))) ? "mapped/${meta.id}/${it}" : null }
         ]
     }
 
-// MARKDUPLICATES
-
-    withName: 'SAMTOOLS_CRAMTOBAM.*'{
-        ext.args         = "-b"
+    withName: 'MERGE_BAM'{
+        ext.prefix = { "${meta.id}.sorted" }
     }
 
+// MARKDUPLICATES
 
-    withName: 'SAMTOOLS_STATS_CRAM'{
-        ext.prefix       = { (params.skip_tools && params.skip_tools.split(',').contains('markduplicates')) ? "${meta.id}.sorted.cram" : "${meta.id}.md.cram" }
+    withName: 'CRAMTOBAM.*'{
+        ext.args         = "-b"
     }
 
-    withName: 'SAMTOOLS_BAMTOCRAM.*' {
-        // BAM provided for step Markduplicates either run through MD or Convert -> then saved as md.cram
+    withName: 'BAMTOCRAM_MARKDUPLICATES' {
+        // BAM provided for step Markduplicates either run through MD or Convert -> then saved as sorted.cram (convert) or md.cram (md directly)
         // BAM files provided for step prepare_recal are converted and run through BQSR -> then saved as md.cram
         // BAM files provided for step recal are converted and run through BQSR II -> then saved as md.cram
         ext.args         = "-C"
-        ext.prefix       = { (params.skip_tools && params.skip_tools.split(',').contains('markduplicates')) ? "${meta.id}.sorted" : "${meta.id}.md" }
+        ext.prefix       = { (params.skip_tools && params.skip_tools.split(',').contains('markduplicates')) ? "${meta.id}.md" : "${meta.id}.sorted" }
         publishDir       = [
             enabled: !params.save_output_as_bam,
             mode: params.publish_dir_mode,
-            path: { (params.skip_tools && params.skip_tools.split(',').contains('markduplicates')) ? "${params.outdir}/preprocessing/mapped/${meta.id}/" : "${params.outdir}/preprocessing/markduplicates/${meta.id}/" },
+            path: { "${params.outdir}/preprocessing" },
+            pattern: "*{cram,crai}",
+            saveAs: { (params.skip_tools && params.skip_tools.split(',').contains('markduplicates')) ? "markduplicates/${meta.id}/${it}" : "mapped/${meta.id}/${it}" },
+        ]
+    }
+
+    withName: 'BAMTOCRAM_MAPPING'{
+        // Run only when mapping should be saved as CRAM or when no MD is done
+        ext.when         = (params.save_mapped && !params.save_output_as_bam)  || (params.skip_tools && params.skip_tools.split(',').contains('markduplicates'))
+        ext.prefix       = { "${meta.id}.sorted" }
+        publishDir       = [
+            // Never publish if BAM only should be published
+            enabled: !params.save_output_as_bam,
+            mode: params.publish_dir_mode,
+            path: { "${params.outdir}/preprocessing/mapped/${meta.id}/" },
             pattern: "*{cram,crai}"
         ]
     }
@@ -446,14 +462,14 @@ process {
 
     withName: 'GATK4_MARKDUPLICATES' {
         ext.args         = '-REMOVE_DUPLICATES false -VALIDATION_STRINGENCY LENIENT --CREATE_INDEX true'
-        ext.prefix       = { "${meta.id}.md" }
+        ext.prefix       = { "${meta.id}.md.cram" }
         ext.when         = { !(params.skip_tools && params.skip_tools.split(',').contains('markduplicates')) }
         publishDir       = [
             [
-                enabled: params.save_output_as_bam,
+                enabled: !params.save_output_as_bam,
                 mode: params.publish_dir_mode,
                 path: { "${params.outdir}/preprocessing/markduplicates/${meta.id}/" },
-                pattern: "*{bam,bai}"
+                pattern: "*{cram,crai}"
             ],
             [
                 mode: params.publish_dir_mode,
@@ -484,17 +500,16 @@ process {
         ]
     }
 
-    if (params.use_gatk_spark && params.use_gatk_spark.split(',').contains('markduplicates')) {
-        withName: 'NFCORE_SAREK:SAREK:MARKDUPLICATES_SPARK:SAMTOOLS_CRAMTOBAM'{
-            ext.prefix       = { "${meta.id}.md" }
-            ext.when         = { params.save_output_as_bam}
-            publishDir       = [
-                enabled: params.save_output_as_bam,
-                mode: params.publish_dir_mode,
-                path: { "${params.outdir}/preprocessing/markduplicates/${meta.id}/" },
-                pattern: "*{md.bam,md.bam.bai}"
-            ]
-        }
+
+    withName: 'CRAMTOBAM_MARKDUPLICATES'{
+        ext.prefix       = { "${meta.id}.md" }
+        ext.when         = { params.save_output_as_bam}
+        publishDir       = [
+            enabled: params.save_output_as_bam,
+            mode: params.publish_dir_mode,
+            path: { "${params.outdir}/preprocessing/markduplicates/${meta.id}/" },
+            pattern: "*{md.bam,md.bam.bai}"
+        ]
     }
 
 // PREPARE_RECALIBRATION
@@ -555,7 +570,7 @@ process {
         }
     }
 
-    withName: 'SAMTOOLS_CRAMTOBAM_RECAL'{
+    withName: 'CRAMTOBAM_RECAL'{
         ext.prefix       = { "${meta.id}.recal" }
         ext.when         = { params.save_output_as_bam}
         publishDir       = [
@@ -566,7 +581,7 @@ process {
         ]
     }
 
-    withName: 'SAMTOOLS_BAMTOCRAM_VARIANTCALLING' {
+    withName: 'BAMTOCRAM_VARIANTCALLING' {
         ext.prefix       = { "${meta.id}.recal" }
         publishDir       = [
             mode: params.publish_dir_mode,
@@ -590,21 +605,22 @@ process {
         ]
     }
 
-    if (params.step == 'mapping' || params.step == 'markduplicates') {
-        withName: 'NFCORE_SAREK:SAREK:CRAM_QC:SAMTOOLS_STATS' {
-            ext.prefix       = { "${meta.id}.cram" }
-        }
-    }
-
     withName: 'SAMTOOLS_STATS' {
         ext.when         = { !(params.skip_tools && params.skip_tools.split(',').contains('samtools')) }
+        ext.prefix       = { "${meta.id}.md.cram" }
         publishDir       = [
             mode: params.publish_dir_mode,
             path: { "${params.outdir}/reports/samtools/${meta.id}" },
             saveAs: { filename -> filename.equals('versions.yml') ? null : filename }
         ]
     }
 
+    if(params.skip_tools && params.skip_tools.split(',').contains('markduplicates')){
+        withName: 'NFCORE_SAREK:SAREK:CRAM_QC_NO_MARKDUPLICATES:SAMTOOLS_STATS'{
+            ext.prefix       = { "${meta.id}.sorted.cram" }
+        }
+    }
+
     withName: 'MOSDEPTH' {
         ext.args = { !params.wes ? "-n --fast-mode --by 500" : ""}
         ext.prefix       = { params.skip_tools && params.skip_tools.split(',').contains('markduplicates') ? "${meta.id}.sorted" : "${meta.id}.md" }
@@ -617,11 +633,11 @@ process {
     }
 
     if ((params.step == 'mapping' || params.step == 'markduplicates'|| params.step == 'prepare_recalibration'|| params.step == 'recalibrate') && (!(params.skip_tools && params.skip_tools.split(',').contains('baserecalibrator')))) {
-        withName: 'NFCORE_SAREK:SAREK:CRAM_QC:MOSDEPTH' {
+        withName: 'NFCORE_SAREK:SAREK:CRAM_QC_RECALIBRATION:MOSDEPTH' {
             ext.prefix       = { "${meta.id}.recal" }
         }
 
-        withName: 'NFCORE_SAREK:SAREK:CRAM_QC:SAMTOOLS_STATS' {
+        withName: 'NFCORE_SAREK:SAREK:CRAM_QC_RECALIBRATION:SAMTOOLS_STATS' {
             ext.prefix       = { "${meta.id}.recal.cram" }
             ext.when         = { !(params.skip_tools && params.skip_tools.split(',').contains('samtools')) }
             publishDir       = [

conf/test.config

@@ -96,8 +96,8 @@ profiles {
         params.step                   = 'recalibrate'
         params.tools                  = null
     }
-    save_bam_mapped {
-        params.save_bam_mapped        = true
+    save_mapped {
+        params.save_mapped        = true
         params.tools                  = null
     }
     skip_bqsr {

docs/output.md

@@ -150,26 +150,26 @@ These files are intermediate and by default not kept in the final files delivere
 
 [BWA](https://github.com/lh3/bwa) is a software package for mapping low-divergent sequences against a large reference genome. The aligned reads are then coordinate-sorted (or name-sorted if [`GATK MarkDuplicatesSpark`](https://gatk.broadinstitute.org/hc/en-us/articles/5358833264411-MarkDuplicatesSpark) is used for duplicate marking) with [samtools](https://www.htslib.org/doc/samtools.html).
 
-These files are intermediate and by default not kept in the final files delivered to users. Set `--save_bam_mapped` to enable publishing.
+These files are intermediate and by default not kept in the final files delivered to users. Set `--save_mapped` to enable publishing.
 
 #### BWA-mem2
 
 [BWA-mem2](https://github.com/bwa-mem2/bwa-mem2) is a software package for mapping low-divergent sequences against a large reference genome.The aligned reads are then coordinate-sorted (or name-sorted if [`GATK MarkDuplicatesSpark`](https://gatk.broadinstitute.org/hc/en-us/articles/5358833264411-MarkDuplicatesSpark) is used for duplicate marking) with [samtools](https://www.htslib.org/doc/samtools.html).
 
-These files are intermediate and by default not kept in the final files delivered to users. Set `--save_bam_mapped` to enable publishing.
+These files are intermediate and by default not kept in the final files delivered to users. Set `--save_mapped` to enable publishing.
 
 #### DragMap
 
 [DragMap](https://github.com/Illumina/dragmap) is an open-source software implementation of the DRAGEN mapper, which the Illumina team created so that we would have an open-source way to produce the same results as their proprietary DRAGEN hardware. The aligned reads are then coordinate-sorted (or name-sorted if [`GATK MarkDuplicatesSpark`](https://gatk.broadinstitute.org/hc/en-us/articles/5358833264411-MarkDuplicatesSpark) is used for duplicate marking) with [samtools](https://www.htslib.org/doc/samtools.html).
 
-These files are intermediate and by default not kept in the final files delivered to users. Set `--save_bam_mapped` to enable publishing.
+These files are intermediate and by default not kept in the final files delivered to users. Set `--save_mapped` to enable publishing.
 
 <details markdown="1">
 <summary>Output files for all mappers and samples</summary>
 
 **Output directory: `{outdir}/preprocessing/mapped/<sample>/`**
 
-- if `--save_bam_mapped`: `<sample>.bam` and `<sample>.bam.bai`
+- if `--save_mapped`: `<sample>.bam` and `<sample>.bam.bai`
   - BAM file and index
 
 </details>
@@ -253,7 +253,7 @@ See the [`--input`](usage.md#--input) section in the usage documentation for fur
 **Output directory: `{outdir}/preprocessing/csv`**
 
 - `mapped.csv`
-  - if `--save_bam_mapped`
+  - if `--save_mapped`
   - CSV containing an entry for each sample with the columns `patient,sample,sex,status,bam,bai`
 - `markduplicates_no_table.csv`
   - CSV containing an entry for each sample with the columns `patient,sample,sex,status,cram,crai`
@@ -824,7 +824,7 @@ The plot will show:
 
 **Output directory: `{outdir}/reports/markduplicates/<sample>`**
 
-- `<sample>.md.metrics`
+- `<sample>.md.cram.metrics`
   - file used by [MultiQC](https://multiqc.info/)
   </details>
 

docs/usage.md

@@ -148,7 +148,7 @@ patient,sample,cram,crai
 patient1,test_sample,test_mapped.cram,test_mapped.cram.crai
 ```
 
-The Sarek-generated CSV file is stored under `results/csv/mapped.csv` if in a previous run `--save_bam_mapped` was set and will automatically be used as an input when specifying the parameter `--step markduplicates`. Otherwise this file will need to be manually generated.
+The Sarek-generated CSV file is stored under `results/csv/mapped.csv` if in a previous run `--save_mapped` was set and will automatically be used as an input when specifying the parameter `--step markduplicates`. Otherwise this file will need to be manually generated.
 
 #### Full samplesheet
 
@@ -526,9 +526,10 @@ Expected run output:
 [-        ] process > NFCORE_SAREK:SAREK:GATK4_MAPPING:BWAMEM2_MEM                                          -
 [-        ] process > NFCORE_SAREK:SAREK:GATK4_MAPPING:DRAGMAP_ALIGN                                        -
 [46/35a640] process > NFCORE_SAREK:SAREK:MARKDUPLICATES:GATK4_MARKDUPLICATES (test)                         [100%] 1 of 1 ✔
-[e0/525bb3] process > NFCORE_SAREK:SAREK:MARKDUPLICATES:BAM_TO_CRAM:SAMTOOLS_BAMTOCRAM (test)               [100%] 1 of 1 ✔
-[46/9fe93a] process > NFCORE_SAREK:SAREK:MARKDUPLICATES:BAM_TO_CRAM:SAMTOOLS_STATS_CRAM (test)              [100%] 1 of 1 ✔
-[77/2c8b1b] process > NFCORE_SAREK:SAREK:MARKDUPLICATES:BAM_TO_CRAM:MOSDEPTH (test)                         [100%] 1 of 1 ✔
+[9a/76cef7] process > NFCORE_SAREK:SAREK:MARKDUPLICATES:INDEX_MARKDUPLICATES (test)                         [100%] 1 of 1 ✔
+[46/9fe93a] process > NFCORE_SAREK:SAREK:MARKDUPLICATES:CRAM_QC:SAMTOOLS_STATS (test)                       [100%] 1 of 1 ✔
+[77/2c8b1b] process > NFCORE_SAREK:SAREK:MARKDUPLICATES:CRAM_QC:MOSDEPTH (test)                             [100%] 1 of 1 ✔
+[-        ] process > NFCORE_SAREK:SAREK:SAMTOOLS_CRAMTOBAM_MARKDUPLICATES                                  -
 [f7/499800] process > NFCORE_SAREK:SAREK:PREPARE_RECALIBRATION:BASERECALIBRATOR (test)                      [100%] 1 of 1 ✔
 [-        ] process > NFCORE_SAREK:SAREK:PREPARE_RECALIBRATION:GATHERBQSRREPORTS                            -
 [9d/3d1fff] process > NFCORE_SAREK:SAREK:RECALIBRATE:APPLYBQSR (test)                                       [100%] 1 of 1 ✔

modules.json

@@ -195,7 +195,7 @@
                 },
                 "gatk4/markduplicates": {
                     "branch": "master",
-                    "git_sha": "df2620cfc7e4c21b14ed03c1c928f09fbabf83c4"
+                    "git_sha": "e726b1730dff525bde4a6839e544dabfea4cd7fd"
                 },
                 "gatk4/markduplicatesspark": {
                     "branch": "master",

nextflow.config

@@ -43,7 +43,7 @@ params {
     // Preprocessing
     aligner            = 'bwa-mem'  // Default is bwa-mem, bwa-mem2 and dragmap can be used too
     use_gatk_spark     = null       // GATK Spark implementation of their tools in local mode not used by default
-    save_bam_mapped    = false      // Mapped BAMs not saved
+    save_mapped    = false      // Mapped BAMs not saved
     save_output_as_bam = false      //Output files from preprocessing are saved as bam and not as cram files
     seq_center         = null       // No sequencing center to be written in read group CN field by aligner
     seq_platform       = 'ILLUMINA' // Default platform written in read group PL field by aligner