7 Differential equations
Models in this family specify a rate law and estimate its parameters. They are presented below in order of increasing parameter count, with the minimum design required for each.
7.1 D1. Uncoupled parametric rate law
Each region is modeled independently. The exponential law
\[ \frac{dv_r}{dt} = -k_r \, v_r \qquad \Longrightarrow \qquad \log v_r(t) = \log v_r(0) - k_r t \]
is the natural baseline: it is scale-free, guarantees positivity, and is linear on the log scale, so it can be fitted as a linear mixed model. Saturating alternatives — logistic, Gompertz — add an asymptote and an inflection point.
Estimation is by nonlinear mixed effects, with fixed effects at the population level and random effects on \(k\) and, where applicable, onset. Parameters: \(O(1)\) per region for fixed effects, plus a random-effect covariance.
Minimum design: \(T = 2\) for the exponential form with \(N\) large; \(T = 3\) or latent time for any form with curvature.
7.2 D2. Coupled linear dynamical system
Regions are modeled jointly:
\[ \frac{d\mathbf{x}}{dt} = \mathbf{A}\mathbf{x} + \mathbf{C}\mathbf{u}, \]
where \(\mathbf{A} \in \mathbb{R}^{R \times R}\) contains self- and cross-region couplings and \(\mathbf{u}\) are exogenous inputs. This is the form used by Ziegler et al. (2017) for gray matter volume in 26 bilateral regions with two to three scans per subject, with inputs representing a manifest pubertal factor and a hidden growth factor, and inference by Bayesian model comparison. That study reported predominantly negative inter-regional couplings and higher model evidence for within-network than between-network connectivity.
The difficulty is dimensional (Chapter 5). \(\mathbf{A}\) has \(R^2\) entries: 676 for 26 regions, 4{,}624 for 68, and 160{,}000 for a 400-region parcellation. Feasible constraints are sparsity penalties, low-rank factorization \(\mathbf{A} = \mathbf{U}\mathbf{V}^\top\) with \(\mathrm{rank} \ll R\), restriction to a predefined network partition, or an empirical Bayes prior shrinking toward a structural connectivity matrix.
Minimum design: \(T = 2\) with large \(N\) and heavy regularization for a constrained \(\mathbf{A}\); \(T \geq 4\) before a general \(\mathbf{A}\) is meaningfully constrained by within-subject information.
7.3 D3. Network diffusion
A structural prior replaces free coupling parameters:
\[ \frac{d\mathbf{x}}{dt} = -\beta \, \mathbf{L} \mathbf{x}, \]
where \(\mathbf{L}\) is the Laplacian of a connectome graph derived independently of the longitudinal data. The entire \(R\)-region system has one or two free parameters (Raj et al. 2012, 2015). Related network-based formulations discriminate among competing spread mechanisms rather than assuming diffusion (Zhou et al. 2012); see Chapter 6.
This is the most parsimonious member of the family and the one least affected by small \(T\). Because its parameter count does not grow with \(R\), it is also the arm most likely to outperform flexible models on extrapolation. Its limitation is that the prior is fixed: if the true coupling is not proportional to structural connectivity, the model cannot represent it.
Minimum design: \(T = 2\). Estimable from cross-sectional data with an assumed initial condition.
7.4 D4. Latent time reparametrization
Rather than a new rate law, this is a reparametrization applicable to D1–D3: each subject is assigned a time shift \(\tau_i\) and rate \(\alpha_i\) mapping observed time to population time (Oxtoby 2023). The population trajectory can then span a far longer interval than any individual’s follow-up.
This is the principal technique that makes nonlinear shape estimation feasible at \(T = 2\)–\(3\), and it should be included in any sparse-design comparison. Its identifiability conditions and anchoring requirements are discussed in Chapter 3.
Minimum design: \(T = 2\), with \(N\) large and good coverage of the stage range.
7.5 Methods that are not available at small \(T\)
The following are frequently proposed and are not identifiable from two or three observations per subject:
- Per-subject sparse regression on a library of candidate terms (SINDy and relatives). Requires dense sampling to estimate derivatives; at \(T = 3\) there is one noisy curvature estimate against a library of many candidate functions.
- Per-subject neural ODEs. The vector field is a neural network with more parameters than the subject has observations.
- Free-form nonparametric rate estimation. Requires enough points per subject to estimate \(dv/dt\) as a function of \(v\).
These become available in the regimes described in Chapter 13. Proposing them for a two- or three-visit design is a parameter-counting error and will be identified as such in review.
7.6 Estimation notes
- Collocation and gradient matching avoid repeated numerical integration by fitting a smooth interpolant and matching derivatives. They require enough points per subject to smooth, and are therefore unavailable at \(T \leq 3\); full likelihood with integration is required there.
- Multi-start optimization is not optional. In the closest systematic study of hybrid ODE estimation, single-start optimization failed to recover parameters even under favorable conditions (Philipps et al. 2025). The same applies to the mechanistic models here whenever the likelihood is non-convex.