Skip to content

Commit ccb13a7

Browse files
pinin4fjordsclaude
andauthored
feat(custom): orfnormalise + orfmerge modules + orftable_fasta_gtf_buildorfcatalogue subworkflow (#11740)
* feat(custom/orfnormalise,orfmerge,orftable_fasta_gtf_buildorfcatalogue): scaffold ORF catalogue chain Add a three-component upstream chain for building cross-caller Ribo-seq ORF catalogues: - custom/orfnormalise: dispatches on meta.caller in {ribocode, ribotish, ribotricer, rpbp, price}, emits unified BED12 + sidecar TSV with the harmonised orf_class vocabulary (canonical_cds / uORF / dORF / novel_u / smORF / other). - custom/orfmerge: class-aware clustering across callers and samples, records called_by_<caller> / score_<caller> provenance with direction-aware best-score aggregation. - orftable_fasta_gtf_buildorfcatalogue: composes normaliser x N callers + merger + bedtools/getfasta + seqkit/translate, emits BED12, catalogue TSV, orf_to_gene TSV, AA FASTA and an MQC sidecar. Templates use the same python/pandas/pyyaml Wave container as custom/bed12codonpositions. Stub tests scaffolded; real-fixture tests follow once the paired test-datasets PR lands. Co-Authored-By: Claude Opus 4.7 (1M context) <noreply@anthropic.com> * feat(orf-catalogue chain): add real-fixture tests, fix ribotricer ORF_ID fallback - orfnormalise: 5 per-caller tests + 1 stub (all green). Fall back to parsing the ribotricer ORF_ID when the optional 'coordinate' column is absent (detect-orfs output does not carry it). Tolerate empty GTF inputs by requiring a real file before trying to read. - orfmerge: chain via setup{} from orfnormalise outputs (ribotish + ribocode on chr20); disambiguate per-caller prefix via a bundled setup_prefix.config so the merger doesn't see name collisions. - subworkflow: bundle nextflow.config that scopes the same caller-aware prefix to ORFTABLE_FASTA_GTF_BUILDORFCATALOGUE:CUSTOM_ORFNORMALISE. Test sets up GUNZIP for the chr20 FASTA since bedtools/getfasta needs uncompressed input. Fixtures live in nf-core/test-datasets#2070; will switch the test URLs from the pinin4fjords fork back to params.modules_testdata_base_path once that PR merges. Co-Authored-By: Claude Opus 4.7 (1M context) <noreply@anthropic.com> * fix(custom/orfnormalise): correct PRICE Location regex + drop empty-md5 stubs PRICE encodes strand as a single char appended to the chrom rather than a separate colon-bracketed field (`19+:start-end`), so the previous regex never matched and the parser returned an empty BED12 on real data. Switch to a non-greedy `chrom` capture so the trailing strand char is correctly extracted. Drop the per-module stub tests: real-fixture tests already cover every caller + the merge + the full subworkflow chain, and the stub snapshots were tripping the test_snap_md5sum lint rule (empty-file md5). Co-Authored-By: Claude Opus 4.7 (1M context) <noreply@anthropic.com> * fix(custom/orfnormalise + subworkflow): self-review pass on inherited parsers + empty-cohort handling Inherited parser bugs (carried over from the local riboseq port): - RiboCode: aa_length always 0 because the parser looked up an AA_length column that does not exist in predicted_orfs.txt. Derive from ORF_length (nt) instead. - Rp-Bp: bf_mean read from column 15 (p_translated_var) instead of column 18 (bayes_factor_mean), and orf_type read from column 14 (a metric, not a category). RPBP's predicted-orfs BED has no orf_type column; default to canonical_cds for the post-selectfinalpredictionset curated set. - Ribo-TISH: combined p-value never read because the parser looked for Pvalcombined/Pvalue/Pvalcom; actual columns are FisherPvalue / RiboPvalue / TISPvalue, with "None" string sentinels. Walk the preference list and skip None strings. Subworkflow gaps: - ch_orf_tables empty case crashed CUSTOM_ORFMERGE (arity '1..*' violated by .collect() emitting an empty list). Filter out the empty case so output channels are simply empty. - BEDTOOLS_GETFASTA and SEQKIT_TRANSLATE got generic ${meta.id}.{fa,fasta} filenames; pin to .catalogue.nt.fa and .catalogue.aa.fasta via bundled ext.prefix. Add ext.args = '-split -s -nameOnly' to getfasta (splice-aware extraction of BED12 blocks) and '--trim' to seqkit translate (drop trailing stop codons). - meta.yml: cite all five caller tools with verified DOIs (RiboCode Xiao 2018, Ribo-TISH Zhang 2017, Ribotricer Choudhary 2020, Rp-Bp Malone 2017, PRICE Erhard 2018); document that ch_orf_tables empty short-circuits; clarify GTF is required for ribocode/ribotish and optional for the rest; note that cohort-level output meta is hardcoded to [id:'cohort']. Co-Authored-By: Claude Opus 4.7 (1M context) <noreply@anthropic.com> * fix(orf-catalogue): second-review pass on parsers, clustering, fixtures, tests Surfaced + addressed via independent third-party assessment of branch 6597190 (parsers vs. canonical upstream output formats; merger clustering logic; test depth). Parsers (orfnormalise): - Ribotricer: detect-orfs output does NOT have the `coordinate` column (that lives in the prepare-orfs index file). The previous fallback treated the ORF_ID span as a single block, which for any multi-exon host transcript emits a BED12 record that spans introns - biologically nonsense and contaminates the catalogue AA FASTA produced via bedtools/getfasta + seqkit/translate. New _ribotricer_blocks_from_id takes the transcript map and intersects the ORF span with the host transcript's exon structure (same pattern ribotish's fallback already uses), recovering proper multi-exon blocks. - PRICE: orf_id_raw.split("__", 1)[0] used double-underscore; actual ID format is `<tid>_<type>_<index>` (single underscores). Was dead code for fixtures with Gene populated, but still wrong; switch to single underscore so the transcript_id lookup against the GTF works. Merger (orfmerge): - cluster_by_transcript collapsed every uORF / dORF / other ORF sharing a transcript_id into one cluster. A transcript can host multiple distinct uORFs (biologically common) so the catalogue under-reported them. Split into cluster_by_transcript (canonical_cds only) + cluster_by_transcript_position (uORF/dORF/other; keyed on outer span as well so distinct positional ORFs stay separate). - Document order-dependence of cluster_by_reciprocal_overlap at the threshold boundary. - Drop dead `_bed_key` helper. Tests: - Snapshot-only assertions caught changes but not correctness. Add content-level checks per caller: aa_length > 0 for all rows, score column populated, orf_class distribution not collapsed to a single bucket, multi-block BED12 records for ribotricer and price. These would have caught every parser bug in the assessment. - Merger test asserts called_by_*/score_* columns are coherent (catches the case where a parser bug populates the indicator but the score column is empty). Documentation (meta.yml): - Cite all five caller tools with verified DOIs (already in previous commit) - add note that Rp-Bp's predicted-orfs BED has no orf_type column so this module defaults Rp-Bp calls to canonical_cds. - Note PRICE iORF/intronic/orphan -> 'other' collapse (lossy). - Note ribotricer multi-exon block recovery via GTF + that GTF is strongly recommended for ribotricer (not just ribocode/ribotish). - Subworkflow: AA fasta now lands as `${cohort}.catalogue.aa.fasta` (from the previous commit). Paired test-datasets PR nf-core/test-datasets#2070 also gets a small fixture update: the ribotish/predict fixture's source run lacked a TIS BAM so every p-value column came out as the literal string "None", which silently masked the broken pvalue parser. Substitute the real RiboPvalue value into the FisherPvalue column via awk so the parser path is exercised; documented in the fixture README. Co-Authored-By: Claude Opus 4.7 (1M context) <noreply@anthropic.com> * feat(custom/orfnormalise): per-caller field overrides via ext.args + provenance header The source columns the parser reads for each derived TSV field (`score`, `orf_type`, `aa_length`, length-for-aa-derivation) were previously hardcoded in each parser. Some callers expose multiple meaningful choices - Ribo-TISH has TISPvalue / RiboPvalue / FisherPvalue (plus Q-values), Rp-Bp has bayes_factor_mean / chi_square_p / p_translated_mean - and there was no way to override the default chain. This change: - Plumbs `task.ext.args` through the module to the template. - Adds `--score-field`, `--orf-type-field`, `--length-field`, `--aa-length-field` CLI options. When set, an override replaces the per-caller default chain (so the user gets exactly what they asked for). - Refactors the rpbp parser from positional column indices to csv.DictReader + named RPBP_COLUMNS so users can address columns by name (e.g. `--score-field bayes_factor_var`). - Centralises field-name preference chains into a DEFAULT_FIELDS table next to the parsers. - Writes a `# parser_columns: caller=<X> score=<col> orf_type=<col> ...` provenance header at the top of the normalised TSV, so consumers (and reviewers) can verify which source column was actually read for each derived field. Standard csv.DictReader and pandas `comment='#'` skip the line automatically. Also updates custom/orfmerge's load_normalised to skip `#`-prefixed comment lines so the new provenance header doesn't leak into the cohort catalogue. Test changes: - Existing per-caller tests gain a provenance-line assertion (verifies `# parser_columns:` is present and reports the expected default column for at least the ribocode case). - New `homo_sapiens [chr20] - ribotish - score-field override` test runs with `ext.args = '--score-field RiboPvalue'` and asserts the provenance line reports `score=RiboPvalue` rather than the default `FisherPvalue`. - Tests now filter `#`-prefixed lines before parsing the TSV header. Verified end-to-end on chr19/chr20 fixtures: default + override paths emit the right provenance line, and all parser-specific assertions (non-zero aa_length, populated score column, non-collapsed class distribution, multi-exon BED12 blocks) continue to pass. Co-Authored-By: Claude Opus 4.7 (1M context) <noreply@anthropic.com> * refactor(custom/orfnormalise): caller as val input (not meta.caller); subworkflow channel now carries it as val too nf-core/modules has a hardcoded `permitted_meta_keys = {"id", "single_end"}` allow-list - any other meta.<key> reference in main.nf fails the `main_nf_meta_key` lint check. The previous design embedded the per-record caller id in meta.caller, which violated that policy and tripped CI lint. Lift the caller out of meta and into a proper val input - the pattern used by raxmlng (val(model)), last/mafconvert (val(format)), clair3 (val(platform)), and similar nf-core modules: - Module input: `tuple val(meta), path(orfs_table), val(caller)` plus the unchanged `tuple val(meta2), path(gtf)`. The `caller` value goes into a per-caller enum in meta.yml (one of ribocode / ribotish / ribotricer / rpbp / price). - Module main.nf no longer references meta.caller; the script binds `caller` directly from the input tuple. - Subworkflow `ch_orf_tables` API also moves caller out of meta: `[ val(meta), path(orf_table), val(caller) ]`. Subworkflow appends caller to meta.id locally for the normalise call so the merger's `beds/*` staging gets unique per-caller filenames; the merger then rebuilds meta as `[id: 'cohort']` so no caller leakage downstream. - Subworkflow's bundled `nextflow.config` drops the `withName: CUSTOM_ORFNORMALISE` prefix override (the meta.id append in main.nf does the work). Pre-commit fixes: - ruff E741: rename `l` to `line` in write_outputs (loop variable). - ruff UP015: drop the redundant `"r"` mode arg from open() in open_text. - ruff format: minor reformatting of long literal lists. Tests: - All 6 orfnormalise tests + 1 merger test + 1 subworkflow test pass module input as the new 3-tuple (meta carries only `id`; caller is the trailing val). Merger setup config still constructs `meta + [caller: ...]` for prefix disambiguation - that meta key is consumed by the *test* setup_prefix.config (not by the module), which is fine since lint only scans modules' main.nf. Local `nf-core modules lint custom/orfnormalise` now passes cleanly (only the cosmetic Wave container-version warnings remain, shared with the reference custom/bed12codonpositions module). Co-Authored-By: Claude Opus 4.7 (1M context) <noreply@anthropic.com> * fix(custom/orfnormalise+orfmerge): prefix convention, args, ribotricer score [skip ci] - orfnormalise prefix: follow #11955 convention — descriptor into default prefix (`${meta.id}.normalised`) and drop it from output glob (`${prefix}.bed12`, `${prefix}.tsv`). Default output filenames are unchanged; users who override ext.prefix still get a correctly-named file. - Both templates: replace `"${args}".split() if...` with `shlex.split("${args}")` (consistent with rpbp chain modules; handles quoted arg values; empty string returns [] without special-casing). - orfnormalise/parse_ribotricer: fix bug where phase_score was read for the BED score column but `raw_score` was always stored as "". All other parsers correctly store the native score; ribotricer rows will now populate `score_ribotricer` in the merged catalogue. Ribotricer test snapshot dropped for regeneration; score assertion added to the ribotricer test to match the other caller tests. - Fix meta.caller -> caller val wording in orfnormalise meta.yml tool description, orfnormalise.py error message, and subworkflow meta.yml + main.nf comment. Co-Authored-By: Claude Sonnet 4.6 <noreply@anthropic.com> * simplify Python: replace x==x NaN idiom and csv.DictReader comment-skip [skip ci] Co-Authored-By: Claude Sonnet 4.6 <noreply@anthropic.com> * refactor(custom/orf): drop pandas, simplify orfnormalise/orfmerge templates Strip pandas from both bundled scripts: it was used only as a typed list-of-dicts container (build DataFrame then iterrows / read_csv then to_dict). Parsers now return plain lists of dicts and orfmerge reads TSVs with csv.DictReader. Container drops to python + pyyaml. orfnormalise: replace five near-identical classifier functions with a keyword-rule table, factor shared aa-length and column-resolution helpers, unify ribotish/ribotricer interval parsing, add a clamp helper, and remove the unused SCORE_DIRECTIONS map. orfmerge: collapse the two transcript-grouping functions into one group_by helper. Stubs write versions.yml via a one-line python -c using pyyaml. Module output is byte-identical: all bed12/tsv/catalogue/fasta snapshot md5s are unchanged; only versions.yml md5s differ (pandas removed). Adds a previously-missing ribotricer snapshot entry. Co-Authored-By: Claude Opus 4.8 (1M context) <noreply@anthropic.com> * feat(custom/orfmerge): record cross-sample recurrence in the catalogue The normalised TSV's sample_id column was previously written by orfnormalise but never consumed. orfmerge now aggregates it per cluster into two catalogue columns, mirroring the existing per-caller consensus: - n_samples: number of distinct samples contributing to the cluster. - samples: sorted, comma-separated list of those sample ids. This gives a cohort recurrence metric for filtering the catalogue. The orfmerge test now merges a second sample (reusing the ribotish fixture) so recurrence (n_samples == 2) is exercised and asserted. Co-Authored-By: Claude Opus 4.8 (1M context) <noreply@anthropic.com> --------- Co-authored-by: Claude Opus 4.7 (1M context) <noreply@anthropic.com>
1 parent 457a7a8 commit ccb13a7

19 files changed

Lines changed: 2692 additions & 0 deletions

File tree

Lines changed: 8 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,8 @@
1+
---
2+
# yaml-language-server: $schema=https://raw.githubusercontent.com/nf-core/modules/master/modules/environment-schema.json
3+
channels:
4+
- conda-forge
5+
- bioconda
6+
dependencies:
7+
- conda-forge::python=3.12.11
8+
- conda-forge::pyyaml=6.0.2
Lines changed: 38 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,38 @@
1+
process CUSTOM_ORFMERGE {
2+
tag "$meta.id"
3+
label 'process_medium'
4+
5+
conda "${moduleDir}/environment.yml"
6+
container "${ workflow.containerEngine in ['singularity', 'apptainer'] && !task.ext.singularity_pull_docker_container ?
7+
'https://community-cr-prod.seqera.io/docker/registry/v2/blobs/sha256/15/15f690593e9d2bd7ebeb72fff4068eeae9cdfc909103de457b81a25d29cfbef3/data' :
8+
'community.wave.seqera.io/library/python_pyyaml:25b0f22c7e7bf847' }"
9+
10+
input:
11+
tuple val(meta), path(bed12s, arity: '1..*', stageAs: 'beds/*'), path(tsvs, arity: '1..*', stageAs: 'tsvs/*')
12+
13+
output:
14+
tuple val(meta), path("${prefix}.catalogue.bed12") , emit: bed12
15+
tuple val(meta), path("${prefix}.catalogue.tsv") , emit: catalogue_tsv
16+
tuple val(meta), path("${prefix}.orf_to_gene.tsv") , emit: orf_to_gene_tsv
17+
tuple val(meta), path("${prefix}.catalogue.mqc.tsv") , emit: multiqc
18+
path "versions.yml" , emit: versions, topic: versions
19+
20+
when:
21+
task.ext.when == null || task.ext.when
22+
23+
script:
24+
prefix = task.ext.prefix ?: "${meta.id}"
25+
args = task.ext.args ?: ''
26+
template 'orfmerge.py'
27+
28+
stub:
29+
prefix = task.ext.prefix ?: "${meta.id}"
30+
"""
31+
touch ${prefix}.catalogue.bed12
32+
touch ${prefix}.catalogue.tsv
33+
touch ${prefix}.orf_to_gene.tsv
34+
touch ${prefix}.catalogue.mqc.tsv
35+
36+
python -c "import platform, yaml; yaml.safe_dump({'${task.process}': {'python': platform.python_version()}}, open('versions.yml', 'w'), default_flow_style=False, sort_keys=False)"
37+
"""
38+
}
Lines changed: 154 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,154 @@
1+
# yaml-language-server: $schema=https://raw.githubusercontent.com/nf-core/modules/master/modules/meta-schema.json
2+
name: "custom_orfmerge"
3+
description: |
4+
Cluster normalised per-sample, per-caller ORF predictions into a single
5+
cohort-level catalogue. Pair with `custom/orfnormalise` upstream and
6+
(typically) `bedtools/getfasta` + `seqkit/translate` downstream to obtain
7+
the AA FASTA.
8+
9+
Strategy is class-aware (operating on the harmonised `orf_class` written
10+
by `custom/orfnormalise`):
11+
12+
- canonical_cds: collapse by (transcript_id, strand). One
13+
canonical CDS per transcript by definition.
14+
- uORF, dORF, other: collapse by (transcript_id, strand, start,
15+
end). A single transcript can host multiple
16+
distinct uORFs / dORFs / internal ORFs, so
17+
keying on the outer span keeps them in
18+
separate clusters while still merging
19+
cross-caller calls that agree on coordinates.
20+
- novel_u, smORF: greedy reciprocal-overlap clustering on the
21+
outer genomic span at `--reciprocal-overlap`
22+
(default 0.8). Catches fuzzy cross-caller
23+
matches and exact-coordinate collapses in
24+
one pass. Order-dependent at the boundary:
25+
a chain A-B-C where A-B and B-C overlap at
26+
~0.85 but A-C only at ~0.75 may cluster as
27+
{A,B,C} or {A,B}+{C} depending on iteration
28+
order. Rare in practice at 0.8.
29+
30+
Cross-caller consensus is recorded in two column families on the
31+
catalogue TSV:
32+
33+
- `called_by_<caller>`: 0/1 indicator per supported caller
34+
(ribotish, ribocode, ribotricer, rpbp, price).
35+
- `score_<caller>`: best score from that caller within the cluster.
36+
Score direction is per-caller (p-values are
37+
minimised; Bayes factors / phase scores are
38+
maximised).
39+
40+
Cross-sample recurrence is recorded in two further columns:
41+
42+
- `n_samples`: number of distinct samples contributing to the
43+
cluster (a cohort recurrence metric).
44+
- `samples`: sorted, comma-separated list of those sample ids.
45+
46+
Emits a small MultiQC custom-content TSV (per-class counts) for
47+
inclusion in downstream MultiQC reports.
48+
keywords:
49+
- orf
50+
- ribo-seq
51+
- catalogue
52+
- merge
53+
- clustering
54+
tools:
55+
- "orfmerge":
56+
description: |
57+
Python helper that clusters normalised ORF BED12+TSV pairs across
58+
callers and samples into one unified catalogue, recording per-caller
59+
provenance and best score in the output table.
60+
tool_dev_url: "https://github.com/nf-core/modules/blob/master/modules/nf-core/custom/orfmerge/main.nf"
61+
licence: ["MIT"]
62+
identifier: ""
63+
input:
64+
- - meta:
65+
type: map
66+
description: |
67+
Groovy Map. Typically a cohort-level `[ id: 'allsamples' ]`.
68+
- bed12s:
69+
type: file
70+
description: |
71+
One or more normalised BED12 files (output of `custom/orfnormalise`).
72+
Files are co-indexed with `tsvs` only via shared `orf_id`s; per-row
73+
alignment is not required.
74+
pattern: "*.bed12"
75+
ontologies:
76+
- edam: http://edamontology.org/format_3586 # BED
77+
- tsvs:
78+
type: file
79+
description: |
80+
One or more normalised sidecar TSV files (output of
81+
`custom/orfnormalise`).
82+
pattern: "*.tsv"
83+
ontologies:
84+
- edam: http://edamontology.org/format_3475 # TSV
85+
output:
86+
bed12:
87+
- - meta:
88+
type: map
89+
description: |
90+
Groovy Map matching the input meta.
91+
- ${prefix}.catalogue.bed12:
92+
type: file
93+
description: Merged ORF catalogue as BED12 (genomic blocks).
94+
pattern: "*.catalogue.bed12"
95+
ontologies:
96+
- edam: http://edamontology.org/format_3586 # BED
97+
catalogue_tsv:
98+
- - meta:
99+
type: map
100+
description: |
101+
Groovy Map matching the input meta.
102+
- ${prefix}.catalogue.tsv:
103+
type: file
104+
description: |
105+
Per-ORF table with `called_by_<caller>` and `score_<caller>`
106+
columns for each supported caller, `n_samples` / `samples`
107+
cross-sample recurrence columns, plus orf_class, aa_length,
108+
and host gene_id / transcript_id.
109+
pattern: "*.catalogue.tsv"
110+
ontologies:
111+
- edam: http://edamontology.org/format_3475 # TSV
112+
orf_to_gene_tsv:
113+
- - meta:
114+
type: map
115+
description: |
116+
Groovy Map matching the input meta.
117+
- ${prefix}.orf_to_gene.tsv:
118+
type: file
119+
description: |
120+
One row per (orf_id, gene_id, transcript_id) mapping. An ORF can
121+
map to multiple host transcripts/genes when callers picked
122+
different annotated isoforms; downstream gene-level aggregation
123+
collapses these.
124+
pattern: "*.orf_to_gene.tsv"
125+
ontologies:
126+
- edam: http://edamontology.org/format_3475 # TSV
127+
multiqc:
128+
- - meta:
129+
type: map
130+
description: |
131+
Groovy Map matching the input meta.
132+
- ${prefix}.catalogue.mqc.tsv:
133+
type: file
134+
description: |
135+
MultiQC custom-content TSV (per-class ORF counts).
136+
pattern: "*.catalogue.mqc.tsv"
137+
ontologies:
138+
- edam: http://edamontology.org/format_3475 # TSV
139+
versions:
140+
- versions.yml:
141+
type: file
142+
description: File containing software versions
143+
pattern: "versions.yml"
144+
ontologies:
145+
- edam: http://edamontology.org/format_3750 # YAML
146+
topics:
147+
versions:
148+
- versions.yml:
149+
type: string
150+
description: The name of the process
151+
authors:
152+
- "@pinin4fjords"
153+
maintainers:
154+
- "@pinin4fjords"

0 commit comments

Comments
 (0)