Alignment process in Nextflow not picking up the index files from the index process

Question

I'm using bowtie2 to build a small index. After building the index, I want to pass it to the bowtie2 alignment process. The issue is bowtie2-build outputs multiple index files with the .bt2 extension. When I try to provide these index files as an input to the Alignmnet process, it gives the following error:

Error executing process > 'SAM_FOR_STRAND (1)'

Caused by: Process SAM_FOR_STRAND (1) terminated with an error exit status (255)

Command executed:

bowtie2 -t -x index.1.bt2 index.2.bt2 index.3.bt2 index.4.bt2 index.rev.1.bt2 index.rev.2.bt2 -1 SRR5655560_1.atria.fastq -2 SRR5655560_2.atria.fastq -S "SRR5655560.sam"

Command exit status: 255

Command output: (empty)

Command error: (ERR): "index.1.bt2" does not exist or is not a Bowtie 2 index Exiting now ...

executor > local (4) [71/c095e4] process > INDEX [100%] 1 of 1 ✔ [1f/35a1e6] process > GTF2BED [100%] 1 of 1 ✔ [fa/b1e933] process > TRIM (SRR5655560) [100%] 1 of 1 ✔ [d6/d47ead] process > SAM_FOR_STRAND (1) [100%] 1 of 1, failed: 1 ✘ Error executing process > 'SAM_FOR_STRAND (1)'

Caused by: Process SAM_FOR_STRAND (1) terminated with an error exit status (255)

Command executed:

bowtie2 -t -x index.1.bt2 index.2.bt2 index.3.bt2 index.4.bt2 index.rev.1.bt2 index.rev.2.bt2 -1 SRR5655560_1.atria.fastq -2 SRR5655560_2.atria.fastq -S "SRR5655560.sam"

The following is simplifies version of the pipeline:

// Declare syntax version
nextflow.enable.dsl=2

/* 
 * pipeline input parameters 
 */


params.reads = "$projectDir/data/*_{1,2}.fastq"
params.ref_genome = "$projectDir/data/Glycine_max.Glycine_max_v2.1.dna.toplevel.fa"
params.gtf_file = "$projectDir/data/Glycine_max.Glycine_max_v2.1.55.gtf"
params.outdir = "results"



log.info """\

        R N A S E Q - N F   P I P E L I N E    
        ===================================
        reads            : ${params.reads}
        reference genome : ${params.ref_genome}
        gtf              : ${params.gtf_file}
        outdir           : ${params.outdir}
        """
        .stripIndent()

/*
 * define the `index` process that creates a binary index
 * given the reference genome file
 */

process INDEX {

  input:
  path ref_genome

  output:
  path "index.*"

  script:
  """
  bowtie2-build $ref_genome index
  """
}

/*
 * define the `GTF2BED` process that creates a BED file
 * given the params.gtf_file as input. 
 */
 
process GTF2BED {

    input:
    path gtf

    output:
    path "${gtf.baseName}.bed"

    script:
    """
    gtf2bed --gtf "${gtf}" --bed "${gtf.baseName}.bed"
    """
}


/*
 * define the `TRIM` process that creates trimmed fastq files
 * given the reads fastq file. 
 */
 
process TRIM {

    tag { pair_id }

    input:
    tuple val(pair_id), path(reads)

    output:
    tuple val(pair_id), path(trimmed_reads), emit: trimmed_fastqs
    path "${pair_id}/*.atria.log*", emit: logs

    script:
    def (r1, r2) = reads

    trimmed_reads = reads.collect {
        "${pair_id}/${it.baseName}.atria.fastq"
    }

    """
    atria \\
        -r "${r1}" \\
        -R "${r2}" \\
        -o "${pair_id}" \\
        -t 1 \\
        -l 12 \\
        -q 30
    """
}

/*
 * define the `SAM_FOR_STRAND` process that creates a SAM file
 * to be used by `STRAND_CHECK` process to determine 
 * the strandedness of each fastq file pair. 
 */
 
process SAM_FOR_STRAND {

    input:
    path index
    tuple val(pair_id), path(reads)

    output:
    path "${pair_id}.sam"

    script:
    """
    bowtie2 -x "${index.baseName}" -1 ${reads[0]} -2 ${reads[1]} -S ${pair_id}.sam
    """
}


workflow {

    index_ch = INDEX(params.ref_genome)

    gtf2bed_ch=GTF2BED(params.gtf_file)

    read_pairs_ch = Channel.fromFilePairs(params.reads, checkIfExists:true)

    trim_out = TRIM(read_pairs_ch)

    sam_for_strand_ch = SAM_FOR_STRAND(index_ch, TRIM.out.trimmed_fastqs)

}

I've tried passing the index argument directly to the script in the SAM_FOR_STRAND Alignmnet process, it gave the "index.1.bt2" does not exist error. Then I tried giving the basename of the index files (as shown in the below). It still does not work, even though the index files are being created properly.

 /*
     * define the `SAM_FOR_STRAND` process that creates a SAM file
     * to be used by `STRAND_CHECK` process to determine 
     * the strandedness of each fastq file pair. 
     */
     
    process SAM_FOR_STRAND {
    
        input:
        path index
        tuple val(pair_id), path(reads)
    
        output:
        path "${pair_id}.sam"
    
        script:
        """
        bowtie2 -x "${index.baseName}" -1 ${reads[0]} -2 ${reads[1]} -S ${pair_id}.sam
        """
    }

Best wishes.

Answer 1

When we run Bowtie2, we just need to make sure we supply the index basename using -x. This is the index filename prefix minus the .X.bt2 extension. A safe assumption is that the first file in the collection of output files (sorted in lexicographical order) has the .1.bt2 suffix. We can then just strip this off to get the required prefix. For example:

process SAM_FOR_STRAND {

    input:
    path index_files, stageAs: "bt2_index/*"
    tuple val(pair_id), path(reads)

    output:
    path "${pair_id}.sam"

    script:
    def idx = index_files[0].getBaseName(2)    

    """
    bowtie2 \\
        -x "${idx}" \\
        -1 "${reads[0]}" \\ 
        -2 "${reads[1]}" \\ 
        -S "${pair_id}.sam"
    """
}