From d01162025ad530b7d8ce450cd7e44c32ab46bc2c Mon Sep 17 00:00:00 2001 From: Neil Horner Date: Tue, 5 Dec 2023 18:42:47 +0000 Subject: [PATCH] Resolve CW-2369 --- CHANGELOG.md | 1 + main.nf | 46 +++++++++++++++++++----------- nextflow.config | 1 - subworkflows/reference_assembly.nf | 15 ++++++---- 4 files changed, 40 insertions(+), 23 deletions(-) diff --git a/CHANGELOG.md b/CHANGELOG.md index bb3f33f..8a9fb0e 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -6,6 +6,7 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0 ## [unreleased] ### Added +- Published minimap2 and pychopper results to output directory. - Two extra pychopper parameters `--cdna_kit` and `--pychopper_backend`. `--pychopper_options` is still available to define any other options. ## [v0.4.2] diff --git a/main.nf b/main.nf index 7d6bcf7..c5e71cd 100644 --- a/main.nf +++ b/main.nf @@ -142,21 +142,32 @@ process preprocess_reads { cpus 4 memory "2 GB" input: - tuple val(meta), path(input_reads) + tuple val(meta), path('seqs.fastq.gz') output: - tuple val("${meta.alias}"), path("${meta.alias}_full_length_reads.fastq"), emit: full_len_reads - path '*.tsv', emit: report + tuple val("${meta.alias}"), + path("${meta.alias}_pychopper_output/${meta.alias}_full_length_reads.fastq"), + emit: full_len_reads + tuple val("${meta.alias}"), + path("${meta.alias}_pychopper_output/"), + emit: pychopper_output + path "${meta.alias}_pychopper_output/${meta.alias}_pychopper.tsv", + emit: report script: def cdna_kit = params.cdna_kit.split("-")[-1] - def extra_params = params.pychopper_opts ?: '' + def extra_params = params.pychopper_opts ?: '' """ - pychopper -t ${params.threads} -k ${cdna_kit} -m ${params.pychopper_backend} ${extra_params} ${input_reads} ${meta.alias}_full_length_reads.fastq + pychopper -t ${params.threads} -k ${cdna_kit} -m ${params.pychopper_backend} ${extra_params} 'seqs.fastq.gz' ${meta.alias}_full_length_reads.fastq mv pychopper.tsv ${meta.alias}_pychopper.tsv workflow-glue generate_pychopper_stats --data ${meta.alias}_pychopper.tsv --output . - # Add sample id column + # Add sample id colum sed "1s/\$/\tsample_id/; 1 ! s/\$/\t${meta.alias}/" ${meta.alias}_pychopper.tsv > tmp mv tmp ${meta.alias}_pychopper.tsv + + + mkdir "${meta.alias}_pychopper_output/" + shopt -s extglob # Allow extended pattern matching so we can exclude files from the mv + mv !("${meta.alias}_pychopper_output"|seqs.fastq.gz) "${meta.alias}_pychopper_output/" """ } @@ -541,7 +552,7 @@ workflow pipeline { } - + results = Channel.empty() software_versions = getVersions() workflow_params = getParams() input_reads = reads.map{ meta, samples, stats -> [meta, samples]} @@ -552,6 +563,8 @@ workflow pipeline { preprocess_reads(input_reads) full_len_reads = preprocess_reads.out.full_len_reads pychopper_report = preprocess_reads.out.report.collectFile(keepHeader: true) + pychopper_results_dir = preprocess_reads.out.pychopper_output.map{ it -> it[1]} + results = results.concat(pychopper_results_dir) } else{ full_len_reads = input_reads.map{ meta, reads -> [meta.alias, reads]} @@ -564,7 +577,7 @@ workflow pipeline { assembly_stats = assembly.stats.map{ it -> it[1]}.collect() - split_bam(assembly.bam) + split_bam(assembly.bam.map {sample_id, bam, bai -> [sample_id, bam]}) assemble_transcripts(split_bam.out.bundles.flatMap(map_sample_ids_cls).combine(ref_annotation),use_ref_ann) @@ -582,14 +595,14 @@ workflow pipeline { gff_compare = run_gffcompare.out.gffcmp_dir.map{ it -> it[1]}.collect() merge_gff = merge_gff_bundles.out.gff.map{ it -> it[1]}.collect() - results = Channel.empty() - }else - { + results = results.concat(assembly.bam.map {sample_id, bam, bai -> [bam, bai]}.flatten()) + } + else{ gff_compare = file("$projectDir/data/OPTIONAL_FILE") merge_gff = file("$projectDir/data/OPTIONAL_FILE") assembly_stats = file("$projectDir/data/OPTIONAL_FILE") use_ref_ann = false - results = Channel.empty() + } if (jaffal_refBase){ gene_fusions(full_len_reads, jaffal_refBase, jaffal_genome, jaffal_annotation) @@ -638,11 +651,9 @@ workflow pipeline { de_report, count_transcripts_file) - report = makeReport.out.report - - - - results = results.concat(makeReport.out.report) + report = makeReport.out.report + + results = results.concat(makeReport.out.report) if (use_ref_ann){ results = run_gffcompare.output.gffcmp_dir.concat( @@ -675,6 +686,7 @@ workflow pipeline { } results.concat(workflow_params.map{ [it, null]}) + emit: results diff --git a/nextflow.config b/nextflow.config index d38d5f8..9911804 100644 --- a/nextflow.config +++ b/nextflow.config @@ -115,7 +115,6 @@ manifest { version = 'v0.4.2' } - epi2melabs { tags = "isoforms, transcriptomics" } diff --git a/subworkflows/reference_assembly.nf b/subworkflows/reference_assembly.nf index 22fc83d..f0bd299 100644 --- a/subworkflows/reference_assembly.nf +++ b/subworkflows/reference_assembly.nf @@ -12,18 +12,23 @@ process map_reads{ tuple val(sample_id), path (fastq_reads), path(index), path(reference) output: - tuple val(sample_id), path("${sample_id}_reads_aln_sorted.bam"), emit: bam + tuple val(sample_id), + path("${sample_id}_reads_aln_sorted.bam"), + path("${sample_id}_reads_aln_sorted.bam.bai"), + emit: bam tuple val(sample_id), path("${sample_id}_read_aln_stats.tsv"), emit: stats script: def ContextFilter = """AlnContext: { Ref: "${reference}", LeftShift: -${params.poly_context}, RightShift: ${params.poly_context}, RegexEnd: "[Aa]{${params.max_poly_run},}", Stranded: True,Invert: True, Tsv: "internal_priming_fail.tsv"} """ + + def mm2_threads = Math.min(task.cpus - 3, 1) """ - minimap2 -t ${params.threads} -ax splice ${params.minimap2_opts} ${index} ${fastq_reads}\ + minimap2 -t ${mm2_threads} -ax splice ${params.minimap2_opts} ${index} ${fastq_reads}\ | samtools view -q ${params.minimum_mapping_quality} -F 2304 -Sb -\ - | seqkit bam -j ${params.threads} -x -T '${ContextFilter}' -\ - | samtools sort -@ ${params.threads} -o "${sample_id}_reads_aln_sorted.bam" - ; - ((cat "${sample_id}_reads_aln_sorted.bam" | seqkit bam -s -j ${params.threads} - 2>&1) | tee ${sample_id}_read_aln_stats.tsv ) || true + | seqkit bam -j 1 -x -T '${ContextFilter}' -\ + | samtools sort --write-index -@ 1 -o "${sample_id}_reads_aln_sorted.bam##idx##${sample_id}_reads_aln_sorted.bam.bai" - ; + ((cat "${sample_id}_reads_aln_sorted.bam" | seqkit bam -s -j 1 - 2>&1) | tee ${sample_id}_read_aln_stats.tsv ) || true if [[ -s "internal_priming_fail.tsv" ]]; then