Reporting faithfulness: audit items B1-B11 (plan B)

CHANGELOG.md

@@ -6,8 +6,20 @@ All notable changes to Tessera are recorded here. The format follows
 
 ## [Unreleased]
 
+Fixes from the post-1.2.0 audit
+(`docs/superpowers/specs/2026-10-01-post-1.2.0-audit-design.md`, items A1-A10 and B1-B11).
+For a given alignment the callers find the same regions, with the same coordinates and
+p-values, as before. Three things around them do change: the `donor_undercovered` flag and
+donor-absent rows where a coverage gap is now recognised as a breakpoint artefact; runs
+that select the barcode caller on a panel where it cannot run (it no longer counts toward
+`--min-methods`, and a barcode-only run is refused); and what reaches the scan -- an
+alignment built with the MAFFT or MAF-based backends, or a panel built with `--curate`,
+may now be a different (corrected) one.
+
 ### Fixed
 
+Data safety and alignment:
+
 - **The MAFFT backend ignored strand.** A genome, or one contig of a draft assembly, on the
   opposite strand to the backbone was aligned as given and came out at chance-level identity
   (about 0.40 against 0.97 for the same genome in forward orientation), which the scan then
@@ -63,6 +75,77 @@ All notable changes to Tessera are recorded here. The format follows
   reference and an empty MAFFT result are reported as input/output errors that name the
   file, instead of "Unexpected error".
 
+Reporting:
+
+- **The run provenance named the wrong callers.** MaxChi, Bootscan, GENECONV and the barcode
+  caller were each described as `heuristic (min ... / margin ... / merge ...)` in
+  `run_provenance.json` and the report, so a default run recorded
+  `hmm + 3seq + heuristic + heuristic`. Each caller is now described under its own name, and
+  the record gains the settings that change what is reported: the agreement gate
+  (`--min-methods`), sibling exclusion, lineage clustering and donor re-attribution.
+- **A clean recombinant between divergent parents was reported as a possible missing
+  reference.** A window straddling a breakpoint matches neither parent well on its own, so
+  its best similarity fell below the coverage threshold, the stretch was called a `divergent`
+  coverage gap, and the region was marked `donor_undercovered` -- on the shipped
+  `divergent` example the headline read "low confidence" for a donor identical to the query.
+  A gap within one window of a called region boundary is now labelled `breakpoint` in
+  `coverage_gaps.tsv` and the report when every under-threshold window in it is matched,
+  to the coverage threshold, by a single switch between the region's two parents. Such a
+  gap does not caveat the region, is not turned into a donor-absent region, and is left out
+  of the headline. A stretch from a source outside the panel stays `divergent` even when it
+  is short and sits beside a breakpoint. **`donor_undercovered` and the confidence wording
+  change for regions whose only gaps were breakpoint artefacts.** One consequence is not yet
+  measured on the hybrid harness: a region whose donor is a close stand-in for an absent
+  lineage (above the coverage threshold) was previously caveated only through its
+  breakpoint gaps, and is now reported without a caveat.
+  Reference recruitment (`fill-references`, `find-references`) is unaffected.
+- The report's "covering N kb (P %) of the query" added the lengths of overlapping regions
+  that name different donors; it now reports their union.
+- After `--reattribute-donors`, `recombination_methods.tsv` and the report's method table
+  kept the donor from before re-attribution.
+- `--lineage-map` pointing at a file that does not exist was ignored (exit 0, untyped
+  report) by `recomb`, `type-lineages`, `detect`, `fill-references` and `build-panel`. It
+  is now an error.
+- **The barcode caller on an untyped panel read as a negative.** It named the first record
+  of the alignment as major parent and its column in the method table said `no`. It now
+  reports no major parent; in an ensemble it is logged and shown as `not run`, the
+  agreement gate counts only callers that ran (a gate lowered for that reason is logged and
+  recorded as, for example, `1 (requested 2; 1 caller ran)`), and a run that selected only
+  `barcode` is refused. The report and provenance give the actual reason, which on a typed
+  panel is that fewer than two clades carry enough markers.
+- The report judged the PHI p-value at alpha 0.05 whatever `--alpha` was.
+- **A PHI test that could not reject was reported as "no signal".** With no more
+  informative sites than the window can hold (`--phi-window`, default 100) the permutation
+  p-value is 1 for any data. This is now reported as `not testable` (`NA` in the
+  `recombination_profile.tsv` header, `phi_p = None` in the API).
+- Plots labelled a donor-absent region "recombinant: <backbone>"; `similarity_pair` showed
+  the two leading window winners rather than the major parent and the leading donor; with
+  `--top-n 1` the donor was drawn grey.
+- The report footer listed `.pdf` plots under `--plot-format png` and files that were not
+  written; the methods text and references described a two-caller ensemble; `--method`
+  help omitted `geneconv`.
+- `sibeliaz` was probed with `-v`, which it rejects, and the error line was recorded as the
+  aligner version. A failed probe is no longer recorded as a version.
+- Input checks: `find-references --msa <missing>` and `recomb -o <existing file>` failed
+  with "Unexpected error"; `--max-rounds 0` exited 0 having built nothing; `reassort`
+  accepted `--ani-floor 500`, `--margin -3` and a `--dataset` key matching no segment;
+  `type-lineages` rejected `.fasta.gz` / `.fas` collections (it now reads what the other
+  commands read, and like them rejects a file in the collection that is not FASTA, where it
+  used to skip it); a multi-line aligner error broke `segment_scan.tsv`.
+
+### Changed
+
+- `RecombinationSignal.phi_p` is `float | None`.
+- `seaborn` is no longer a dependency; nothing imported it.
+
+### Documentation
+
+- `docs/detection-methods.md` describes lineage clustering (including its limit on panels
+  below about 1.5 % divergence), the `breakpoint` coverage kind and the untestable-PHI case.
+- `validation/README.md` states that the specificity harness defaults to `--min-methods 2`
+  while the CLI defaults to 1, gives the measured rate at each, and records that
+  `hcv_clonal_1b` currently fails.
+
 ## [1.2.0] - 2026-10-01
 
 ### Fixed

docs/detection-methods.md

@@ -78,6 +78,23 @@ Recombination regions (major parent: cowpox_KC813504):
   variola       cowpox_KC813504  66268        147150     80882       0.999      0.977      0.97     2e-300   66768
 ```
 
+### Near-duplicate references (lineage clustering)
+
+A recruited panel often holds several near-identical genomes of one lineage. Competed
+individually they tie in every window and fragment the call, so the HMM caller first
+pools them (`--cluster-lineages`, on by default; hmm only). Two references are pooled
+when their identity stays at or above 98.5 % in every window, with no region-sized run
+below it. The pooled lineage competes as one state under the label of its best-covering
+member, and a region names that member. Clustering is skipped for panels of fewer than 4
+or more than 200 references.
+
+One limitation follows from the absolute threshold. On a panel in which *every* pair of
+references is at least 98.5 % identical in every window -- mpox, VZV and other sets below
+roughly 1.5 % divergence -- all references pool into a single lineage and the HMM has
+nothing left to compete, so it calls no region whatever the data. The site-based callers
+(3SEQ, MaxChi) are unaffected. On such a panel, run with `--no-cluster-lineages` to keep
+the HMM's vote. `run_provenance.json` records whether clustering was on.
+
 ### Low-divergence panels (intra-species sets, DNA viruses)
 
 When the references are nearly identical -- e.g. mpox clades (~0.5 %), VZV (~0.2 %),
@@ -181,8 +198,10 @@ and scans the query window by window for the clade whose local markers it carrie
 where a non-major clade dominates is a region attributed to that lineage. A marker is
 denoised across the clade's members and the comparison is multi-way, so the call is robust
 to a single near-identical adjacent-clade genome winning by chance -- the failure mode of
-the genome-level callers at low divergence. It is silent on untyped panels, and opt-in
-(needs typed references). It composes with `--pool-consensus` but needs only a typed panel.
+the genome-level callers at low divergence. It is opt-in and needs typed references. On an
+untyped panel it cannot run, which is not the same as finding nothing: in an ensemble it is
+logged and shown as `not run` in the method comparison, and a run that selected only
+`barcode` is refused. It composes with `--pool-consensus` but needs only a typed panel.
 
 ## Parent-free recombination signal (PHI + Rmin)
 
@@ -232,6 +251,17 @@ localizes the signal rather than averaging it away -- are the more informative
 parent-free outputs. The diagnostic runs for every `--method`; disable with
 `--no-phi`, or widen its window with `--phi-window`.
 
+The PHI p-value is judged at the run's `--alpha`, in the report and for the per-region
+flag alike. The test is **not testable** when the alignment has too few informative sites
+for the window (`--phi-window`, default 100 site ranks; at least window + 2 sites are
+needed): every pair of sites then falls inside one window, so reordering the sites cannot
+change the statistic and the permutation p-value would be 1 whatever the data. Tessera
+reports this as `not testable` (`NA` in the `recombination_profile.tsv` header) rather
+than as a non-significant result, and no region is flagged `parent_free_support`; lower
+`--phi-window` to test such an alignment. Just above that limit the test is defined but
+has little power: the shipped `divergent` example (159 informative sites) gives p = 1 at
+the default window and p = 0.001 at `--phi-window 20`.
+
 This all remains an **indicative screen**: the built-in HMM and 3SEQ tests are fast
 triplet/segmentation screens, not a full tree-based analysis (such as GARD). Treat
 regions as candidates to confirm.
@@ -276,16 +306,16 @@ The practical readings:
 |---|---|
 | `report.html` | Self-contained report: run provenance, the region table, the per-dataset stats, and an embedded interactive plot |
 | `recombination_regions.tsv` | Called regions: minor/major parent, start/end in **both MSA columns and query bases**, `length_bp` / `length_msa`, `support`, `pvalue` / `qvalue` with the `test` and `statistic` that produced them, mean similarities, the calling `methods`, and `parent_free_support` |
-| `recombination_methods.tsv` | Ensemble breakdown (only when several methods run): one row per region with a Y/n per method and the parent-free flag |
-| `recombination_profile.tsv` | Parent-free signal: header with the PHI p-value and Rmin, then per-informative-site local incompatibility (the PHI profile) |
+| `recombination_methods.tsv` | Ensemble breakdown (only when several methods run): one row per region with `yes` / `no` per method (`not run` for a selected caller that could not run) and the parent-free flag |
+| `recombination_profile.tsv` | Parent-free signal: header with the PHI p-value (`NA` when not testable) and Rmin, then per-informative-site local incompatibility (the PHI profile) |
 | `similarity_windows.tsv` | Full per-window matrix: `msa_position`, `query_position`, `winner`, and one similarity column per dataset. Always base-pair windows (identity over all comparable columns) |
 | `informative_site_windows.tsv` | Only under informative-site windowing: the windows the HMM segmented -- `msa_position`, `query_position`, the window's `msa_start` / `msa_end` (end-exclusive), `winner`, and per dataset the query's identity **at polymorphic columns only**. Not comparable with similarity or ANI |
 | `similarity_stats.tsv` | Per-dataset similarity statistics (median, windows above identity thresholds) |
 | `window_winners.tsv` | Per-dataset count of windows won (ties included) |
-| `coverage_gaps.tsv` | Stretches where even the closest reference is a poor match -- possible missing references |
+| `coverage_gaps.tsv` | Stretches where even the closest reference is below the best-similarity threshold, with a `kind`: `divergent` (the query is far from every reference -- a possible missing reference), `low_information` (too few comparable bases to judge), or `breakpoint` (windows straddling a called breakpoint: each under-threshold window is matched, to the threshold, by a single switch between the region's two parents; not a missing reference, and it does not caveat the region) |
 | `similarity_top{N}.{fmt}` | Static plot of the nearest `--top-n` datasets, called regions shaded |
 | `informative_sites_top{N}.{fmt}` | Only under informative-site windowing: the same plot for identity at informative sites |
-| `similarity_pair.{fmt}` | Static plot of the major vs leading minor parent, region shaded |
+| `similarity_pair.{fmt}` | Static plot of the major parent against the donor of the longest called region, regions shaded (the two leading window winners when no donor was called) |
 | `run_provenance.json` | Machine-readable record of the run: Tessera version, parameters, caller description, and -- when the alignment came from `tessera msa` -- the aligner, its version and its arguments |
 
 ### Reading `pvalue`, `support` and the lengths

docs/reference-panels.md

@@ -244,9 +244,10 @@ masking sibling (e.g. SARS-CoV-2 sublineages), add `--seed-keep-siblings`. NCBI
 By default (`--auto-diversify`), BLAST seeding will **switch to the `ncbi-virus`
 diversity path automatically** when it finds only siblings -- i.e. when the query's
 lineage saturates `nt` and no parental lineage can be recruited by similarity. A broad
-fetch is capped (`--fetch-limit`, default 2000) and dereplicated; for a heavily
-sequenced taxon the capped sample may miss lineages, so a curated `--candidate-pool`
-is recommended (and the run says so). Whether the diversity panel actually contains the
+fetch is not truncated: the whole set is downloaded and dereplicated locally, and the run
+logs a notice when it exceeds `--fetch-limit` (default 2000) because that step can take
+a few minutes. For a heavily sequenced taxon a curated `--candidate-pool` is the faster
+route. Whether the diversity panel actually contains the
 parents depends on the taxon: it works when they are genotype/lineage representatives,
 less so for fine genotype-specific recombinants. Disable with `--no-auto-diversify`.
 This complements the caller-side defence: `tessera recomb` excludes whole-genome

docs/superpowers/specs/2026-10-01-post-1.2.0-audit-design.md

@@ -327,6 +327,22 @@ this section disagree with the text above, the plan and this section hold.
   7/12, positive control 3/3, 55 bp); `run_validation.py` 6 PASS, 1 FAIL
   (`hcv_clonal_1b`), 1 SKIP; the divergent example reports high confidence.
 
+**Plan B, after its whole-branch review** (fix commit `0c9daf2`)
+
+- *B2 tightened again:* "the two parents together explain the gap" was measured over the
+  whole gap, where the well-matching flanks dilute a short stretch from a source outside
+  the panel (150 bp at 15 % divergence beside a breakpoint was relabelled). Each
+  under-threshold window in the gap must now be matched, to the threshold, by the best
+  single-switch mosaic of the two parents.
+- *B2, unmeasured:* a region whose donor is a close stand-in for an absent lineage (above
+  the coverage threshold) used to be caveated only through its breakpoint gaps and is now
+  uncaveated. `run_hybrids.py` scores `panel_donor_absent` cases on that caveat and was
+  not run; it is a gate before merge.
+- *B6:* the reason a caller did not run is carried into the provenance and report; a
+  `--min-methods` gate lowered because a selected caller could not run is logged and
+  recorded with the requested value.
+- The changelog's "no region call changes" now lists its exceptions.
+
 **Plan C**
 
 - *C2:* the sign test stays on the HMM segment and only the reported span is trimmed.

example_data/README.md

@@ -16,9 +16,10 @@ tessera recomb --msa divergent.msa.fasta --query query --output out_divergent \
     --window-size 300 --window-step 30
 ```
 
-Both callers call `parent_B` over the insert (q-value ~1e-29) with a sharp breakpoint,
+The four default callers all call `parent_B` over the insert with a sharp breakpoint,
 so the region is flagged as agreeing (high confidence); the similarity plot shows an
-obvious crossover.
+obvious crossover. The two short stretches listed under reference coverage are the
+windows straddling the breakpoints (kind `breakpoint`), not missing references.
 
 ## `cryptic_insert.msa.fasta` -- why the ensemble exists
 
@@ -31,12 +32,14 @@ tessera recomb --msa cryptic_insert.msa.fasta --query query --output out_hmm \
     --window-size 1000 --window-step 100 --method hmm     # 0 regions
 ```
 
-The default ensemble also runs 3SEQ, which pools the discriminating sites into an exact
-triplet test and recovers the event (q-value ~1e-12, `methods` = 3seq):
+The default ensemble also runs the site-based callers, which pool the discriminating
+sites into triplet tests and recover the event (q-value ~1e-12, `methods` = 3seq,maxchi):
 
 ```
 tessera recomb --msa cryptic_insert.msa.fasta --query query --output out_cryptic \
     --window-size 1000 --window-step 100                  # finds parent_B insert
 ```
 
-Both runs also report the parent-free PHI / Rmin signal in `recombination_profile.tsv`.
+Both runs also write the parent-free PHI / Rmin signal to `recombination_profile.tsv`.
+The cryptic example has 54 informative sites, fewer than the default PHI window of 100, so
+its PHI test is reported as `not testable`; add `--phi-window 5` to test it.

pyproject.toml

@@ -35,7 +35,6 @@ dependencies = [
     "numpy>=1.26,<3.0",
     "pandas>=2.0,<4.0",
     "matplotlib>=3.8,<4.0",
-    "seaborn>=0.13,<1.0",
     "plotly>=5.18,<7.0",
     "biopython>=1.83,<2.0",
 ]

src/tessera/aligners/sibeliaz.py

@@ -47,7 +47,9 @@ class SibeliazAligner(Aligner):
     capabilities = ToolCapabilities(
         name="sibeliaz",
         conda=("bioconda::sibeliaz",),
-        required_binaries=(BinarySpec("sibeliaz", version_args=("-v",)),),
+        # The sibeliaz wrapper has no version option ("-v" is rejected as an illegal
+        # option), so there is nothing to probe; its version is recorded as unknown.
+        required_binaries=(BinarySpec("sibeliaz", version_args=None),),
         recommended_max_genomes=2000,
         threads_param="-t",
     )

src/tessera/cli/cmd_build_panel.py

@@ -13,7 +13,16 @@
 
 import typer
 
-from .main import _require_choice, _require_file, app, get_logger, stage_errors
+from .main import (
+    _require_choice,
+    _require_file,
+    _require_lineage_map,
+    _require_range,
+    _require_scan_windows,
+    app,
+    get_logger,
+    stage_errors,
+)
 
 
 def _seed_source(candidate_pool, nextclade: bool, nextclade_dataset: str | None) -> str:
@@ -84,6 +93,12 @@ def build_panel(
     logger = get_logger(output)
     with stage_errors(logger):
         _require_file(query, "Query file")
+        _require_lineage_map(lineage_map)
+        # Checked here, before any network or aligner work: a round count of zero
+        # builds nothing and still exits 0, and a bad window otherwise surfaces only
+        # after the panel has been recruited and aligned.
+        _require_range(max_rounds, "--max-rounds", lo=1)
+        _require_scan_windows(window_size, window_step)
         _require_choice(aligner, set(aligner_registry.names()), "--aligner")
         params = FillParams(
             query=query, collection=collection, output=output,