Alignment to fastq fix and add tests

.github/workflows/ci.yml

@@ -32,6 +32,7 @@ jobs:
             NXF_EDGE: "1"
         test:
           - "aligner"
+          - "alignment_to_fastq"
           - "annotation"
           - "cnvkit"
           - "controlfreec"

CHANGELOG.md

@@ -103,6 +103,7 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0
 - [#599](https://github.com/nf-core/sarek/pull/599) - Add checks for correct data type for `params.step`
 - [#599](https://github.com/nf-core/sarek/pull/599) - Add checks for no empty `--tools` with `--step variant_calling` or `--step annotation`
 - [#600](https://github.com/nf-core/sarek/pull/600) - Remove `nf-core lint` warnings
+- [#602](https://github.com/nf-core/sarek/pull/602/) - Fixed bug in `alignment_to_fastq` and added tests
 
 ### Deprecated
 

conf/test.config

@@ -26,7 +26,7 @@ params {
     max_time   = '8.h'
 
     // Input data
-    input = "${baseDir}/tests/csv/3.0/fastq_single.csv"
+    input = "${projectDir}/tests/csv/3.0/fastq_single.csv"
 
     // Small reference genome
     genome            = null
@@ -48,37 +48,37 @@ params {
 
 profiles {
     annotation {
-        params.input               = "${baseDir}/tests/csv/3.0/vcf_single.csv"
+        params.input               = "${projectDir}/tests/csv/3.0/vcf_single.csv"
         params.step                = 'annotate'
     }
     no_intervals {
         params.no_intervals        = true
     }
     pair {
-        params.input               = "${baseDir}/tests/csv/3.0/fastq_pair.csv"
+        params.input               = "${projectDir}/tests/csv/3.0/fastq_pair.csv"
     }
     markduplicates_bam {
-        params.input               = "${baseDir}/tests/csv/3.0/mapped_single_bam.csv"
+        params.input               = "${projectDir}/tests/csv/3.0/mapped_single_bam.csv"
         params.step                = 'markduplicates'
     }
     markduplicates_cram {
-        params.input               = "${baseDir}/tests/csv/3.0/mapped_single_cram.csv"
+        params.input               = "${projectDir}/tests/csv/3.0/mapped_single_cram.csv"
         params.step                = 'markduplicates'
     }
     prepare_recalibration_bam {
-        params.input               = "${baseDir}/tests/csv/3.0/mapped_single_bam.csv"
+        params.input               = "${projectDir}/tests/csv/3.0/mapped_single_bam.csv"
         params.step                = 'prepare_recalibration'
     }
     prepare_recalibration_cram {
-        params.input               = "${baseDir}/tests/csv/3.0/mapped_single_cram.csv"
+        params.input               = "${projectDir}/tests/csv/3.0/mapped_single_cram.csv"
         params.step                = 'prepare_recalibration'
     }
     recalibrate_bam {
-        params.input               = "${baseDir}/tests/csv/3.0/prepare_recalibration_single_bam.csv"
+        params.input               = "${projectDir}/tests/csv/3.0/prepare_recalibration_single_bam.csv"
         params.step                = 'recalibrate'
     }
     recalibrate_cram {
-        params.input               = "${baseDir}/tests/csv/3.0/prepare_recalibration_single_cram.csv"
+        params.input               = "${projectDir}/tests/csv/3.0/prepare_recalibration_single_cram.csv"
         params.step                = 'recalibrate'
     }
     save_bam_mapped {
@@ -100,7 +100,7 @@ profiles {
         params.nucleotides_per_second  = 20
     }
     tools {
-        params.input               = "${baseDir}/tests/csv/3.0/recalibrated.csv"
+        params.input               = "${projectDir}/tests/csv/3.0/recalibrated.csv"
         params.dbsnp               = params.test_data['homo_sapiens']['genome']['dbsnp_138_hg38_21_vcf_gz']
         params.fasta               = params.test_data['homo_sapiens']['genome']['genome_21_fasta']
         params.germline_resource   = params.test_data['homo_sapiens']['genome']['gnomad_r2_1_1_21_vcf_gz']
@@ -113,7 +113,7 @@ profiles {
         params.nucleotides_per_second = 20
     }
     tools_germline {
-        params.input               = "${baseDir}/tests/csv/3.0/recalibrated_germline.csv"
+        params.input               = "${projectDir}/tests/csv/3.0/recalibrated_germline.csv"
         params.dbsnp               = params.test_data['homo_sapiens']['genome']['dbsnp_138_hg38_21_vcf_gz']
         params.known_indels        = params.test_data['homo_sapiens']['genome']['mills_and_1000g_indels_21_vcf_gz']
         params.fasta               = params.test_data['homo_sapiens']['genome']['genome_21_fasta']
@@ -124,7 +124,7 @@ profiles {
         params.nucleotides_per_second = 20
     }
     tools_tumoronly {
-        params.input               = "${baseDir}/tests/csv/3.0/recalibrated_tumoronly.csv"
+        params.input               = "${projectDir}/tests/csv/3.0/recalibrated_tumoronly.csv"
         params.dbsnp               = params.test_data['homo_sapiens']['genome']['dbsnp_138_hg38_21_vcf_gz']
         params.fasta               = params.test_data['homo_sapiens']['genome']['genome_21_fasta']
         params.germline_resource   = params.test_data['homo_sapiens']['genome']['gnomad_r2_1_1_21_vcf_gz']
@@ -137,7 +137,7 @@ profiles {
         params.nucleotides_per_second = 20
     }
     tools_somatic {
-        params.input               = "${baseDir}/tests/csv/3.0/recalibrated_somatic.csv"
+        params.input               = "${projectDir}/tests/csv/3.0/recalibrated_somatic.csv"
         params.chr_dir             = params.test_data['homo_sapiens']['genome']['genome_21_chromosomes_dir']
         params.dbsnp               = params.test_data['homo_sapiens']['genome']['dbsnp_138_hg38_21_vcf_gz']
         params.fasta               = params.test_data['homo_sapiens']['genome']['genome_21_fasta']
@@ -158,19 +158,22 @@ profiles {
         params.trim_fastq          = true
     }
     umi {
-        params.input               = "${baseDir}/tests/csv/3.0/fastq_umi.csv"
+        params.input               = "${projectDir}/tests/csv/3.0/fastq_umi.csv"
         params.umi_read_structure  = '7M1S+T'
     }
     use_gatk_spark {
         params.use_gatk_spark      = 'baserecalibrator,markduplicates'
     }
     variantcalling_channels {
-        params.input               = "${baseDir}/tests/csv/3.0/recalibrated.csv"
+        params.input               = "${projectDir}/tests/csv/3.0/recalibrated.csv"
         params.fasta               = params.test_data['homo_sapiens']['genome']['genome_21_fasta']
         params.intervals           = params.test_data['homo_sapiens']['genome']['genome_21_multi_interval_bed']
         params.wes                 = true
         params.step                = 'variant_calling'
 
         params.nucleotides_per_second = 20
     }
+    alignment_to_fastq {
+        params.input               = "${projectDir}/tests/csv/3.0/bam_for_remapping.csv"
+    }
 }

conf/test_full.config

@@ -15,7 +15,7 @@ params {
     config_profile_description = 'Full test dataset to check pipeline function'
 
     // Input data for full size test
-    input = "${baseDir}/tests/csv/3.0/test_full_data.csv"
+    input = "${projectDir}/tests/csv/3.0/test_full_data.csv"
 
     wes   = true
 

tests/csv/3.0/bam_for_remapping.csv

@@ -0,0 +1,2 @@
+patient,gender,status,sample,lane,bam,bai
+test,XX,0,test,1,https://raw.githubusercontent.com/nf-core/test-datasets/modules/data/genomics/homo_sapiens/illumina/bam/test.paired_end.sorted.bam,https://raw.githubusercontent.com/nf-core/test-datasets/modules/data/genomics/homo_sapiens/illumina/bam/test.paired_end.sorted.bam.bai

tests/test_bam_remap.yml

@@ -0,0 +1,17 @@
+- name: Run alignment to fastq and then remap on bam files
+  command: nextflow run main.nf -profile test,alignment_to_fastq,docker -c ./tests/nextflow.config
+  tags:
+    - alignment_to_fastq
+
+  files:
+    - path: results/cat/test-1_1.merged.fastq.gz
+    - path: results/cat/test-1_2.merged.fastq.gz
+    - path: results/csv/markduplicates.csv
+    - path: results/csv/markduplicates_no_table.csv
+    - path: results/csv/recalibrated.csv
+    - path: results/multiqc/multiqc_report.html
+    - path: results/pipeline_info
+    - path: results/preprocessing/test
+    - path: results/reports
+    - path: results/samtools
+    - path: results/collate

workflows/sarek.nf

@@ -1020,18 +1020,22 @@ def extract_csv(csv_file) {
 
         // start from BAM
         } else if (row.lane && row.bam) {
+            if (!row.bai) {
+                log.error "BAM index (bai) should be provided."
+            }
             meta.id         = "${row.sample}-${row.lane}".toString()
             def bam         = file(row.bam,   checkIfExists: true)
+            def bai         = file(row.bai,   checkIfExists: true)
             def CN          = params.seq_center ? "CN:${params.seq_center}\\t" : ''
-            def read_group  = "\"@RG\\tID:${row_sample}_${row.lane}\\t${CN}PU:${row.lane}\\tSM:${row.sample}\\tLB:${row.sample}\\tPL:${params.seq_platform}\""
+            def read_group  = "\"@RG\\tID:${row.sample}_${row.lane}\\t${CN}PU:${row.lane}\\tSM:${row.sample}\\tLB:${row.sample}\\tPL:${params.seq_platform}\""
 
             meta.numLanes   = numLanes.toInteger()
             meta.read_group = read_group.toString()
             meta.data_type  = "bam"
 
             meta.size       = 1 // default number of splitted fastq
 
-            if (params.step == 'mapping') return [meta, bam]
+            if (params.step == 'mapping') return [meta, bam, bai]
             else {
                 log.error "Samplesheet contains ubam files but step is `$params.step`. Please check your samplesheet or adjust the step parameter.\nhttps://nf-co.re/sarek/usage#input-samplesheet-configurations"
                 System.exit(1)