Feat/analytic optimal tc

[rneher]

No description provided.

[ivan-aksamentov]

⚠️ AI-generated content below. Verify all claims.

Review of PR #855 "Feat/analytic optimal tc" (neherlab/treetime).

• Head: 723f86c7780d757601604ed1c3b2624aa0a95e99
• Base: rust

Overview

Click to expand

The branch replaces the numerical optimization of the constant coalescent timescale ¹ $T_c$ with its closed-form maximum-likelihood estimate. v0 finds $T_c$ by minimizing the negative total log-likelihood with Brent's method; this branch computes it directly. The change and its rationale are recorded in kb/decisions/coalescent-analytic-tc-optimization.md, and the same estimand is preserved (see Notes).

Scope: the coalescent optimizer and its per-edge accumulation, the merger-rate integration helper, and knowledge-base documentation. No CLI surface is touched by this branch.

Background

The analytic coalescent MLE for a constant timescale [click to expand]

Under the Kingman coalescent Kingman 1982 [1], the genealogy likelihood factorizes over intervals of constant lineage count. Each merger contributes a hazard factor and each interval a survival factor, so for a single constant $T_c$ the likelihood collapses to

$$\mathcal{L}(T_c) \propto T_c^{-M}, \exp!\left(-,\frac{I}{T_c}\right)$$

where $M$ is the number of merger events and $I = \int \binom{k}{2},dt$ is the pairwise-merger-rate integral over the tree with $T_c$ factored out, and $k$ is the lineage count. Maximizing $\log \mathcal{L} = -M\ln T_c - I/T_c$ gives

$$\frac{d\log\mathcal{L}}{dT_c} = \frac{I}{T_c^2} - \frac{M}{T_c} = 0 \quad\Longrightarrow\quad \hat{T_c} = \frac{I}{M},$$

the unique interior maximum. This is the standard scale-parameter MLE for the variable-population-size coalescent Griffiths and Tavaré 1994 [2]: total coalescent opportunity divided by the number of mergers, closed form and needing no iteration.

Blocking issues

Correctness concerns worth resolving before merge.

None. The analytic estimator is sound (see Notes).

Non-blocking issues

Documentation drift. Fix if time allows.

🔵 N1. `--coalescent-opt` help still describes Brent's method and an initial guess [click to expand]

The help for --coalescent-opt reads "find the optimal Tc using Brent's method ... If --coalescent is also provided, that value is used as the initial guess; otherwise defaults to 1.0" [src]. The decision this branch records replaces the Brent search with the closed-form $\hat{T_c} = I/M$ solve [src], so there is no iteration and no initial guess; the help was not updated to match.

Suggestions:

• Update the help to describe the closed-form solve and drop the "initial guess" wording.

Notes

Positive observations about the change.

Click to expand

• The closed form $\hat{T_c} = I/M$ is the correct maximum-likelihood estimate and is derived in the decision record and the docstring [src]. Removing Brent's [-20, 2] log-space bracket also removes a bounds artifact that silently capped extreme timescales.
• Accumulating $I$ and $M$ per edge is consistent and general: an $m$-child node contributes $m-1$ mergers through the sum of $(\texttt{n_siblings}-1)/\texttt{n_siblings}$ over its edges [src], so hard polytomies contribute their full merger multiplicity and unary (degree-two) nodes contribute zero, matching the theory. The per-edge form respects the same bad-branch exclusion as the likelihood and carries over to the skyline case.
• The success guard rejects non-finite or non-positive optima and falls back visibly rather than returning a silently wrong $T_c$ [src].
• The reported likelihood is evaluated through the shared coalescent model, so it matches fn compute_coalescent_total_lh() used elsewhere [src].

Glossary

Click to expand

1. Coalescent timescale ($T_c$). The generation-scaled parameter proportional to effective population size that sets the pairwise coalescence rate $1/T_c$ in the Kingman coalescent (Kingman 1982 [1]). ↩

References

Click to expand

1. Kingman, John Frank Charles. 1982. "The coalescent." Stochastic Processes and their Applications 13:235-248. https://doi.org/10.1016/0304-4149(82)90011-4 ↩
2. Griffiths, Robert C., and Simon Tavaré. 1994. "Sampling theory for neutral alleles in a varying environment." Philosophical Transactions of the Royal Society B 344:403-410. https://doi.org/10.1098/rstb.1994.0079 ↩