profiling.nf

//
// Run profiling
//

include { MALT_RUN                                      } from '../../modules/nf-core/malt/run/main'
include { MEGAN_RMA2INFO as MEGAN_RMA2INFO_TSV          } from '../../modules/nf-core/megan/rma2info/main'
include { KRAKEN2_KRAKEN2                               } from '../../modules/nf-core/kraken2/kraken2/main'
include { KRAKEN2_STANDARD_REPORT                       } from '../../modules/local/kraken2_standard_report'
include { BRACKEN_BRACKEN                               } from '../../modules/nf-core/bracken/bracken/main'
include { CENTRIFUGE_CENTRIFUGE                         } from '../../modules/nf-core/centrifuge/centrifuge/main'
include { CENTRIFUGE_KREPORT                            } from '../../modules/nf-core/centrifuge/kreport/main'
include { METAPHLAN_METAPHLAN                           } from '../../modules/nf-core/metaphlan/metaphlan/main'
include { KAIJU_KAIJU                                   } from '../../modules/nf-core/kaiju/kaiju/main'
include { KAIJU_KAIJU2TABLE as KAIJU_KAIJU2TABLE_SINGLE } from '../../modules/nf-core/kaiju/kaiju2table/main'
include { DIAMOND_BLASTX                                } from '../../modules/nf-core/diamond/blastx/main'
include { MOTUS_PROFILE                                 } from '../../modules/nf-core/motus/profile/main'
include { MOTUS_PREPLONG                                } from '../../modules/nf-core/motus/preplong/main'
include { KRAKENUNIQ_PRELOADEDKRAKENUNIQ                } from '../../modules/nf-core/krakenuniq/preloadedkrakenuniq/main'
include { KMCP_SEARCH                                   } from '../../modules/nf-core/kmcp/search/main'
include { KMCP_PROFILE                                  } from '../../modules/nf-core/kmcp/profile/main'
include { GANON_CLASSIFY                                } from '../../modules/nf-core/ganon/classify/main'
include { GANON_REPORT                                  } from '../../modules/nf-core/ganon/report/main'
include { SYLPH_PROFILE                                 } from '../../modules/nf-core/sylph/profile/main'
include { SYLPHTAX_TAXPROF                              } from '../../modules/nf-core/sylphtax/taxprof/main'
include { MELON                                         } from '../../modules/nf-core/melon/main'
include { SINGLEM_PIPE                                  } from '../../modules/nf-core/singlem/singlem/main'

workflow PROFILING {
    take:
    reads // [ [ meta ], [ reads ] ]
    databases // [ [ meta ], path ]

    main:
    ch_versions = channel.empty()
    ch_multiqc_files = channel.empty()
    ch_raw_classifications = channel.empty()
    // These per-read ID taxonomic assingment
    ch_raw_profiles = channel.empty()
    // These are count tables

    /*
        COMBINE READS WITH POSSIBLE DATABASES
    */

    // Separate default 'short;long' (when necessary) databases when short/long specified in database sheet
    ch_dbs = databases
        .map { meta_db, db ->
            [[meta_db.db_type.split(";")].flatten(), meta_db, db]
        }
        .transpose(by: 0)
        .map { type, meta_db, db ->
            [[type: type], meta_db.subMap(meta_db.keySet() - 'db_type') + [type: type], db]
        }

    // Join short and long reads with their corresponding short/long database
    // Note that for not-specified `short;long`, it will match with the database.
    // E.g. if there is no 'long' reads the above generated 'long' database channel element
    //  will have nothing to join to and will be discarded
    // Final output [DUMP: reads_plus_db] [['id':'2612', 'run_accession':'combined', 'instrument_platform':'ILLUMINA', 'single_end':false, 'is_fasta':false, 'type':'short'], <reads_path>/2612.merged.fastq.gz, ['tool':'malt', 'db_name':'malt95', 'db_params':'"-id 90"', 'type':'short'], <db_path>/malt95]

    ch_input_for_profiling = reads
        .map { meta, input_reads ->
            [[type: meta.type], meta, input_reads]
        }
        .combine(ch_dbs, by: 0)
        .map { _db_type, meta, input_reads, db_meta, db ->
            [meta, input_reads, db_meta, db]
        }
        .branch { _meta, _input_reads, db_meta, _db ->
            centrifuge: db_meta.tool == 'centrifuge'
            diamond: db_meta.tool == 'diamond'
            kaiju: db_meta.tool == 'kaiju'
            kraken2: db_meta.tool == 'kraken2' || db_meta.tool == 'bracken'
            krakenuniq: db_meta.tool == 'krakenuniq'
            malt: db_meta.tool == 'malt'
            metaphlan: db_meta.tool == 'metaphlan'
            motus: db_meta.tool == 'motus'
            kmcp: db_meta.tool == 'kmcp'
            ganon: db_meta.tool == 'ganon'
            sylph: db_meta.tool == 'sylph'
            melon: db_meta.tool == 'melon'
            singlem: db_meta.tool == 'singlem'
            unknown: true
        }

    /*
        PREPARE PROFILER INPUT CHANNELS & RUN PROFILING
    */

    // Each tool as a slightly different input structure and generally separate
    // input channels for reads vs databases. We restructure the channel tuple
    // for each tool and make liberal use of multiMap to keep reads/databases
    // channel element order in sync with each other

    if (params.run_malt) {

        // MALT: We groupTuple to have all samples in one channel for MALT as database
        // loading takes a long time, so we only want to run it once per database
        ch_input_for_malt = ch_input_for_profiling.malt
            .map { _meta, input_reads, db_meta, db ->

                // Reset entire input meta for MALT to just database name,
                // as we don't run run on a per-sample basis due to huge datbaases
                // so all samples are in one run and so sample-specific metadata
                // unnecessary. Set as database name to prevent `null` job ID and prefix.
                def new_meta = db_meta + [id: db_meta.db_name]

                // Extend database parameters to specify whether to save alignments or not
                def sam_format = params.malt_save_reads ? ' --alignments ./ -za false' : ""
                new_meta.db_params = db_meta.db_params + sam_format

                [new_meta, input_reads, db]
            }
            .groupTuple(by: [0, 2])
            .multiMap { meta, input_reads, db ->
                reads: [meta, input_reads.flatten()]
                db: db
            }

        MALT_RUN(ch_input_for_malt.reads, ch_input_for_malt.db)

        ch_maltrun_for_megan = MALT_RUN.out.rma6
            .transpose()
            .map { meta, rma ->
                // re-extract meta from file names, use filename without rma to
                // ensure we keep paired-end information in downstream filenames
                // when no pair-merging
                def meta_new = meta + [db_name: meta.id, id: rma.baseName]

                [meta_new, rma]
            }

        MEGAN_RMA2INFO_TSV(ch_maltrun_for_megan, params.malt_generate_megansummary)
        ch_multiqc_files = ch_multiqc_files.mix(MALT_RUN.out.log)
        ch_versions = ch_versions.mix(MALT_RUN.out.versions.first(), MEGAN_RMA2INFO_TSV.out.versions.first())
        ch_raw_classifications = ch_raw_classifications.mix(ch_maltrun_for_megan)
        ch_raw_profiles = ch_raw_profiles.mix(MEGAN_RMA2INFO_TSV.out.txt)
    }

    if (params.run_kraken2 || params.run_bracken) {
        // Have to pick first element of db_params if using bracken,
        // as db sheet for bracken must have ; sep list to
        // distinguish between kraken and bracken parameters
        def ch_prepare_for_kraken2 = ch_input_for_profiling.kraken2
            .map { meta, input_reads, db_meta, db ->

                // Only take first element if one exists
                def parsed_params = db_meta['db_params'].split(";")
                def db_meta_new = [:]
                if (parsed_params.size() == 2) {
                    db_meta_new = db_meta + [db_params: parsed_params[0]]
                }
                else if (parsed_params.size() == 0) {
                    db_meta_new = db_meta + [db_params: ""]
                }
                else {
                    db_meta_new = db_meta + [db_params: parsed_params[0]]
                }

                [meta, input_reads, db_meta_new, db]
            }
            .filter { meta, _input_reads, db_meta_new, _db ->
                if (db_meta_new.tool == 'bracken' && (meta.instrument_platform == 'OXFORD_NANOPORE' || meta.instrument_platform == 'PACBIO_SMRT')) {
                    log.warn("[nf-core/taxprofiler] Bracken has not been evaluated for long-read datasets. Skipping Bracken for sample ${meta.id}.")
                }
                db_meta_new.tool == 'kraken2' || (db_meta_new.tool == 'bracken' && (meta.instrument_platform != 'OXFORD_NANOPORE' && meta.instrument_platform != 'PACBIO_SMRT'))
            }

        ch_input_for_kraken2 = ch_prepare_for_kraken2.multiMap { it ->
            reads: [it[0] + it[2], it[1]]
            db: it[3]
        }

        KRAKEN2_KRAKEN2(ch_input_for_kraken2.reads, ch_input_for_kraken2.db, params.kraken2_save_reads, params.kraken2_save_readclassifications)
        ch_multiqc_files = ch_multiqc_files.mix(KRAKEN2_KRAKEN2.out.report)
        ch_versions = ch_versions.mix(KRAKEN2_KRAKEN2.out.versions.first())
        ch_raw_classifications = ch_raw_classifications.mix(KRAKEN2_KRAKEN2.out.classified_reads_assignment)
        ch_raw_profiles = ch_raw_profiles.mix(
            KRAKEN2_KRAKEN2.out.report.map { meta, report ->
                [meta + [tool: meta.tool == 'bracken' ? 'kraken2-bracken' : meta.tool], report]
            }
        )
    }

    if (params.run_kraken2 && params.run_bracken) {
        // Remove files from 'pure' kraken2 runs, so only those aligned against Bracken & kraken2 database are used.
        ch_kraken2_output = KRAKEN2_KRAKEN2.out.report

        // If necessary, convert the eight column output to six column output.
        if (params.kraken2_save_minimizers) {
            ch_kraken2_output = KRAKEN2_STANDARD_REPORT(ch_kraken2_output).report
        }

        // Extract the database name to combine by.
        ch_bracken_databases = databases
            .filter { meta, _db -> meta.tool == 'bracken' }
            .map { meta, db -> [meta.db_name, meta, db] }

        // Combine back with the reads
        ch_input_for_bracken = ch_kraken2_output
            .map { meta, report -> [meta.db_name, meta, report] }
            .combine(ch_bracken_databases, by: 0)
            .map { key, meta, input_reads, db_meta, db ->

                // // Have to make a completely fresh copy here as otherwise
                // // was getting db_param loss due to upstream meta parsing at
                // // kraken2 input channel manipulation step
                def db_meta_keys = db_meta.keySet()
                def db_meta_new = db_meta.subMap(db_meta_keys)

                // Have to pick second element if using bracken, as first element
                // contains kraken parameters
                if (db_meta.tool == 'bracken') {

                    // Only take second element if one exists
                    def parsed_params = db_meta['db_params'].split(";")

                    if (parsed_params.size() == 2) {
                        db_meta_new = db_meta + [db_params: parsed_params[1]]
                    }
                    else {
                        db_meta_new = db_meta + [db_params: ""]
                    }
                }
                else {
                    db_meta_new['db_params']
                }

                [key, meta, input_reads, db_meta_new, db]
            }
            .multiMap { _key, meta, report, db_meta, db ->
                report: [meta + db_meta, report]
                db: db
            }

        BRACKEN_BRACKEN(ch_input_for_bracken.report, ch_input_for_bracken.db)
        ch_versions = ch_versions.mix(BRACKEN_BRACKEN.out.versions.first())
        ch_raw_profiles = ch_raw_profiles.mix(BRACKEN_BRACKEN.out.reports)
    }

    if (params.run_centrifuge) {

        ch_input_for_centrifuge = ch_input_for_profiling.centrifuge
            .filter {
                if (it[0].is_fasta) {
                    log.warn("[nf-core/taxprofiler] Centrifuge currently does not accept FASTA files as input. Skipping Centrifuge for sample ${it[0].id}.")
                }
                !it[0].is_fasta
            }
            .multiMap { it ->
                reads: [it[0] + it[2], it[1]]
                db: it[3]
            }

        CENTRIFUGE_CENTRIFUGE(ch_input_for_centrifuge.reads, ch_input_for_centrifuge.db, params.centrifuge_save_reads, params.centrifuge_save_reads)
        ch_versions = ch_versions.mix(CENTRIFUGE_CENTRIFUGE.out.versions.first())
        ch_raw_classifications = ch_raw_classifications.mix(CENTRIFUGE_CENTRIFUGE.out.results)

        // Ensure the correct database goes with the generated report for KREPORT
        ch_database_for_centrifugekreport = databases
            .filter { meta, _db -> meta.tool == 'centrifuge' }
            .map { meta, db -> [meta.db_name, meta, db] }

        ch_input_for_centrifuge_kreport = combineProfilesWithDatabase(CENTRIFUGE_CENTRIFUGE.out.results, ch_database_for_centrifugekreport)

        // Generate profile
        CENTRIFUGE_KREPORT(ch_input_for_centrifuge_kreport.profile, ch_input_for_centrifuge_kreport.db)
        ch_versions = ch_versions.mix(CENTRIFUGE_KREPORT.out.versions.first())
        ch_raw_profiles = ch_raw_profiles.mix(CENTRIFUGE_KREPORT.out.kreport)
        ch_multiqc_files = ch_multiqc_files.mix(CENTRIFUGE_KREPORT.out.kreport)
    }

    if (params.run_metaphlan) {

        ch_input_for_metaphlan = ch_input_for_profiling.metaphlan.multiMap { it ->
            reads: [it[0] + it[2], it[1]]
            db: it[3]
        }

        METAPHLAN_METAPHLAN(ch_input_for_metaphlan.reads, ch_input_for_metaphlan.db, params.metaphlan_save_samfiles)
        ch_versions = ch_versions.mix(METAPHLAN_METAPHLAN.out.versions.first())
        ch_raw_profiles = ch_raw_profiles.mix(METAPHLAN_METAPHLAN.out.profile)
        ch_multiqc_files = ch_multiqc_files.mix(METAPHLAN_METAPHLAN.out.profile)
    }

    if (params.run_kaiju) {

        ch_input_for_kaiju = ch_input_for_profiling.kaiju.multiMap { it ->
            reads: [it[0] + it[2], it[1]]
            db: it[3]
        }

        KAIJU_KAIJU(ch_input_for_kaiju.reads, ch_input_for_kaiju.db)
        ch_versions = ch_versions.mix(KAIJU_KAIJU.out.versions.first())
        ch_raw_classifications = ch_raw_classifications.mix(KAIJU_KAIJU.out.results)

        // Ensure the correct database goes with the generated report for KAIJU2TABLE
        ch_database_for_kaiju2table = databases
            .filter { meta, _db -> meta.tool == 'kaiju' }
            .map { meta, db -> [meta.db_name, meta, db] }

        ch_input_for_kaiju2table = combineProfilesWithDatabase(KAIJU_KAIJU.out.results, ch_database_for_kaiju2table)
        // Generate profile
        KAIJU_KAIJU2TABLE_SINGLE(ch_input_for_kaiju2table.profile, ch_input_for_kaiju2table.db, params.kaiju_taxon_rank)
        ch_versions = ch_versions.mix(KAIJU_KAIJU2TABLE_SINGLE.out.versions)
        ch_multiqc_files = ch_multiqc_files.mix(KAIJU_KAIJU2TABLE_SINGLE.out.summary)
        ch_raw_profiles = ch_raw_profiles.mix(KAIJU_KAIJU2TABLE_SINGLE.out.summary)
    }

    if (params.run_diamond) {
        ch_input_for_diamond = ch_input_for_profiling.diamond.multiMap { meta, input_reads, meta_db, db ->
            if (!meta.single_end) {
                log.warn("[nf-core/taxprofiler] DIAMOND does not accept paired-end files as input. Only read 1 will be used for profiling. Running DIAMOND for sample ${meta.id} using only read 1.")
            }
            reads: [meta + meta_db, meta.single_end ? input_reads : input_reads[0]]
            db: [meta_db, db]
        }

        // diamond only accepts single output file specification, therefore
        // this will replace output file!
        ch_diamond_reads_format = params.diamond_save_reads ? 'sam' : params.diamond_output_format

        DIAMOND_BLASTX(ch_input_for_diamond.reads, ch_input_for_diamond.db, ch_diamond_reads_format, [])
        ch_versions = ch_versions.mix(DIAMOND_BLASTX.out.versions.first())
        ch_raw_profiles = ch_raw_profiles.mix(DIAMOND_BLASTX.out.tsv)
        ch_multiqc_files = ch_multiqc_files.mix(DIAMOND_BLASTX.out.log)
    }

    if (params.run_motus) {

        ch_input_for_motus = ch_input_for_profiling.motus
            .filter {
                if (it[0].is_fasta) {
                    log.warn("[nf-core/taxprofiler] mOTUs currently does not accept FASTA files as input. Skipping mOTUs for sample ${it[0].id}.")
                }
                !it[0].is_fasta
            }
            .branch {
                longread: it[0].instrument_platform == 'OXFORD_NANOPORE' || it[0].instrument_platform == 'PACBIO_SMRT'
                shortread: it[0].instrument_platform != 'OXFORD_NANOPORE' && it[0].instrument_platform != 'PACBIO_SMRT'
            }
        ch_input_for_motus_longread = ch_input_for_motus.longread.multiMap { it ->
            reads: [it[0] + it[2], it[1]]
            db: it[3]
        }

        MOTUS_PREPLONG(ch_input_for_motus_longread.reads, ch_input_for_motus_longread.db)

        ch_database_for_motus = databases
            .filter { meta, _db -> meta.tool == 'motus' }
            .map { meta, db -> [meta.db_name, meta, db] }


        ch_prepped_to_motus = MOTUS_PREPLONG.out.out
            .map { meta, in_reads -> [meta.db_name, [meta, in_reads]] }
            .combine(ch_database_for_motus, by: 0)
            .map { item ->
                def _db_name = item[0]
                def meta = item[1][0]
                def in_reads = item[1][1]
                def db_meta = item[2]
                def db = item[3]
                def meta_keys = ['id', 'run_accession', 'instrument_platform', 'single_end', 'is_fasta', 'type']
                def new_meta = meta.subMap(meta_keys)

                def db_meta_keys = db_meta.keySet()
                def new_db_meta = db_meta.subMap(db_meta_keys) + [type: meta.type]

                [new_meta, [in_reads], new_db_meta, db]
            }
            .mix(ch_input_for_motus.shortread)
            .multiMap { it ->
                reads: [it[0] + it[2], it[1]]
                db: it[3]
            }


        MOTUS_PROFILE(ch_prepped_to_motus.reads, ch_prepped_to_motus.db)

        ch_versions = ch_versions.mix(MOTUS_PREPLONG.out.versions.first())
        ch_versions = ch_versions.mix(MOTUS_PROFILE.out.versions.first())
        ch_raw_profiles = ch_raw_profiles.mix(MOTUS_PROFILE.out.out)
        ch_multiqc_files = ch_multiqc_files.mix(MOTUS_PROFILE.out.log)
    }

    if (params.run_krakenuniq) {

        ch_input_for_krakenuniq = ch_input_for_profiling.krakenuniq
            .map { meta, input_reads, db_meta, db ->
                def seqtype = input_reads[0].name ==~ /.+?\.f\w{0,3}a(\.gz)?$/ ? 'fasta' : 'fastq'
                // We bundle the sample identifier with the sequencing files to undergo batching.
                def prefix = params.perform_runmerging ? meta.id : "${meta.id}_${meta.run_accession}"
                prefix = meta.single_end ? "${prefix}.se" : "${prefix}.pe"
                [[id: db_meta.db_name, single_end: meta.single_end, seqtype: seqtype], [input_reads].flatten() + [prefix], db_meta, db]
            }
            .groupTuple(by: [0, 2, 3])
            .flatMap { single_meta, input_reads, db_meta, db ->
                // Sort reads array by comparing last element, prefix. This will ensure batch membership remains
                // constant across runs, enabling retrieval of cached tasks.
                input_reads.sort { a, b -> a[-1] <=> b[-1] }
                def batches = input_reads.collate(params.krakenuniq_batch_size)
                return batches.collect { batch ->
                    // We split the sample identifier from the reads again after batching.
                    def reads_batch = batch.collect { elements -> elements.take(elements.size() - 1) }.flatten()
                    def prefixes = batch.collect { elements -> elements[-1] }
                    return [single_meta + db_meta, reads_batch, prefixes, db]
                }
            }
            .multiMap { meta, input_reads, prefixes, db ->
                reads: [meta, input_reads, prefixes]
                db: db
                seqtype: meta.seqtype
            }
        // Hardcode to _always_ produce the report file (which is our basic output, and goes into)
        KRAKENUNIQ_PRELOADEDKRAKENUNIQ(ch_input_for_krakenuniq.reads, ch_input_for_krakenuniq.seqtype, ch_input_for_krakenuniq.db, params.krakenuniq_save_reads, true, params.krakenuniq_save_readclassifications)
        ch_multiqc_files = ch_multiqc_files.mix(KRAKENUNIQ_PRELOADEDKRAKENUNIQ.out.report)
        ch_versions = ch_versions.mix(KRAKENUNIQ_PRELOADEDKRAKENUNIQ.out.versions.first())
        ch_raw_classifications = ch_raw_classifications.mix(KRAKENUNIQ_PRELOADEDKRAKENUNIQ.out.classified_assignment.map { meta, profiles -> [meta - meta.subMap('seqtype'), profiles] })
        ch_raw_profiles = ch_raw_profiles.mix(KRAKENUNIQ_PRELOADEDKRAKENUNIQ.out.report.map { meta, profiles -> [meta - meta.subMap('seqtype'), profiles] })
    }

    if (params.run_kmcp) {

        ch_input_for_kmcp = ch_input_for_profiling.kmcp
            .filter { meta, _input_reads, meta_db, _db ->
                if (meta.instrument_platform == 'OXFORD_NANOPORE' || meta.instrument_platform == 'PACBIO_SMRT') {
                    log.warn("[nf-core/taxprofiler] KMCP is only suitable for short-read metagenomic profiling, with much lower sensitivity on long-read datasets. Skipping KMCP for sample ${meta.id}.")
                }
                meta_db.tool == 'kmcp' && (meta.instrument_platform != 'OXFORD_NANOPORE' && meta.instrument_platform != 'PACBIO_SMRT')
            }
            .map { meta, input_reads, db_meta, db ->
                def db_meta_keys = db_meta.keySet()
                def db_meta_new = db_meta.subMap(db_meta_keys)

                // Split the string, the arguments before semicolon should be parsed into kmcp search
                def parsed_params = db_meta_new.db_params.split(";")
                if (parsed_params.size() == 2) {
                    db_meta_new.db_params = parsed_params[0]
                }
                else if (parsed_params.size() == 0) {
                    db_meta_new.db_params = ""
                }
                else {
                    db_meta_new.db_params = parsed_params[0]
                }

                [meta, input_reads, db_meta_new, db]
            }
            .multiMap { it ->
                reads: [it[0] + it[2], it[1]]
                db: it[3]
            }

        KMCP_SEARCH(ch_input_for_kmcp.reads, ch_input_for_kmcp.db)

        ch_versions = ch_versions.mix(KMCP_SEARCH.out.versions.first())
        ch_raw_classifications = ch_raw_classifications.mix(KMCP_SEARCH.out.result)

        ch_database_for_kmcp_profile = databases
            .filter { meta, _db -> meta.tool == 'kmcp' }
            .map { meta, db -> [meta.db_name, meta, db] }

        ch_input_for_kmcp_profile = KMCP_SEARCH.out.result
            .map { meta, report -> [meta.db_name, meta, report] }
            .combine(ch_database_for_kmcp_profile, by: 0)
            .map { key, meta, input_reads, db_meta, db ->

                // Same as kraken2/bracken logic here. Arguments after semicolon are going into KMCP_PROFILE
                def db_meta_keys = db_meta.keySet()
                def db_meta_new = db_meta.subMap(db_meta_keys)

                def parsed_params = db_meta.db_params.split(";")

                if (parsed_params.size() == 2) {
                    db_meta_new = db_meta + [db_params: parsed_params[1]]
                }
                else {
                    db_meta_new = db_meta + [db_params: ""]
                }

                [key, meta, input_reads, db_meta_new, db]
            }
            .multiMap { _key, meta, report, db_meta, db ->
                report: [meta + db_meta, report]
                db: db
            }

        //Generate kmcp profile
        KMCP_PROFILE(ch_input_for_kmcp_profile.report, ch_input_for_kmcp.db)
        ch_versions = ch_versions.mix(KMCP_PROFILE.out.versions.first())
        ch_raw_profiles = ch_raw_profiles.mix(KMCP_PROFILE.out.profile)
        ch_multiqc_files = ch_multiqc_files.mix(KMCP_PROFILE.out.profile)
    }


    if (params.run_ganon) {

        ch_input_for_ganonclassify = ch_input_for_profiling.ganon
            .filter { meta, _input_reads, meta_db, _db ->
                if (meta.instrument_platform == 'OXFORD_NANOPORE' || meta.instrument_platform == 'PACBIO_SMRT') {
                    log.warn("[nf-core/taxprofiler] Ganon has not been evaluated for long-read datasets. Skipping Ganon for sample ${meta.id}.")
                }
                meta_db.tool == 'ganon' && (meta.instrument_platform != 'OXFORD_NANOPORE' && meta.instrument_platform != 'PACBIO_SMRT')
            }
            .multiMap { it ->
                reads: [it[0] + it[2], it[1]]
                db: it[3]
            }

        ch_input_for_ganonclassify.reads

        GANON_CLASSIFY(ch_input_for_ganonclassify.reads, ch_input_for_ganonclassify.db)
        ch_versions = ch_versions.mix(GANON_CLASSIFY.out.versions.first())

        ch_database_for_ganonreport = databases
            .filter { meta, _db -> meta.tool == "ganon" }
            .map { meta, db -> [meta.db_name, meta, db] }

        ch_report_for_ganonreport = combineProfilesWithDatabase(GANON_CLASSIFY.out.report, ch_database_for_ganonreport)

        GANON_REPORT(ch_report_for_ganonreport.profile, ch_report_for_ganonreport.db)
        ch_versions = ch_versions.mix(GANON_REPORT.out.versions.first())

        // Might be flipped - check/define what is a profile vs raw classification
        ch_raw_profiles = ch_raw_profiles.mix(GANON_REPORT.out.tre)
        ch_raw_classifications = ch_raw_classifications.mix(GANON_CLASSIFY.out.all)
    }

    if (params.run_sylph) {
        ch_input_for_sylph = ch_input_for_profiling.sylph
            .filter {
                if (it[0].is_fasta) {
                    log.warn("[nf-core/taxprofiler] sylph currently does not accept FASTA files as input. Skipping sylph for sample ${it[0].id}.")
                }
                !it[0].is_fasta
            }
            .map { meta, in_reads, db_meta, db ->
                def db_meta_keys = db_meta.keySet()
                def db_meta_new = db_meta.subMap(db_meta_keys)

                // Split the string, the arguments before semicolon should be parsed into sylph profile
                def parsed_params = db_meta_new['db_params'].split(";")
                if (parsed_params.size() == 2) {
                    db_meta_new['db_params'] = parsed_params[0]
                }
                else if (parsed_params.size() == 0) {
                    db_meta_new['db_params'] = ""
                }
                else {
                    db_meta_new['db_params'] = parsed_params[0]
                }

                [meta, in_reads, db_meta_new, db]
            }
            .multiMap { it ->
                reads: [it[0] + it[2], it[1]]
                db: it[3]
            }

        SYLPH_PROFILE(ch_input_for_sylph.reads, ch_input_for_sylph.db)
        ch_versions = ch_versions.mix(SYLPH_PROFILE.out.versions.first())

        ch_database_for_sylph_profile = databases
            .filter { meta, _db -> meta.tool == 'sylph' }
            .map { meta, db -> [meta.db_name, meta, db] }

        ch_input_for_sylphtax = SYLPH_PROFILE.out.profile_out
            .map { meta, report -> [meta.db_name, meta, report] }
            .combine(ch_database_for_sylph_profile, by: 0)
            .map { key, meta, in_reads, db_meta, db ->

                // Same as kraken2/bracken logic here. Arguments after semicolon are going into sylph-tax taxprof
                def db_meta_keys = db_meta.keySet()
                def db_meta_new = db_meta.subMap(db_meta_keys)

                def parsed_params = db_meta['db_params'].split(";")

                if (parsed_params.size() == 2) {
                    db_meta_new = db_meta + [db_params: parsed_params[1]]
                }
                else {
                    db_meta_new = db_meta + [db_params: ""]
                }

                [key, meta, in_reads, db_meta_new, db]
            }
            .multiMap { _key, meta, report, db_meta, db ->
                report: [meta + db_meta, report]
                db: db
            }

        SYLPHTAX_TAXPROF(ch_input_for_sylphtax.report, file(params.sylph_taxonomy, checkExists: true))
        ch_versions = ch_versions.mix(SYLPHTAX_TAXPROF.out.versions.first())
        ch_raw_profiles = ch_raw_profiles.mix(SYLPHTAX_TAXPROF.out.taxprof_output)
    }

    if (params.run_melon) {

        ch_input_for_melon = ch_input_for_profiling.melon
            .filter { meta, _in_reads, meta_db, _db ->
                if (meta['type'] != 'long') {
                    log.warn("[nf-core/taxprofiler] Melon is only suitable for long-read metagenomic profiling. Skipping Melon for sample ${meta.id}.")
                }
                meta_db['tool'] == 'melon' && meta['type'] == 'long'
            }
            .multiMap { it ->
                reads: [it[0] + it[2], it[1]]
                db: it[3]
            }

        MELON(ch_input_for_melon.reads, ch_input_for_melon.db, [])
        ch_versions = ch_versions.mix(MELON.out.versions.first())
        ch_raw_classifications = ch_raw_classifications.mix(MELON.out.json_output)
        ch_raw_profiles = ch_raw_profiles.mix(MELON.out.tsv_output)
    }

if (params.run_singlem) {

    ch_input_for_singlem = ch_input_for_profiling.singlem
        .map { sample_meta, input_reads, db_meta, db ->
            def meta = sample_meta + [
                tool    : 'singlem',           // profiler name expected by Taxpasta
                db_name : db_meta.db_name,     // for grouping/ID
                db_type : db_meta.type ?: null // optional, if you use it elsewhere
            ]
            tuple(meta, input_reads, db)
        }

    SINGLEM_PIPE(ch_input_for_singlem)

    ch_raw_profiles = ch_raw_profiles.mix(SINGLEM_PIPE.out.results)
    ch_versions     = ch_versions.mix(SINGLEM_PIPE.out.versions)
    }
    
    emit:
    classifications = ch_raw_classifications
    profiles        = ch_raw_profiles // channel: [ val(meta), [ reads ] ] - should be text files or biom
    versions        = ch_versions // channel: [ versions.yml ]
    motus_version   = params.run_motus ? MOTUS_PROFILE.out.versions.first() : channel.empty()
    mqc             = ch_multiqc_files
}


// Custom Functions

/**
* Combine profiles with their original database, then separate into two channels.
*
* The channel elements are assumed to be tuples one of [ meta, profile ], and the
* database to be of [db_key, meta, database_file].
*
* @param ch_profile A channel containing a meta and the profilign report of a given profiler
* @param ch_database A channel containing a key, the database meta, and the database file/folders itself
* @return A multiMap'ed output channel with two sub channels, one with the profile and the other with the db
*/
def combineProfilesWithDatabase(ch_profile, ch_database) {

    return ch_profile
        .map { meta, profile -> [meta.db_name, meta, profile] }
        .combine(ch_database, by: 0)
        .multiMap { _key, meta, profile, _db_meta, db ->
            profile: [meta, profile]
            db: db
        }
}