diff --git a/.github/workflows/workflowhub_update.yml b/.github/workflows/workflowhub_update.yml new file mode 100644 index 00000000..9b760491 --- /dev/null +++ b/.github/workflows/workflowhub_update.yml @@ -0,0 +1,18 @@ +name: Publish workflows on WorkflowHub + +on: + release: + types: [published] + +jobs: + wfh-submit: + name: WorkflowHub submission + runs-on: ubuntu-latest + steps: + - uses: actions/checkout@v4 + - name: Submit workflows + uses: workflowhub-eu/submission-action@v0 + env: + API_TOKEN: ${{ secrets.WORKFLOWHUB_API_TOKEN }} + with: + team_id: 204 diff --git a/CHANGELOG.md b/CHANGELOG.md index 874bd446..d9ca573d 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -14,6 +14,9 @@ Release 11 of sanger-tol/ascc - `KMER_COUNTER` has been replaced with `COBIONTID_KMERCOUNTER` to increase efficiency. - This required the addition of `REFORMAT_NPY_2_CSV` to generate the kmer table. - `FCSGX_RUNGX` has been updated to not depend on `modulecmd` instead, production profiles will instead default to a local installation of fcs_gx. Avoiding containerised options provided in the module. +- Samtools modules have been updated to `1.23.1`. +- `ORGANELLAR` and `GENOMIC` subworkflows have been merged into a single `ASCC` workflow. + - This was an artifact from when the two would have been doing significantly different processes. ### `Dependencies` @@ -24,6 +27,7 @@ Release 11 of sanger-tol/ascc | `REFORMAT_NPY_2_CSV` | npy_2_csv.py | NA | 1.0.0 | | `SAMTOOLS_DICT` | samtools | 1.22.1 | 1.23.1 | | `SAMTOOLS_FAIDX` | samtools | 1.22.1 | 1.23.1 | +| `SAMTOOLS_SORT` | samtools | 1.22.1 | 1.23.1 | | `MINIMAP2_ALIGN2` | minimap2 + samtools | 2.29 + 1.21 | 2.30 + 1.23.1 | ## [0.6.0] - Red Notebook [28/01/2025] diff --git a/assets/farm_testing/samplesheet.csv b/assets/farm_testing/samplesheet.csv index 7c704184..d0c1af05 100644 --- a/assets/farm_testing/samplesheet.csv +++ b/assets/farm_testing/samplesheet.csv @@ -1,4 +1,2 @@ sample,assembly_type,assembly_file -odCymConc1,PRIMARY,/lustre/scratch122/tol/data/d/8/a/6/e/b/Cymbastela_concentrica/assembly/draft/odCymConc1.20241101/odCymConc1.20241101.fa.gz -odCymConc1,HAPLO,/lustre/scratch122/tol/data/d/8/a/6/e/b/Cymbastela_concentrica/assembly/draft/odCymConc1.20241101/odCymConc1.20241101.haplotigs.decontaminated.fa.gz -odCymConc1,MITO,/lustre/scratch122/tol/data/d/8/a/6/e/b/Cymbastela_concentrica/assembly/draft/odCymConc1.20241101/odCymConc1-MT.fa.gz +odCymConc1,PRIMARY,/nfs/treeoflife-01/resources/nextflow/test-data/Laetiporus_sulphureus/assembly/release/gfLaeSulp1.1/insdc/GCA_927399515.1.fasta.gz diff --git a/bin/abnormal_contamination_check.py b/bin/abnormal_contamination_check.py index c3ed9143..6d39a02b 100755 --- a/bin/abnormal_contamination_check.py +++ b/bin/abnormal_contamination_check.py @@ -1,12 +1,13 @@ #!/usr/bin/env python -import general_purpose_functions as gpf -import sys +import argparse import os.path import pathlib -import argparse +import sys import textwrap +import general_purpose_functions as gpf + VERSION = "V1.2.0" DESCRIPTION = """ @@ -38,11 +39,33 @@ def parse_args(): parser.add_argument("assembly", type=str, help="Path to the fasta assembly file") parser.add_argument("summary_path", type=str, help="Path to the tiara summary file") parser.add_argument("-q", "--out_prefix", type=str, help="Output file prefix for the report") - parser.add_argument("-o", "--output", type=str, help="Path to output file", default="fcs-gx_alarm_indicator_file.txt") - parser.add_argument("-p", "--alarm_percentage", type=int, help="Percentage of putative contaminant sequence in genomic assembly that will trip the alarm", default=3) - parser.add_argument("-l", "--alarm_length_removed", type=int, help="Length of removed sequence is greater than default, greater than this will trip the alarm.", default=1e7) - parser.add_argument("-s", "--alarm_scaff_length", type=int, help="Length of largest scaffold removed to trip alarm.", default=1.8e6) - parser.add_argument("-t", "--alarm_scaff_percent_removed", type=float, help="Percentage of Scaffolds set for removal from assembly to trip the alarm.", default=10.0) + parser.add_argument( + "-o", "--output", type=str, help="Path to output file", default="fcs-gx_alarm_indicator_file.txt" + ) + parser.add_argument( + "-p", + "--alarm_percentage", + type=int, + help="Percentage of putative contaminant sequence in genomic assembly that will trip the alarm", + default=3, + ) + parser.add_argument( + "-l", + "--alarm_length_removed", + type=int, + help="Length of removed sequence is greater than default, greater than this will trip the alarm.", + default=1e7, + ) + parser.add_argument( + "-s", "--alarm_scaff_length", type=int, help="Length of largest scaffold removed to trip alarm.", default=1.8e6 + ) + parser.add_argument( + "-t", + "--alarm_scaff_percent_removed", + type=float, + help="Percentage of Scaffolds set for removal from assembly to trip the alarm.", + default=10.0, + ) parser.add_argument("-r", "--review_info", type=int, help="Number of REVIEW/INFO to the trigger alarm", default=0) parser.add_argument("-v", "--version", action="version", version=VERSION) return parser.parse_args() @@ -66,9 +89,7 @@ def load_fcs_gx_results(seq_dict, fcs_gx_and_tiara_summary_path): Loads FCS-GX actions from the FCS-GX and Tiara results summary file, adds them to the dictionary that contains sequence lengths """ fcs_gx_and_tiara_summary_data = gpf.l(fcs_gx_and_tiara_summary_path) - fcs_gx_and_tiara_summary_data = fcs_gx_and_tiara_summary_data[ - 1 : len(fcs_gx_and_tiara_summary_data) - ] + fcs_gx_and_tiara_summary_data = fcs_gx_and_tiara_summary_data[1 : len(fcs_gx_and_tiara_summary_data)] for line in fcs_gx_and_tiara_summary_data: split_line = line.split(",") assert len(split_line) == 5 @@ -87,9 +108,7 @@ def main(): sys.exit(1) if os.path.isfile(args.assembly) is False: - sys.stderr.write( - f"The assembly FASTA file was not found at the expected location ({args.assembly})\n" - ) + sys.stderr.write(f"The assembly FASTA file was not found at the expected location ({args.assembly})\n") sys.exit(1) seq_dict = get_sequence_lengths(args.assembly) @@ -115,7 +134,7 @@ def main(): "PERCENTAGE_LENGTH_REMOVED": args.alarm_percentage, "LARGEST_SCAFFOLD_REMOVED": args.alarm_scaff_length, "PERCENTAGE_SCAFFOLDS_REMOVED": args.alarm_scaff_percent_removed, - "REVIEW_OR_INFO": args.review_info + "REVIEW_OR_INFO": args.review_info, } report_dict = { @@ -124,7 +143,7 @@ def main(): "LARGEST_SCAFFOLD_REMOVED": max(lengths_removed, default=0), "SCAFFOLDS_REMOVED": scaffolds_removed, "PERCENTAGE_SCAFFOLDS_REMOVED": 100 * scaffolds_removed / scaffold_count, - "REVIEW_OR_INFO": review_info + "REVIEW_OR_INFO": review_info, } # Seperated out to ensure that the file is written in one go and doesn't confuse Nextflow diff --git a/bin/general_purpose_functions.py b/bin/general_purpose_functions.py index f53ef221..2a1de0d1 100755 --- a/bin/general_purpose_functions.py +++ b/bin/general_purpose_functions.py @@ -7,8 +7,7 @@ # # Copyright (c) 2020-2021 Genome Research Ltd. # -# Author: Eerik Aunin (eeaunin@gmail.com) -# +# Author: Eerik Aunin (eeaunin@gmail.com # This file is a part of the Genome Decomposition Analysis (GDA) pipeline. # # Permission is hereby granted, free of charge, to any person obtaining a copy @@ -29,13 +28,13 @@ # OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE # SOFTWARE. +import argparse import os -from os.path import isfile -import sys -import subprocess import signal +import subprocess +import sys from datetime import datetime -import argparse +from os.path import isfile def l(path): @@ -172,17 +171,15 @@ def string_to_chunks(line, n): return [line[i : i + n] for i in range(0, len(line), n)] -def run_system_command( - system_command, verbose=True, dry_run=False, tries=1, expected_exit_code=0 -): +def run_system_command(system_command, verbose=True, dry_run=False, tries=1, expected_exit_code=0): """ Executes a system command and checks its exit code """ triggering_script_name = sys.argv[0].split("/")[-1] try_counter_string = "" - if dry_run == False: + if not dry_run: for i in range(0, tries): - if verbose == True: + if verbose: time_now = datetime.now().strftime("%Y-%m-%d %H:%M:%S") if i > 0: try_counter_string = ", try {}".format(i + 1) @@ -201,11 +198,7 @@ def run_system_command( break except subprocess.CalledProcessError as exc: out_errormessage = ( - "<" - + triggering_script_name - + "> " - + " exited with error code " - + str(exc.returncode) + "<" + triggering_script_name + "> " + " exited with error code " + str(exc.returncode) ) if exc.output.isspace() == False: out_errormessage += ". Error message: " + exc.output @@ -234,11 +227,7 @@ def get_file_paths(in_folder_path, extension): onlyfiles = list() selected_file_paths = list() if os.path.isdir(in_folder_path): - onlyfiles = [ - f - for f in os.listdir(in_folder_path) - if os.path.isfile(os.path.join(in_folder_path, f)) - ] + onlyfiles = [f for f in os.listdir(in_folder_path) if os.path.isfile(os.path.join(in_folder_path, f))] for file_item in onlyfiles: if "." + extension in file_item: file_item_split = file_item.split(".") diff --git a/conf/modules.config b/conf/modules.config index 61991313..34e0f0b1 100644 --- a/conf/modules.config +++ b/conf/modules.config @@ -155,8 +155,8 @@ process { ] } - withName: ".*:ASCC_ORGANELLAR:RUN_DECONTAMINATE_FASTA:DECONTAMINATE_GENERATE_BED" { - ext.args = "--is_organelle True" + withName: ".*:ASCC_ASSEMBLY:RUN_DECONTAMINATE_FASTA:DECONTAMINATE_GENERATE_BED" { + ext.args = { meta.assembly_type == "organellar" ? "--is_organelle True" : "" } } withName: "ASCC_MERGE_TABLES|DECONTAMINATE_GENERATE_BED|GZIP" { diff --git a/conf/production.config b/conf/production.config index 1efe4d91..bbdd1da9 100644 --- a/conf/production.config +++ b/conf/production.config @@ -20,8 +20,6 @@ process { // Allowing it to use default resources, stops it showing up in LSF? withName: FCSGX_RUNGX { container = "" - module = "" - //module = "fcs-gx/farm/0.5.5" cpus = { 32 } memory = { 520.GB } time = { 12.h } @@ -67,13 +65,13 @@ params { // for sanger decontamination. run_essentials = "both" - run_kmers = "off" //genomic + run_kmers = "genomic" //genomic run_tiara = "both" run_coverage = "both" - run_nt_blast = "off" //both - run_nr_diamond = "off" //both - run_uniprot_diamond = "off" //both - run_kraken = "off" //both + run_nt_blast = "both" //both + run_nr_diamond = "both" //both + run_uniprot_diamond = "both" //both + run_kraken = "both" //both run_fcsgx = "both" //both run_fcs_adaptor = "both" //both run_vecscreen = "both" //both diff --git a/conf/test.config b/conf/test.config index 6d98331b..a190ff5e 100644 --- a/conf/test.config +++ b/conf/test.config @@ -10,14 +10,33 @@ ---------------------------------------------------------------------------------------- */ + process { - resourceLimits = [ - cpus: 4, - memory: '15.GB', - time: '1.h' - ] + // Adding the module info here + // stopped the resources from base.config + // from being used, and causes the job + // to only use the defaulted process_* tag + // which is just not enough for FCS_GX + // TODO + // Allowing it to use default resources, stops it showing up in LSF? + withName: FCSGX_RUNGX { + container = "" + cpus = { 32 } + memory = { 520.GB } + time = { 12.h } + } + + withName: SANGER_TOL_BTK { + clusterOptions = {"-J NF_ASCC::BLOBTOOLKIT(${meta.id})"} + queue = "oversubscribed" + cpus = { 2 } + memory = { 1200.MB * task.attempt } + time = { 96.h * task.attempt } + } } + + params { config_profile_name = 'Test profile' config_profile_description = 'Minimal test dataset to check pipeline function' @@ -82,7 +101,7 @@ params { run_organellar_blast = "genomic" run_autofilter_assembly = "both" run_create_btk_dataset = "both" - run_merge_datasets = "genomic" + run_merge_datasets = "both" run_decontaminate_fasta = "both" run_html_report = "both" } diff --git a/conf/test_full.config b/conf/test_full.config index 4f0b9408..a73207cb 100644 --- a/conf/test_full.config +++ b/conf/test_full.config @@ -64,21 +64,21 @@ params { run_essentials = "both" // both run_kmers = "genomic" // genomic run_tiara = "both" // both - run_coverage = "both" // both - run_nt_blast = "both" // both - run_nr_diamond = "both" // both - run_uniprot_diamond = "both" // both - run_kraken = "both" // both + run_coverage = "off" // both + run_nt_blast = "off" // both + run_nr_diamond = "off" // both + run_uniprot_diamond = "off" // both + run_kraken = "off" // both run_fcsgx = "both" // both - run_fcs_adaptor = "both" // both - run_vecscreen = "both" // both - run_btk_busco = "genomic" // genomic - run_pacbio_barcodes = "both" // both - run_organellar_blast = "genomic" // genomic + run_fcs_adaptor = "off" // both + run_vecscreen = "off" // both + run_btk_busco = "off" // genomic + run_pacbio_barcodes = "off" // both + run_organellar_blast = "off" // genomic run_autofilter_assembly = "both" // both - run_create_btk_dataset = "both" // both - run_merge_datasets = "genomic" // genomic + run_create_btk_dataset = "off" // both + run_merge_datasets = "off" // genomic run_decontaminate_fasta = "both" // both - run_html_report = "both" // both + run_html_report = "off" // both } diff --git a/modules.json b/modules.json index 2e224b12..71830a1b 100644 --- a/modules.json +++ b/modules.json @@ -37,12 +37,6 @@ "installed_by": ["modules"], "patch": "modules/nf-core/fcs/fcsadaptor/fcs-fcsadaptor.diff" }, - "fcsgx/rungx": { - "branch": "master", - "git_sha": "6d46786420b4d7bc88eba026eb389c0c5535d120", - "installed_by": ["modules"], - "patch": "modules/nf-core/fcsgx/rungx/fcsgx-rungx.diff" - }, "gnu/sort": { "branch": "master", "git_sha": "7ad318b7b2c0ef6a101c01b6083b8acd6a9a63de", @@ -151,6 +145,16 @@ "git_sha": "f919028603ca42cb01a59e45475c19106840372d", "installed_by": ["fastx_map_long_reads"] }, + "fcsgx/parseresults": { + "branch": "main", + "git_sha": "17d7d065459873b60bb795714fac19137304b793", + "installed_by": ["modules"] + }, + "fcsgx/rungx": { + "branch": "main", + "git_sha": "864e75dc4bf5b4fc50f7a5e568a4cb6236bfb528", + "installed_by": ["modules"] + }, "samtools/mergedup": { "branch": "main", "git_sha": "729e9c8cfa83bc64c95ea6024dd5477c936e305e", diff --git a/modules/local/ascc/merge_tables/main.nf b/modules/local/ascc/merge_tables/main.nf index 5aac0459..e0e291e9 100644 --- a/modules/local/ascc/merge_tables/main.nf +++ b/modules/local/ascc/merge_tables/main.nf @@ -34,7 +34,7 @@ process ASCC_MERGE_TABLES { script: def args = task.ext.args ?: "" - def empty_file_size = 80 + def empty_file_size = 50 def coverage_data = coverage.size() > empty_file_size ? "-c ${coverage}" : "" def tiara_data = tiara.size() > empty_file_size ? "-t ${tiara}" : "" def nt_kraken_data = nt_kraken.size() > empty_file_size ? "-nk ${nt_kraken}" : "" diff --git a/modules/nf-core/fcsgx/rungx/fcsgx-rungx.diff b/modules/nf-core/fcsgx/rungx/fcsgx-rungx.diff deleted file mode 100644 index 1c6f7567..00000000 --- a/modules/nf-core/fcsgx/rungx/fcsgx-rungx.diff +++ /dev/null @@ -1,106 +0,0 @@ -Changes in component 'nf-core/fcsgx/rungx' -'modules/nf-core/fcsgx/rungx/meta.yml' is unchanged -'modules/nf-core/fcsgx/rungx/environment.yml' is unchanged -Changes in 'fcsgx/rungx/main.nf': ---- modules/nf-core/fcsgx/rungx/main.nf -+++ modules/nf-core/fcsgx/rungx/main.nf -@@ -1,5 +1,5 @@ - process FCSGX_RUNGX { -- tag "$meta.id" -+ tag "${meta.id}" - label 'process_high' - - conda "${moduleDir}/environment.yml" -@@ -8,9 +8,10 @@ - 'quay.io/biocontainers/ncbi-fcs-gx:0.5.5--h9948957_0' }" - - input: -- tuple val(meta), val(taxid), path(fasta) -+ tuple val(meta), path(fasta) - path gxdb -- val ramdisk_path -+ path ramdisk_path -+ val production_mode - - output: - tuple val(meta), path("*.fcs_gx_report.txt"), emit: fcsgx_report -@@ -25,32 +26,44 @@ - script: - def args = task.ext.args ?: '' - def prefix = task.ext.prefix ?: "${meta.id}" -+ -+ // At Sanger we have a permenant home for the DB on NVME storage -+ // def mv_database_to_ram = ramdisk_path ? "rclone copy $gxdb $ramdisk_path/$task.index/" : '' -+ // def database = ramdisk_path ? "$ramdisk_path/$task.index/" : gxdb // Use task.index to make memory location unique - def database = ramdisk_path ?: gxdb -- ( ramdisk_path ? -- """ -- if [ -d "${database}" ]; then -- echo "ERROR: Database exists in memory, and may be in use by another process" >&2 -- ls -l ${database} -- exit 1 -- fi -- # Clean up shared memory on exit -- trap "rm -rf ${database}" EXIT -- # Copy DB to RAM-disk when supplied. Otherwise, rungx is very slow. -- rclone copy ${gxdb} ${database} - -- """: "") -- << -- """ -- export GX_NUM_CORES=${task.cpus} -- run_gx.py \\ -- --fasta ${fasta} \\ -- --gx-db ${database} \\ -- --tax-id ${taxid} \\ -- --generate-logfile true \\ -- --out-basename ${prefix} \\ -- --out-dir . \\ -- ${args} -- """ -+ if ( production_mode ) { -+ """ -+ echo "Using Production FCSGX with local installation" -+ -+ export GX_NUM_CORES=${task.cpus} -+ export GX_INSTANTIATE_FASTA=1 -+ -+ run_gx \\ -+ --fasta ${fasta} \\ -+ --gx-db ${database} \\ -+ --tax-id ${meta.taxid} \\ -+ --generate-logfile true \\ -+ --out-basename ${prefix} \\ -+ --out-dir . \\ -+ ${args} -+ -+ """ -+ } else { -+ """ -+ echo "Using Standard FCSGX with container" -+ -+ run_gx.py \\ -+ --fasta ${fasta} \\ -+ --gx-db ${database} \\ -+ --tax-id ${meta.taxid} \\ -+ --generate-logfile true \\ -+ --out-basename ${prefix} \\ -+ --out-dir . \\ -+ ${args} -+ -+ """ -+ } - - stub: - // def args = task.ext.args ?: '' -@@ -60,5 +73,6 @@ - touch ${prefix}.taxonomy.rpt - touch ${prefix}.summary.txt - echo "" | gzip > ${prefix}.hits.tsv.gz -+ - """ - } - -'modules/nf-core/fcsgx/rungx/tests/main.nf.test' is unchanged -'modules/nf-core/fcsgx/rungx/tests/main.nf.test.snap' is unchanged -************************************************************ diff --git a/modules/sanger-tol/fcsgx/parseresults/environment.yml b/modules/sanger-tol/fcsgx/parseresults/environment.yml new file mode 100644 index 00000000..380e2d27 --- /dev/null +++ b/modules/sanger-tol/fcsgx/parseresults/environment.yml @@ -0,0 +1,7 @@ +--- +# yaml-language-server: $schema=https://raw.githubusercontent.com/nf-core/modules/master/modules/environment-schema.json +channels: + - conda-forge + - bioconda +dependencies: + - conda-forge::python=3.14.3 diff --git a/modules/sanger-tol/fcsgx/parseresults/main.nf b/modules/sanger-tol/fcsgx/parseresults/main.nf new file mode 100644 index 00000000..1d122355 --- /dev/null +++ b/modules/sanger-tol/fcsgx/parseresults/main.nf @@ -0,0 +1,37 @@ +process FCSGX_PARSERESULTS { + tag "${meta.id}" + label 'process_low' + + conda "${moduleDir}/environment.yml" + container "${ workflow.containerEngine == 'singularity' && !task.ext.singularity_pull_docker_container ? + 'https://depot.galaxyproject.org/singularity/python:3.14' : + 'quay.io/biocontainers/python:3.14' }" + + input: + tuple val(meta), path(taxonomy_report) + tuple val(meta2), path(fcs_report) + path ncbi_rankedlineage_path + + output: + tuple val(meta), path( "*_parsed_fcsgx.csv" ), emit: fcsgxresult + tuple val("${task.process}"), val('python'), eval('python --version | sed "s/Python //"'), emit: versions_python, topic: versions + tuple val("${task.process}"), val('parse_fcsgx_result'), eval("parse_fcsgx_result.py --version"), emit: versions_parsefcsgx, topic: versions + + when: + task.ext.when == null || task.ext.when + + script: + def prefix = task.ext.prefix ?: "${meta.id}" + """ + parse_fcsgx_result.py \\ + ${taxonomy_report} \\ + ${fcs_report} \\ + ${ncbi_rankedlineage_path} > ${prefix}_parsed_fcsgx.csv + """ + + stub: + def prefix = task.ext.prefix ?: "${meta.id}" + """ + touch ${prefix}_parsed_fcsgx.csv + """ +} diff --git a/modules/sanger-tol/fcsgx/parseresults/meta.yml b/modules/sanger-tol/fcsgx/parseresults/meta.yml new file mode 100644 index 00000000..fe43e9d4 --- /dev/null +++ b/modules/sanger-tol/fcsgx/parseresults/meta.yml @@ -0,0 +1,102 @@ +name: fcsgx_parseresults +description: | + Parses FCS-GX (Fungal Contamination Screening - GX) result files and augments them + with NCBI taxonomy lineage information to produce a comprehensive CSV output containing + contamination screening results with taxonomic classification. + +keywords: + - fcs-gx + - contamination + - taxonomy + - parsing +tools: + - parse_fcsgx_result: + description: Parse FCS-GX result files and augment with NCBI taxonomy lineage information + homepage: https://github.com/sanger-tol/treeval + documentation: https://github.com/sanger-tol/treeval + tool_dev_url: https://github.com/sanger-tol/treeval + doi: no DOI available + licence: + - MIT + identifier: "" +input: + - - meta: + type: map + description: | + Groovy Map containing sample information + e.g. [ id:'sample1' ] + - taxonomy_report: + type: file + description: FCS-GX taxonomy report file + pattern: "*.taxonomy.rpt" + ontologies: [] + - - meta2: + type: map + description: | + Groovy Map containing sample information + e.g. [ id:'sample1' ] + - fcs_report: + type: file + description: FCS-GX report file + pattern: "*.fcs_gx_report.txt" + ontologies: + - edam: http://edamontology.org/format_2330 # TXT + - ncbi_rankedlineage_path: + type: file + description: NCBI rankedlineage.dmp file + pattern: "*rankedlineage*" + ontologies: [] +output: + fcsgxresult: + - - meta: + type: map + description: | + Groovy Map containing sample information + e.g. `[ id:'sample1' ]` + - "*_parsed_fcsgx.csv": + type: file + description: Parsed FCS-GX results with NCBI taxonomy lineage information in CSV format + ontologies: + - edam: http://edamontology.org/format_3752 # CSV + versions_python: + - - ${task.process}: + type: string + description: The name of the process + - python: + type: string + description: The name of the tool + - python --version | sed "s/Python //": + type: eval + description: The expression to obtain the version of the tool + versions_parsefcsgx: + - - ${task.process}: + type: string + description: The name of the process + - parse_fcsgx_result: + type: string + description: The name of the tool + - parse_fcsgx_result.py --version: + type: eval + description: The expression to obtain the version of the tool +topics: + versions: + - - ${task.process}: + type: string + description: The name of the process + - python: + type: string + description: The name of the tool + - python --version | sed "s/Python //": + type: eval + description: The expression to obtain the version of the tool + - - ${task.process}: + type: string + description: The name of the process + - parse_fcsgx_result: + type: string + description: The name of the tool + - parse_fcsgx_result.py --version: + type: eval + description: The expression to obtain the version of the tool +authors: + - "@DLBPointon" diff --git a/modules/sanger-tol/fcsgx/parseresults/nextflow.config b/modules/sanger-tol/fcsgx/parseresults/nextflow.config new file mode 100644 index 00000000..651f0b86 --- /dev/null +++ b/modules/sanger-tol/fcsgx/parseresults/nextflow.config @@ -0,0 +1 @@ +nextflow.enable.moduleBinaries = true diff --git a/modules/sanger-tol/fcsgx/parseresults/resources/usr/bin/parse_fcsgx_result.py b/modules/sanger-tol/fcsgx/parseresults/resources/usr/bin/parse_fcsgx_result.py new file mode 100755 index 00000000..1edcd87d --- /dev/null +++ b/modules/sanger-tol/fcsgx/parseresults/resources/usr/bin/parse_fcsgx_result.py @@ -0,0 +1,224 @@ +#!/usr/bin/env python3 +""" +Script for parsing the result files of FCS-GX, originally written by Eerik Aunin (ea10) +further refactoring/modifications by Yumi Sims (yy5) +Updated by Damon-Lee Pointon (dp24) +""" + +import argparse +import os +import sys +from pathlib import Path + + +def file_to_generator(input_file: str): + """ + Brought in from Eeriks GPF library + + Originally as gpf.l(input_file) + """ + file_path = Path(input_file) + + if file_path.exists(): + with open(file_path) as f: + for line in f: + line = line.rstrip() + yield line + + +def file_to_list(path: str) -> list: + """ + Function for loading text file as a list and removing line breaks from line ends + """ + try: + return list(file_to_generator(path)) + except FileNotFoundError: + sys.stderr.write("Error: file not found (" + path + ")\n") + sys.exit(1) + + +def load_taxids_data(taxonomy_file: str) -> dict: + """ + Parses the *.taxonomy.rpt to find taxids that correspond to species names. Returns this as a dictionary + """ + taxonomy_data: list = file_to_list(taxonomy_file)[2:] + assert taxonomy_data, "Taxonomy data is empty" + + collection_dict = dict() + + for line in taxonomy_data: + split_line = line.split("\t") + assert len(split_line) == 34 + + scaff = split_line[0].split("~")[0] + tax_name_1, tax_id_1, div_1, cvg_by_div_1, cvg_by_tax_1, score_1 = ( + split_line[5], + split_line[6], + split_line[7], + int(split_line[8]) if split_line[8] else None, + int(split_line[9]) if split_line[9] else None, + int(split_line[10]) if split_line[10] else 0, + ) + + if scaff in collection_dict and int(collection_dict[scaff]["fcs_gx_score"] or 0) < score_1: + row_dict = { + "fcs_gx_top_tax_name": tax_name_1, + "fcs_gx_top_taxid": tax_id_1, + "fcs_gx_div": div_1, + "fcs_gx_coverage_by_div": cvg_by_div_1, + "fcs_gx_coverage_by_tax": cvg_by_tax_1, + "fcs_gx_score": score_1, + "fcs_gx_multiple_divs_per_scaff": True, + "fcs_gx_action": "NA", + } + collection_dict[scaff] = row_dict + elif scaff not in collection_dict: + row_dict = { + "fcs_gx_top_tax_name": tax_name_1, + "fcs_gx_top_taxid": tax_id_1, + "fcs_gx_div": div_1, + "fcs_gx_coverage_by_div": cvg_by_div_1, + "fcs_gx_coverage_by_tax": cvg_by_tax_1, + "fcs_gx_score": score_1, + "fcs_gx_multiple_divs_per_scaff": False, + "fcs_gx_action": "NA", + } + collection_dict[scaff] = row_dict + + return collection_dict + + +def load_report_data(report_file: str, collection_dict: dict) -> dict: + """ + Parses the *.fcs_gx_report.txt to add entries from the 'action' column to the collection of entries per scaffold that is stored in collection_dict + """ + report_data = file_to_list(report_file) + if len(report_data) > 2: + report_data = report_data[2 : len(report_data)] + for line in report_data: + split_line = line.split("\t") + assert len(split_line) == 8 + scaff: str = split_line[0] + fcs_gx_action: str = split_line[4] + collection_dict[scaff]["fcs_gx_action"] = fcs_gx_action + return collection_dict + + +def get_taxids_list(fcs_gx_taxonomy_file_path: str) -> list: + """ + Goes through FCS-GX taxonomy output file and returns a list of unique taxIDs found in the file + """ + if not os.path.isfile(fcs_gx_taxonomy_file_path): + sys.stderr.write( + f"The FCS-GX taxonomy file was not found at the expected location ({fcs_gx_taxonomy_file_path})\n" + ) + sys.exit(1) + taxids_list = list() + for line in file_to_list(fcs_gx_taxonomy_file_path): + if not line.startswith("#"): + split_line = line.split("\t") + assert len(split_line) == 34 + taxid = split_line[6] + if taxid not in taxids_list: + taxids_list.append(taxid) + return taxids_list + + +def get_lineages_by_taxid(taxids_list: list, rankedlineage_path: str) -> dict: + """ + Takes a list of taxIDs and the path to the NCBI rankedlineage.dmp file as the input. Returns the lineage corresponding to each taxID + """ + lineages_dict = dict() + rankedlineage_col_names = ( + "taxid", + "fcs_gx_name", + "fcs_gx_species", + "fcs_gx_genus", + "fcs_gx_family", + "fcs_gx_order", + "fcs_gx_class", + "fcs_gx_phylum", + "fcs_gx_kingdom", + "fcs_gx_domain", + ) + + for line in file_to_list(rankedlineage_path): + split_line: list = line.split("|") + split_line: list = [n.strip() for n in split_line] + assert len(split_line) >= 11, ( + f"Expected at least 11 columns in rankedlineage.dmp, got {len(split_line)}" + ) # This should now handle both new and old formats (as of April 2025) + taxid = split_line[0] + if taxid in taxids_list: + current_lineage_dict = dict() + for i in range(1, 10): + current_lineage_dict[rankedlineage_col_names[i]] = split_line[i] + lineages_dict[taxid] = current_lineage_dict + return lineages_dict + + +def main(taxonomy_report: str, fcs_report: str, ncbi_rankedlineage_path: str) -> None: + collection_dict = load_taxids_data(taxonomy_report) + collection_dict = load_report_data(fcs_report, collection_dict) + taxids_list = get_taxids_list(taxonomy_report) + lineages_dict = get_lineages_by_taxid(taxids_list, ncbi_rankedlineage_path) + rankedlineage_col_names = ( + "taxid", + "fcs_gx_name", + "fcs_gx_species", + "fcs_gx_genus", + "fcs_gx_family", + "fcs_gx_order", + "fcs_gx_class", + "fcs_gx_phylum", + "fcs_gx_kingdom", + "fcs_gx_domain", + ) + + out_header = "scaff,fcs_gx_top_tax_name,fcs_gx_top_taxid,fcs_gx_div,fcs_gx_coverage_by_div,fcs_gx_coverage_by_tax,fcs_gx_score,fcs_gx_multiple_divs_per_scaff,fcs_gx_action" + out_header += "," + ",".join(rankedlineage_col_names) + print(out_header) + for scaff, row_dict in collection_dict.items(): + out_line = [ + scaff, + row_dict["fcs_gx_top_tax_name"], + row_dict["fcs_gx_top_taxid"], + row_dict["fcs_gx_div"], + row_dict["fcs_gx_coverage_by_div"], + row_dict["fcs_gx_coverage_by_tax"], + row_dict["fcs_gx_score"], + row_dict["fcs_gx_multiple_divs_per_scaff"], + row_dict["fcs_gx_action"], + ] + row_top_taxid: str = row_dict["fcs_gx_top_taxid"] + if row_top_taxid in lineages_dict: + current_lineage_dict: dict = lineages_dict[row_dict["fcs_gx_top_taxid"]] + for i in range(1, 10): + out_line.append(current_lineage_dict[rankedlineage_col_names[i]]) + else: + for _ in range(1, 10): + out_line.append("") + print(*out_line, sep=",") + + +if __name__ == "__main__": + parser = argparse.ArgumentParser(description=__doc__) + parser.add_argument( + "taxonomy_report", + type=str, + help="Path to directory with FCS-GX output files *.taxonomy.rpt", + ) + parser.add_argument( + "fcs_report", + type=str, + help="Path to directory with FCS-GX output files *.fcs_gx_report.txt", + ) + parser.add_argument( + "ncbi_rankedlineage_path", + type=str, + help="Path to the rankedlineage.dmp of NCBI taxonomy", + ) + parser.add_argument("--version", action="version", version="1.1.0") + args = parser.parse_args() + + main(args.taxonomy_report, args.fcs_report, args.ncbi_rankedlineage_path) diff --git a/modules/sanger-tol/fcsgx/parseresults/tests/main.nf.test b/modules/sanger-tol/fcsgx/parseresults/tests/main.nf.test new file mode 100644 index 00000000..38186b90 --- /dev/null +++ b/modules/sanger-tol/fcsgx/parseresults/tests/main.nf.test @@ -0,0 +1,101 @@ +nextflow_process { + + name "Test Process FCSGX_PARSERESULTS" + script "../main.nf" + process "FCSGX_PARSERESULTS" + + tag "modules" + tag "modules_sangertol" + tag "fcsgx" + tag "fcsgx/parseresults" + + setup { + println "\nDownloading the new_taxdump..." + def new_taxdump_dir = new File("${launchDir}/new_taxdump/") + new_taxdump_dir.mkdirs() + + def new_taxdump_url = "https://ftp.ncbi.nlm.nih.gov/pub/taxonomy/new_taxdump/new_taxdump.tar.gz" + + def command = ['bash', '-c', "curl ${new_taxdump_url} -o ${launchDir}/new_taxdump.tar.gz"] + def process = command.execute() + process.waitFor() + + def command2 = ['bash', '-c', "gunzip -c ${launchDir}/new_taxdump.tar.gz | tar xf - -C ${new_taxdump_dir}/"] + def process2 = command2.execute() + process2.waitFor() + + if (process2.exitValue() != 0) { + throw new RuntimeException("Error - failed to download new_taxdump: ${process2.err.text}") + } + } + + test("Laetiporus_sulphureus - csv") { + + config "../nextflow.config" + + when { + process { + """ + input[0] = channel.of( + [ + [ id: "Laetiporus_sulphureus" ], + file(params.modules_testdata_base_path + "Laetiporus_sulphureus/analysis/gfLaeSulp1.1/fcsgx/gfLaeSulp1.1_PRIMARY.taxonomy.rpt", checkIfExists: true) + ] + ) + input[1] = channel.of( + [ + [ id: "Laetiporus_sulphureus" ], + file(params.modules_testdata_base_path + "Laetiporus_sulphureus/analysis/gfLaeSulp1.1/fcsgx/gfLaeSulp1.1_PRIMARY.fcs_gx_report.txt", checkIfExists: true) + ] + ) + input[2] = channel.of( + file("${launchDir}/new_taxdump/rankedlineage.dmp", checkIfExists: true) + ) + """ + } + } + + then { + assert process.success + assertAll( + { assert snapshot(process.out).match() } + ) + } + } + + test("Laetiporus_sulphureus - fasta - stub") { + + options "-stub" + + config "../nextflow.config" + + when { + process { + """ + input[0] = channel.of( + [ + [ id: "Laetiporus_sulphureus" ], + file(params.modules_testdata_base_path + "Laetiporus_sulphureus/analysis/gfLaeSulp1.1/fcsgx/gfLaeSulp1.1_PRIMARY.taxonomy.rpt", checkIfExists: true) + ] + ) + input[1] = channel.of( + [ + [ id: "Laetiporus_sulphureus" ], + file(params.modules_testdata_base_path + "Laetiporus_sulphureus/analysis/gfLaeSulp1.1/fcsgx/gfLaeSulp1.1_PRIMARY.fcs_gx_report.txt", checkIfExists: true) + ] + ) + input[2] = channel.of( + file("${launchDir}/new_taxdump/rankedlineage.dmp", checkIfExists: true) + ) + """ + } + } + + then { + assert process.success + assertAll( + { assert snapshot(process.out).match() } + ) + } + } +} diff --git a/modules/sanger-tol/fcsgx/parseresults/tests/main.nf.test.snap b/modules/sanger-tol/fcsgx/parseresults/tests/main.nf.test.snap new file mode 100644 index 00000000..c9ec4fa1 --- /dev/null +++ b/modules/sanger-tol/fcsgx/parseresults/tests/main.nf.test.snap @@ -0,0 +1,112 @@ +{ + "Laetiporus_sulphureus - fasta - stub": { + "content": [ + { + "0": [ + [ + { + "id": "Laetiporus_sulphureus" + }, + "Laetiporus_sulphureus_parsed_fcsgx.csv:md5,d41d8cd98f00b204e9800998ecf8427e" + ] + ], + "1": [ + [ + "FCSGX_PARSERESULTS", + "python", + "3.14.3" + ] + ], + "2": [ + [ + "FCSGX_PARSERESULTS", + "parse_fcsgx_result", + "1.1.0" + ] + ], + "fcsgxresult": [ + [ + { + "id": "Laetiporus_sulphureus" + }, + "Laetiporus_sulphureus_parsed_fcsgx.csv:md5,d41d8cd98f00b204e9800998ecf8427e" + ] + ], + "versions_parsefcsgx": [ + [ + "FCSGX_PARSERESULTS", + "parse_fcsgx_result", + "1.1.0" + ] + ], + "versions_python": [ + [ + "FCSGX_PARSERESULTS", + "python", + "3.14.3" + ] + ] + } + ], + "timestamp": "2026-07-09T14:48:37.325284408", + "meta": { + "nf-test": "0.9.5", + "nextflow": "25.10.4" + } + }, + "Laetiporus_sulphureus - csv": { + "content": [ + { + "0": [ + [ + { + "id": "Laetiporus_sulphureus" + }, + "Laetiporus_sulphureus_parsed_fcsgx.csv:md5,2e0a628cce100843cc36d184a1d92653" + ] + ], + "1": [ + [ + "FCSGX_PARSERESULTS", + "python", + "3.14.3" + ] + ], + "2": [ + [ + "FCSGX_PARSERESULTS", + "parse_fcsgx_result", + "1.1.0" + ] + ], + "fcsgxresult": [ + [ + { + "id": "Laetiporus_sulphureus" + }, + "Laetiporus_sulphureus_parsed_fcsgx.csv:md5,2e0a628cce100843cc36d184a1d92653" + ] + ], + "versions_parsefcsgx": [ + [ + "FCSGX_PARSERESULTS", + "parse_fcsgx_result", + "1.1.0" + ] + ], + "versions_python": [ + [ + "FCSGX_PARSERESULTS", + "python", + "3.14.3" + ] + ] + } + ], + "timestamp": "2026-07-09T14:48:17.817628452", + "meta": { + "nf-test": "0.9.5", + "nextflow": "25.10.4" + } + } +} \ No newline at end of file diff --git a/modules/nf-core/fcsgx/rungx/environment.yml b/modules/sanger-tol/fcsgx/rungx/environment.yml similarity index 100% rename from modules/nf-core/fcsgx/rungx/environment.yml rename to modules/sanger-tol/fcsgx/rungx/environment.yml diff --git a/modules/nf-core/fcsgx/rungx/main.nf b/modules/sanger-tol/fcsgx/rungx/main.nf similarity index 52% rename from modules/nf-core/fcsgx/rungx/main.nf rename to modules/sanger-tol/fcsgx/rungx/main.nf index 29fd415b..c5ad8420 100644 --- a/modules/nf-core/fcsgx/rungx/main.nf +++ b/modules/sanger-tol/fcsgx/rungx/main.nf @@ -8,9 +8,9 @@ process FCSGX_RUNGX { 'quay.io/biocontainers/ncbi-fcs-gx:0.5.5--h9948957_0' }" input: - tuple val(meta), path(fasta) + tuple val(meta), val(taxid), path(fasta) path gxdb - path ramdisk_path + val ramdisk_path val production_mode output: @@ -27,52 +27,35 @@ process FCSGX_RUNGX { def args = task.ext.args ?: '' def prefix = task.ext.prefix ?: "${meta.id}" - // At Sanger we have a permenant home for the DB on NVME storage - // def mv_database_to_ram = ramdisk_path ? "rclone copy $gxdb $ramdisk_path/$task.index/" : '' - // def database = ramdisk_path ? "$ramdisk_path/$task.index/" : gxdb // Use task.index to make memory location unique def database = ramdisk_path ?: gxdb - if ( production_mode ) { + ( ramdisk_path ? """ - echo "Using Production FCSGX with local installation" + # Copy DB to RAM-disk when supplied. Otherwise, rungx is very slow. + rclone copy --ignore-existing ${gxdb} ${database} + """ : "" + ) - export GX_NUM_CORES=${task.cpus} - export GX_INSTANTIATE_FASTA=1 - - run_gx \\ - --fasta ${fasta} \\ - --gx-db ${database} \\ - --tax-id ${meta.taxid} \\ - --generate-logfile true \\ - --out-basename ${prefix} \\ - --out-dir . \\ - ${args} - - """ - } else { - """ - echo "Using Standard FCSGX with container" - - run_gx.py \\ - --fasta ${fasta} \\ - --gx-db ${database} \\ - --tax-id ${meta.taxid} \\ - --generate-logfile true \\ - --out-basename ${prefix} \\ - --out-dir . \\ - ${args} - - """ - } + """ + export GX_NUM_CORES=${task.cpus} + export GX_INSTANTIATE_FASTA=1 + + run_gx.py \\ + --fasta ${fasta} \\ + --gx-db ${database} \\ + --tax-id ${taxid} \\ + --generate-logfile true \\ + --out-basename ${prefix} \\ + --out-dir . \\ + ${args} + """ stub: - // def args = task.ext.args ?: '' def prefix = task.ext.prefix ?: "${meta.id}" """ touch ${prefix}.fcs_gx_report.txt touch ${prefix}.taxonomy.rpt touch ${prefix}.summary.txt echo "" | gzip > ${prefix}.hits.tsv.gz - """ } diff --git a/modules/nf-core/fcsgx/rungx/meta.yml b/modules/sanger-tol/fcsgx/rungx/meta.yml similarity index 97% rename from modules/nf-core/fcsgx/rungx/meta.yml rename to modules/sanger-tol/fcsgx/rungx/meta.yml index b7a23508..4530c450 100644 --- a/modules/nf-core/fcsgx/rungx/meta.yml +++ b/modules/sanger-tol/fcsgx/rungx/meta.yml @@ -39,6 +39,9 @@ input: - ramdisk_path: type: string description: Path to RAM disk for improved performance (optional) + - production_mode: + type: boolean + description: Whether to run in production mode (optional) output: fcsgx_report: - - meta: diff --git a/modules/nf-core/fcsgx/rungx/tests/main.nf.test b/modules/sanger-tol/fcsgx/rungx/tests/main.nf.test similarity index 76% rename from modules/nf-core/fcsgx/rungx/tests/main.nf.test rename to modules/sanger-tol/fcsgx/rungx/tests/main.nf.test index 9224b201..b5b4c0b9 100644 --- a/modules/nf-core/fcsgx/rungx/tests/main.nf.test +++ b/modules/sanger-tol/fcsgx/rungx/tests/main.nf.test @@ -5,14 +5,19 @@ nextflow_process { process "FCSGX_RUNGX" tag "modules" - tag "modules_nfcore" + tag "modules_sangertol" tag "fcsgx" - tag "fcsgx/fetchdb" + tag "nf-core/fcsgx/fetchdb" + tag "nf-core/fetchdb" tag "fcsgx/rungx" setup { + nfcoreInitialise("${launchDir}/library/") + nfcoreInstall("${launchDir}/library/", ["fcsgx/fetchdb"]) + nfcoreLink("${launchDir}/library/", "${baseDir}/modules") + run("FCSGX_FETCHDB"){ - script "../../fetchdb/main.nf" + script "../../../../nf-core/fcsgx/fetchdb/main.nf" process { """ input[0] = file('https://ftp.ncbi.nlm.nih.gov/genomes/TOOLS/FCS/database/test-only/test-only.manifest', checkIfExists: true) @@ -29,10 +34,11 @@ nextflow_process { input[0] = [ [ id:'test', single_end:false ], // meta map '2697049', // taxid for SARS-CoV-2 - file(params.modules_testdata_base_path + 'genomics/sarscov2/genome/genome.fasta', checkIfExists: true), + file(params.nfcore_modules_testdata_base_path + 'genomics/sarscov2/genome/genome.fasta', checkIfExists: true), ] input[1] = FCSGX_FETCHDB.out.database input[2] = [] + input[3] = channel.of(false) """ } } @@ -65,10 +71,11 @@ nextflow_process { input[0] = [ [ id:'test', single_end:false ], // meta map '2697049', // taxid for SARS-CoV-2 - file(params.modules_testdata_base_path + 'genomics/sarscov2/genome/genome.fasta', checkIfExists: true), + file(params.nfcore_modules_testdata_base_path + 'genomics/sarscov2/genome/genome.fasta', checkIfExists: true), ] input[1] = FCSGX_FETCHDB.out.database input[2] = [] + input[3] = channel.of(false) """ } } diff --git a/modules/nf-core/fcsgx/rungx/tests/main.nf.test.snap b/modules/sanger-tol/fcsgx/rungx/tests/main.nf.test.snap similarity index 100% rename from modules/nf-core/fcsgx/rungx/tests/main.nf.test.snap rename to modules/sanger-tol/fcsgx/rungx/tests/main.nf.test.snap diff --git a/nextflow.config b/nextflow.config index 250d4125..97847bcd 100644 --- a/nextflow.config +++ b/nextflow.config @@ -108,6 +108,15 @@ process { includeConfig 'conf/base.config' profiles { + sanger { + process { + withName: 'FCSGX_RUNGX' { + module = 'fcs-gx/farm/0.5.5' + container = null + conda = null + } + } + } debug { dumpHashes = true process.beforeScript = 'echo $HOSTNAME' diff --git a/subworkflows/local/generate_html_report/main.nf b/subworkflows/local/generate_html_report/main.nf index 2bf04477..fc89fc74 100644 --- a/subworkflows/local/generate_html_report/main.nf +++ b/subworkflows/local/generate_html_report/main.nf @@ -45,7 +45,11 @@ workflow GENERATE_HTML_REPORT_WORKFLOW { .join(fcs_adaptor .map { meta, files -> // LOGIC: SORT INTO PREDICTABLE ORDER - def sorted_files = files.sort { file -> + // NOTE: `sort(false)` returns a NEW list. The default `sort{}` + // mutates in place, and the list emitted by `groupTuple()` + // is the SAME instance shared with every other consumer of + // this channel -> `ConcurrentModificationException`. + def sorted_files = (files ?: []).sort(false) { file -> file.toString().contains('_euk') ? 0 : file.toString().contains('_prok') ? 1 : 2 } diff --git a/subworkflows/local/run_fcsadaptor/main.nf b/subworkflows/local/run_fcsadaptor/main.nf index d4f815cf..83abfb2a 100644 --- a/subworkflows/local/run_fcsadaptor/main.nf +++ b/subworkflows/local/run_fcsadaptor/main.nf @@ -36,7 +36,12 @@ workflow RUN_FCSADAPTOR { .map{ meta, file -> [meta.id, file] } ) .groupTuple() - .map { id, files -> [[id: id], files] } + // NOTE: `groupTuple()` emits its internal ArrayBag, and that SAME instance is + // handed to every consumer of this channel. Any in-place operation + // (sort/unique/add) by one consumer corrupts it for the others, which + // surfaces as `Unexpected error [ConcurrentModificationException]`. + // Emit a defensive copy so each downstream branch owns its own list. + .map { id, files -> [[id: id], new ArrayList(files)] } .set { ch_fcsadapt } emit: diff --git a/subworkflows/local/run_fcsgx/main.nf b/subworkflows/local/run_fcsgx/main.nf index 68cf9a07..07eeb74c 100644 --- a/subworkflows/local/run_fcsgx/main.nf +++ b/subworkflows/local/run_fcsgx/main.nf @@ -1,11 +1,8 @@ -// -// MODULE IMPORT BLOCK -// include { SAMTOOLS_DICT } from '../../../modules/nf-core/samtools/dict/main' -include { FCSGX_RUNGX } from '../../../modules/nf-core/fcsgx/rungx/main' -include { PARSE_FCSGX_RESULT } from '../../../modules/local/fcsgx/parse_results/main' +include { FCSGX_RUNGX } from '../../../modules/sanger-tol/fcsgx/rungx/main' +include { FCSGX_PARSERESULTS } from '../../../modules/sanger-tol/fcsgx/parseresults/main' -workflow RUN_FCSGX { +workflow FCSGX_PARSECSV { take: reference // channel [ val(meta), path(file) ] @@ -13,10 +10,11 @@ workflow RUN_FCSGX { ncbi_rankedlineage_path // channel path(file) main: - ch_versions = channel.empty() // - // MODULE: Use SAMTOOLS_DICT to get origin file and md5sum of each sequence + // MODULE: USE SAMTOOLS_DICT TO GET THE ORIGIN FILE OF EACH SEQUENCE + // ITS NOT NEEDED BY ANYTHING IN THE SUBWORKFLOW BUT IS A NICE + // VERIFICATION OF SEQ ORIGIN // SAMTOOLS_DICT( reference @@ -24,17 +22,10 @@ workflow RUN_FCSGX { // - // MODULE: FCSGX_RUNGX RUN ON ASSEMBLY FASTA TUPLE WITH THE TAXID AGAINST THE FCSGXDB - // PRODUCTION_MODE WILL CHANGE HOW THE FCSGX MODULE IS RUN E.G IT IS SPECIFIC FOR `module` + // MODULE: RUN FCSGX FOR CLASSIFICATION OF SEQUENCES IN ASSEMBLY // - SAMTOOLS_DICT.out.dict - .map { meta, ref, _dict -> - [meta, ref] - } - .set { samtools_reference } - FCSGX_RUNGX ( - samtools_reference, + reference.map { meta, ref -> tuple(meta, meta.taxid, ref) }, fcsgxpath, [], "production" in workflow.profile.tokenize(',') @@ -42,44 +33,32 @@ workflow RUN_FCSGX { fcsgx_report_txt = FCSGX_RUNGX.out.fcsgx_report .map { meta, file -> - file ? [[ id: meta.id ], file] : [[:], []] + file ? tuple([ id: meta.id ], file) : [[:], []] } fcsgx_taxonomy_rpt = FCSGX_RUNGX.out.taxonomy_report .map { meta, file -> - file ? [[ id: meta.id ], file] : [[:], []] + file ? tuple([ id: meta.id ], file) : [[:], []] } - // - // MODULE: CREATE INPUT CHANNEL FOR PARSING RESULT MODULE - // - fcsgx_report_txt - .map{ meta, file -> - [meta, file.getParent()] - } - .set { report_path } - // - // MODULE: PARSE_FCSGX_RESULT to parse the FCSGX_RUNGX result output in csv format. + // MODULE: CONVER FCSGX_RUNGX RESULTS INTO A SINGLE CSV FILE // - PARSE_FCSGX_RESULT ( - report_path, + FCSGX_PARSERESULTS ( + fcsgx_taxonomy_rpt, + fcsgx_report_txt, ncbi_rankedlineage_path ) - ch_versions = ch_versions.mix( PARSE_FCSGX_RESULT.out.versions ) - fcsgxresult = PARSE_FCSGX_RESULT.out.fcsgxresult + fcsgxresult = FCSGX_PARSERESULTS.out.fcsgxresult .map { meta, file -> - file ? [[id: meta.id], file] : [[:], []] + file ? tuple(meta, file) : [[:], []] } emit: - fcsgxresult - genomedict = samtools_reference + genomedict = SAMTOOLS_DICT.out.dict fcsgx_report_txt fcsgx_taxonomy_rpt - versions = ch_versions - } diff --git a/tests/default.nf.test b/tests/default.nf.test index 72c4e7cb..06d74409 100644 --- a/tests/default.nf.test +++ b/tests/default.nf.test @@ -243,7 +243,7 @@ nextflow_pipeline { // sizes have been written to show that there are X numbers across the Y relevant input assemblies // ascc_main_output // file output has minute changes per run - ascc_main_output.size() == 2, + ascc_main_output.size() == 3, ascc_decontaminated_fa.size() == 3, // kmer_data_files, // ML method outputs, changes per run, included in stable name diff --git a/tests/default.nf.test.snap b/tests/default.nf.test.snap index 3588e2c2..7b2b11d2 100644 --- a/tests/default.nf.test.snap +++ b/tests/default.nf.test.snap @@ -1,7 +1,7 @@ { "-profile test": { "content": [ - 139, + 140, { "ASCC_MERGE_TABLES": { "python": "3.11.0", @@ -330,6 +330,7 @@ "asccTinyTest_V2_MITO/KRAKEN2/asccTinyTest_V2_MITO_nt_kraken_lineage_file.txt", "asccTinyTest_V2_MITO/KRAKEN2/versions.yml", "asccTinyTest_V2_MITO/ascc_main_output", + "asccTinyTest_V2_MITO/ascc_main_output/asccTinyTest_V2_MITO_contamination_check_merged_table.csv", "asccTinyTest_V2_MITO/ascc_main_output/asccTinyTest_V2_MITO_filtered.fasta", "asccTinyTest_V2_MITO/ascc_main_output/asccTinyTest_V2_MITO_filtered.fasta.abnormal_details.txt", "asccTinyTest_V2_MITO/ascc_main_output/asccTinyTest_V2_MITO_filtered.fasta.contamination.bed", diff --git a/workflows/ascc.nf b/workflows/ascc.nf index e380e2a8..6466a7be 100644 --- a/workflows/ascc.nf +++ b/workflows/ascc.nf @@ -3,8 +3,7 @@ IMPORT MODULES / SUBWORKFLOWS / FUNCTIONS ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ */ -include { ASCC_GENOMIC as GENOMIC } from './ascc_genomic' -include { ASCC_ORGANELLAR as ORGANELLAR } from './ascc_organellar' +include { ASCC_ASSEMBLY as ASSEMBLY } from './ascc_assembly' include { softwareVersionsToYAML } from '../subworkflows/nf-core/utils_nfcore_pipeline' include { methodsDescriptionText } from '../subworkflows/local/utils_nfcore_ascc_pipeline' @@ -48,13 +47,15 @@ workflow ASCC { ch_versions = channel.empty() // - // WORKFLOW: Run main workflow for GENOMIC samples + // WORKFLOW: Run main workflow for all assemblies (genomic + organellar) // - GENOMIC ( - genomic_genomes, + ASSEMBLY ( + params.genomic_only + ? genomic_genomes + : genomic_genomes.mix(organellar_genomes), organellar_genomes, fcs_override, - genomic_fcs_samplesheet, + genomic_fcs_samplesheet.mix(organellar_fcs_samplesheet), fcs_db, collected_reads, scientific_name, @@ -75,36 +76,7 @@ workflow ASCC { barcodes, val_reads_per_chunk ) - ch_versions = ch_versions.mix(GENOMIC.out.versions) - - - // - // WORKFLOW: Run main workflow for ORGANELLAR samples - // - if ( !params.genomic_only ) { - ORGANELLAR ( - organellar_genomes, - fcs_override, - organellar_fcs_samplesheet, - fcs_db, - collected_reads, - scientific_name, - pacbio_database, - ncbi_taxonomy_path, - ncbi_ranked_lineage_path, - nt_database_path, - diamond_nr_db_path, - diamond_uniprot_db_path, - taxid, - nt_kraken_db_path, - vecscreen_database_path, - reads_path, - reads_type, - barcodes, - val_reads_per_chunk - ) - ch_versions = ch_versions.mix(ORGANELLAR.out.versions) - } + ch_versions = ch_versions.mix(ASSEMBLY.out.versions) // diff --git a/workflows/ascc_genomic.nf b/workflows/ascc_assembly.nf similarity index 69% rename from workflows/ascc_genomic.nf rename to workflows/ascc_assembly.nf index b133d86f..7dff4a1f 100644 --- a/workflows/ascc_genomic.nf +++ b/workflows/ascc_assembly.nf @@ -22,7 +22,7 @@ include { PACBIO_BARCODE_CHECK } from '../subworkflows/ include { RUN_READ_COVERAGE } from '../subworkflows/local/run_read_coverage/main' include { RUN_VECSCREEN } from '../subworkflows/local/run_vecscreen/main' include { RUN_NT_KRAKEN } from '../subworkflows/local/run_nt_kraken/main' -include { RUN_FCSGX } from '../subworkflows/local/run_fcsgx/main' +include { FCSGX_PARSECSV as RUN_FCSGX } from '../subworkflows/local/run_fcsgx/main' include { RUN_FCSADAPTOR } from '../subworkflows/local/run_fcsadaptor/main' include { RUN_DIAMOND as NR_DIAMOND } from '../subworkflows/local/run_diamond/main' include { RUN_DIAMOND as UP_DIAMOND } from '../subworkflows/local/run_diamond/main' @@ -34,20 +34,19 @@ include { GENERATE_HTML_REPORT_WORKFLOW } from '../subworkflows/ include { getEmptyPlaceholder } from '../functions/local/ascc_utils' - /* ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ RUN MAIN WORKFLOW ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ */ -workflow ASCC_GENOMIC { +workflow ASCC_ASSEMBLY { take: - ch_samplesheet // channel: samplesheet read in from --input - organellar_genomes // channel: tuple(meta, reference) - _fcs_ov // params.fcs_override - fcs_samplesheet // The FCS override samplesheet for override + ch_samplesheet // channel: combined genomic + organellar assemblies; meta.assembly_type identifies each + organellar_genomes // channel: tuple(meta, reference) – organellar only, for ORGANELLAR_BLAST against genomic + _fcs_ov // params.fcs_override + fcs_samplesheet // The FCS override samplesheet (combined genomic + organellar entries) fcs_db // [path(path)] _reads scientific_name // val(name) @@ -72,9 +71,13 @@ workflow ASCC_GENOMIC { ch_versions = channel.empty() // - // LOGIC: CREATE run_conditional LIST + // LOGIC: HELPER CLOSURE AND CONDITIONAL LISTS + // isOrganellar resolves the correct run_conditionals list per assembly item, + // preserving the original per-workflow gating logic in a single unified workflow. // - run_conditionals = ["both", "genomic"] + def genomicConditionals = ["both", "genomic"] + def organellarConditionals = ["both", "organellar"] + def isOrganellar = { meta -> meta.assembly_type in ["MITO", "PLASTID"] } // @@ -82,10 +85,21 @@ workflow ASCC_GENOMIC { // ch_samplesheet .map { meta, sample -> - log.info "[ASCC INFO]: GENOMIC WORKFLOW:\n\t-- $meta\n\t-- $sample\n" + def type = isOrganellar(meta) ? "ORGANELLAR" : "GENOMIC" + log.info "[ASCC INFO]: ${type} WORKFLOW:\n\t-- $meta\n\t-- $sample\n" } + // + // LOGIC: BUILD ASSEMBLY TYPE LOOKUP FROM THE ORIGINAL SAMPLESHEET BEFORE ESSENTIAL_JOBS. + // ESSENTIAL_JOBS rebuilds meta internally with only id/sliding/window/taxid, + // stripping assembly_type. This lookup is joined back onto the output channels + // so that all isOrganellar() checks and the type branch works correctly. + // + ch_assembly_type_lookup = ch_samplesheet + .map { meta, _file -> [[id: meta.id], meta.assembly_type] } + + //------------------------------------------------------------------------- // // SUBWORKFLOW: RUNS FILTER_FASTA, GENERATE .GENOME, CALCS GC_CONTENT AND FINDS RUNS OF N's @@ -95,7 +109,6 @@ workflow ASCC_GENOMIC { ch_samplesheet ) ch_versions = ch_versions.mix(ESSENTIAL_JOBS.out.versions) - ej_reference_tuple = ESSENTIAL_JOBS.out.reference_tuple_from_GG ej_dot_genome = ESSENTIAL_JOBS.out.dot_genome ej_gc_coverage = ESSENTIAL_JOBS.out.gc_content_txt ej_trailing_ns = ESSENTIAL_JOBS.out.trailing_ns_report @@ -103,21 +116,55 @@ workflow ASCC_GENOMIC { ej_fasta_filter_log = ESSENTIAL_JOBS.out.filter_fasta_length_filtering_log + // + // LOGIC: RESTORE assembly_type TO ESSENTIAL_JOBS OUTPUT CHANNELS. + // ESSENTIAL_JOBS strips assembly_type when it rebuilds meta to inject seqkit + // sliding/window params. Rejoin against ch_assembly_type_lookup so that + // isOrganellar() and the type branch below correctly route each assembly. + // + ej_reference_tuple = ESSENTIAL_JOBS.out.reference_tuple_from_GG + .map { meta, file -> [[id: meta.id], meta, file] } + .join(ch_assembly_type_lookup) + .map { _id_meta, meta, file, assembly_type -> + [meta + [assembly_type: assembly_type], file] + } + + ej_seqkit_reference = ESSENTIAL_JOBS.out.reference_with_seqkit + .map { meta, file -> [[id: meta.id], meta, file] } + .join(ch_assembly_type_lookup) + .map { _id_meta, meta, file, assembly_type -> + [meta + [assembly_type: assembly_type], file] + } + + //------------------------------------------------------------------------- // - // LOGIC: CONVERT THE CHANNEL I AN EPOCH COUNT FOR THE GET_KMER_PROFILE + // LOGIC: BRANCH THE REFERENCE CHANNEL BY ASSEMBLY TYPE FOR TYPE-SPECIFIC PROCESS GATING // ej_reference_tuple + .branch { meta, _f -> + genomic: !isOrganellar(meta) + organellar: true + } + .set { ch_type_branch } + + + //------------------------------------------------------------------------- + // + // LOGIC: CONVERT THE CHANNEL INTO AN EPOCH COUNT FOR GET_KMER_PROFILE (GENOMIC ONLY) + // + ch_type_branch.genomic .map { _meta, file -> file.countFasta() * 3 } - .set {autoencoder_epochs_count} + .set { autoencoder_epochs_count } + // - // SUBWORKFLOW: COUNT KMERS, THEN REDUCE DIMENSIONS USING SELECTED METHODS + // SUBWORKFLOW: COUNT KMERS, THEN REDUCE DIMENSIONS USING SELECTED METHODS (GENOMIC ONLY) // GET_KMERS_PROFILE ( - ej_reference_tuple.filter{ _meta, _file -> params.run_kmers in run_conditionals }, + ch_type_branch.genomic.filter{ _meta, _file -> params.run_kmers in genomicConditionals }, params.kmer_length, params.dimensionality_reduction_methods, autoencoder_epochs_count @@ -137,7 +184,12 @@ workflow ASCC_GENOMIC { // SUBWORKFLOW: EXTRACT RESULTS HITS FROM TIARA // TIARA_TIARA ( - ej_reference_tuple.filter{ _meta, _file -> params.run_tiara in run_conditionals } + ch_type_branch.genomic + .filter{ _meta, _file -> params.run_tiara in genomicConditionals } + .mix( + ch_type_branch.organellar + .filter{ _meta, _file -> params.run_tiara in organellarConditionals } + ) ) ch_versions = ch_versions.mix( TIARA_TIARA.out.versions ) ch_tiara = TIARA_TIARA.out.classifications @@ -149,11 +201,59 @@ workflow ASCC_GENOMIC { // // SUBWORKFLOW: EXTRACT RESULTS HITS FROM NT-BLAST // + // Genomic assemblies: use ej_reference_tuple directly (no length gate) + // Organellar assemblies: apply seqkit length filter first (valid_length_fasta) + // + + // + // LOGIC: FOR ORGANELLAR ASSEMBLIES, WE NEED TO MAKE SURE THAT THE INPUT SEQUENCE + // IS OF AT LEAST LENGTH OF params.seqkit_window BEFORE RUNNING BLAST/DIAMOND + // + valid_length_fasta = ej_seqkit_reference + .filter { meta, _f -> isOrganellar(meta) } + // + // NOTE: Here we are using the un-filtered genome, any filtering may (accidently) cause an empty fasta + // + .map{ meta, file -> + def total_length = 0 + file.eachLine { line -> + if (line && !line.startsWith('>')) { + total_length += line.length() + } + } + + def meta2 = [ + id: meta.id, + sliding: meta.sliding, + window: meta.window, + seq_count: total_length + ] + + [meta2, file] + } + .filter { meta, _file -> + meta.seq_count >= params.seqkit_window + } + + valid_length_fasta + .map{ meta, _file -> + log.info "[ASCC INFO]: Running BLAST (NT, DIAMOND, NR) on VALID ORGANELLE: \n\t-- ${meta.id}'s sequence ($meta.seq_count bases) is >= seqkit_window $params.seqkit_window\n" + } + + assemblies_to_blast = ch_type_branch.genomic + .filter{ _meta, _file -> params.run_nt_blast in genomicConditionals } + .mix( + valid_length_fasta + .filter{ _meta, _file -> params.run_nt_blast in organellarConditionals } + ) EXTRACT_NT_BLAST ( - ej_reference_tuple.filter{ _meta, _file -> - params.run_nt_blast in run_conditionals - }, + ch_type_branch.genomic + .filter{ _meta, _file -> params.run_nt_blast in genomicConditionals } + .mix( + valid_length_fasta + .filter{ _meta, _file -> params.run_nt_blast in organellarConditionals } + ), nt_database_path.first(), ncbi_ranked_lineage_path.first() ) @@ -169,9 +269,12 @@ workflow ASCC_GENOMIC { // NR_DIAMOND ( - ej_reference_tuple.filter{ _meta, _file -> - params.run_nr_diamond in run_conditionals - }, + ch_type_branch.genomic + .filter{ _meta, _file -> params.run_nr_diamond in genomicConditionals } + .mix( + valid_length_fasta + .filter{ _meta, _file -> params.run_nr_diamond in organellarConditionals } + ), diamond_nr_db_path.first() ) ch_versions = ch_versions.mix(NR_DIAMOND.out.versions) @@ -192,9 +295,12 @@ workflow ASCC_GENOMIC { // qseqid sseqid pident length mismatch gapopen qstart qend sstart send // evalue bitscore staxids sscinames sskingdoms sphylums salltitles UP_DIAMOND ( - ej_reference_tuple.filter{ _meta, _file -> - params.run_uniprot_diamond in run_conditionals - }, + ch_type_branch.genomic + .filter{ _meta, _file -> params.run_uniprot_diamond in genomicConditionals } + .mix( + valid_length_fasta + .filter{ _meta, _file -> params.run_uniprot_diamond in organellarConditionals } + ), diamond_uniprot_db_path.first() ) ch_versions = ch_versions.mix(UP_DIAMOND.out.versions) @@ -210,11 +316,11 @@ workflow ASCC_GENOMIC { //------------------------------------------------------------------------- // - // LOGIC: CHECK WHETHER THERE IS A MITO AND BRANCH + // LOGIC: CHECK WHETHER THERE IS A MITO AND BRANCH (GENOMIC ONLY) // organellar_check = organellar_genomes .filter{ _meta, _file -> - params.run_organellar_blast in run_conditionals + params.run_organellar_blast in genomicConditionals } .branch { meta, _assembly -> mito: meta.assembly_type == "MITO" @@ -224,20 +330,20 @@ workflow ASCC_GENOMIC { // - // SUBWORKFLOW: BLASTING FOR MITO ASSEMBLIES IN GENOME + // SUBWORKFLOW: BLASTING FOR MITO ASSEMBLIES IN GENOME (GENOMIC ONLY) // MITO_ORGANELLAR_BLAST ( - ej_reference_tuple, + ch_type_branch.genomic, organellar_check.mito ) ch_versions = ch_versions.mix(MITO_ORGANELLAR_BLAST.out.versions) // - // SUBWORKFLOW: BLASTING FOR PLASTID ASSEMBLIES IN GENOME + // SUBWORKFLOW: BLASTING FOR PLASTID ASSEMBLIES IN GENOME (GENOMIC ONLY) // PLASTID_ORGANELLAR_BLAST ( - ej_reference_tuple, + ch_type_branch.genomic, organellar_check.plastid ) ch_versions = ch_versions.mix(PLASTID_ORGANELLAR_BLAST.out.versions) @@ -278,8 +384,9 @@ workflow ASCC_GENOMIC { .set { duplicated_db } PACBIO_BARCODE_CHECK ( - duplicated_db.reference.filter{ _meta, _file -> - params.run_pacbio_barcodes in run_conditionals + duplicated_db.reference.filter{ meta, _file -> + def conds = isOrganellar(meta) ? organellarConditionals : genomicConditionals + params.run_pacbio_barcodes in conds }, ch_barcodes, duplicated_db.pacbio_db @@ -293,29 +400,35 @@ workflow ASCC_GENOMIC { // SUBWORKFLOW: RUN FCS-ADAPTOR TO IDENTIDY ADAPTOR AND VECTORR CONTAMINATION // RUN_FCSADAPTOR ( - ej_reference_tuple.filter{ _meta, _file -> - params.run_fcs_adaptor in run_conditionals + ej_reference_tuple.filter{ meta, _file -> + def conds = isOrganellar(meta) ? organellarConditionals : genomicConditionals + params.run_fcs_adaptor in conds } ) - ch_versions = ch_versions.mix(RUN_FCSADAPTOR.out.versions) - ch_fcsadapt = RUN_FCSADAPTOR.out.ch_joint_report.ifEmpty{ [[:],[]] } + ch_versions = ch_versions.mix(RUN_FCSADAPTOR.out.versions) + ch_fcsadapt = RUN_FCSADAPTOR.out.ch_joint_report.ifEmpty{ [[:],[]] } //------------------------------------------------------------------------- // // SUBWORKFLOW: RUN FCS-GX TO IDENTIFY CONTAMINATION IN THE ASSEMBLY // - if ( params.run_fcsgx in run_conditionals && !params.fcs_override ) { + + if ( params.run_fcsgx != "off" && !params.fcs_override ) { joint_channel = ej_reference_tuple + .filter { meta, _f -> + def conds = isOrganellar(meta) ? organellarConditionals : genomicConditionals + params.run_fcsgx in conds + } .combine(fcs_db) .combine(taxid) .combine(ncbi_ranked_lineage_path) .multiMap { meta, ref, db, _tax_id, tax_path -> def new_meta = [id: meta.id, taxid: meta.taxid] - reference: [new_meta, ref] - fcs_db_path: db - ncbi_tax_path: tax_path + reference: [new_meta, ref] + fcs_db_path: db + ncbi_tax_path: tax_path } RUN_FCSGX ( @@ -323,7 +436,6 @@ workflow ASCC_GENOMIC { joint_channel.fcs_db_path, joint_channel.ncbi_tax_path ) - ch_versions = ch_versions.mix(RUN_FCSGX.out.versions) ch_fcsgx = RUN_FCSGX.out.fcsgxresult ch_fcsgx_report = RUN_FCSGX.out.fcsgx_report_txt @@ -353,8 +465,9 @@ workflow ASCC_GENOMIC { // SUBWORKFLOW: CALCULATE AVERAGE READ COVERAGE // RUN_READ_COVERAGE ( - ej_reference_tuple.filter{ _meta, _file -> - params.run_coverage in run_conditionals + ej_reference_tuple.filter{ meta, _file -> + def conds = isOrganellar(meta) ? organellarConditionals : genomicConditionals + params.run_coverage in conds }, reads_path, reads_type, //Subworkflow uses the param, not this value... as soon as it's in a channel it can't be used for a comparator. @@ -370,8 +483,9 @@ workflow ASCC_GENOMIC { // SUBWORKFLOW: SCREENING FOR VECTOR SEQUENCE // RUN_VECSCREEN ( - ej_reference_tuple.filter{ _meta, _file -> - params.run_vecscreen in run_conditionals + ej_reference_tuple.filter{ meta, _file -> + def conds = isOrganellar(meta) ? organellarConditionals : genomicConditionals + params.run_vecscreen in conds }, vecscreen_database_path.first() ) @@ -384,8 +498,9 @@ workflow ASCC_GENOMIC { // SUBWORKFLOW: RUN THE KRAKEN CLASSIFIER // RUN_NT_KRAKEN( - ej_reference_tuple.filter{ _meta, _file -> - params.run_kraken in run_conditionals + ej_reference_tuple.filter{ meta, _file -> + def conds = isOrganellar(meta) ? organellarConditionals : genomicConditionals + params.run_kraken in conds }, nt_kraken_db_path.first(), ncbi_ranked_lineage_path.first() @@ -397,28 +512,47 @@ workflow ASCC_GENOMIC { //------------------------------------------------------------------------- - if ( params.run_create_btk_dataset in run_conditionals ) { + if ( params.run_create_btk_dataset != "off" ) { + + // + // LOGIC: FILTER fcsgx TO GENOMIC ITEMS ONLY FOR BTK DATASET INPUT + // The organellar workflow intentionally excludes fcsgx from CREATE_BTK_DATASET; + // organellar items won't find a match here and will get null → placeholder via join remainder. + // + ch_fcsgx_for_btk = ch_fcsgx + .map { meta, f -> [[id: meta.id], f] } + .join( + ch_type_branch.genomic + .map { meta, _f -> [[id: meta.id], true] } + ) + .map { meta, f, _flag -> [meta, f] } // // LOGIC: FOUND RACE CONDITION EFFECTING LONG RUNNING JOBS // AND INPUT TO HERE ARE NOW MERGED AND MAPPED // EMPTY CHANNELS ARE CHECKED AND DEFAULTED TO [[:],[]] // + // ch_fcsgx_for_btk : genomic items only → organellar items → null → placeholder + // ch_kmers : genomic items only → organellar items → null → placeholder // ej_reference_tuple + .filter { meta, _f -> + def conds = isOrganellar(meta) ? organellarConditionals : genomicConditionals + params.run_create_btk_dataset in conds + } .map{meta, file -> [[id: meta.id], file]} - .join(ch_nt_blast, remainder: true) - .join(ch_tiara, remainder: true) - .join(ej_dot_genome,remainder: true) - .join(ch_fcsgx, remainder: true) - .join(ch_bam, remainder: true) - .join(ch_coverage, remainder: true) - .join(ch_kmers, remainder: true) - .join(ch_kraken1, remainder: true) - .join(ch_kraken2, remainder: true) - .join(ch_kraken3, remainder: true) - .join(nr_full, remainder: true) - .join(un_full, remainder: true) + .join(ch_nt_blast, remainder: true) + .join(ch_tiara, remainder: true) + .join(ej_dot_genome, remainder: true) + .join(ch_fcsgx_for_btk, remainder: true) + .join(ch_bam, remainder: true) + .join(ch_coverage, remainder: true) + .join(ch_kmers, remainder: true) + .join(ch_kraken1, remainder: true) + .join(ch_kraken2, remainder: true) + .join(ch_kraken3, remainder: true) + .join(nr_full, remainder: true) + .join(un_full, remainder: true) .filter { items -> def meta = items[0] meta != null && @@ -469,28 +603,29 @@ workflow ASCC_GENOMIC { // OR BY include_steps CONTAINING ALL AND EXCLUDE NOT CONTAINING autofilter_assembly. // if ( - ( params.run_tiara in run_conditionals ) && - ( params.run_fcsgx in run_conditionals ) && - ( params.run_autofilter_assembly in run_conditionals ) + params.run_tiara != "off" && + params.run_fcsgx != "off" && + params.run_autofilter_assembly != "off" ) { // - // LOGIC: FILTER THE INPUT FOR THE AUTOFILTER STEP - // - We can't just combine on meta.id as some of the channel. have other data + // LOGIC: FILTER THE INPUT FOR THE AUTOFILTER STEP PER ASSEMBLY TYPE + // - We can't just combine on meta.id as some of the Channels have other data // in there too so we just sanitise, and _then_ combine on 0, and // _then_ add back in the taxid as we need that for this process. // Thankfully taxid is a param so easy enough to add back in. // Actually, it just makes more sense to passs in as its own channel. // - ej_reference_tuple + .filter { meta, _f -> + def conds = isOrganellar(meta) ? organellarConditionals : genomicConditionals + params.run_tiara in conds && params.run_fcsgx in conds && params.run_autofilter_assembly in conds + } .map{ meta, file -> [[id: meta.id], file] } .combine( - ch_tiara.map{ meta, file -> [[id: meta.id], file] }, - by: 0 + ch_tiara.map{ meta, file -> [[id: meta.id], file] }, by: 0 ) .combine( - ch_fcsgx.map{ meta, file -> [[id: meta.id], file] }, - by: 0 + ch_fcsgx.map{ meta, file -> [[id: meta.id], file] }, by: 0 ) .combine( ncbi_ranked_lineage_path @@ -508,6 +643,7 @@ workflow ASCC_GENOMIC { } .set { autofilter_input_formatted } + // // MODULE: AUTOFILTER ASSEMBLY BY TIARA AND FCSGX RESULTS // @@ -517,14 +653,30 @@ workflow ASCC_GENOMIC { autofilter_input_formatted.fcs_file, autofilter_input_formatted.ncbi_rank ) - ch_autofilt_assem = AUTOFILTER_AND_CHECK_ASSEMBLY.out.decontaminated_assembly.map{_meta, file -> file} + ch_versions = ch_versions.mix(AUTOFILTER_AND_CHECK_ASSEMBLY.out.versions) + ch_autofilt_assem = AUTOFILTER_AND_CHECK_ASSEMBLY.out.decontaminated_assembly ch_autofilt_indicator = AUTOFILTER_AND_CHECK_ASSEMBLY.out.indicator_file + ch_autofilt_removed_seqs= AUTOFILTER_AND_CHECK_ASSEMBLY.out.removed_seqs + ch_autofilt_raw_report = AUTOFILTER_AND_CHECK_ASSEMBLY.out.raw_report + + ch_autofilt_alarm_file = AUTOFILTER_AND_CHECK_ASSEMBLY.out.alarm_file + .map{ meta, _file -> [[ id: meta.id ], _file] } + + ch_autofilt_fcs_tiara = AUTOFILTER_AND_CHECK_ASSEMBLY.out.fcs_tiara_summary + .map{ meta, _file -> [[ id: meta.id ], _file] } // - // LOGIC: BRANCH THE CHANNEL ON WHETHER OR NOT THERE IS ABNORMAL CONTAMINATION IN THE - // OUTPUT FILE. - // CHANGE OUTPUT NAME TO BE REFERENCE NAME AND THEN ALARM FILE + // LOGIC: BRANCH THE ALARM FILE ON ABNORMAL CONTAMINATION FOR THE GENOMIC BTK PIPELINE. + // Only genomic items feed the BTK conditional; organellar items are excluded + // by joining against ch_type_branch.genomic before branching. + // btk_bool = AUTOFILTER_AND_CHECK_ASSEMBLY.out.alarm_file + .map { meta, file -> [[id: meta.id], file] } + .join( + ch_type_branch.genomic + .map { meta, _f -> [[id: meta.id], true] } + ) + .map { meta, file, _flag -> [meta, file] } .map { meta, file -> [meta, file.text.trim()] } .branch { meta, data -> log.info("[ASCC INFO]: Run for ${meta.id} has:\n${data}\n") @@ -532,16 +684,8 @@ workflow ASCC_GENOMIC { run_btk : data.contains("YES_ABNORMAL_CONTAMINATION") dont_run : true // only other lines to be produced are "NO_ABNORMAL_CONTAMINATION" } - ch_versions = ch_versions.mix(AUTOFILTER_AND_CHECK_ASSEMBLY.out.versions) - btk_bool_run_btk = btk_bool.run_btk - ch_autofilt_removed_seqs= AUTOFILTER_AND_CHECK_ASSEMBLY.out.removed_seqs - ch_autofilt_raw_report = AUTOFILTER_AND_CHECK_ASSEMBLY.out.raw_report - - ch_autofilt_alarm_file = AUTOFILTER_AND_CHECK_ASSEMBLY.out.alarm_file - .map{ meta, file -> [[id: meta.id], file ] } - ch_autofilt_fcs_tiara = AUTOFILTER_AND_CHECK_ASSEMBLY.out.fcs_tiara_summary - .map{ meta, _file -> [[id: meta.id], _file] } + btk_bool_run_btk = btk_bool.run_btk } else { btk_bool_run_btk = channel.of([[id: "NA"], "false"]) @@ -556,13 +700,11 @@ workflow ASCC_GENOMIC { //------------------------------------------------------------------------- // - // LOGIC: DETERMINE WHETHER BLOBTOOLKIT SHOULD RUN BASED ON CONDITIONALS + // LOGIC: DETERMINE WHETHER BLOBTOOLKIT SHOULD RUN BASED ON CONDITIONALS (GENOMIC ONLY) // - ALWAYS RUN IF params.btk_busco_run_mode == "mandatory" AND BTK - run_btk_conditional = ej_reference_tuple - .map { meta, file -> - [[id: meta.id, taxid: meta.taxid], file] - } + run_btk_conditional = ch_type_branch.genomic + .map { meta, file -> [[id: meta.id], file] } // below is combined into the tuple to enforce the block to only run when channel is present. .combine ( btk_bool_run_btk .map{ meta, data -> @@ -571,7 +713,7 @@ workflow ASCC_GENOMIC { .replaceAll(/\s+/, "-") // Replace remaining spaces with "-" .replaceAll(/_+/, "_") // Keep underscores as they are .replaceAll(/-+/, "-") // Clean up multiple dashes - [[id: meta.id, taxid: meta.taxid], joined_content] + [[id: meta.id], joined_content] }, by: [0] ) @@ -593,7 +735,6 @@ workflow ASCC_GENOMIC { log.info "\t- You can verify here: $file" return [meta, file] } - //.set { skipped_btk_ch } if (params.run_autofilter_assembly == "off" && params.run_btk_busco != "off") { log.warn "[ASCC WARN]: run_autofilter_assembly is off, but run_btk_busco != off" @@ -601,8 +742,7 @@ workflow ASCC_GENOMIC { log.warn " `--btk_busco_run_mode mandatory`" } - // Noticed a race condition, this should fix that. - // + // NOTE: Noticed a race condition, this should fix that. run_btk_conditional.run_btk .map { meta, file, _data -> [meta.id, meta, file] } .join( @@ -632,6 +772,7 @@ workflow ASCC_GENOMIC { } .set { ch_meta_reads } + BLOBTOOLKIT_GENERATECSV ( ch_meta_reads, [[],[]], @@ -666,6 +807,7 @@ workflow ASCC_GENOMIC { } .combine(btk_samplesheet, by: 0) + combined_input .map{ meta, ref, samplesheet -> log.info("[ASCC INFO]: BTK will run for $meta\n\t| REF: ${ref}\n\t| SST: ${samplesheet}\n") @@ -695,37 +837,34 @@ workflow ASCC_GENOMIC { //------------------------------------------------------------------------- if ( - ( params.run_merge_datasets in run_conditionals ) && - ( params.run_btk_busco in run_conditionals ) - ) { - // - // MODULE: MERGE THE TWO BTK FORMATTED DATASETS INTO ONE DATASET FOR EASIER USE - // - - merged_channel = ch_create_btk_dataset - .map { meta, file -> [meta.id, meta, file] } - .combine( - SANGER_TOL_BTK.out.dataset - .map { meta, file -> - [meta.id, meta, file] - }, by: 0) - .map { _id, ref_meta, ref_file, _btk_meta, btk_file -> - [ref_meta, ref_file, btk_file] - } - - - MERGE_BTK_DATASETS ( - merged_channel + ( params.run_merge_datasets in genomicConditionals ) && + ( params.run_btk_busco in genomicConditionals ) + ) { + // + // MODULE: MERGE THE TWO BTK FORMATTED DATASETS INTO ONE DATASET FOR EASIER USE + // + merged_channel = ch_create_btk_dataset + .map { meta, file -> tuple(meta.id, [meta, file]) } + .join( + SANGER_TOL_BTK.out.dataset + .map { meta, file -> tuple(meta.id, [meta, file]) } ) - ch_versions = ch_versions.mix(MERGE_BTK_DATASETS.out.versions) - busco_merge_btk = MERGE_BTK_DATASETS.out.busco_summary_tsv - .map{ meta, _tsv -> [[id: meta.id], _tsv] } - merged_ds = MERGE_BTK_DATASETS.out.merged_datasets - } else { - busco_merge_btk = channel.of( [[:],[]] ) - merged_ds = channel.of( [[:],[]] ) + .map { + _id, ref, btk -> tuple(ref[0], ref[1], btk[1]) + } - } + MERGE_BTK_DATASETS ( + merged_channel + ) + ch_versions = ch_versions.mix(MERGE_BTK_DATASETS.out.versions) + busco_merge_btk = MERGE_BTK_DATASETS.out.busco_summary_tsv + .map{ meta, _tsv -> [[id: meta.id], _tsv] } + merged_ds = MERGE_BTK_DATASETS.out.merged_datasets + } else { + busco_merge_btk = channel.of( [[:],[]] ) + merged_ds = channel.of( [[:],[]] ) + + } //------------------------------------------------------------------------- @@ -733,15 +872,22 @@ workflow ASCC_GENOMIC { // LOGIC: EACH SUBWORKFLOW OUTPUTS EITHER AN EMPTY CHANNEL OR A FILE CHANNEL DEPENDING ON THE RUN RULES // SO THE RULES FOR THIS ONLY NEED TO BE A SIMPLE "DO YOU WANT IT OR NOT" // + // ch_kmers and busco_merge_btk are genomic-only; organellar items will find no match + // in the join and receive null → replaced by placeholder via getEmptyPlaceholder. + // if ( - ( params.run_essentials in run_conditionals ) && - ( params.run_merge_datasets in run_conditionals ) + params.run_essentials != "off" && + params.run_merge_datasets != "off" ) { // // LOGIC: JOIN CHANNELS INTO ONE BASED ON META.ID WHILST RETAINING EMPTY CHANNELS // ej_reference_tuple + .filter { meta, _f -> + def conds = isOrganellar(meta) ? organellarConditionals : genomicConditionals + params.run_essentials in conds && params.run_merge_datasets in conds + } .map{meta, file -> [[id: meta.id], file]} .join(ej_gc_coverage .map{meta, file -> [[id: meta.id], file]}, remainder: true) @@ -797,22 +943,24 @@ workflow ASCC_GENOMIC { // SUBWORKFLOW: GENERATE DECONTAMINATION FILES AND POTENTIALLY A DECONTAMINATED FASTA // THIS SHOULD ONLY RUN IF STANDARD CONDITIONALS ARE MET // AND ABNORMAL CONTAMINATION IS FOUND - // + // AUTOFILTERING THE ASSEMBLY IS ESSENTIAL FOR DECON TO RUN - // We only want the EUKARYOTIC report - // Not using the collection will result in a `Unexpected error [ConcurrentModificationException]` - // `ch_fcsadapt` because it is a mix channel, is technically still mutable + // NOTE: We only want the EUKARYOTIC report. euk_fcsadapt = ch_fcsadapt.map{ meta, files -> - def filesCopy = (files ?: []).collect() // defensive copy - [meta, filesCopy.find{ file -> file.name.endsWith('_euk.fcs_adaptor_report.txt') }] + [meta, (files ?: []).find{ file -> file.name.endsWith('_euk.fcs_adaptor_report.txt') }] } ej_reference_tuple_filtered = ej_reference_tuple - .filter{ meta, file -> - params.run_decontaminate_fasta in run_conditionals && params.run_autofilter_assembly in run_conditionals + .filter{ meta, _file -> + def conds = isOrganellar(meta) ? organellarConditionals : genomicConditionals + params.run_decontaminate_fasta in conds && params.run_autofilter_assembly in conds } .map{ meta, file -> [[id: meta.id], file] } + // + // ch_mito_full and ch_chloro_full are genomic-only; organellar items will find no match + // in the join inside RUN_DECONTAMINATE_FASTA and receive null → replaced by placeholder. + // RUN_DECONTAMINATE_FASTA( ej_reference_tuple_filtered, ch_fcsgx, @@ -831,6 +979,9 @@ workflow ASCC_GENOMIC { // SUBWORKFLOW: GENERATE HTML REPORT (minimal wiring, opt-in) // Gate with params.run_html_report to avoid altering default behavior. // + // ch_kmers_results is genomic-only; organellar items will find no match + // in the join inside GENERATE_HTML_REPORT_WORKFLOW and receive a placeholder. + // // Params file ch_params_file = params.params_file ? channel.fromPath(params.params_file) : channel.value([]) @@ -844,8 +995,9 @@ workflow ASCC_GENOMIC { merged_table, merged_phylum_count, ch_kmers_results, - ej_reference_tuple.filter{ _meta, _file -> - params.run_html_report in run_conditionals + ej_reference_tuple.filter{ meta, _file -> + def conds = isOrganellar(meta) ? organellarConditionals : genomicConditionals + params.run_html_report in conds }, ej_fasta_sanitation_log, ej_fasta_filter_log, diff --git a/workflows/ascc_organellar.nf b/workflows/ascc_organellar.nf deleted file mode 100644 index 7e55e3ad..00000000 --- a/workflows/ascc_organellar.nf +++ /dev/null @@ -1,654 +0,0 @@ -/* -~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ - IMPORT MODULES / SUBWORKFLOWS / FUNCTIONS -~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ -*/ - -include { CREATE_BTK_DATASET } from '../modules/local/blobtoolkit/create_dataset/main' -include { AUTOFILTER_AND_CHECK_ASSEMBLY } from '../modules/local/autofilter/autofilter/main' - -include { TIARA_TIARA } from '../modules/nf-core/tiara/tiara/main' - -include { ESSENTIAL_JOBS } from '../subworkflows/local/essential_jobs/main' -include { EXTRACT_NT_BLAST } from '../subworkflows/local/extract_nt_blast/main' -include { PACBIO_BARCODE_CHECK } from '../subworkflows/local/pacbio_barcode_check/main' -include { RUN_READ_COVERAGE } from '../subworkflows/local/run_read_coverage/main' -include { RUN_VECSCREEN } from '../subworkflows/local/run_vecscreen/main' -include { RUN_NT_KRAKEN } from '../subworkflows/local/run_nt_kraken/main' -include { RUN_FCSGX } from '../subworkflows/local/run_fcsgx/main' -include { RUN_FCSADAPTOR } from '../subworkflows/local/run_fcsadaptor/main' -include { RUN_DIAMOND as NR_DIAMOND } from '../subworkflows/local/run_diamond/main' -include { RUN_DIAMOND as UP_DIAMOND } from '../subworkflows/local/run_diamond/main' -include { ASCC_MERGE_TABLES } from '../modules/local/ascc/merge_tables/main' -include { RUN_DECONTAMINATE_FASTA } from '../subworkflows/local/run_decontaminate_fasta' -include { GENERATE_HTML_REPORT_WORKFLOW } from '../subworkflows/local/generate_html_report/main' - -// FUNCTION IMPORTS -// NOTE: IN FUTURE SHOULD ALSO CONTAIN DATA-MAPPER FUNCTIONS -include { getEmptyPlaceholder } from '../functions/local/ascc_utils' - - -/* -~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ - RUN MAIN WORKFLOW -~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ -*/ - -workflow ASCC_ORGANELLAR { - - take: - ch_samplesheet // channel: samplesheet read in from --input - _fcs_ov // params.fcs_override - fcs_samplesheet // The FCS override samplesheet for override - fcs_db // [path(path)] - _reads - scientific_name // val(name) - pacbio_database // tuple [[meta.id], pacbio_database] - ncbi_taxonomy_path - ncbi_ranked_lineage_path - nt_database_path - diamond_nr_db_path - diamond_uniprot_db_path - taxid - nt_kraken_db_path - vecscreen_database_path - reads_path - reads_type - ch_barcodes - val_reads_per_chunk - - main: - ch_versions = channel.empty() - - // - // LOGIC: CREATE run_conditional LIST - // - run_conditionals = ["both", "organellar"] - - - // - // LOGIC: PRETTY NOTIFICATION OF FILES AT STAGE - // - ch_samplesheet - .map { meta, sample -> - log.info "[ASCC INFO]: ORGANELLAR WORKFLOW:\n\t-- $meta\n\t-- $sample\n" - } - - - //------------------------------------------------------------------------- - // - // SUBWORKFLOW: RUNS FILTER_FASTA, GENERATE .GENOME, CALCS GC_CONTENT AND FINDS RUNS OF N's - // THIS SHOULD NOT RUN ONLY WHEN SPECIFICALLY REQUESTED - // - ESSENTIAL_JOBS( - ch_samplesheet - ) - ch_versions = ch_versions.mix(ESSENTIAL_JOBS.out.versions) - ej_reference_tuple = ESSENTIAL_JOBS.out.reference_tuple_from_GG - ej_seqkit_reference = ESSENTIAL_JOBS.out.reference_with_seqkit - ej_dot_genome = ESSENTIAL_JOBS.out.dot_genome - ej_gc_coverage = ESSENTIAL_JOBS.out.gc_content_txt - ej_trailing_ns = ESSENTIAL_JOBS.out.trailing_ns_report - ej_fasta_sanitation_log = ESSENTIAL_JOBS.out.filter_fasta_sanitation_log - ej_fasta_filter_log = ESSENTIAL_JOBS.out.filter_fasta_length_filtering_log - - - //------------------------------------------------------------------------- - // - // SUBWORKFLOW: EXTRACT RESULTS HITS FROM TIARA - // - TIARA_TIARA ( - ej_reference_tuple.filter{ _meta, _file -> params.run_tiara in run_conditionals } - ) - ch_versions = ch_versions.mix( TIARA_TIARA.out.versions ) - ch_tiara = TIARA_TIARA.out.classifications - .map { meta, file -> [[id: meta.id ], file] } - .ifEmpty { [[:],[]] } - - - //------------------------------------------------------------------------- - // - // SUBWORKFLOW: IDENTITY PACBIO BARCODES IN INPUT DATA - // - ej_reference_tuple - .combine(pacbio_database) - .multiMap{ - ref_meta, ref_data, pdb_meta, pdb_data -> - reference: [ref_meta, ref_data] - pacbio_db: [pdb_meta, pdb_data] - } - .set { duplicated_db } - - PACBIO_BARCODE_CHECK ( - duplicated_db.reference.filter{ _meta, _file -> - params.run_pacbio_barcodes in run_conditionals - }, - ch_barcodes, - duplicated_db.pacbio_db - ) - ch_versions = ch_versions.mix(PACBIO_BARCODE_CHECK.out.versions) - ch_barcode_check = PACBIO_BARCODE_CHECK.out.filtered.ifEmpty{ [[:],[]] } - - - //------------------------------------------------------------------------- - // - // SUBWORKFLOW: RUN FCS-ADAPTOR TO IDENTIDY ADAPTOR AND VECTORR CONTAMINATION - // - RUN_FCSADAPTOR ( - ej_reference_tuple.filter{ _meta, _file -> - params.run_fcs_adaptor in run_conditionals - } - ) - ch_versions = ch_versions.mix(RUN_FCSADAPTOR.out.versions) - ch_fcsadapt = RUN_FCSADAPTOR.out.ch_joint_report - - - //------------------------------------------------------------------------- - // - // SUBWORKFLOW: RUN FCS-GX TO IDENTIFY CONTAMINATION IN THE ASSEMBLY - // - - if ( params.run_fcsgx in run_conditionals && !params.fcs_override) { - - joint_channel = ej_reference_tuple - .combine(fcs_db) - .combine(taxid) - .combine(ncbi_ranked_lineage_path) - .multiMap { meta, ref, db, _tax_id, tax_path -> - def new_meta = [id: meta.id, taxid: meta.taxid] - reference: [new_meta, ref] - fcs_db_path: db - ncbi_tax_path: tax_path - } - - RUN_FCSGX ( - joint_channel.reference, - joint_channel.fcs_db_path, - joint_channel.ncbi_tax_path - ) - ch_versions = ch_versions.mix(RUN_FCSGX.out.versions) - - ch_fcsgx = RUN_FCSGX.out.fcsgxresult - ch_fcsgx_report = RUN_FCSGX.out.fcsgx_report_txt - ch_fcsgx_taxonomy = RUN_FCSGX.out.fcsgx_taxonomy_rpt - - } else if ( params.fcs_override ) { - - fcs_samplesheet.map{ meta, file -> - log.info("[ASCC INFO]: Overriding Internal FCSGX with ${file}") - [[id: meta.id], file] - - } - .set { ch_fcsgx } - - ch_fcsgx_report = channel.of( [[:],[]] ) - ch_fcsgx_taxonomy = channel.of( [[:],[]] ) - - } else { - ch_fcsgx = channel.of( [[:],[]] ) - ch_fcsgx_report = channel.of( [[:],[]] ) - ch_fcsgx_taxonomy = channel.of( [[:],[]] ) - } - - - //------------------------------------------------------------------------- - // - // SUBWORKFLOW: CALCULATE AVERAGE READ COVERAGE - // - RUN_READ_COVERAGE ( - ej_reference_tuple.filter{ _meta, _file -> - params.run_coverage in run_conditionals - }, - reads_path, - reads_type, //Subworkflow uses the param, not this value... as soon as it's in a channel it can't be used for a comparator. - val_reads_per_chunk - ) - ch_versions = ch_versions.mix(RUN_READ_COVERAGE.out.versions) - ch_coverage = RUN_READ_COVERAGE.out.tsv_ch.ifEmpty{ [[:], []] } - ch_bam = RUN_READ_COVERAGE.out.bam_ch.ifEmpty{ [[:], []] } - - - //------------------------------------------------------------------------- - // - // SUBWORKFLOW: SCREENING FOR VECTOR SEQUENCE - // - RUN_VECSCREEN ( - ej_reference_tuple.filter{ _meta, _file -> - params.run_vecscreen in run_conditionals - }, - vecscreen_database_path.first() - ) - ch_versions = ch_versions.mix(RUN_VECSCREEN.out.versions) - ch_vecscreen = RUN_VECSCREEN.out.vecscreen_contam.ifEmpty{ [[:],[]] } - - - //------------------------------------------------------------------------- - // - // SUBWORKFLOW: RUN THE KRAKEN CLASSIFIER - // - RUN_NT_KRAKEN( - ej_reference_tuple.filter{ _meta, _file -> - params.run_kraken in run_conditionals - }, - nt_kraken_db_path.first(), - ncbi_ranked_lineage_path.first() - ) - ch_versions = ch_versions.mix(RUN_NT_KRAKEN.out.versions) - ch_kraken1 = RUN_NT_KRAKEN.out.classified.ifEmpty{ [[:], []] } - ch_kraken2 = RUN_NT_KRAKEN.out.report.ifEmpty{ [[:], []] } - ch_kraken3 = RUN_NT_KRAKEN.out.lineage.ifEmpty{ [[:], []] } - - - // - // LOGIC: WE NEED TO MAKE SURE THAT THE INPUT SEQUENCE IS OF AT LEAST LENGTH OF params.seqkit_window - // - valid_length_fasta = ej_seqkit_reference - // - // NOTE: Here we are using the un-filtered genome, any filtering may (accidently) cause an empty fasta - // - .map{ meta, file -> - def total_length = 0 - file.eachLine { line -> - if (line && !line.startsWith('>')) { - total_length += line.length() - } - } - - def meta2 = [ - id: meta.id, - sliding: meta.sliding, - window: meta.window, - seq_count: total_length - ] - - [meta2, file] - } - .filter { meta, _file -> - meta.seq_count >= params.seqkit_window - } - - valid_length_fasta - .map{ meta, _file -> - log.info "[ASCC INFO]: Running BLAST (NT, DIAMOND, NR) on VALID ORGANELLE: \n\t-- ${meta.id}'s sequence ($meta.seq_count bases) is >= seqkit_window $params.seqkit_window\n" - } - - - //------------------------------------------------------------------------- - // - // SUBWORKFLOW: EXTRACT RESULTS HITS FROM NT-BLAST - // - EXTRACT_NT_BLAST ( - valid_length_fasta.filter{ _meta, _file -> - params.run_nt_blast in run_conditionals - }, - nt_database_path.first(), - ncbi_ranked_lineage_path.first() - ) - ch_versions = ch_versions.mix(EXTRACT_NT_BLAST.out.versions) - ch_nt_blast = EXTRACT_NT_BLAST.out.ch_blast_hits.ifEmpty { [[:],[]] } - ch_blast_lineage = EXTRACT_NT_BLAST.out.ch_top_lineages.ifEmpty { [[:],[]] } - ch_btk_format = EXTRACT_NT_BLAST.out.ch_btk_format.ifEmpty { [[:],[]] } - - - //------------------------------------------------------------------------- - // - // SUBWORKFLOW: DIAMOND BLAST FOR INPUT ASSEMBLY - // - NR_DIAMOND ( - valid_length_fasta.filter{ _meta, _file -> - params.run_nr_diamond in run_conditionals - }, - diamond_nr_db_path.first() - ) - ch_versions = ch_versions.mix(NR_DIAMOND.out.versions) - nr_full = NR_DIAMOND.out.reformed - .map { meta, file -> [[id: meta.id ], file] } - .ifEmpty { [[:],[]] } - - nr_hits = NR_DIAMOND.out.hits_file - .map { meta, file -> [[id: meta.id ], file] } - .ifEmpty { [[:],[]] } - - - //------------------------------------------------------------------------- - // - // SUBWORKFLOW: DIAMOND BLAST FOR INPUT ASSEMBLY - // - // NOTE: HEADER FORMAT WILL BE - - // qseqid sseqid pident length mismatch gapopen qstart qend sstart send - // evalue bitscore staxids sscinames sskingdoms sphylums salltitles - UP_DIAMOND ( - valid_length_fasta.filter{ _meta, _file -> - params.run_uniprot_diamond in run_conditionals - }, - diamond_uniprot_db_path.first() - ) - ch_versions = ch_versions.mix(UP_DIAMOND.out.versions) - un_full = UP_DIAMOND.out.reformed - .map { meta, file -> [[id: meta.id], file ] } - .ifEmpty { [[:],[]] } - - un_hits = UP_DIAMOND.out.hits_file - .map { meta, file -> [[id: meta.id ], file ] } - .ifEmpty { [[:],[]] } - - - //------------------------------------------------------------------------- - if ( params.run_create_btk_dataset in run_conditionals ) { - - // - // LOGIC: FOUND RACE CONDITION EFFECTING LONG RUNNING JOBS - // AND INPUT TO HERE ARE NOW MERGED AND MAPPED - // EMPTY CHANNELS ARE CHECKED AND DEFAULTED TO [[:],[]] - // - // - ej_reference_tuple - .map{meta, file -> [[id: meta.id], file]} - .join(ch_nt_blast, remainder: true) - .join(ch_tiara, remainder: true) - .join(ej_dot_genome,remainder: true) - .join(channel.of([[:],[]]), remainder: true) //ch_fcsgx - .join(ch_bam, remainder: true) - .join(ch_coverage, remainder: true) - .join(channel.of([[:],[]]), remainder: true) //ch_kmers - .join(ch_kraken1, remainder: true) - .join(ch_kraken2, remainder: true) - .join(ch_kraken3, remainder: true) - .join(nr_full, remainder: true) - .join(un_full, remainder: true) - .filter { items -> - def meta = items[0] - meta != null && - meta != [] && - !(meta instanceof Map && (meta.id == null || meta.isEmpty())) - } - .map { items -> - // Replace null values with placeholder file - items.withIndex().collect { item, index -> - if (item == null) { - getEmptyPlaceholder(index) - } else if (item instanceof List && item.isEmpty()) { - getEmptyPlaceholder(index) - } else { - item - } - } - } - .set{ create_input_channel} - - - // - // MODULE: CREATE A BTK COMPATIBLE DATASET FOR NEW DATA - // - CREATE_BTK_DATASET ( - create_input_channel, - params.taxid, - ncbi_taxonomy_path.first(), - scientific_name - - ) - ch_versions = ch_versions.mix(CREATE_BTK_DATASET.out.versions) - - ch_create_summary = CREATE_BTK_DATASET.out.create_summary - .map{ meta, _file -> [[ id: meta.id ], _file] } - - ch_create_btk_dataset = CREATE_BTK_DATASET.out.btk_datasets - .map{ meta, _file -> [[ id: meta.id ], _file] } - } else { - ch_create_summary = channel.of( [[:],[]] ) - ch_create_btk_dataset = channel.of( [[:],[]] ) - } - - - //------------------------------------------------------------------------- - // - // LOGIC: AUTOFILTER ASSEMBLY BY TIARA AND FCSGX RESULTS SO THE SUBWORKLOW CAN EITHER BE TRIGGERED BY THE VALUES tiara, fcs-gx, autofilter_assemlby AND EXCLUDE STEPS NOT CONTAINING autofilter_assembly - // OR BY include_steps CONTAINING ALL AND EXCLUDE NOT CONTAINING autofilter_assembly. - // - if ( - ( params.run_tiara in run_conditionals ) && - ( params.run_fcsgx in run_conditionals ) && - ( params.run_autofilter_assembly in run_conditionals ) - ) { - // - // LOGIC: FILTER THE INPUT FOR THE AUTOFILTER STEP - // - We can't just combine on meta.id as some of the Channels have other data - // in there too so we just sanitise, and _then_ combine on 0, and - // _then_ add back in the taxid as we need that for this process. - // Thankfully taxid is a param so easy enough to add back in. - // Actually, it just makes more sense to passs in as its own channel. - // - ej_reference_tuple - .map{ meta, file -> [[id: meta.id], file] } - .combine( - ch_tiara.map{ meta, file -> [[id: meta.id], file] }, by: 0 - ) - .combine( - ch_fcsgx.map{ meta, file -> [[id: meta.id], file] }, by: 0 - ) - .combine( - ncbi_ranked_lineage_path - ) - .combine( - taxid - ) - .multiMap{ - meta, ref, tiara, fcs, ncbi, thetaxid -> - def new_meta = [id: meta.id, taxid: thetaxid] - reference: [new_meta, ref] - tiara_file: [new_meta, tiara] - fcs_file: [new_meta, fcs] - ncbi_rank: ncbi - } - .set { autofilter_input_formatted } - - - // - // MODULE: AUTOFILTER ASSEMBLY BY TIARA AND FCSGX RESULTS - // - AUTOFILTER_AND_CHECK_ASSEMBLY ( - autofilter_input_formatted.reference, - autofilter_input_formatted.tiara_file, - autofilter_input_formatted.fcs_file, - autofilter_input_formatted.ncbi_rank - ) - ch_versions = ch_versions.mix(AUTOFILTER_AND_CHECK_ASSEMBLY.out.versions) - ch_autofilt_assem = AUTOFILTER_AND_CHECK_ASSEMBLY.out.decontaminated_assembly - ch_autofilt_indicator = AUTOFILTER_AND_CHECK_ASSEMBLY.out.indicator_file - ch_autofilt_removed_seqs= AUTOFILTER_AND_CHECK_ASSEMBLY.out.removed_seqs - ch_autofilt_raw_report = AUTOFILTER_AND_CHECK_ASSEMBLY.out.raw_report - - ch_autofilt_alarm_file = AUTOFILTER_AND_CHECK_ASSEMBLY.out.alarm_file - .map{ meta, _file -> [[ id: meta.id ], _file] } - - ch_autofilt_fcs_tiara = AUTOFILTER_AND_CHECK_ASSEMBLY.out.fcs_tiara_summary - .map{ meta, _file -> [[ id: meta.id ], _file] } - - } else { - ch_autofilt_alarm_file = channel.of( [[:],[]] ) - ch_autofilt_removed_seqs= channel.of( [[:],[]] ) - ch_autofilt_assem = channel.of( [[:],[]] ) - ch_autofilt_indicator = channel.of( [[:],[]] ) - ch_autofilt_fcs_tiara = channel.of( [[:],[]] ) - ch_autofilt_raw_report = channel.of( [[:],[]] ) - } - - - // - // LOGIC: EACH SUBWORKFLOW OUTPUTS EITHER AN EMPTY CHANNEL OR A FILE CHANNEL DEPENDING ON THE RUN RULES - // SO THE RULES FOR THIS ONLY NEED TO BE A SIMPLE "DO YOU WANT IT OR NOT" - // - if ( - ( params.run_essentials in run_conditionals ) && - ( params.run_merge_datasets in run_conditionals ) - ) { - - // - // LOGIC: FOUND RACE CONDITION EFFECTING LONG RUNNING JOBS - // AND INPUT TO HERE ARE NOW MERGED AND MAPPED - // EMPTY CHANNELS ARE CHECKED AND DEFAULTED TO [[:],[]] - // - ej_reference_tuple - .map{meta, file -> [[id: meta.id], file]} - .join(ej_gc_coverage - .map{meta, file -> [[id: meta.id], file]}, remainder: true) - .join(ch_coverage, remainder: true) - .join(ch_tiara, remainder: true) - .join(ch_kraken3, remainder: true) - .join(ch_blast_lineage, remainder: true) - .join(channel.of([[:],[]]), remainder: true) //ch_kmers - not in organellar - .join(nr_hits, remainder: true) - .join(un_hits, remainder: true) - .join(ch_create_summary,remainder: true) - .join(channel.of([[:],[]]), remainder: true) //busco_merge_btk - not in organellar - .join(ch_fcsgx, remainder: true) - .filter { items -> - def meta = items[0] - meta != null && - meta != [] && - !(meta instanceof Map && (meta.id == null || meta.isEmpty())) - } - .map { items -> - // Replace null values with placeholder file - items.withIndex().collect { item, index -> - if (item == null) { - getEmptyPlaceholder(index) - } else if (item instanceof List && item.isEmpty()) { - getEmptyPlaceholder(index) - } else { - item - } - } - } - .set{ merge_input_channel} - - ASCC_MERGE_TABLES ( - merge_input_channel - ) - ch_versions = ch_versions.mix(ASCC_MERGE_TABLES.out.versions) - org_merged_table = ASCC_MERGE_TABLES.out.merged_table - .map{ meta, _file -> [[id:meta.id ], _file] } - - org_merged_phylum_count = ASCC_MERGE_TABLES.out.phylum_counts - .map{ meta, _file -> [[id:meta.id], _file] } - - } else { - org_merged_table = channel.of( [[:],[]] ) - //merged_extended_table = channel.empty() - org_merged_phylum_count = channel.of( [[:],[]] ) - } - - - //------------------------------------------------------------------------- - // - // SUBWORKFLOW: GENERATE HTML REPORT (minimal wiring, opt-in) - // Gate with params.run_html_report to avoid altering default behavior. - // - - // Params file - ch_params_file = params.params_file ? channel.fromPath(params.params_file) : channel.value([]) - - GENERATE_HTML_REPORT_WORKFLOW ( - ch_barcode_check, - ch_fcsadapt, - ej_trailing_ns, - ch_vecscreen, - ch_autofilt_fcs_tiara, - org_merged_table, - org_merged_phylum_count, - channel.of( [[:],[]] ), - ej_reference_tuple.filter{ _meta, _file -> - params.run_html_report in run_conditionals - }, - ej_fasta_sanitation_log, - ej_fasta_filter_log, - ch_params_file, - ch_fcsgx_report, - ch_fcsgx_taxonomy, - ch_create_btk_dataset - ) - ch_versions = ch_versions.mix(GENERATE_HTML_REPORT_WORKFLOW.out.versions) - - - // - // SUBWORKFLOW: GENERATE DECONTAMINATION FILES AND POTENTIALLY A DECONTAMINATED FASTA - // THIS SHOULD ONLY RUN IF STANDARD CONDITIONALS ARE MET - // AND ABNORMAL CONTAMINATION IS FOUND - // AUTOFILTERING THE ASSEMBLY IS ESSENTIAL FOR DECON TO RUN - - // We only want the EUKARYOTIC report - // Not using the collection will result in a `Unexpected error [ConcurrentModificationException]` - // `ch_fcsadapt` because it is a mix channel, is technically still mutable - euk_fcsadapt = ch_fcsadapt.map{ meta, files -> - def filesCopy = (files ?: []).collect() // defensive copy - [meta, filesCopy.find{ file -> file.name.endsWith('_euk.fcs_adaptor_report.txt') }] - } - - ej_reference_tuple_filtered = ej_reference_tuple - .filter{ meta, file -> - params.run_decontaminate_fasta in run_conditionals && params.run_autofilter_assembly in run_conditionals - } - .map{ meta, file -> [[id: meta.id], file] } - - RUN_DECONTAMINATE_FASTA( - ej_reference_tuple_filtered, - ch_fcsgx, - ch_autofilt_fcs_tiara, - euk_fcsadapt, - ej_trailing_ns, - ch_barcode_check, - channel.of( [[:],[]] ), - channel.of( [[:],[]] ) - ) - ch_versions = ch_versions.mix(RUN_DECONTAMINATE_FASTA.out.versions) - - - emit: - - essential_reference = ej_reference_tuple - essential_genome_file = ej_dot_genome - essential_gc_cov = ej_gc_coverage - - blast_output = ch_nt_blast - blast_lineage = ch_blast_lineage - blast_btk_formatted = ch_btk_format - - diamond_nr_blast_full = nr_full - diamond_nr_blast_hits = nr_hits - - diamond_un_blast_full = un_full - diamond_un_blast_hits = un_hits - - read_coverage_output = ch_coverage - read_coverage_bam = ch_bam - - fcsadaptor_prok_euk = ch_fcsadapt - fcsgx_output = ch_fcsgx - - autofilter_deconned_assm = ch_autofilt_assem - autofilter_fcs_tiar_smry = ch_autofilt_fcs_tiara - autofilter_removed_seqs = ch_autofilt_removed_seqs - autofilter_alarm_file = ch_autofilt_alarm_file - autofilter_indicator_file = ch_autofilt_indicator - autofilter_raw_report = ch_autofilt_raw_report - - create_btk_ds_dataset = ch_create_btk_dataset - create_btk_ds_create_smry = ch_create_summary - - kraken2_classified = ch_kraken1 - kraken2_report = ch_kraken2 - kraken2_lineage = ch_kraken3 - - vecscreen_contam = ch_vecscreen - - tiara_output = ch_tiara - - versions = ch_versions - -} - -/* -~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ - THE END -~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ -*/