Propagating Cross-Species Parameter Uncertainty

Why Propagate Cross-Species Uncertainty?

A human PK/PD prediction combines quantities that are uncertain for different reasons: an animal parameter is estimated from finite data, a scaling relationship has uncertain coefficients, and the animal and human may differ in ways that are not fully captured by the scaling model. Reporting only a single projected concentration or response hides how those uncertainties combine. In nonlinear PK/PD models, even the average input values may not produce the average predicted response.

The introductory material on mechanistic PBPK-PD integration for human projection establishes how a mechanistic model maps inputs to human concentration and response. Here the focus is on carrying uncertainty through that mapping. The result is a distribution of predictions under explicit assumptions—not a guarantee that the model includes every source of uncertainty.

Definition: Uncertainty propagation is the process of representing uncertainty in model inputs by a joint distribution, evaluating the prediction for draws from that distribution, and summarizing the resulting distribution of predicted quantities. A predicted response profile is the set of predicted responses across specified times or conditions, with uncertainty summarized at each point or for a defined profile-level quantity.

Keep three sources conceptually distinct. Parameter uncertainty concerns unknown values in the chosen model, such as an animal clearance estimate or an allometric exponent. Individual variability concerns real differences among people, such as differences in clearance between patients. Model discrepancy concerns systematic inadequacy of the model, such as an unrepresented species-specific pathway. A prediction distribution includes only the sources represented in its inputs or model. Adding a distribution for patient variability does not, by itself, account for model discrepancy.

Specify the Joint Uncertainty Before Simulating

Write the human prediction as \(R(t;\boldsymbol{\theta})\), where \(t\) is time and \(\boldsymbol{\theta}\) contains uncertain inputs: animal parameter estimates, scaling coefficients, species-difference terms, and any other quantities treated as uncertain. For each draw of the entire input vector, run the same specified PK/PD model and retain the predicted response at the times of interest.

Assumptions for a simulation-based uncertainty summary: The model equations and input definitions are fixed; the joint distribution assigned to \(\boldsymbol{\theta}\) reasonably represents the uncertainty being considered; dependence among inputs is represented; and model evaluation is numerically reliable for the simulated draws. The resulting summaries are conditional on these choices. They do not automatically include uncertainty in the model structure or unrepresented biological processes.

The inputs must be drawn jointly, not independently by default. For example, an animal clearance estimate and an allometric exponent may be statistically related because both were estimated using overlapping data. Species-specific terms may also be correlated with scaled parameters. Ignoring such dependence changes the uncertainty in the prediction, sometimes substantially. Use a joint sampling distribution, bootstrap draws, or another justified joint representation rather than combining marginal distributions without checking their dependence.

For \(M\) simulated input vectors, compute \(R_m(t)=R(t;\boldsymbol{\theta}_m)\), for \(m=1,\ldots,M\). At a specified time, the sample mean of the predictions estimates the distribution's mean; the median and chosen quantiles describe its center and spread in ways that may be more suitable for skewed predictions. Pointwise quantiles across time form an uncertainty band for the response profile. They summarize uncertainty at each time separately; they are not necessarily a simultaneous band that covers the entire trajectory at a stated rate.

\boldsymbol{\theta}_m \sim F_{\boldsymbol{\theta}},\qquad R_m(t)=R(t;\boldsymbol{\theta}_m),\qquad \widehat{q}_{\gamma}(t)=\text{empirical } \gamma\text{ quantile of } \{R_1(t),\ldots,R_M(t)\}.

Here \(F_{\boldsymbol{\theta}}\) is the specified joint uncertainty distribution. If its quantiles represent uncertainty about unknown parameters, a central interval between its 2.5th and 97.5th percentiles is a 95% interval under that distribution and model. It is not automatically a confidence interval with repeated-sampling coverage, nor does it include uncertainties that were omitted. State which inputs were varied and how the interval was constructed.

Worked Examples

Worked Example: Uncertainty in an Allometrically Scaled Clearance

Suppose a hypothetical animal-to-human calculation gives a log-scale human clearance prediction with mean \(4.919\), where the modeled quantity is \(\log(CL_H/(1\ \mathrm{L/h}))\). Assume three independent, normally distributed uncertainty contributions: the animal log-clearance estimate has standard deviation \(0.20\), the allometric exponent has standard deviation \(0.080\), and a log-scale species-difference term has standard deviation \(0.250\). The animal-to-human body-weight ratio is 280. The normality and independence assumptions are stipulated for this illustration; in an application they require justification.

For \(Y=\log(CL_H/(1\ \mathrm{L/h}))\), the exponent uncertainty is multiplied by the log weight ratio. Independence therefore gives

\operatorname{Var}(Y)=0.20^2+\left[0.080\log(280)\right]^2+0.250^2.

Since \(\log(280)\approx5.63479\), the exponent contribution to the standard deviation is \(0.080(5.63479)=0.450783\). Thus

\operatorname{Var}(Y)=0.040000+0.203205+0.062500=0.305705,\qquad SD(Y)=\sqrt{0.305705}=0.552906.

A central 90% interval under the assumed normal log-scale distribution uses \(1.645\) standard deviations on either side of the mean. Back-transforming its endpoints gives

\left[ \exp\{4.919-1.645(0.552907)\}, \exp\{4.919+1.645(0.552907)\} \right] = [\exp(4.00947),\exp(5.82853)] \approx[55.1,\ 340]\ \mathrm{L/h}.

This interval is wide and asymmetric on the clearance scale because exponentiation maps a symmetric log-scale interval to an asymmetric interval. Its interpretation is conditional: under the specified independent normal uncertainties and scaling model, the central 90% of the modeled clearance distribution lies approximately between 55.1 and 340 L/h. The calculation does not establish that the assumptions are realistic.

Worked Example: Simulating an Uncertain Response Profile

Consider a hypothetical one-compartment concentration profile \(C(t)=C_0e^{-kt}\), with fixed \(k=0.5\ \mathrm{h}^{-1}\). Let \(C_0\) be equally likely to be 1 or 2 mg/L, and let \(EC_{50}\) independently be equally likely to be 1 or 3 mg/L. At each time, predict response magnitude using \(R(t)=80C(t)/(EC_{50}+C(t))\), measured in percentage points. This deliberately small discrete joint distribution makes every simulated outcome inspectable.

The four equally likely joint draws are the four combinations of \(C_0\) and \(EC_{50}\). At \(t=0\), the concentrations are 1 or 2 mg/L. The four responses are

\frac{80(1)}{1+1}=40,\quad \frac{80(1)}{3+1}=20,\quad \frac{80(2)}{1+2}=53.33,\quad \frac{80(2)}{3+2}=32 \quad\text{percentage points}.

The equally weighted mean response at time zero is \((40+20+53.33+32)/4=36.33\) percentage points, and the simulated range is 20 to 53.33 percentage points.

At \(t=1\) h, the two possible concentrations are \(1e^{-0.5}=0.6065\) mg/L and \(2e^{-0.5}=1.2131\) mg/L. Substitution for the same four joint draws gives approximately \(30.20,\ 13.45,\ 43.85,\) and \(23.03\) percentage points. The mean, calculated from the unrounded responses, is approximately 27.64 percentage points, with a range of about 13.45 to 43.85 percentage points. Repeating this calculation at more times produces a pointwise summary of the profile. With only four draws, these values are an exact summary of the stipulated discrete example, not a stable approximation to a continuous uncertainty distribution.

Worked Example: Why Correlated Inputs Must Be Drawn Together

Suppose a hypothetical log-clearance prediction has the form \(Y=\log(CL_A)+5b\), where \(\log(CL_A)\) is an animal log-clearance estimate and \(b\) is an uncertain scaling exponent. Let their standard deviations be 0.20 and 0.10, respectively, with correlation \(0.60\). Assume no other uncertain terms. The covariance is \(0.60(0.20)(0.10)=0.0120\).

Using the variance formula for a sum of correlated random variables,

\operatorname{Var}(Y)=\operatorname{Var}\{\log(CL_A)\}+25\operatorname{Var}(b) +2(5)\operatorname{Cov}\{\log(CL_A),b\}.

Substitution gives \(\operatorname{Var}(Y)=0.20^2+25(0.10^2)+10(0.0120)=0.0400+0.250+0.120=0.410\), so \(SD(Y)=\sqrt{0.410}=0.6403\). If the two inputs were incorrectly treated as independent, the covariance term would be omitted, giving variance \(0.290\) and standard deviation \(0.5385\). Here the positive dependence increases the uncertainty in predicted log clearance. In simulation, the correction is to generate paired values from the joint distribution, preserving the correlation, and apply the prediction equation to each pair.

Derivation, Diagnostics, and Alternatives

The covariance calculation used above follows by expanding the variance. For random variables \(A\) and \(B\), and constants \(c\) and \(d\),

\begin{aligned} \operatorname{Var}(cA+dB) &=E\left[\{c(A-E[A])+d(B-E[B])\}^2\right]\\ &=c^2\operatorname{Var}(A)+d^2\operatorname{Var}(B) +2cd\operatorname{Cov}(A,B). \end{aligned}

The cross-term remains unless the covariance is zero. This result explains why independently sampling correlated inputs is not a harmless computational shortcut: it changes the distribution being propagated. More generally, simulation is valuable when the mapping from inputs to predictions is nonlinear, involves several linked PK/PD components, or produces skewed or bounded predictions for which a simple normal approximation is questionable.

An analytic approximation such as the delta method is a useful alternative when the prediction is smooth and input uncertainty is modest. It approximates the prediction locally using derivatives and the input covariance matrix. It can be faster and helps identify which local input contributions affect prediction variance. Its accuracy may deteriorate for strong nonlinearity, broad uncertainty, boundaries, or multimodal input distributions. Simulation instead evaluates the stated model repeatedly, but it cannot repair a misspecified input distribution or model.

Before interpreting simulated summaries, check that the computation is behaving as intended:

Inspect input draws. Confirm plausible ranges, units, transformations, and correlations. Unexpected negative clearances or impossible probabilities can indicate that an unconstrained distribution was used on an inappropriate scale.
Check the model output for every draw. Look for solver failures, non-finite values, discontinuities, or predictions outside biological or mathematical bounds. Do not silently discard failed draws; investigate their cause and report any justified exclusion rule.
Check Monte Carlo stability. Increase the number of draws and compare means and quantiles. Tail quantiles generally need more draws for stable estimation than central summaries. A smooth-looking profile alone does not demonstrate numerical convergence.
Check dependence and scale. Confirm that joint draws preserve intended dependence and that operations are on the correct scale. For example, a normal distribution for log clearance implies a lognormal clearance distribution after exponentiation.

A frequent mistake is to plug the mean animal estimate and mean scaling exponent into the model and call that the mean human prediction. In general, \(E\{R(t;\boldsymbol{\theta})\}\ne R(t;E[\boldsymbol{\theta}])\) when the prediction is nonlinear. The correction is to propagate the joint uncertainty through the model and summarize the resulting predictions. A second pitfall is to label a pointwise 95% interval at each time a 95% band for the whole curve; pointwise coverage does not imply simultaneous coverage of an entire profile.

Finally, distinguish an uncertainty interval for the model's predicted mean response from a predictive distribution for future individuals. The latter may additionally require a model for between-person variability and residual error. Neither type automatically covers uncertainty from omitted mechanisms or an incorrect species-translation assumption. The reported distribution should therefore identify its included uncertainty sources, its conditioning assumptions, and whether its summaries are pointwise or profile-level.

Key takeaway: Human translation uncertainty is determined by both the uncertainty in each input and the joint structure linking them. Preserve that joint structure, propagate it through the full prediction model, and interpret the resulting response distribution only for the uncertainties and assumptions actually represented.

Check Your Understanding

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

Why is it generally incorrect to vary an animal estimate and scaling exponent independently when their uncertainty is correlated?
For \(Y=A+3B\), suppose \(SD(A)=0.10\), \(SD(B)=0.20\), and \(\operatorname{Cov}(A,B)=0.005\). Find \(\operatorname{Var}(Y)\).
A simulation reports pointwise 95% intervals at every time. Does that establish 95% coverage for the complete response trajectory?
What does a simulation interval represent if it varies only parameter estimates and does not include between-person variability?

Solutions. (1) Independent draws remove the covariance present in the joint uncertainty distribution, changing the variance and potentially the shape of the predicted distribution. Draw paired values from a joint representation. (2) \(\operatorname{Var}(Y)=0.10^2+9(0.20^2)+2(3)(0.005)=0.010+0.360+0.030=0.400\). (3) No. The intervals describe pointwise uncertainty; simultaneous coverage of the whole trajectory requires a profile-level summary and a corresponding method. (4) It summarizes uncertainty due to the varied parameter estimates conditional on the model and its other fixed inputs. It is not a predictive distribution for future individuals unless individual variability and any required residual error are also represented.

Propagating Cross-Species Parameter Uncertainty

What You'll Learn