Skip to content

hmmsearch -T prefilter is global (shared HMM module) — loosens VOG/DRAM-v and sulfur/metals; make it per-DB #13

Description

@tpall

Summary

PR #12 switched the HMM prefilter from -E ${e_value} to a fixed bit-score floor -T 10 --domT 10. That change lives in the shared HMM_SEARCH process (modules/local/annotate/hmmsearch.nf), so it applies to every HMM search — kofam, sulfur, vog, canthyd, fegenie, metals, camper-hmm — not only the pooled kofam/sulfur searches the PR targeted. The per-DB e_value inputs (e.g. vog_e_value = 1e-10) are now silently ignored.

-T is required for the pooled path (bit scores are search-space/Z-independent, so pooled == per-genome), so a global revert to -E isn't the answer — the prefilter needs to be per-DB.

Why it matters (benign for some DBs, not others)

Whether the swap is harmless depends on the downstream filtering in bin/hmm_parser.py:

  • kofam (plus canthyd via A/B rank, fegenie via its cutoffs) carry per-profile bit-score thresholds, so hmm_parser keeps score >= threshold and the -T 10 prefilter is a near no-op — the "+0.08%" the PR reported. ✅
  • vog, sulfur, metals pass the dummy sheet (bin/assets/forms/distill_sheets/dummy.tsv, contents NULL), so hmm_parser finds no threshold and falls through to only perc_cov >= 0.35. Their real significance gate was the hmmsearch prefilter itself. For VOG that was -E 1e-10 (stringent); -T 10 (bit score ≥ 10) is far looser. ⚠️

Measured impact (DRAM-v mode)

Ran DRAM-v on a 161-contig viral voTU catalog (geNomad gene input, catalog mode, --max_auxiliary_score 3, --use_dramv --use_kofam --use_merops) on two refs differing only by this prefilter line:

metric -E 1e-10 -T 10 change
VOG-annotated genes 2871 5872 +105% (doubles)
V-flag genes 999 1737 +74%
auxiliary_score histogram {1: 6848, 4: 1416, 5: 1212} identical 0
AMG calls (M & aux ≤ 3) 294 294 0

Raw VOG annotations and the V flag roughly double, but the AMG endpoint is unchangedauxiliary_score is geNomad/VirSorter-derived (independent of VOG), and the AMG gate keys off amg_flags (M, lacking A/P/T/N) + auxiliary_score, none of which depend on VOG.

Net: AMG-identification output is safe, but any downstream use of raw VOG annotations / V-flag counts / "% viral genes" is now ~2× inflated with low-confidence hits.

Proposed fix

Make the prefilter per-DB instead of hardcoded in the shared module. HMM_SEARCH already reads def args = task.ext.args ?: "" but never uses it in the command — wire the prefilter through that (or an explicit process input), so:

  • pooled kofam/sulfur searches use -T 10 --domT 10 (needed for Z-independent pooling; provably benign there),
  • per-genome vog/sulfur/metals keep -E ${e_value} (e.g. VOG -E 1e-10).

This preserves PR #12's pooling-correctness win while restoring the intended significance filter for the threshold-less DBs.

Refs: modules/local/annotate/hmmsearch.nf, bin/hmm_parser.py, docs/dev/search-bundling.md.


🤖 Filed with Claude Code after an -E-vs--T DRAM-v validation run.

Metadata

Metadata

Assignees

No one assigned

    Labels

    No labels
    No labels

    Projects

    No projects

    Milestone

    No milestone

    Relationships

    None yet

    Development

    No branches or pull requests

    Issue actions