diff --git a/.github/workflows/nf-test.yml b/.github/workflows/nf-test.yml index b53498b..9b512f5 100644 --- a/.github/workflows/nf-test.yml +++ b/.github/workflows/nf-test.yml @@ -36,6 +36,7 @@ jobs: - "only_trimming/skip_trimming" - "only_alignment/alignment" - "only_moveumitoheader" + - "only_cleanheaders" steps: - name: Check out pipeline code uses: actions/checkout@v3 diff --git a/Documentation.md b/Documentation.md index db3ffe1..0a5c45f 100644 --- a/Documentation.md +++ b/Documentation.md @@ -37,6 +37,12 @@ MultiQC report of important quality control metrics. Produced by MULTIQC process ## Commonly used parameters +**Cleaning spaces from fastq headers** + +Some fastq files carry spaces in their header lines (for example `@READID 1:N:0:INDEX`). Tools further down the pipeline truncate read names at the first space, which can break UMI extraction and read-name matching after alignment. + +To enable this option set `clean_fastq_headers = true`. Every space in a header line is replaced with an underscore, and the cleaned reads are used for all subsequent steps. This runs before any other analysis, including UMI moving and trimming. + **Moving UMI from fastq reads to read header** If you are analysing demultiplexed sample files, then depending on where you have sourced your fastq from, the UMI and experimental barcode might still be present at the 5' end of reads, which will cause errors in downstream analysis. Note this is common for historical public data downloaded from ArrayExpress, for example. diff --git a/conf/logic.config b/conf/logic.config index 44c0311..5748255 100644 --- a/conf/logic.config +++ b/conf/logic.config @@ -7,6 +7,7 @@ // Defaults params { + run_clean_fastq_headers = false run_genome_prep = true run_input_check = true run_move_umi_to_header = false @@ -20,6 +21,7 @@ params { } // Set other logic +if(params.clean_fastq_headers) { params.run_clean_fastq_headers = true } if(params.move_umi_to_header) { params.run_move_umi_to_header = true } if(params.skip_umi_dedupe) { params.run_umi_dedup = false } @@ -36,6 +38,7 @@ if(params.only_input) { } if(params.only_genome) { + params.run_clean_fastq_headers = false params.run_input_check = false params.run_trim_galore_fastqc = false params.run_alignment = false diff --git a/conf/modules.config b/conf/modules.config index 62e5729..0ae451e 100644 --- a/conf/modules.config +++ b/conf/modules.config @@ -196,6 +196,18 @@ if(params.run_genome_prep) { FASTQ + TRIMMING ======================================================================================== */ +if(params.run_clean_fastq_headers) { + process { + withName: 'CLIPSEQ:CLEAN_FASTQ_HEADERS' { + publishDir = [ + path: { "${params.outdir}/01_prealign/clean_fastq_headers" }, + mode: "${params.publish_dir_mode}", + enabled: false + ] + } + } +} + if(params.run_move_umi_to_header) { process { withName: 'CLIPSEQ:UMITOOLS_EXTRACT' { diff --git a/main.nf b/main.nf index 42d943a..d6350d2 100644 --- a/main.nf +++ b/main.nf @@ -94,6 +94,7 @@ ch_multiqc_config = file("$projectDir/assets/multiqc_config.yml", checkIfExists: // MODULEs // +include { CLEAN_FASTQ_HEADERS } from './modules/local/clean_fastq_headers/main' include { MULTIQC } from './modules/local/multiqc' include { GET_CROSSLINKS as CALC_SMRNA_K1_CROSSLINKS } from './modules/local/get_crosslinks' include { GET_CROSSLINKS as CALC_GENOME_CROSSLINKS } from './modules/local/get_crosslinks' @@ -279,6 +280,17 @@ workflow CLIPSEQ { } //EXAMPLE CHANNEL STRUCT: [[id:h3k27me3_R1, group:h3k27me3, replicate:1, single_end:false], [FASTQ]] //ch_fastq | view + if(params.run_clean_fastq_headers) { + /* + * MODULE: Replace spaces in FASTQ headers with underscores before any analysis + */ + CLEAN_FASTQ_HEADERS ( + ch_fastq + ) + ch_versions = ch_versions.mix(CLEAN_FASTQ_HEADERS.out.versions) + ch_fastq = CLEAN_FASTQ_HEADERS.out.reads + } + if(params.encode_eclip){ ENCODE_MOVEUMI ( ch_fastq diff --git a/modules/local/clean_fastq_headers/main.nf b/modules/local/clean_fastq_headers/main.nf new file mode 100644 index 0000000..6554649 --- /dev/null +++ b/modules/local/clean_fastq_headers/main.nf @@ -0,0 +1,22 @@ +process CLEAN_FASTQ_HEADERS { + tag "$meta.id" + label "process_single" + + conda "bioconda::biopython=1.78 pigz=2.6" + container "quay.io/biocontainers/mulled-v2-877c4e5a8fad685ea5bde487e04924ac447923b9:b7daa641364165419b9a87d9988bc803f913c5b6-0" + + input: + tuple val(meta), path(reads) + + output: + tuple val(meta), path("*.clean.fastq.gz"), emit: reads + path "versions.yml" , emit: versions + + when: + task.ext.when == null || task.ext.when + + shell: + prefix = task.ext.prefix ?: "${meta.id}" + process_name = task.process + template 'clean_fastq_headers.py' +} diff --git a/modules/local/clean_fastq_headers/meta.yml b/modules/local/clean_fastq_headers/meta.yml new file mode 100644 index 0000000..623cb54 --- /dev/null +++ b/modules/local/clean_fastq_headers/meta.yml @@ -0,0 +1,37 @@ +name: clean_fastq_headers +description: Replace spaces in FASTQ header lines with underscores before any analysis +tools: + - pigz: + description: | + A parallel implementation of gzip for modern multi-processor, + multi-core machines. + homepage: https://zlib.net/pigz/ + documentation: https://zlib.net/pigz/pigz.pdf + licence: ["Zlib"] +input: + - meta: + type: map + description: | + Groovy Map containing sample information + e.g. [ id:'test', single_end:false ] + - reads: + type: file + description: A single FASTQ file (this pipeline is single-end only) + pattern: "*.{fastq,fastq.gz}" + +output: + - meta: + type: map + description: | + Groovy Map containing sample information + e.g. [ id:'test', single_end:false ] + - reads: + type: file + description: FASTQ file with spaces in the header lines replaced by underscores + pattern: "*.clean.fastq.gz" + - versions: + type: file + description: File containing software versions + pattern: "versions.yml" +authors: + - "@Chromojones" diff --git a/modules/local/clean_fastq_headers/templates/clean_fastq_headers.py b/modules/local/clean_fastq_headers/templates/clean_fastq_headers.py new file mode 100644 index 0000000..5b0696b --- /dev/null +++ b/modules/local/clean_fastq_headers/templates/clean_fastq_headers.py @@ -0,0 +1,76 @@ +#!/usr/bin/env python3 + +import platform +import subprocess + + +def process_fastq(input_file, output_file, threads): + """Process FASTQ file to remove spaces from headers using pigz.""" + is_gzipped = input_file.endswith('.gz') + + in_file = None + if is_gzipped: + decompress = subprocess.Popen(['pigz', '-p', str(threads), '-dc', input_file], + stdout=subprocess.PIPE) + else: + in_file = open(input_file, 'rb') + decompress = subprocess.Popen(['cat'], stdin=in_file, stdout=subprocess.PIPE) + + with open(output_file, 'wb') as out_handle: + compress = subprocess.Popen(['pigz', '-p', str(threads), '-c'], + stdin=subprocess.PIPE, stdout=out_handle) + try: + line_count = 0 + for line in decompress.stdout: + if line_count % 4 == 0: # Header line + line = line.decode().strip().replace(' ', '_').encode() + b'\n' + compress.stdin.write(line) + line_count += 1 + finally: + decompress.stdout.close() + decompress.wait() + if in_file is not None: + in_file.close() + compress.stdin.close() + compress.wait() + + # Both ends of the pipe must be checked. pigz reports a corrupt or + # truncated input by exiting non-zero after emitting a partial stream, so + # without this the task "succeeds" and hands a silently truncated FASTQ to + # the rest of the pipeline. + if decompress.returncode != 0: + raise RuntimeError( + "reading %s failed (exit %d) - the file is likely truncated or corrupt" + % (input_file, decompress.returncode)) + if compress.returncode != 0: + raise RuntimeError( + "writing %s failed (exit %d)" % (output_file, compress.returncode)) + + if line_count % 4 != 0: + raise RuntimeError( + "%s ended mid-record after %d lines - the file is likely truncated" + % (input_file, line_count)) + + +def pigz_version(): + """pigz reports its version on stderr on some builds, stdout on others.""" + proc = subprocess.run(['pigz', '--version'], stdout=subprocess.PIPE, + stderr=subprocess.STDOUT) + return proc.stdout.decode().strip().replace('pigz ', '') + + +# This pipeline is single-end only: single_end is parsed into meta but never +# read downstream, and each sample is one FASTQ. Fail loudly rather than +# quietly emit files the rest of the pipeline cannot use. +reads = "!{reads}".replace('[', '').replace(']', '').replace(',', ' ').split() +if len(reads) != 1: + raise SystemExit( + "CLEAN_FASTQ_HEADERS expects exactly one FASTQ per sample, got %d (%s)" + % (len(reads), ", ".join(reads))) + +process_fastq(reads[0], "!{prefix}.clean.fastq.gz", !{task.cpus}) + +with open("versions.yml", "w") as out_f: + out_f.write("!{process_name}" + ":\n") + out_f.write(" python: " + platform.python_version() + "\n") + out_f.write(" pigz: " + pigz_version() + "\n") diff --git a/nextflow.config b/nextflow.config index e857fac..7db3e4e 100644 --- a/nextflow.config +++ b/nextflow.config @@ -14,7 +14,7 @@ params { publish_dir_mode = 'symlink' monochrome_logs = false debug = false - ignore_params = "run_genome_prep,run_input_check,run_trim_galore_fastqc,run_alignment,run_read_filter,run_umi_dedup,run_calc_crosslinks,run_peak_calling,run_reporting,only_input,only_genome,only_trimming,only_alignment,only_filtering,only_dedup,only_crosslinks,only_peakcalling" + ignore_params = "run_clean_fastq_headers,run_genome_prep,run_input_check,run_trim_galore_fastqc,run_alignment,run_read_filter,run_umi_dedup,run_calc_crosslinks,run_peak_calling,run_reporting,only_input,only_genome,only_trimming,only_alignment,only_filtering,only_dedup,only_crosslinks,only_peakcalling" multiqc_title = null // Max resource options @@ -70,6 +70,7 @@ params { // Pipeline params crosslink_position = "start" + clean_fastq_headers = false encode_eclip = false move_umi_to_header = false umi_header_format = null diff --git a/schema/clipseq.json b/schema/clipseq.json index 55c0c9b..59e4187 100644 --- a/schema/clipseq.json +++ b/schema/clipseq.json @@ -303,6 +303,11 @@ "description": "Format of the UMI header to extract.", "type": "string" }, + "clean_fastq_headers": { + "name": "Clean FASTQ Headers", + "description": "Replace spaces in FASTQ header lines with underscores before any analysis.", + "type": "boolean" + }, "skip_umi_dedupe": { "name": "Skip UMI Deduplication", "description": "Skip the UMI deduplication of BAMs.", diff --git a/tests/nf_test/only_cleanheaders/main.nf.test b/tests/nf_test/only_cleanheaders/main.nf.test new file mode 100644 index 0000000..b05a62a --- /dev/null +++ b/tests/nf_test/only_cleanheaders/main.nf.test @@ -0,0 +1,28 @@ +nextflow_pipeline { + + name "only_cleanheaders" + script "main.nf" + + test("run_test") { + when { + params { + outdir = "$outputDir" + only_trimming = true + save_trimmed = true + clean_fastq_headers = true + } + } + + then { + assert workflow.success + + // The cleaned fastqs are intermediates and are not published, so + // assert the process ran. TrimGalore re-links its input to + // "${prefix}.fastq.gz", so downstream file names are unchanged. + assert workflow.trace.succeeded().any { it.process ==~ /.*CLEAN_FASTQ_HEADERS.*/ } + + assert new File("$outputDir/01_prealign/trimgalore/sample_1_R1.fastq.gz_trimming_report.txt").exists() + assert new File("$outputDir/01_prealign/trimgalore/sample_1_R1_trimmed.fq.gz").exists() + } + } +}