nextflow: merge lane files from samples forward and reverse read

Question

I'd like to merge files from different lanes for forward and reverse reads. For e.g., if a sample A has three lanes in forward, and three reverse, ultimately one reverse and one forward is to be generated. I have normal and tumor for each sample. N as normal and T as tumor.
The lanes files are for e.g. as:

MM-265-DNA-N-01-01_L003_R1_001.fastq.gz
MM-265-DNA-N-01-01_L003_R2_001.fastq.gz
MM-265-DNA-N-01-01_L004_R1_001.fastq.gz
MM-265-DNA-N-01-01_L004_R2_001.fastq.gz

MM-265-DNA-T-01-01_L003_R1_001.fastq.gz
MM-265-DNA-T-01-01_L003_R2_001.fastq.gz
MM-265-DNA-T-01-01_L004_R1_001.fastq.gz
MM-265-DNA-T-01-01_L004_R2_001.fastq.gz

MM-465-DNA-N-01-01_L003_R1_001.fastq.gz
MM-465-DNA-N-01-01_L003_R2_001.fastq.gz
MM-465-DNA-N-01-01_L004_R1_001.fastq.gz
MM-465-DNA-N-01-01_L004_R2_001.fastq.gz

MM-465-DNA-T-01-01_L003_R1_001.fastq.gz
MM-465-DNA-T-01-01_L003_R2_001.fastq.gz
MM-465-DNA-T-01-01_L004_R1_001.fastq.gz
MM-465-DNA-T-01-01_L004_R2_001.fastq.gz  

The final output needs to be as:

MM-465_normal_R1.fastq MM-465_normal_R2.fastq
MM-465_tumor_R1.fastq MM-465_tumor_R2.fastq
MM-265_normal_R1.fastq MM-265_normal_R2.fastq
MM-265_tumor_R1.fastq MM-265_tumor_R2.fastq

I've code below broken in three parts:

1. regular expression
2. process that iterates over files
3. workflow that calls the above two code blocks

regex:

def get_sample_id( file ) {
def matcher = (file =~ /(MM-\d+)-DNA-N-.+/)
    if (matcher) {
        return matcher[0][1]
    } else {
        return ''
    }
}

process:

process file_print {

debug true

input:

tuple val(sample_id), path(datasetFile)
val(sample_status) 
val(read_value)

output:
publishDir '/arion/learn_new_concepts/', mode:'copy'

tuple val (sample_id), path("${output_file_name}")

script:

output_file_name = "${sample_id}_${sample_status}_${read_value}.fastq"

for (g in datasetFile){
println(g)

"zcat ${g} >> ${output_file_name} "

}

}

Workflow:

workflow {

normal_reads=Channel.fromFilePairs("/sc/arion/MM-*-N-*_L00*_R{1,2}_*.fastq.gz")
//MM-2658-DNA-N-01-01_L004_R2_001.fastq.gz  //

normal_forward_reads=normal_reads
        .map{ file -> return tuple( get_sample_id(file), file[1][0] )}.groupTuple(sort:true ) 

file_print(normal_forward_reads,"normal","R1").view() /// similar for reverse [1][1], R2 .. for tumor pass reads, "tumor", "R1"/R2"

}

I have not put code for reverse, but I get error as:

Caused by:   Process `file_print (1)` terminated with an error exit status (127) Command executed:   null

Command exit status:   127

I tried many things such as ${ output_file_name }, { output_file_name }, $output_file_name but in vain.

I'd like to improve few things in code:

a) The iteration/loop over file in datasetFile I think something like each could help, but I need order to be maintained

b) Order in which loop of datasetFile is created. For instance, I'd like to keep lanes L003, L004, L005 in this order. Or, L001, L002, L003 order, I hope this is clear.

I cannot move ahead or able to sort this by myself.

I do not quite understand how and when to use .map{}

Edited: after output added tuple for val and path. I still get same error.

Answer 1

I think you are on the right track using groupTuple(sort: true). Another way to get what you want would be to instead parse the fromFilePairs grouping key. Unless you're looking to do more than just concatenate (and decompress) FASTQ files, I think all the process really requires is an output prefix (in addition to the FASTQ files of course), for example:

params.reads = "/sc/arion/MM-*_R{1,2}_001.fastq.gz"


process concat_reads {

    tag { prefix }

    input:
    tuple val(prefix), path(reads, stageAs: 'reads/*')

    output:
    tuple val(prefix), path("${prefix}.fastq")

    """
    zcat -f ${reads} > "${prefix}.fastq"
    """
}

workflow {

    Channel.fromFilePairs( params.reads ) \
        | branch { key, fastqs ->

            def matcher = key =~ /^(MM-\w+)-DNA-(\w+)-.*_L(\d+).*$/

            matches: matcher
                def (_, sample_id, sample_type, lane) = matcher[0].findAll()

                def new_sample_id = "${sample_id}_${sample_type}"

                return tuple( new_sample_id, lane as int, *fastqs )

            other: true
                return tuple( key, *fastqs )
        } \
        | set { samples }

    samples.matches \
        | groupTuple( sort: true ) \
        | flatMap { sample_id, lanes, fq1s, fq2s ->
            [
                tuple( "${sample_id}_R1", fq1s ),
                tuple( "${sample_id}_R2", fq2s ),
            ]
        } \
        | concat_reads \
        | branch { prefix, fastq ->

            def idx = prefix.lastIndexOf('_')

            def sample_id = prefix.substring(0, idx)
            def suffix = prefix.substring(idx + 1)

            read1: suffix == 'R1'
                return tuple( sample_id, fastq )
            read2: suffix == 'R2'
                return tuple( sample_id, fastq )
        } \
        | set { result }

    result.read1 \
        | join( result.read2 ) \
        | map { sample_id, fastq1, fastq2 ->

            tuple( sample_id, [ fastq1, fastq2 ] )
        } \
        | view()
}

Results:

$ nextflow run main.nf
N E X T F L O W  ~  version 23.04.1
Launching `main.nf` [jolly_noether] DSL2 - revision: 4a81f7a875
executor >  local (8)
[b8/a0b534] process > concat_reads (MM-265_T_R2) [100%] 8 of 8 ✔
[MM-465_T, [/path/to/work/9c/8b6fd328c848379423a0aba098cb74/MM-465_T_R1.fastq, /path/to/work/34/913780b1da902b5035f3bf56aac038/MM-465_T_R2.fastq]]
[MM-465_N, [/path/to/work/69/595120a1ee0a42217d06855ca9d1bc/MM-465_N_R1.fastq, /path/to/work/8f/31c406c5180073de577502fdb88fbf/MM-465_N_R2.fastq]]
[MM-265_N, [/path/to/work/1b/9c17cefb5aa51b3ee69832dc65d915/MM-265_N_R1.fastq, /path/to/work/d0/aac9f46e4a7f7165693b2537ad5630/MM-265_N_R2.fastq]]
[MM-265_T, [/path/to/work/d5/6295ff109df3ccb829a6f0a97b6e45/MM-265_T_R1.fastq, /path/to/work/b8/a0b53499ca13a523d811ecf772449b/MM-265_T_R2.fastq]]

Note that unless the FASTQ files are small, I would try to avoid decompressing the FASTQ files if possible. Most bioinformatics tools that deal with large FASTQ files will likely have gzip support, meaning you can replace the zcat -f with a cat in the example above (and remember to add the .gz extension to the output filename).