Deriving the Likelihood for Branching Trajectories

From a Tree Model to an Observed-Data Likelihood

A tree-structured trajectory model describes cells using two quantities that are usually not observed directly: a cell’s progress along the trajectory and its identity as a member of a particular branch. The observed data are expression measurements, not the latent coordinates themselves. A likelihood for those data must therefore account for every possible latent progress and branch identity that could have generated each cell.

The central operation is marginalization: multiply the probability or density of an observation under a latent state by the probability of that state, then sum over discrete latent states and integrate over continuous ones. A branch contributes through both its probability and its observation distribution. Omitting either part generally gives a different model.

Definition: For cell $i$, let $\boldsymbol y_i$ be its observed expression vector, $T_i$ its latent progress, and $B_i\in\{1,\ldots,K\}$ its latent branch identity. Let $g(t)$ be the density of progress, $q_b(t)=P(B_i=b\mid T_i=t)$ the branch probability, and $f_b(\boldsymbol y_i\mid t)$ the observation density on branch $b$ at progress $t$. The marginal density for one cell is $$ p(\boldsymbol y_i)=\sum_{b=1}^{K}\int g(t)q_b(t)f_b(\boldsymbol y_i\mid t)\,dt. $$

The density $g(t)$ and probabilities $q_b(t)$ describe the distribution of latent states among the cells represented by the sampling process. They need not describe the speed of biological development or the fraction of all cells in each state unless the sampling design supports that interpretation. At each progress value, $q_b(t)\geq 0$ and $\sum_{b=1}^{K}q_b(t)=1$. The observation distributions must also be valid probability distributions: for continuous measurements, each density integrates to one over possible observations; for discrete counts, the corresponding probability mass function sums to one.

In a tree with a shared trunk, the branch labels can represent eventual branch identity even before the paths physically diverge. Before the branchpoint, the observation distribution can be shared, so $f_b(\boldsymbol y\mid t)=f_0(\boldsymbol y\mid t)$ for every $b$ in that region. Summing over the branch probabilities then gives $\sum_b q_b(t)f_0(\boldsymbol y\mid t)=f_0(\boldsymbol y\mid t)$. After divergence, branch-specific means, variances, or other features of the observation distributions can differ.

Assumptions Behind the Likelihood

The formula is a generative model, not merely a bookkeeping identity. Its interpretation relies on specifying how a cell’s latent state is drawn and how its observed data arise from that state. For a sample of $n$ cells, a common working assumption is that cells are independent draws from the same modeled sampling distribution. Conditional on the latent states and model parameters, their observations are also independent. Under these assumptions, the observed-data likelihood is

L(\theta)=\prod_{i=1}^{n}\left\{\sum_{b=1}^{K}\int g(t;\theta)q_b(t;\theta)f_b(\boldsymbol y_i\mid t;\theta)\,dt\right\}, \qquad \ell(\theta)=\sum_{i=1}^{n}\log p(\boldsymbol y_i;\theta).

Here $\theta$ collects the model quantities that are treated as unknown or variable, such as parameters of the progress distribution, branch probabilities, and observation distributions. The likelihood is the joint density of the observed data as a function of $\theta$; it is not a probability distribution over $\theta$. If cell sampling makes observations dependent, or if cells are drawn from different groups with different distributions, this simple product does not represent the full joint model without further structure.

The product assumes the observations are modeled on the scale actually used in the likelihood. For example, an expression vector of normalized continuous measurements could be assigned a continuous density, while raw count data would require a discrete probability model. A continuous density evaluated at an observed value is not the probability of observing that exact value; probabilities concern intervals or regions and are obtained by integrating the density.

Deriving the Marginal Contribution

The marginal expression follows directly from conditioning on latent variables. For a small progress interval around $t$, the joint contribution of progress near $t$, branch $b$, and observation $\boldsymbol y_i$ is represented by $g(t)q_b(t)f_b(\boldsymbol y_i\mid t)\,dt$. Integrating adds contributions over progress; summing adds contributions over branches. More formally, the law of total probability gives the following justification.

Condition on progress and branch.
For fixed $t$ and $b$, the observation has density $f_b(\boldsymbol y_i\mid t)$.

Weight by the latent-state distribution.
The joint density for progress $t$ and branch $b$ is $g(t)q_b(t)$, since $P(B_i=b\mid T_i=t)=q_b(t)$.

Remove the unobserved state.
Sum the weighted density over $b$, then integrate over all allowed progress values. The result is $p(\boldsymbol y_i)=\sum_b\int g(t)q_b(t)f_b(\boldsymbol y_i\mid t)\,dt$.

An equivalent parameterization assigns each branch a prior weight $\pi_b$ and a branch-specific progress density $g_b(t)$. In that case the joint density of branch and progress is $\pi_b g_b(t)$, and

p(\boldsymbol y_i)=\sum_{b=1}^{K}\pi_b\int g_b(t)f_b(\boldsymbol y_i\mid t)\,dt, \qquad \sum_{b=1}^{K}\pi_b=1,\quad \int g_b(t)\,dt=1.

These forms are consistent when the joint latent-state distribution is factored as $g(t)q_b(t)$. The choice of factorization should match the scientific model: branch probabilities can depend on progress in the first form, whereas the second form makes explicit a branch-specific progress distribution. Neither factorization removes the need to specify the observation model $f_b$.

Worked Examples

Worked Example: A Finite Progress-and-Branch Mixture

Suppose progress has two possible values, $t_1$ and $t_2$, each with probability $0.5$. At either value, branch A has probability $0.6$ and branch B has probability $0.4$. A single expression measurement $y$ has a normal observation distribution with standard deviation 1. At $t_1$, the means are 0 on A and 1 on B; at $t_2$, they are 2 on A and 3 on B. Calculate the marginal density for a hypothetical cell with $y=1$.

For a normal density with mean $\mu$ and standard deviation 1, $f(y\mid\mu)=\frac{1}{\sqrt{2\pi}}\exp[-(y-\mu)^2/2]$. At $y=1$, the densities for means 0 and 2 are both approximately $0.2420$, while those for means 1 and 3 are approximately $0.3989$ and $0.0540$, respectively. Thus,

\begin{aligned} p(y=1) &=0.5\{0.6(0.3989)+0.4(0.3989)\}\\ &\quad+0.5\{0.6(0.2420)+0.4(0.0540)\}\\ &=0.19945+0.08340\\ &=0.28285. \end{aligned}

The first pair of terms corresponds to progress $t_1$; the second pair corresponds to $t_2$. Within each pair, the branch probability weights the branch-specific density. The result is a density, not the probability of the exact measurement $y=1$.

Worked Example: Branch-Specific Observation Variance

Consider a cell whose progress is fixed at a known value for this calculation. There are two branches, with probabilities $q_A=0.7$ and $q_B=0.3$. The observed scalar measurement is $y=0.5$. Suppose branch A uses $N(0,1)$, and branch B uses $N(1,0.25)$, where the second argument is the variance. Calculate the mixture density, retaining the branch-specific variances.

For a normal distribution with mean $\mu$ and standard deviation $\sigma$, the density is $\phi(y;\mu,\sigma^2)=\frac{1}{\sqrt{2\pi}\sigma}\exp[-(y-\mu)^2/(2\sigma^2)]$. Branch A has standard deviation 1, giving $f_A(0.5)=\frac{1}{\sqrt{2\pi}}\exp(-0.125)\approx0.3521$. Branch B has standard deviation $0.5$, giving $f_B(0.5)=\frac{1}{\sqrt{2\pi}(0.5)}\exp(-0.5)\approx0.4839$. Therefore,

$$ p(y=0.5)=0.7(0.3521)+0.3(0.4839)=0.24647+0.14517=0.39164. $$

The narrower branch-B distribution assigns a relatively high density at this observation, but its contribution is still weighted by its branch probability. Replacing both distributions with one shared observation density would discard a specified part of the model and generally change the result.

Worked Example: From Cell Contributions to a Sample Likelihood

Use the two-progress, two-branch model from the first example for two hypothetical independent cells, with measurements $y_1=1$ and $y_2=2$. The first marginal density was $p(y_1)=0.28285$. For $y_2=2$, the densities at means 0, 1, 2, and 3 are approximately $0.0540$, $0.2420$, $0.3989$, and $0.2420$, respectively. The marginal contribution for the second cell is

\begin{aligned} p(y_2=2) &=0.5\{0.6(0.0540)+0.4(0.2420)\}\\ &\quad+0.5\{0.6(0.3989)+0.4(0.2420)\}\\ &=0.5(0.1292)+0.5(0.33614)\\ &=0.23267. \end{aligned}

Under the stated independent, common-mixture sampling assumption, the two-cell likelihood is the product of the marginal contributions: $L=0.28285(0.23267)\approx0.06581$. It is not obtained by adding the two densities, nor by summing over one shared branch label for both cells. The model assigns each cell its own latent progress and branch identity.

Continuous Progress and Numerical Evaluation

When progress is continuous, the integral may have no closed-form solution. A standard numerical approach is quadrature: choose progress points $t_m$ and nonnegative weights $w_m$ that approximate integration with respect to $g(t)$. The weights incorporate the progress distribution and are normalized to sum to one. The cell contribution is then approximated by

p(\boldsymbol y_i)\approx\sum_{m=1}^{M}w_m\sum_{b=1}^{K}q_b(t_m)f_b(\boldsymbol y_i\mid t_m).

This is a weighted average of observation densities across plausible latent states. Increasing the number or changing the placement of quadrature points can improve accuracy, especially when the integrand is sharply peaked or changes quickly near the branchpoint. A useful numerical check is to recompute contributions with a finer grid and verify that the log-likelihood changes negligibly at the precision needed for the analysis. All branch weights at a given point must sum to one, and the quadrature weights must represent the intended progress distribution.

For a multivariate normal observation model, one possible branch-specific specification is

f_b(\boldsymbol y\mid t)= \frac{\exp\left[-\frac{1}{2}(\boldsymbol y-\boldsymbol\mu_b(t))^\mathsf{T} \boldsymbol\Sigma_b(t)^{-1}(\boldsymbol y-\boldsymbol\mu_b(t))\right]} {(2\pi)^{d/2}|\boldsymbol\Sigma_b(t)|^{1/2}},

where $d$ is the number of measured features, $\boldsymbol\mu_b(t)$ is the branch-specific mean vector, and $\boldsymbol\Sigma_b(t)$ is a positive-definite covariance matrix. Branch-specific observation distributions need not differ only in their means: they may also differ in dispersion or dependence among features. The chosen family must be appropriate for the measurement scale and must have valid parameters at every modeled progress value.

Diagnostics, Alternatives, and Common Pitfalls

A likelihood can be calculated exactly and still describe the data poorly. Before relying on it, check that the branch probabilities are valid, the progress distribution is normalized, and each branch-specific observation distribution is defined over the data’s support. Inspect numerical stability: very small contributions can underflow when many cells are multiplied, so calculations are typically performed as a sum of log marginal densities. This uses $\log L=\sum_i\log p(\boldsymbol y_i)$ and avoids direct multiplication of many small numbers.

The marginal likelihood is also sensitive to the assumed sampling distribution of latent progress. If captured cells are enriched for particular states, a progress distribution intended to represent all cells may not describe the sampled cells. Likewise, a branch probability estimated or specified for a different population cannot automatically be treated as the branch probability among observed cells. These are modeling and sampling issues, not numerical integration errors.

A common alternative is a complete-data likelihood in which each cell’s $T_i$ and $B_i$ are treated as observed. Its contribution would be $g(T_i)q_{B_i}(T_i)f_{B_i}(\boldsymbol y_i\mid T_i)$. That expression is appropriate only when the latent states are known or when the analysis explicitly conditions on assigned states; substituting uncertain assignments as if they were known can understate uncertainty. The observed-data likelihood instead sums and integrates over the states. A hard assignment to one branch is another approximation, not the same calculation: it removes the contributions from other plausible branches.

Common mistake: “A cell is most likely on branch A, so I will use only branch A’s density and ignore branch B.” The correction is to retain all branch contributions in the marginal density: $\sum_b\int g(t)q_b(t)f_b(\boldsymbol y_i\mid t)\,dt$. A branch with a smaller prior probability can still contribute substantially if its observation density at the measured value is comparatively large. Conversely, the branch probability alone is not the likelihood contribution; it must be combined with the observation density and progress distribution.

Key takeaway: For each cell, the observed-data likelihood averages its branch-specific observation densities over the joint distribution of latent progress and branch identity. For independent cells, multiply those marginal contributions, or equivalently add their logarithms.

Check Your Understanding

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

Write the marginal density for one cell when progress has density $g(t)$, branch probabilities are $q_b(t)$, and branch-specific observation densities are $f_b(\boldsymbol y\mid t)$.
Why must branch-specific observation densities be included in addition to branch probabilities?
If two cells are independent draws from the same model and have marginal densities $0.20$ and $0.35$, what is their joint likelihood?
What changes in the likelihood calculation if expression measurements are discrete counts rather than continuous measurements?

Solutions. (1) The marginal density is $\sum_b\int g(t)q_b(t)f_b(\boldsymbol y\mid t)\,dt$. (2) Branch probabilities describe how likely the latent branch is, while $f_b$ describes how compatible the observed measurement is with that branch and progress; both are required by marginalization. (3) Under the stated independence assumption, the likelihood is $0.20(0.35)=0.070$. (4) Use a probability mass function for the count observation in place of a continuous density, and sum over discrete outcomes where marginalization over outcomes is required.