Conversation
Smooths over every bipartition of a clade rather than only those in the CCD graph, with the pseudocount set by how many novel clades the bipartition introduces: a per-split alpha on the CCD0 split set (classes 0 and 1 pooled), and class totals alpha1, alpha2 for the one- and two-novel-clade classes. The exponentially large class is obtained by subtraction from 2^(m-1)-1, so it is never enumerated and there is no #P-hard counting. Per-class totals rather than per-split constants are essential: with a constant pseudocount the two-novel class swamps the data, leaving observed splits 6e-9 of the probability at a 40-taxon root clade. Properties, all by construction: exactly normalised by the chain rule (no partition function, no truncation); full support; regCCD nested exactly at alpha1 = alpha2 = 0; an observed split always outranks an expanded one at the same clade; and the mass held back shrinks as f(C) grows. Implements scoring, exact sampling, MAP and entropy. Sampling draws a class then a member, using rejection for the two-novel class, and scores each draw with the same routine the scorer uses, so the sampled and scored distributions coincide. MAP runs a DP over the observed-clade DAG with an optional wider search admitting one-novel-clade splits, plus a bound that can certify global optimality. Entropy is unbiased by Monte Carlo; a deterministic recursion is also provided, which approximates novel subclades as structureless and is therefore optimistic. Verified in CRegCCDTest and CRegCCDMapEntropyTest: total mass 1.000000000000 by enumeration on 4-7 taxa; agreement with RegCCD to 1.8e-15 at alpha1 = alpha2 = 0; class sizes partitioning all 2^(m-1)-1 bipartitions; sampled frequencies matching scored probabilities across all 105 five-taxon topologies; and the MAP search matching brute force in twelve configurations with the certificate firing in each. Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
… real data KRegNormalisationTest enumerates every rooted topology and asserts that KRegCCD is exactly normalised on four taxa but sub-normalised beyond it, with the deficit scaling as mu^2 and never exceeding one. This is the Theta(mu^2) maximality deficit, which is distinct from, and larger than, the O(mu^(k+1)) reserve truncation. SplitClassSizeAnalysis reports the four split-class sizes at a root clade, showing that a constant per-split pseudocount leaves the observed splits 6.3e-9 of the probability by 40 taxa. ClassUsageAnalysis attributes a held-out tree's log probability to the split classes and measures how much observed structure a two-novel-clade split destroys. RealDataHeadToHeadTest compares CCD1, regCCD, KRegCCD, MRegCCD and CRegCCD on a real posterior using the manuscript's RSV2 protocol, selecting each model's hyperparameters on a validation split disjoint from the fitted set. Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
MRegCCD enumerated novel resolutions recursively, which blew up at boundary size 4 and silently truncated against its op budget: on RSV2 at depth 4 the score moved from -60.6219 to -60.6443 when the budget was raised 20x, so the default was returning a truncated answer without saying so. The new MRegCCD overrides only countsFor, indexing disjoint pairs of observed subclades by the bitset they cover and counting each boundary order against that index rather than by enumeration. Boundary 4 needs a canonical form -- a 4-part boundary decomposes into two pairs in three ways -- or every such boundary is counted three times; the agreement test caught exactly that. The previous implementation is kept as MRegCCDSlow so the two can be checked against each other, and MRegCCDAgreementTest holds them to identical counts on 2199 clades. Both now default to reserve depth 4, matching KRegCCD's default reserve so the two models look equally far past the CCD graph. Measured 4.3x faster at depth 4 and no longer truncating. BitSet.newBitSet used other.length() rather than other.size() to size the copy. length() is the index of the highest set bit plus one, so copying a set whose top words are empty returned an undersized BitSet, and the bitwise operations index the operand by this.words.length. Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
The comparison was driven by a JUnit test forked once per dataset per model via Maven. Class loading of BEAST and CCD takes 30ms, so there was nothing to amortise by forking; project scan, dependency resolution, plugin loading and the surefire fork cost several seconds each, and a sweep spent over 99% of wall clock on startup -- 11m39s of which 2.4s was modelling. HeadToHeadBenchmark is a plain main that loops all datasets and models in one JVM. The full corpus of 101 datasets now takes about 26 minutes. It resumes: a sweep reads back the rows already in the CSV, skips those datasets and appends. Three separate runs were lost to a crash partway through (heap exhaustion, and a SIGBUS once swap was exhausted), each costing every dataset already scored. Failures are caught per model per dataset and recorded in the row, so a model that fails on one dataset costs that cell rather than the row or the run -- KRegCCD and MRegCCD both exhaust an 8GB heap on the 441-taxon stervander set while CRegCCD completes it. Hyperparameters are selected on a validation quarter disjoint from the fitted quarter. Drawing the validation trees from within the fitted set drove every parameter to its grid floor. The alpha grid is capped at 5: alpha is a per-split pseudocount, so beyond alpha = 1 the prior outweighs a real observation, and the only principled reason to exceed it is MCMC autocorrelation inflating counts by N/ESS -- the cap admits that down to ESS = 200 of 1000 training trees. Parameters are serialised with Double.toString. A "%.4f" format collapsed the alpha2 bands at 1.6e-5, 8e-5 and 4e-4 all to "0.0000", so the selected value could not be recovered from the CSV at all. Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
CRegCCD's four split classes -- observed, expanded, one novel clade, two novel clades -- are a property of the split and the training set rather than of any model, so all three full-support models can be asked what they charge for the same split. This tool classifies every split of every held-out tree directly from the training trees (no model internals) and regresses the per-tree log probability on the class counts crossed with the size of the clade being split. The regression has no intercept: the four class counts sum to the number of internal nodes, which is fixed for a given taxon set, so a design matrix with an intercept is singular and the split between intercept and coefficients is arbitrary. Without it each coefficient is directly the mean log probability contributed by one split of that class. The fitted mean is printed alongside the actual mean, since OLS forces them equal and the check catches a solver that is quietly wrong. On pela-zego this reproduces each model's benchmark mean exactly and explains R^2 = 0.996-0.998 of the per-tree variance, which is what licenses reading the coefficients as per-split costs even though KRegCCD and MRegCCD price maximal novel regions rather than individual splits. The result it was written to test: CRegCCD spreads the class total alpha2 uniformly over |A_2(C)| ~ 2^(m-1) splits, so its two-novel penalty falls linearly in the clade size m with slope -log 2, measured at -0.659 nats per taxon against the predicted -0.693. KRegCCD and MRegCCD raise an escape rate eps ~ mu/O(s^2) to a fixed power, so theirs grows only like log m. At m = 193 that is -129 nats per split against -14 and -6. Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
CRegCCD contributes one conditional split factor per internal node, including at novel clades. KRegCCD and MRegCCD cluster novel nodes into maximal regions and score each region once at its top, so nodes interior to a region never get a conditional distribution at all. On pela-zego that is 377.0 factors per held-out tree against 369.4. An attempt to measure this by regressing the score on class counts crossed with whether the clade was observed produced a coefficient of +43 nats for a single split, which is impossible: the interior-of-region counts are largely determined by the novel splits above them, so those columns are collinear and the coefficients are not identifiable. Counting the nodes is exact and needs no fit, so that is what this reports. Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
The slope of the two-novel penalty against clade size was read off three binned means, which needed an arbitrary "too few splits to report" threshold and then fitted a line through the bins as though each carried equal weight. On pela-zego those bins hold 1, 4, 23, 859, 6 and 438 splits, and the top bin spans clade sizes 74 to 252, so a single coefficient for it averages a predicted penalty ranging from -51 to -174 nats. The size now enters as a covariate instead: per tree, the counts of the first three classes, the count of two-novel splits, and the sum over those splits of (m-1) in one fit and log m in the other. The coefficient on that column is the slope, estimated from every two-novel split at its proper weight, and the two fits can be compared on residuals to test which functional form holds. This sharpens the result rather than softening it. CRegCCD's slope is -0.7002 per taxon against the -log 2 = -0.6931 predicted by dividing a class total across |A_2(C)| ~ 2^(m-1), within 1% where the binned estimate was 5% off, and the linear form beats the logarithmic one on residuals 146316 to 250084. KRegCCD and MRegCCD are flat in m at -0.023 and -0.010 per taxon, and for them the two forms are indistinguishable, as expected when the dependence enters only through log m. Cell counts and clade-size ranges are printed alongside the binned table that remains, so a thin cell is visible rather than merely suppressed. Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
An earlier commit added hs_err_pid*.log by appending to a .gitignore whose last line had no trailing newline, producing ".projecths_err_pid*.log". That matched nothing, so both patterns were silently disabled: .project stopped being ignored and the crash logs it was meant to add never were. Also drops a __pycache__ entry that was appended here by a command run from the wrong working directory; it belongs to the benchmark scripts in the CCD_regularisation repository, not to this one. Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
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
CRegCCD, a full-support CCD that replaces KRegCCD's escape mechanism with additive smoothing over every bipartition of a clade. Paper section is in a PR on the manuscript repo (yangsoph/CCD_regularisation#10).The model
The pseudocount depends only on how many novel clades a bipartition introduces:
alphaper splitalpha1totalalpha2total|A_2|follows by subtraction from2^(m-1)-1, so the exponentially large class is never enumerated. There is no #P-hard counting, no reserve equation and no truncation — the cost is oneO(#clades)scan per clade, comparable to the CCD0 expansion.Per-class totals rather than per-split constants are essential: with a constant pseudocount the two-novel class swamps the data, leaving observed splits
6.3e-9of the probability at a 40-taxon root clade (SplitClassSizeAnalysis).Properties hold by construction rather than by correction:
Theta(mu^2)maximality deficit KRegCCD carries cannot arise, because there are no regionsalpha1 = alpha2 = 0f(C)growsWhat is implemented
Scoring, sampling, MAP and entropy.
Sampling draws a class then a member — classes 0 and 1 by enumeration, class 2 by rejection from uniform bipartitions with an enumeration fallback when acceptance would be poor. Each draw is scored with the same routine the scorer uses, so the sampled and scored distributions coincide exactly; there is no sampling-fidelity choice to make.
MAP runs a DP over the observed-clade DAG, with an optional wider search admitting one-novel-clade splits, plus a bound that can certify global optimality.
Entropy is unbiased by Monte Carlo, since sampler and scorer agree. A deterministic recursion is also provided; it approximates novel subclades as structureless and is therefore optimistic, measured at +2.0% to +10.4% on real data, so Monte Carlo is the estimator of record.
Verification
In
CRegCCDTestandCRegCCDMapEntropyTest:1.000000000000by enumeration on 4–7 taxa, for several training sets and arbitrary pseudocountsRegCCDto1.8e-15atalpha1 = alpha2 = 02^(m-1)-1bipartitions on 135 cladesAlso included
Measurement tests used to characterise the models:
KRegNormalisationTest(KRegCCD is exactly normalised on four taxa, sub-normalised beyond, deficit scaling asmu^2and never exceeding one),SplitClassSizeAnalysis,ClassUsageAnalysis, andRealDataHeadToHeadTest, which compares CCD1, regCCD, KRegCCD, MRegCCD and CRegCCD on a real posterior with hyperparameters selected on a disjoint validation split.Known limitations
alpha1andalpha2are class totals, not per-split pseudocounts, so thealphanaming is arguably misleading —mu1/mu2is under consideration and not yet applied.🤖 Generated with Claude Code