Propagating Random-Effect Variability Into Exposure Metrics

Why Propagate Random-Effect Variability?

A typical parameter set describes a reference individual, not the spread of exposures across individuals. If clearance and volume vary between individuals—and their random effects are correlated—then AUC, peak concentration, and half-life inherit that variability in different ways. AUC depends on clearance, peak concentration after an intravenous bolus depends on volume, and half-life depends on their ratio. In particular, the correlation between clearance and volume can substantially change half-life variability even when the marginal variability of each parameter stays fixed.

The calculations below use a one-compartment model after a single intravenous bolus dose. This deliberately simple setting makes the propagation transparent; the same general strategy applies to other exposure metrics by evaluating them for draws of the individual parameters.

Model and Assumptions

Definition: Let $CL_{\mathrm{typ}}$ and $V_{\mathrm{typ}}$ denote typical clearance and volume. For individual $i$, let $CL_i=CL_{\mathrm{typ}}\exp(\eta_{CL,i})$ and $V_i=V_{\mathrm{typ}}\exp(\eta_{V,i})$, where $(\eta_{CL,i},\eta_{V,i})\sim N(\mathbf{0},\Omega)$ and $$ \Omega= \begin{pmatrix} \omega_{CL}^{2} & \omega_{CL,V}\\ \omega_{CL,V} & \omega_{V}^{2} \end{pmatrix}. $$ Here $\omega_{CL}^{2}$ and $\omega_V^{2}$ are variances on the log-parameter scale, and $\omega_{CL,V}$ is their covariance. The corresponding correlation is $\rho_{CL,V}=\omega_{CL,V}/(\omega_{CL}\omega_V)$. For dose $D$, the exposure metrics are $AUC_i=D/CL_i$, $C_{\max,i}=D/V_i$, and $t_{1/2,i}=\ln(2)V_i/CL_i$.

This model assumes positive clearance and volume, represented by exponentiating normally distributed random effects. It also assumes the stated one-compartment, linear-elimination setting: dose is delivered instantaneously, parameters remain constant during the interval, and AUC is the total exposure from the dose. The formulas for $C_{\max}$ refer to the concentration immediately after an intravenous bolus. They are not formulas for an oral peak, whose timing and value also depend on absorption.

The covariance matrix must be valid: its diagonal entries must be nonnegative, and $\omega_{CL,V}^{2}\leq \omega_{CL}^{2}\omega_V^{2}$. Equivalently, when both marginal variances are positive, $-1\leq\rho_{CL,V}\leq1$. These conditions ensure that the proposed joint normal distribution exists.

Exact Propagation on the Log Scale

Taking logarithms turns each metric into a linear combination of the random effects. That is useful because any linear combination of jointly normal random variables is normal.

\begin{aligned} \ln(AUC_i)&=\ln(D/CL_{\mathrm{typ}})-\eta_{CL,i},\\ \ln(C_{\max,i})&=\ln(D/V_{\mathrm{typ}})-\eta_{V,i},\\ \ln(t_{1/2,i})&=\ln\!\left(\ln(2)V_{\mathrm{typ}}/CL_{\mathrm{typ}}\right) +\eta_{V,i}-\eta_{CL,i}. \end{aligned}

Thus, the metrics are lognormal under the stated assumptions. Their typical values, obtained by setting both random effects to zero, are also their medians. They are not generally their arithmetic means. For a lognormal variable $M$ with $\ln(M)\sim N(m,v)$, its median is $\exp(m)$, its arithmetic mean is $\exp(m+v/2)$, and its coefficient of variation is $\sqrt{\exp(v)-1}$. The mean follows from the normal moment-generating function: $E[\exp(Y)]=\exp(m+v/2)$ for $Y\sim N(m,v)$. Also $E[M^2]=\exp(2m+2v)$, so $\operatorname{Var}(M)=\exp(2m+v)(\exp(v)-1)$; dividing by $E[M]^2$ gives the stated squared coefficient of variation.

Applying the variance rule for linear combinations gives:

\begin{aligned} \operatorname{Var}(\ln AUC_i)&=\omega_{CL}^{2},\\ \operatorname{Var}(\ln C_{\max,i})&=\omega_V^{2},\\ \operatorname{Var}(\ln t_{1/2,i})&=\omega_V^{2}+\omega_{CL}^{2}-2\omega_{CL,V}. \end{aligned}

The minus sign in the half-life expression matters: $t_{1/2}$ increases with volume but decreases with clearance. Positive covariance between the two log random effects therefore reduces the variance of their difference. The marginal distributions of AUC and $C_{\max}$ do not depend on this covariance, although their joint distribution with each other or with half-life does.

Worked Example: Analytic AUC, Peak, and Half-Life Variability

Suppose, hypothetically, that $D=100$ mg, $CL_{\mathrm{typ}}=5$ L/h, and $V_{\mathrm{typ}}=50$ L. Let $\omega_{CL}=0.30$, $\omega_V=0.20$, and $\rho_{CL,V}=0.60$. The covariance is $0.60(0.30)(0.20)=0.0360$, so the half-life log-variance is $0.30^2+0.20^2-2(0.0360)=0.0580$.

For AUC, the typical value is $100/5=20.0$ mg·h/L and the log-variance is $0.30^2=0.0900$. Its population arithmetic mean and coefficient of variation are $20.0\exp(0.0900/2)=20.92$ mg·h/L and $\sqrt{\exp(0.0900)-1}=0.3069$, or 30.69%.

For the immediate post-bolus peak, the typical value is $100/50=2.00$ mg/L and the log-variance is $0.20^2=0.0400$. The arithmetic mean is $2.00\exp(0.0400/2)=2.040$ mg/L, and the coefficient of variation is $\sqrt{\exp(0.0400)-1}=0.2020$, or 20.20%.

The typical half-life is $\ln(2)(50/5)=6.931$ h. Its log-variance is 0.0580, giving an arithmetic mean of $6.931\exp(0.0580/2)=7.135$ h and a coefficient of variation of $\sqrt{\exp(0.0580)-1}=0.2444$, or 24.44%. Each arithmetic mean exceeds its typical value because the corresponding metric is lognormal with positive log-scale variance.

Correlation Can Change Half-Life Variability

Holding the marginal variances fixed, the half-life log-variance can be written as $\omega_{CL}^{2}+\omega_V^{2}-2\rho_{CL,V}\omega_{CL}\omega_V$. This makes the direction of the correlation effect explicit. Positive correlation tends to make clearance and volume move together, reducing variability in their ratio; negative correlation tends to make them move in opposite directions, increasing that variability. The effect is on the distribution of the ratio, not on the typical half-life.

Worked Example: Comparing Correlation Scenarios

Keep the typical parameters and marginal standard deviations from the previous hypothetical example. The typical half-life remains 6.931 h in all three scenarios. For $\rho_{CL,V}=0.60$, the covariance is 0.0360 and the log-variance is $0.0900+0.0400-2(0.0360)=0.0580$; the half-life coefficient of variation is $\sqrt{\exp(0.0580)-1}=0.2444$, or 24.44%.

If $\rho_{CL,V}=0$, the covariance is zero and the log-variance is $0.0900+0.0400=0.1300$. The coefficient of variation is $\sqrt{\exp(0.1300)-1}=0.3726$, or 37.26%. If $\rho_{CL,V}=-0.60$, the covariance is $-0.0360$ and the log-variance is $0.0900+0.0400-2(-0.0360)=0.2020$. The coefficient of variation is $\sqrt{\exp(0.2020)-1}=0.4731$, or 47.31%.

Thus, changing only the correlation from positive to negative nearly doubles the half-life coefficient of variation, while leaving the marginal AUC and peak distributions unchanged. These are model-implied population distributions, not estimates of uncertainty in the typical parameters.

Simulation When Metrics Are More Complicated

For these three simple metrics, the lognormal results are exact and simulation is not needed to obtain their theoretical moments. Simulation is useful when a target metric has no convenient closed form, when a model has more parameters, or when the quantity of interest is a nonlinear function of a full concentration-time profile. The simulation must preserve the estimated dependence structure: drawing clearance and volume random effects independently would impose zero covariance and could give a misleading half-life distribution.

Specify the joint random-effect distribution.
Use the fitted or assumed $\Omega$, and check that it is a valid covariance matrix.

Generate correlated draws.
For each simulated individual, draw $(\eta_{CL},\eta_V)$ jointly from $N(\mathbf{0},\Omega)$.

Transform to individual parameters and metrics.
Calculate $CL$, $V$, AUC, $C_{\max}$, and half-life using the model equations.

Summarize and check stability.
Compute the desired summaries across draws. Increase the number of draws and compare summaries to assess Monte Carlo stability.

Worked Example: A Small Correlated-Draw Simulation

For a transparent, hypothetical illustration, take the same parameters as above. A bivariate-normal draw can be generated using independent standard normal values $z_1,z_2$ as $\eta_{CL}=0.30z_1$ and $\eta_V=0.20(0.60z_1+0.80z_2)$. This construction gives marginal standard deviations 0.30 and 0.20 and correlation 0.60. Consider three illustrative pairs $(z_1,z_2)=(-1,0),(0,1),(1,-1)$. A real analysis would use many more draws.

The resulting random effects are $(-0.30,-0.12),(0,0.16),(0.30,-0.04)$. For the first draw, $AUC=20\exp(0.30)=27.00$ mg·h/L, $C_{\max}=2\exp(0.12)=2.255$ mg/L, and $t_{1/2}=6.931\exp(0.18)=8.298$ h. Applying the same transformations to all three gives:

Draw	$\eta_{CL}$	$\eta_V$	AUC (mg·h/L)	$C_{\max}$ (mg/L)	Half-life (h)
1	-0.30	-0.12	27.00	2.255	8.298
2	0.00	0.16	20.00	1.704	8.135
3	0.30	-0.04	14.82	2.082	4.934

For example, the simulated AUC mean is $(27.00+20.00+14.82)/3=20.61$ mg·h/L. The half-life mean is $(8.298+8.135+4.934)/3=7.122$ h. These three draws are only a demonstration of the calculation; their summaries are not stable estimates of the population moments. With many draws, simulated summaries should approach the analytic benchmarks, such as the 20.92 mg·h/L arithmetic mean AUC in the preceding example.

Diagnostics, Alternatives, and Limitations

A useful check is to compare simulation summaries with the analytic results for any metric whose distribution is known. Material disagreement can signal too few draws, a covariance matrix implemented incorrectly, a sign error in the metric transformation, or a mismatch between the simulated and intended random-effect distribution. Also inspect the generated random effects: their empirical means, variances, and correlation should be reasonably close to the specified values, with discrepancies expected to shrink as the number of draws grows.

The exact lognormal calculation is efficient and gives clear formulas for medians, arithmetic means, and coefficients of variation. A first-order delta-method calculation is an alternative approximation: for a smooth function $g(\boldsymbol{\eta})$, it approximates the variance of $g$ near $\boldsymbol{\eta}=\mathbf{0}$ by $\nabla g(\mathbf{0})^\mathsf{T}\Omega\nabla g(\mathbf{0})$. This can be useful when an exact distribution is unavailable, but approximating variability on the original scale may be inaccurate for large random-effect variances or strongly nonlinear metrics. Direct simulation accommodates more complicated calculations, but its precision depends on the number of draws and it does not correct an incorrect model.

All approaches inherit the assumptions used to define the individual parameters and the structural model. The lognormal results rely on normally distributed log random effects; they need not describe a population with a different random-effect distribution. The bolus peak formula does not apply directly to oral dosing, and the stated AUC formula presumes linear elimination and the specified dose setting. Residual assay error is not included in these population metric calculations: it affects observed concentrations, while the propagation here concerns variability induced by the parameter random effects.

Common mistake: using the marginal variances alone for half-life. The incorrect calculation is $\omega_{CL}^{2}+\omega_V^{2}$ regardless of dependence. The correct expression includes $-2\omega_{CL,V}$, because the log half-life contains the difference $\eta_V-\eta_{CL}$.
Common mistake: calling the typical value the population mean. Setting random effects to zero gives the median of these lognormal metrics. The arithmetic mean is larger and includes the factor $\exp(v/2)$.
Common mistake: simulating the two effects independently by default. Independent draws set their covariance to zero. Use joint draws reflecting the specified covariance whenever dependence is part of the model.

Key takeaway: Write the exposure metric as a function of the individual parameters, then propagate the full random-effect covariance through that function. For lognormal clearance and volume, AUC, bolus peak, and half-life have exact lognormal distributions; their log variances show directly where correlation matters.

Check Your Understanding

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

If the log-clearance variance is 0.09, what is the log-variance of AUC?
With $\omega_{CL}^{2}=0.09$, $\omega_V^{2}=0.04$, and $\omega_{CL,V}=0.02$, find the log-variance of half-life.
For a metric with typical value 10 units and log-variance 0.04, find its arithmetic mean and coefficient of variation.
Why can changing clearance–volume correlation alter half-life variability without changing the marginal AUC distribution?
What should be checked if simulation summaries differ substantially from the analytic AUC mean?

Solutions. (1) AUC has log-variance equal to the log-clearance variance, 0.09. (2) The half-life log-variance is $0.09+0.04-2(0.02)=0.09$. (3) The mean is $10\exp(0.04/2)=10.20$ units, and the coefficient of variation is $\sqrt{\exp(0.04)-1}=0.2020$, or 20.20%. (4) Half-life depends on the difference of the log random effects, so its variance includes their covariance; AUC depends only on clearance, so its marginal distribution uses only the clearance variance. (5) Check the number of draws, the random-effect transformation and covariance, and the AUC formula; also compare the simulated AUC mean with the exact lognormal mean.