Why Multimodality Matters in QSP Calibration
A posterior distribution can assign substantial probability to two or more separated regions of parameter space. In a QSP model, this can occur when different parameter combinations produce similar predicted biomarker profiles, or when distinct mechanistic explanations fit the observations comparably well. A chain that moves readily within one region but rarely crosses to another may yield a stable-looking trace and a precise-looking summary that describe only part of the posterior.
The earlier material on posterior geometry in nonlinear ODE models motivates the possibility of separated regions; the material on rank-normalized R-hat and effective sample size provides complementary chain diagnostics. Here the focus is visual: inspect each chain’s trace, marginal density, and pair plots to assess whether modes are being explored and whether pooled summaries retain a useful interpretation.
A plot can suggest multimodality, but it cannot by itself establish that the chains have sampled the full posterior. The relevant assumption for ordinary posterior summaries is that the retained draws represent the target posterior, including the relative probability assigned to its important regions. In practice, this requires adequate exploration, not merely a large number of saved draws. MCMC draws are dependent, and chains trapped in separate modes do not become representative simply by being run longer without transitions.
What to Look for in the Plots
Start with a trace plot for each parameter, showing iterations on the horizontal axis and sampled values on the vertical axis. Plot chains separately or use clearly distinguishable colours. A chain that stays near one value range while another stays near a distant range is a warning that the chains may occupy different modes. Long stretches in one region separated by rare jumps can indicate difficult movement between modes. Conversely, a trace that moves rapidly within a region does not show that the chain has explored other regions.
Next inspect marginal density estimates for individual parameters, preferably with chain-specific densities as well as a pooled density. Separate peaks in a parameter’s marginal density are evidence that its posterior distribution may be multimodal. However, bandwidth and plotting range affect the appearance: excessive smoothing can merge peaks, while insufficient smoothing can create spurious bumps. Compare reasonable bandwidths and use the trace to check whether apparent peaks correspond to persistent regions of chain occupancy.
A pair plot places draws for two parameters against each other. It can reveal clusters or a curved relationship that is hard to see in either marginal alone. In mechanistic models, distinct combinations of clearance and potency, for example, may form separate clusters even when neither parameter’s marginal plot makes the scientific distinction obvious. Examine pair plots for parameters expected to trade off in the model, and colour points by chain or by time segment. A pooled pair plot can hide the fact that different chains occupy different clusters.
No single display is decisive. A chain-specific trace indicates where each chain travelled; a marginal density summarizes the sampled distribution of one parameter; a pair plot shows how two parameters vary jointly. Their combined evidence is more informative than any one pooled plot.
Worked Example: Chains Occupying Different Regions
Suppose four hypothetical chains each contain 1,000 retained draws of a dimensionless transformed parameter \(\theta\). The trace plots show chains 1 and 2 fluctuating near \(\theta=1\), and chains 3 and 4 fluctuating near \(\theta=5\), with no visible transitions between those regions. Their sample means are 1.02, 0.98, 4.98, and 5.02, respectively. For illustration, suppose the within-region standard deviation is about 0.10.
The pooled mean is the average of the four equally sized chain means:
A value near 3 is not where any chain spends appreciable time. The pooled mean is therefore a poor description of the sampled parameter values. The chain-specific traces and densities show the important feature: two groups of chains occupy separated regions. Equal numbers of chains in each region do not establish that the two posterior modes have equal probability. The chains’ allocation may reflect their starting values and failure to move between regions, rather than the posterior mass of each mode.
A pooled density would display two peaks if the draws are combined without adjustment, but that plot alone still cannot establish their relative posterior masses. Before reporting a pooled mean or interpreting a pooled density as representative of the target posterior, investigate why the chains did not transition and whether the sampler explored both regions adequately.
Why a Mean Can Be Misleading
The effect of mode separation on a summary can be shown by treating the posterior as a mixture of two mode-specific distributions. Let \(w\) be the posterior probability assigned to mode A, let \(1-w\) be the probability assigned to mode B, and let their conditional means and variances for a scalar parameter be \(\mu_A,\mu_B\) and \(\sigma_A^2,\sigma_B^2\). The posterior mean is the weighted average of the mode-specific means:
The total variance includes both variation within the modes and variation between their centres:
To see the decomposition, write the variance as the expected squared deviation from the overall mean and condition on the mode label \(M\). Within each mode, the squared deviation splits into deviation from that mode’s mean plus the squared difference between that mode’s mean and the overall mean. The within-mode terms average to \(w\sigma_A^2+(1-w)\sigma_B^2\). The remaining terms are \(w(\mu_A-E(\theta\mid y))^2+(1-w)(\mu_B-E(\theta\mid y))^2\). Substituting \(E(\theta\mid y)=w\mu_A+(1-w)\mu_B\) gives \(w(1-w)(\mu_A-\mu_B)^2\). Thus the between-mode contribution can be large even when each mode is individually narrow.
Worked Example: A Narrow Mode Pair With a Broad Overall Summary
Suppose a hypothetical scalar posterior has two equally weighted modes centred at 1 and 5, each with within-mode variance \(0.01\) (standard deviation \(0.10\)). The posterior mean is
The variance is
The overall standard deviation is \(\sqrt{4.01}\approx2.002\), despite a standard deviation of just \(0.10\) within each mode. The mean lies midway between the modes, where the posterior may have little density. Reporting only the mean and standard deviation would obscure the two concentrated regions; showing the density or mode-specific summaries is more informative.
Pair Plots and Mechanistic Combinations
A posterior mode can represent a combination of parameter values, not simply a high or low value of one parameter. Pair plots are useful because QSP calibration often involves compensation: changing one parameter can be offset by changing another, yielding similar model predictions. A plot of each parameter separately may not communicate which combinations co-occur.
Worked Example: Two Separated Parameter Combinations
Consider hypothetical draws for two transformed QSP parameters, \(a\) and \(b\). A pair plot shows two similarly sized clusters centred approximately at \((-1,1)\) and \((1,-1)\). Within each cluster, suppose both coordinates have variance \(0.04\). The cluster centres are opposites, so their equally weighted average is \((0,0)\), which lies between the clusters.
For the pooled draws, each coordinate has variance \(1+0.04=1.04\): one unit comes from the separation of the cluster centres and \(0.04\) from within-cluster variation. The covariance is \(-1\), because the centre pairs vary in opposite directions and within-cluster covariance is assumed to be zero. The resulting correlation is
The pair plot makes the two parameter combinations and their negative association visible. A single overall mean for each parameter, both near zero, would describe neither combination. In an actual analysis, inspect the chain identity and iteration order of the points: a pooled plot cannot show whether every chain visited both clusters or different chains remained in different ones.
When Posterior Predictions Also Need Care
Multimodality can matter beyond parameter summaries. If a model prediction is a nonlinear function of a parameter, evaluating the model at the posterior mean need not give the posterior mean prediction. Moreover, predictions from distinct modes may form separate groups, even if a single average prediction appears plausible.
Worked Example: Averaging Parameters Is Not Averaging Predictions
Suppose a hypothetical calibration has two well-explored modes for a dimensionless parameter \(\theta\): mode A is centred at 1 with posterior weight 0.70, and mode B is centred at 5 with weight 0.30. For illustration, define the model’s 24-hour biomarker prediction as \(C_{24}(\theta)=10+2\theta^2\) ng/mL.
The posterior mean parameter is
The predictions at the mode centres are \(C_{24}(1)=12\) ng/mL and \(C_{24}(5)=60\) ng/mL. Their posterior-weighted mean is
By contrast, evaluating the model at the posterior mean parameter gives
These quantities differ because the prediction is nonlinear. If predictions vary little within each mode, the prediction distribution will also have groups near 12 and 60 ng/mL. A single average does not communicate that separation. This calculation is illustrative; in an analysis, compute predictions for each retained posterior draw when the draws adequately represent the posterior.
Diagnostics, Follow-Up, and Limitations
Treat visual evidence as a reason to investigate, not as a verdict. Use the following checks together:
- Compare chains separately. Check whether each chain visits the same regions and whether apparent transitions are sustained or isolated.
- Compare chain-specific densities. Distinct density shapes or peaks can reveal chain separation that a pooled plot hides.
- Inspect relevant pairs. Colour pair plots by chain and, when helpful, by iteration segment to distinguish clusters from gradual movement or a curved relationship.
- Check related diagnostics. Rank-normalized R-hat and effective sample size can provide important evidence, but they do not replace inspection of mode exploration. Their interpretation depends on the chains representing the target posterior.
- Investigate the model and computation. Check whether parameter constraints, transformations, or mechanistic trade-offs explain the clusters. Consider whether dispersed initial values or a different sampler configuration provides evidence of movement between regions.
A common mistake is to pool all draws, observe a stable mean and a large draw count, and conclude that the posterior is well explored. That reasoning confuses the amount of saved output with coverage of the target distribution. If chains are trapped in separate modes, their pooled mean can lie in a low-density region, and their relative counts may not estimate the modes’ posterior masses.
The converse error is also possible: a chain-specific trace with long stretches at different values does not prove that the posterior has multiple modes. Slow movement along a continuous, highly correlated region, poor scaling, or a transient before stationarity can produce similar plots. A density estimate can also change substantially with smoothing choices. Investigate these alternatives before labelling a pattern as multimodality.
Finally, separated modes need not imply distinct scientific conclusions. They may reflect a symmetry or alternative parameter representations that yield nearly identical quantities of interest. Conversely, modes that imply different mechanisms or predictions may be scientifically consequential even when a scalar average appears unremarkable. Report the evidence for separation, whether chains explored the regions, and which parameter combinations or predictions differ. If mode weights are uncertain because transitions are rare, do not present pooled mode proportions as reliable posterior probabilities.
Check Your Understanding
Answer each question, then compare with the solutions.
- Two chains remain near different parameter values for their entire retained runs. What does this suggest, and what does it not establish?
- For two equally weighted modes with means 2 and 8 and within-mode variances \(0.25\), calculate the overall mean and variance.
- Why might a pair plot be more informative than two separate marginal density plots for mechanistic parameters?
- A model prediction is nonlinear in a parameter. Is prediction at the posterior mean generally equal to the posterior mean prediction? Explain.
- Why is the fraction of draws in each region not necessarily a reliable estimate of posterior mode weight when chains do not transition between regions?
Solutions. (1) It suggests chain separation or trapping in different regions. It does not prove that the target posterior is multimodal; slow exploration or nonstationarity are alternatives. (2) The mean is \(0.5(2)+0.5(8)=5\). The variance is \(0.5(0.25)+0.5(0.25)+0.25(8-2)^2=0.25+9=9.25\). (3) A pair plot shows which parameter values occur together and can reveal separated combinations or trade-offs hidden by separate one-dimensional views. (4) No. For a nonlinear function, \(E[f(\theta)\mid y]\) generally differs from \(f(E[\theta\mid y])\); compute the function on posterior draws when sampling is adequate. (5) If chains remain in their starting regions, the draw counts reflect chain allocation and initialization as well as posterior mass, so they need not estimate the modes’ relative probabilities.