Laplace Approximation for Sparse-Data Posteriors

Why Approximate a Sparse-Data Posterior?

Individual pharmacokinetic data are often sparse: a patient may have only a few measured concentrations, while the model contains several individual parameters. Bayesian updating combines the individual likelihood with a population-based prior to obtain a posterior distribution. The posterior mode gives one useful summary, but it does not by itself describe how uncertain the individual parameters remain.

The Laplace approximation summarizes that uncertainty by fitting a Gaussian distribution near the posterior mode. Its key input is local curvature: a sharply curved log posterior implies a narrow approximation, while a flatter one implies a wider approximation. The approximation can be convenient for uncertainty calculations and forecasting, but local curvature is not a substitute for examining the whole posterior. With sparse data, the posterior may be skewed, bounded, or have multiple separated regions of high density.

The Hessian is the matrix of second partial derivatives of the log posterior. Its negative, \(J\), is sometimes called the observed posterior curvature matrix. The approximation is based on a second-order Taylor expansion at the mode. It requires a smooth log posterior near that mode, a finite, isolated interior maximum, and positive curvature in every direction there. Those are local mathematical conditions; for the approximation to summarize the posterior well, the posterior mass must also be concentrated in a region that is reasonably represented by that local quadratic shape.

From the Hessian to a Gaussian

Write \(\boldsymbol{\delta}=\boldsymbol{\eta}-\hat{\boldsymbol{\eta}}\). Because the gradient of \(g\) is zero at an interior mode, the second-order Taylor expansion is

Exponentiating gives a density proportional to \(\exp[-\boldsymbol{\delta}^{T}J\boldsymbol{\delta}/2]\). This is the kernel of a multivariate normal density centered at \(\hat{\boldsymbol{\eta}}\), with precision matrix \(J\) and covariance matrix \(J^{-1}\). In one dimension, if the log posterior has second derivative \(g''(\hat{\eta})<0\), the approximate variance is \(-1/g''(\hat{\eta})\).

The derivation identifies both the strength and the limitation of the method. It uses the posterior’s height and curvature at one point. It does not use the full shape of distant tails, detect a second mode from the curvature at the first, or guarantee that a locally Gaussian approximation gives sensible probabilities far from the mode.

Intervals and Parameter Coordinates

For a scalar parameter in the chosen coordinates, an approximate central 95% posterior interval is \(\hat{\eta}\pm1.96\sqrt{J^{-1}}\). For a vector, the approximate covariance \(J^{-1}\) supplies marginal standard deviations and correlations. A multivariate Gaussian also defines joint ellipsoidal regions, but a marginal interval for one parameter is not the same object as a joint region for the whole vector.

The coordinates matter. Positive parameters are often represented on a log scale, and correlated random effects may be represented in unconstrained coordinates. If a density is transformed, it must include the appropriate change-of-variables Jacobian. The mode of a density generally changes under a nonlinear transformation; transforming a mode in one coordinate system is not necessarily the mode of the density in another.

When Sparse-Data Uncertainty Is Poorly Approximated

Sparse information can leave some parameter directions weakly constrained by the likelihood. The prior may then have substantial influence, and the posterior can retain features that a local quadratic expansion does not capture. Even if the Hessian can be computed and is positive definite, that alone does not establish that the approximation is accurate.

Common warning signs include a visibly skewed marginal posterior, a long tail, a parameter close to a constraint, strong curvature changes across the region of interest, or separated regions of high posterior density. A nearly flat direction can also make the Hessian ill-conditioned: small numerical errors may produce large changes in its inverse. If the Hessian has a non-positive eigenvalue at a purported mode, it is not a valid positive-definite precision matrix for this approximation; investigate the optimization, parameterization, and local posterior shape rather than treating the inverse as a covariance.

Diagnostics, Alternatives, and Limits

A Laplace approximation is most useful when the posterior has one dominant, smooth, interior mode and is close to quadratic over the region containing most of its mass. Before relying on it for an individual forecast or dose decision, check that the optimization found a genuine mode, that the Hessian is symmetric and positive definite up to numerical tolerance, and that its eigenvalues are not so small or ill-conditioned that the estimated covariance is unstable. Inspect the implied intervals for impossible values and compare summaries across sensible parameterizations when transformation is relevant.

The principal alternatives answer the same uncertainty question with different compromises. Direct numerical integration can characterize a low-dimensional posterior without imposing a Gaussian shape, although it becomes difficult as dimension grows. Posterior simulation can represent skewness and multiple modes when exploration is adequate, at additional computational and diagnostic cost. A Gaussian approximation is fast and compact, but can miss tails and modes. None of these methods repairs a misspecified likelihood, an inappropriate prior, or data that do not identify a parameter; they represent the posterior implied by the specified model.