feat(concordance): Nelson-Ladiges "trit" currency and chance correction for QuartetConcordance - #276
Merged
Merged
Conversation
…ance
Add `unit = c("quartet", "trit")` to QuartetConcordance(). "quartet"
(default) is unchanged; "trit" scores in Nelson & Ladiges (1992)
redundancy-corrected currency, where a split of sizes (k, t-k) carries
only (k-1)(t-k-1) independent quartets (4-cycles of K_{k,t-k}; the N&L
entailment is GF(2) cycle addition).
Scoring uses coverage ("option b"): per state-pair, concordant trits
A over the reported unit's own content (Wk for edges, Wc for chars),
pooled by shared information M = min(Wc, Wk). Only a split *displayed*
by the character scores 1; nested/crossing get strict partial credit.
Multistate is handled in the same currency, not deferred: a quartet is
decisive only within a state-pair, those K_{n_i,n_j} blocks are
edge-disjoint, and the entailment never crosses them, so trits add over
pairs (W = sum (n_i-1)(n_j-1), weight 4/(n_i n_j)). Verified by GF(2)
rank of the decisive/concordant sets on 2/3/4-state characters. Missing
/ ambiguous / inapplicable tokens drop per character by construction.
Validation (dev/benchmarks/frac-quart/): cell extraction cross-checked
against the C++ kernel; A <= min(Wc, Wk) everywhere; spec ladder
reproduced; package == an independent re-implementation on binary,
multistate and missing-data inputs across all weight x return combos.
Co-Authored-By: Claude Opus 4.8 <noreply@anthropic.com>
Add a random-expectation baseline to QuartetConcordance(unit = "trit"), mirroring ClusteringConcordance()'s normalize interface: normalize = FALSE (default) no correction; published measure unchanged normalize = TRUE exact expected concordance under a fixed-marginal null normalize = <int> Monte-Carlo tip-shuffle estimate of that expectation Null model: reassign each character's tokens across the leaves at random, holding state counts and split sizes fixed (the multivariate-hypergeometric confusion table; the same null as ClusteringConcordance's miRand and Consistency's ExpectedLength). This is the only coherent null for the trit currency -- it re-zeros each (split, char, state-pair) against its own combinatorial floor and extends unchanged to multistate + missing data. The flat 1/3 SCF baseline is a raw-quartet heuristic that does not survive the trit floors; it is documented as the rejected alternative. Implementation: only A (and, for multistate pairs where a pair's side-A count is random, wk and M = min(wc,wk)) vary under the null, so the exact estimator accumulates E[m], E[m*A/wk], E[m*A/wc] from the hypergeometric pmf (cached on (n_i,n_j,M,t)) and re-zeros the pooled edge/char ratio against the pooled expectation, with the same weight aggregation as the observed score. .Rezero guards the z->1 blow-up (-> NA) and returns values unclamped; .Rezero(1,z) = 1 so a displayed split still scores exactly 1 (ceiling invariant preserved), while conflicting characters go negative. Chance correction for the raw unit="quartet" currency is deferred (it needs the per-state-pair cells the C++ kernel does not expose, and a decision on its null) -- normalize errors there for now. Validated: exact vs MC tip-shuffle agree within MC error at the cell level (dev/benchmarks/frac-quart/chance_baseline.R) and end-to-end (test-Concordance.R). Design settled in dev/plans/frac-quart-weight-agreement.md section 7. Co-Authored-By: Claude Opus 4.8 <noreply@anthropic.com>
…istate
Two gaps in the unit="trit" chance correction, both in the multistate path:
- weight=FALSE re-zero averaged observed and baseline over DIFFERENT
(split, char) cells: observed drops a cell when denM==0 (na.rm), but the
baseline kept it whenever baseDenM>0 (a degenerate multistate pair can have
observed wk=0 yet E[m]>0). Mask both to NA on the union of NA cells so the
edge/char row-means cover the same population before .Rezero. weight=TRUE
(ratio-of-sums) is unaffected.
- The exact-vs-MC end-to-end test used only binary characters, so the
random-wk multistate regime -- the reason E[A] was elaborated into
E[m*A/wk] + E[m] -- had no oracle. Add a 3-state column and cover
edge/char x weight in {TRUE, FALSE}, asserting exact and MC agree on which
cells are NA (guards the cell-matching fix) and within MC error elsewhere.
Co-Authored-By: Claude Opus 4.8 <noreply@anthropic.com>
…tet" Extend the chance-correction option to the raw quartet currency, using the same fixed-marginal (hypergeometric) null as the trit path. Per state-pair the raw concordant/decisive counts are polynomials in the 2x2 cells (no floors), so the exact expectation E[conc], E[dec] is a clean sum over the same trivariate-hypergeometric pmf (.ExpectedQuartet / .QuartetExpect); MC (.QuartetMC) reshuffles tokens and re-scores through the C++ kernel. The pooled conc/dec ratio is re-zeroed against E[conc]/E[dec] (weighted) or with the same NA-cell matching as the trit weight=FALSE path (unweighted). Removes the guard that errored on unit="quartet", normalize != FALSE. .Rezero(1, z) = 1, so a displayed split still scores exactly 1; conflicting characters go negative. normalize = FALSE is byte-identical to the published raw measure. Validated exact vs MC tip-shuffle across return x weight x binary/multistate (test-Concordance.R); all Concordance tests green. Design note in dev/plans/frac-quart-weight-agreement.md section 7 updated (raw unit DONE). This is the R reference for the forthcoming C++ port of the same expectation. Co-Authored-By: Claude Opus 4.8 <noreply@anthropic.com>
Move the exact hypergeometric expectation for QuartetConcordance(normalize = TRUE) from R to C++ (src/concordance_expect.cpp): `quartet_expect` returns E[conc], E[dec]; `trit_expect` returns E[m], E[m*A/wk], E[m*A/wc]. Both sum the per-state-pair trivariate-hypergeometric double sum, hoisting the log-choose terms out of the M/p/r loops (the dominant cost). The R exact helpers (.ExpectedTrit/.CharTritExpect/.ExpectedQuartet and the per-char .QuartetExpect loop) are replaced by thin wrappers / a single all-characters call in .TritConcordance; the Monte-Carlo (normalize = <int>) path stays in R. New Rcpp exports are registered in the manual R_CallMethodDef table (src/TreeSearch-init.c; the package uses R_useDynamicSymbols(FALSE)). Speedup at split-support scale (48 tips, 96 chars): exact normalize=TRUE quartet 2.88s -> 0.14s (~20x), trit 3.48s -> 0.45s (~8x); a full 1000-matrix run drops from ~100 min to ~10 min. Validated bit-for-bit (max |C++ - R| ~1e-14) against a self-contained R reference across binary/multistate/missing data (dev/benchmarks/frac-quart/cpp_expect_parity.R); the exact-vs-MC end-to-end tests in test-Concordance.R now exercise the C++ path and pass. Co-Authored-By: Claude Opus 4.8 <noreply@anthropic.com>
`normalize` becomes `chanceCorrect` in `ClusteringConcordance()`, `QuartetConcordance()` and `ConcordanceTable()`. The argument moves the zero point to the value expected under a fixed-marginal null; the scaling to a maximum of 1 was never governed by it, so the old name described something the argument does not do. `unit = "trit"` becomes `unit = "nrqs"`, after the quantity it counts: non-redundant quartet statements. Of the quartets a split of sizes (k, t - k) resolves, (k - 1)(t - k - 1) are non-redundant under the Nelson-Ladiges entailment. `QuartetConcordance()` chance-corrects by default, as `ClusteringConcordance()` already does, so values may be negative where agreement falls below chance expectation. `chanceCorrect = FALSE` gives the uncorrected measure. Documentation gains the caveat that the NRQS ceiling of 1 bounds each split individually, so it holds for `return = "edge"` but not for `return = "char"`, which averages over every split in the tree. # Conflicts: # NEWS.md
`ms609/MaxMin` is now `ms609/Coreset`; every exported solver keeps its name, but the progress option moved from `MaxMin.progress` to `Coreset.progress`. `WideSample()`'s dispatch, its dependency guard, the `Suggests`/`Remotes` entries and the Shiny app's option are updated accordingly, as is the `skip_if_not_installed()` guard in `test-WideSample.R` -- which had been skipping silently, since `MaxMin` can no longer be installed under that name. Three references stay as they are: `ExactMaxMin()` is a function name, the "MaxMin optimum" is the objective rather than the package, and `inst/REFERENCES.bib` holds the title of the cited paper.
`spell_check_package()` flags four technical terms, all legitimate: the Nelson-Ladiges surname and the NRQS acronym arrive with the concordance rename, while `multistate` and `unclamped` predate it and were already failing the check.
The argument was described in both `@details` and `@param`, and the two disagreed on what an integer samples: the details said `n` uniformly random trees, the parameter entry said `n` reshuffles of the character's tokens. Both are true, of different functions. The one surviving statement sits under `@param`, and names the null each function actually draws from: `ClusteringConcordance()` scores each character against `n` trees from `RandomTree()`, whereas `QuartetConcordance()` permutes each character's tokens across the scored leaves, holding its state frequencies and the split sizes fixed -- the same fixed-marginal null that its `chanceCorrect = TRUE` expectation evaluates in closed form. `man/SiteConcordance.Rd` is updated to match.
Four files needed resolving. `cpp-search` had already followed the MaxMin -> Coreset rename (089b851), and gone further: `FarFirst()` is called with named arguments, and `Remotes:` has become `Additional_repositories:`, which is what lets pak resolve the package on the Windows CI job. Its versions win throughout; the duplicate `Coreset` entry in `Suggests` is dropped. `cpp-search` also renamed `QuartetConcordance()`'s `return` options from `c("character", "site", "default")` to `c("edge", "char")`, matching `ClusteringConcordance()`, and validates them rather than falling back on a default. The NRQS path adopts that vocabulary: `return = "default"` becomes `return = "edge"`, and the validation is hoisted above the `unit == "nrqs"` branch so that path receives a checked value. The chance-correction block keeps its place, and `miRand` is now read through `.ConcSlice()`. Two tests were left stale by the merge, both now updated: - the `return` alias assertions asserted that `"default"`, `"character"` and `"site"` are accepted, which the argument-checking test on the same file now asserts they are not; - the `{0-}` ambiguity test bounded results to [0, 1], which no longer holds now that `QuartetConcordance()` chance-corrects by default. The ceiling of 1 still stands, so the lower bound gives way to a finiteness check. `test-Concordance.R` and `test-WideSample.R` pass in full against a build of the merged tree; `spell_check_package()` is clean.
ms609
enabled auto-merge
August 21, 2026 12:51
This file contains hidden or bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment
Add this suggestion to a batch that can be applied as a single commit.This suggestion is invalid because no changes were made to the code.Suggestions cannot be applied while the pull request is closed.Suggestions cannot be applied while viewing a subset of changes.Only one suggestion per line can be applied in a batch.Add this suggestion to a batch that can be applied as a single commit.Applying suggestions on deleted lines is not supported.You must change the existing code in this line in order to create a valid suggestion.Outdated suggestions cannot be applied.This suggestion has been applied or marked resolved.Suggestions cannot be applied from pending reviews.Suggestions cannot be applied on multi-line comments.Suggestions cannot be applied while the pull request is queued to merge.Suggestion cannot be applied right now. Please check back later.
Adds a redundancy-corrected currency and an exact chance correction to
QuartetConcordance(), via two new arguments:unit = c("quartet", "trit")andnormalize.What this adds
unit = "trit"scores concordance in non-redundant quartet statements. The quartets a bipartition of sizes (k, t−k) resolves are the 4-cycles of K_{k,t−k}; Nelson–Ladiges entailment is GF(2) cycle addition, so only the cyclomatic number (k−1)(t−k−1) are independent. Multistate characters decompose into edge-disjoint state-pair blocks, so trits sum over pairs.normalize = TRUEsubtracts the exact expectation under a fixed-marginal hypergeometric null (both split sizes and the character's own state frequencies held fixed), computed as a closed double sum rather than simulated.normalize = <integer>gives a Monte-Carlo baseline instead.src/concordance_expect.cpp) for the expectation, mirroring the R reference to machine precision.Why chance correction matters here
Under the raw quartet currency the chance level is ~1/3 largely independently of split size, so correcting it is nearly a monotone transform — rank statistics are unchanged (Spearman between corrected and uncorrected = 1.0000000000 for the weighted score over 354,638 splits). Under the trit currency the
(x−1)+floors make the chance level drift with split size, so the correction reorders substantially (Spearman 0.39–0.56). Pooling uncorrected trit scores across split sizes therefore confounds signal with a size artefact; corrected, it does not.Validation
Rcpp::compileAttributes()produces no further changes.quartetandwQuartetreproduce previously cached values from an independent 999-alignment simulation study to 5.55e-16 across all 354,638 splits.Empirical behaviour
Scored against 999 simulated 48-taxon datasets (C-index vs graded similarity to the generating tree; ROC-AUC vs presence in it), on the parsimony-statistics subset (n = 36,761) the correction lifts unweighted quartet concordance from below chance (C 0.489, AUC 0.492) to 0.658/0.659, and the weighted score from 0.572 to 0.678. On the ML subset (n = 44,955) the weighted score goes 0.633/0.643 to 0.706/0.715.
Base
Targets
cpp-searchrather thanmaindeliberately: the character-recoding change inQuartetConcordance()introduced by 78b7414 (which alters whether{0,-}-style tokens count as grouping information) is a prerequisite for this work and is present oncpp-searchbut not onmain.🤖 Generated with Claude Code