Comparing Alternative Lineage Topologies

Why Topology Comparison Needs More Than a Best-Fitting Tree

A lineage topology describes how candidate cell fates are related by branches. Two topologies can fit the observed trajectories differently, but a more elaborate tree can often fit at least as well simply because it has more adjustable quantities. Those quantities may include branchpoint timing, trajectory coefficients, noise levels, and fate probabilities. A fair comparison must therefore account for model complexity and for uncertainty in those continuous parameters.

The comparison also depends on what is held fixed. Candidate topologies should be evaluated against the same observations, with compatible likelihood definitions and comparable assumptions about measurement noise and cell assignments. The likelihood-based criteria below are not interchangeable: AIC targets predictive performance, BIC is a large-sample approximation to a model-evidence comparison under restrictive conditions, and a Bayes factor compares likelihoods averaged over specified parameter priors.

Definition: Let \(T\) denote a candidate lineage topology and let \(\theta_T\) collect its continuous parameters. For observed cell-level data \(Y=(Y_1,\ldots,Y_n)\), the topology-specific likelihood is \(L_T(\theta_T)=P(Y\mid T,\theta_T)\). If a cell's branch label \(Z_i\) is unobserved, its contribution is obtained by summing over possible branches: \(P(Y_i\mid T,\theta_T)=\sum_b P(Z_i=b\mid T,\theta_T)P(Y_i\mid Z_i=b,T,\theta_T)\). The log likelihood is \(\ell_T(\theta_T)=\log L_T(\theta_T)\), and its maximized value is \(\hat{\ell}_T=\ell_T(\hat{\theta}_T)\), where \(\hat{\theta}_T\) maximizes the likelihood for that topology.

The sum over branch labels matters: assigning every cell to its most likely branch and then pretending those assignments were observed generally gives a different likelihood. Clonal labels or other observed lineage constraints can be incorporated into the likelihood, provided the same data and treatment of those data are used for each candidate. For the basic criteria below, suppose the observations used in the likelihood are conditionally independent given the model parameters. If that assumption is inappropriate, the likelihood and the effective number of independent observations must be reconsidered rather than mechanically applying a criterion.

Penalized Likelihood Criteria

AIC and BIC start with the maximized log likelihood and subtract a penalty for the number of estimated parameters. Let \(k_T\) be the number of free parameters estimated under topology \(T\), and let \(n\) be the number of independent observational units represented by the likelihood.

Definition: The Akaike information criterion and Bayesian information criterion are \( \mathrm{AIC}_T=-2\hat{\ell}_T+2k_T \) and \( \mathrm{BIC}_T=-2\hat{\ell}_T+k_T\log(n) \). For either criterion, smaller values indicate a preferred candidate among those compared, according to that criterion. Neither value is a probability that the topology is true.

AIC uses a penalty of two units per estimated parameter. Under its regularity conditions, it estimates relative expected out-of-sample information loss, so it is principally a predictive-selection criterion. BIC's penalty grows with sample size; it approximates a Bayesian model comparison when the models are regular, identifiable, and have fixed dimension, with suitable proper priors and a large sample. These purposes differ, so disagreement between the criteria is possible and is not an arithmetic error.

For BIC, \(n\) refers to the independent units in the likelihood, not the number of genes, branches, or latent labels. In a cell-level likelihood with conditionally independent cells, \(n\) is the number of cells used. The parameter count must likewise reflect what was estimated, including topology-specific noise or branch parameters, rather than only the trajectory coefficients that are easiest to list.

Worked Example: AIC and BIC Favor Different Trees

Suppose two candidate topologies are fitted to the same hypothetical data set of \(n=100\) conditionally independent cells using the same likelihood. Topology A has \(k_A=8\) fitted parameters and \(\hat{\ell}_A=-120\). Topology B has \(k_B=12\) and \(\hat{\ell}_B=-114\). The higher log likelihood for B reflects its closer fit, but it also uses four additional parameters.

For AIC, the calculations are \( \mathrm{AIC}_A=-2(-120)+2(8)=240+16=256 \) and \( \mathrm{AIC}_B=-2(-114)+2(12)=228+24=252 \). AIC favors B by 4 criterion units.

For BIC, \(\log(100)\approx4.605\), so \( \mathrm{BIC}_A=240+8(4.605)=276.84 \) and \( \mathrm{BIC}_B=228+12(4.605)=283.26 \). BIC favors A by about 6.42 units. Both calculations use the same fitted likelihoods and include the topology-specific parameter counts; the disagreement comes from the different penalty strengths. The result is a criterion-dependent preference, not proof that either tree is biologically correct.

Bayes Factors Integrate Over Parameter Uncertainty

A maximized likelihood evaluates each topology at its best-fitting parameter value. A Bayes factor instead compares marginal likelihoods, which average the likelihood over parameter values using a prior distribution. Thus, a topology does not receive credit only for a narrow, high-likelihood fit: the prior-weighted fit across its parameter space also matters.

Definition: For topology \(T\) with proper prior density \(\pi_T(\theta_T)\), its marginal likelihood is \(m_T(Y)=\int L_T(\theta_T)\pi_T(\theta_T)\,d\theta_T\). For two topologies \(T_1\) and \(T_0\), the Bayes factor is \(BF_{10}=m_{T_1}(Y)/m_{T_0}(Y)\). Values above 1 favor \(T_1\) in this comparison; values below 1 favor \(T_0\). Posterior odds equal the Bayes factor multiplied by prior odds.

The priors must be proper for ordinary marginal likelihoods and Bayes factors to be defined. They should also be specified on meaningful scales: a prior that allows implausibly extreme branch times or trajectory differences can change the comparison. A Bayes factor is not the posterior probability of a topology. For two candidates only, that probability also depends on their prior probabilities, and it is conditional on those two candidates and their models.

Worked Example: Integrating Over a Branch Difference

Consider a hypothetical summary \(d=1.2\) of the expression contrast between two descendant branches, with known standard error \(0.5\). Under a no-split topology \(T_0\), suppose \(d\sim N(0,0.25)\), where the second argument is the variance. Under a split topology \(T_1\), let the unknown branch contrast be \(\delta\), with \(d\mid\delta\sim N(\delta,0.25)\) and prior \(\delta\sim N(0,1)\). These assumptions describe a one-dimensional illustration, not a general likelihood for all lineage data.

Under \(T_0\), the marginal distribution of \(d\) is \(N(0,0.25)\). Under \(T_1\), integrating over the normal prior adds the prior variance to the sampling variance, giving \(d\sim N(0,1.25)\). Therefore,

BF_{10} =\frac{f(d=1.2\mid T_1)}{f(d=1.2\mid T_0)} =\sqrt{\frac{0.25}{1.25}} \exp\left\{-\frac{1.2^2}{2(1.25)}+\frac{1.2^2}{2(0.25)}\right\} \approx 0.4472\,e^{2.304} \approx 4.48.

The observed contrast is about 4.48 times as probable under the split model's marginal distribution as under the no-split model's distribution. With equal prior probabilities on these two topologies, the posterior probability of \(T_1\) within this two-model comparison is \(4.48/(1+4.48)\approx0.818\). That is not a probability that the split is true without qualification: it depends on the two-model set, the equal prior odds, and the stated prior for \(\delta\).

For example, widening the prior to \(\delta\sim N(0,4)\) gives marginal variance \(4.25\) under \(T_1\) and a Bayes factor of about 3.65. The change illustrates prior sensitivity: integrating over a wider range can reduce the relative support from data concentrated near a particular contrast.

Derivation: Why BIC Has a Parameter Penalty

The BIC penalty can be motivated by approximating the marginal likelihood near the maximum-likelihood estimate. Assume a regular, identifiable model with \(k_T\) continuous parameters, an interior maximum, a smooth positive prior near that maximum, and enough independent observations for a local quadratic approximation. Near \(\hat{\theta}_T\), the log likelihood is approximately its maximum minus a quadratic term whose curvature grows proportionally to \(n\). The likelihood is therefore concentrated in a region with width proportional to \(n^{-1/2}\) in each parameter direction.

Integrating over \(k_T\) such directions contributes a volume factor proportional to \(n^{-k_T/2}\). Up to terms that do not grow with \(n\), the Laplace approximation gives

\log m_T(Y) \approx \hat{\ell}_T-\frac{k_T}{2}\log(n)+O(1), \qquad -2\log m_T(Y)\approx -2\hat{\ell}_T+k_T\log(n)+O(1).

Dropping the bounded terms yields the BIC expression. The derivation explains both its complexity penalty and its limits: if parameters are unidentifiable, lie on boundaries, or have non-quadratic likelihood shapes, the local-volume argument can fail. Branching models can have precisely these irregularities. BIC may still be reported as a convenient approximation, but its usual evidence interpretation should then be qualified.

Likelihood-Ratio Comparisons and Irregular Topologies

For nested models, one can compare the maximized likelihoods using the likelihood-ratio statistic \(D=2(\hat{\ell}_1-\hat{\ell}_0)\), where model 0 is the restricted model and model 1 contains it as a special case. Under regular conditions, including identifiable parameters and null values in the interior, \(D\) has approximately a chi-squared distribution with degrees of freedom equal to the difference in parameter counts under \(H_0\). This reference distribution does not automatically apply just because one tree can be described as a limiting case of another.

For instance, when a branch disappears under the null, its branchpoint time or branch-specific trajectory may be undefined under that null. A branch probability may also be on a boundary. These features violate regularity conditions, so the usual chi-squared calibration can be inaccurate. For nonnested topologies, the standard nested-model likelihood-ratio test is not applicable in the first place.

Worked Example: Calibrating a Nonregular Split Comparison

Suppose a hypothetical no-split topology is nested in a bifurcating topology when the branch contrast is zero. The observed fitted log likelihoods give \(D_{\mathrm{obs}}=6.4\). Because the branchpoint timing is not identified under the no-split model, assume the ordinary chi-squared reference is not justified. Instead, use a parametric bootstrap under the fitted no-split model: simulate data under that null, refit both topologies to every simulated data set, and calculate the likelihood-ratio statistic each time. This procedure must repeat the fitting steps, including estimating continuous parameters, so that the simulated statistics reflect their estimation uncertainty.

If 9 of 199 simulated statistics are at least 6.4, a common finite-simulation estimate is \( \hat{p}=(9+1)/(199+1)=10/200=0.05 \). The added one in numerator and denominator avoids reporting a zero estimate when no simulated statistic exceeds the observed value. Under the null model and the specified simulation procedure, this is borderline evidence against the no-split topology at \(\alpha=0.05\); Monte Carlo error is substantial with only 199 simulations. It does not establish that the bifurcating topology is true, and a bootstrap result cannot repair a misspecified null model.

Practical Checks, Alternatives, and Limitations

Before interpreting any topology ranking, check that the fitted comparisons are numerically and scientifically meaningful. A topology with a larger maximized likelihood should not be assumed better overall; assess its penalty or integrated evidence as appropriate. Useful checks include:

Optimization: refit from multiple starting values. A lower local maximum can make a topology appear artificially poor.
Parameter uncertainty: inspect profile likelihoods or posterior distributions for branchpoint times and trajectory parameters. Very broad or boundary-concentrated estimates signal weak identification.
Label and branch behavior: verify that latent branch assignments are summed or integrated as specified, and that equivalent label permutations are not counted as distinct biological topologies.
Model compatibility: use the same observations, response representation, and likelihood assumptions across candidates. A numerical criterion cannot correct a poor expression or noise model.
Stability: examine whether the ranking changes under plausible prior scales, alternative optimization starts, or reasonable changes in the fitted features.

The main alternatives answer different questions. AIC is useful when relative predictive performance is the priority and its regularity assumptions are reasonable. BIC provides a simpler, increasingly strong complexity penalty and an approximate evidence comparison in regular large-sample settings. Bayes factors explicitly integrate parameter uncertainty but require defensible priors and can be prior-sensitive. Likelihood-ratio tests provide calibrated evidence for nested comparisons when their reference distribution is valid; simulation-based calibration is preferable when it is not. Predictive validation can assess performance on held-out observations, but it does not by itself produce posterior odds or establish a biological branching mechanism.

No single criterion resolves biological uncertainty. If candidate trees have similar support, report that ambiguity and the parameter uncertainty rather than presenting a unique tree as settled. A topology comparison is conditional on the candidate set: a preferred tree can still be wrong if all candidates omit an important branching structure or share a misspecified trajectory model.

Common mistake: “The tree with the highest maximized likelihood must be the best topology.” A larger model can achieve a higher maximum simply by having more adjustable parameters, and that maximum ignores how uncertain those parameters are. Compare models using a stated criterion or an explicitly specified marginal-likelihood analysis, verify the assumptions behind it, and describe the result as relative support among the candidates considered.

Check Your Understanding

Answer each question, then compare your reasoning with the solutions below.

Two topologies have the same maximized log likelihood, but one estimates more parameters. What happens to their AIC values, and why?
In a cell-level likelihood with 250 conditionally independent cells, what value of \(n\) enters BIC? Explain why the number of genes is not automatically the answer.
A Bayes factor of 5 favors \(T_1\) over \(T_0\). Does this mean \(P(T_1\mid Y)=0.83\) in every analysis? Explain.
A split branch collapses to zero length under a null topology, and its timing is then undefined. Why might the usual chi-squared likelihood-ratio calibration fail?

Solutions. (1) The more complex topology has the larger AIC because its \(2k_T\) penalty is larger while the likelihood terms are equal. (2) Use \(n=250\), given the stated conditional-independence assumption and cell-level likelihood; genes are response dimensions, not automatically independent observational units. (3) No. Posterior odds equal the Bayes factor times prior odds. With equal prior odds and only these two candidates, the posterior probability would be \(5/(1+5)=0.833\); other prior odds or candidate models change it. (4) The branch timing is not identified under the null, so regularity conditions for the standard chi-squared reference distribution fail. A null-based simulation that refits both models can calibrate the statistic under the specified model.