Why Pair Scenario Simulations?
A simulation comparison often asks for a difference between two dosing or treatment scenarios: for example, the change in mean exposure, or the change in the probability of attaining a target. If each scenario is simulated using unrelated virtual patients and random variability, the two estimated outcomes fluctuate independently. Those fluctuations both contribute noise to the estimated contrast.
A paired comparison instead evaluates both scenarios using the same underlying random draws for each replicate. A virtual patient with a particular covariate and parameter draw in scenario A is matched to that same draw in scenario B. When the two scenario outputs respond similarly to those shared features, their difference can be less variable than a difference formed from unrelated simulations. This variance reduction improves Monte Carlo precision; it does not change the scenario-specific targets or make the pharmacometric model itself more certain.
Constructing a Valid Pair
A pair is meaningful only if its two outputs are generated from aligned underlying inputs. A simulation replicate may include a parameter draw, covariates, random effects, and residual variability. For a paired comparison, use the same draw for each input that represents the same simulated patient or uncertainty realization across scenarios. Apply each scenario’s own dosing regimen, model equations, and scenario-specific transformations to those shared inputs.
Sharing inputs is a computational coupling, not a claim that two different treatments occur simultaneously in one real patient. The aim is to preserve the correct distribution of each scenario’s output while making the two outputs move together enough to stabilize their contrast. In particular, the marginal distribution of \(Y_{Ai}\) must be the intended distribution under A, and the marginal distribution of \(Y_{Bi}\) must be the intended distribution under B. Shared draws must not replace a scenario-specific random mechanism with an inappropriate one.
For example, if residual error is generated by applying scenario-specific transformations to a standardized random draw, sharing that standardized draw may preserve both residual distributions while inducing a useful dependence. If the scenarios require distinct residual mechanisms, those mechanisms must still be represented correctly; sharing a draw is not a reason to force identical errors. The earlier material on conditional versus marginal random-effects sampling and on residual-error simulation is relevant to choosing the distributions being drawn. The earlier material on reproducible random-number streams supports repeatable implementation, but repeatability alone does not guarantee correctly aligned pairs.
State the outcome and contrast, such as the difference in expected AUC or in target-attainment probabilities, and specify the population and uncertainty draws represented.
For each \(i\), generate the parameter, covariate, and applicable variability draws intended to be common to scenarios A and B.
Use the same replicate’s shared inputs in both scenario calculations, while changing the scenario-specific regimen or other defined intervention.
Calculate \(D_i=Y_{Bi}-Y_{Ai}\) for each aligned pair, then estimate the mean contrast and its Monte Carlo standard error from the \(D_i\) values.
Why Pairing Can Reduce Monte Carlo Error
Suppose the replicate pairs are independent and identically distributed, with finite output variances. Within each pair, \(Y_A\) and \(Y_B\) may be dependent; that dependence is precisely what the coupling is designed to create. By the variance rule for a difference,
The sample standard deviation of the paired differences is \(s_D\). The Monte Carlo standard error of their mean is \(SE(\bar{D})=s_D/\sqrt{N}\). If the outputs are positively correlated within pairs, the covariance term reduces the variance of each difference. If their marginals are the same as in independent simulation, an unpaired difference of sample means has variance \(\{\operatorname{Var}(Y_A)+\operatorname{Var}(Y_B)\}/N\). Pairing is more precise when its covariance is positive; it is less precise if the induced covariance is negative.
The result follows directly by expanding the variance of the difference. For independent replicate pairs, the variance of a mean is the variance of one replicate divided by \(N\). Therefore:
This derivation also shows why the pairing structure matters: the covariance must be calculated between outputs from the same replicate. Covariance between outputs matched arbitrarily after simulation is not the covariance induced by the CRN design. In practice, the Monte Carlo uncertainty is estimated from the observed paired differences, \(SE(\bar{D})=s_D/\sqrt{N}\), rather than by treating the two scenario means as independent.
Worked Example: Comparing Scenario Means From Paired Outputs
Suppose a hypothetical simulation produces the following exposure metric, in mg·h/L, for five matched replicates. The values are illustrative only.
| Replicate | Scenario A | Scenario B | \(D_i=B-A\) |
|---|---|---|---|
| 1 | 8 | 6 | -2 |
| 2 | 10 | 9 | -1 |
| 3 | 9 | 7 | -2 |
| 4 | 11 | 10 | -1 |
| 5 | 12 | 9 | -3 |
The scenario means are \(\bar{Y}_A=50/5=10.0\) mg·h/L and \(\bar{Y}_B=41/5=8.2\) mg·h/L. Thus the estimated contrast is \(\bar{D}=8.2-10.0=-1.8\) mg·h/L. The paired differences are \(-2,-1,-2,-1,-3\); their mean is \(-1.8\).
The squared deviations of the differences from \(-1.8\) sum to \(2.8\), so their sample variance is \(s_D^2=2.8/(5-1)=0.700\) (mg·h/L)\(^2\), and \(s_D=\sqrt{0.700}=0.837\) mg·h/L. The Monte Carlo standard error is
For comparison, the sample variances in A and B are \(2.50\) and \(2.70\) (mg·h/L)\(^2\). Treating the scenario means as independent would give \(SE=\sqrt{2.50/5+2.70/5}=1.02\) mg·h/L. The paired estimate is much less noisy for these outputs because large and small exposure realizations tend to occur together across the scenarios. The contrast estimates a reduction of 1.8 mg·h/L; its standard error describes simulation error in estimating that contrast, not patient-level variability or uncertainty in the model being simulated.
Contrasts Beyond Continuous Outcomes
The paired method applies to any replicate-level outcome for which a contrast is defined. For a threshold-attainment probability, let \(I_{Ai}\) equal 1 if replicate \(i\) attains the threshold under A and 0 otherwise, and define \(I_{Bi}\) similarly. The estimated difference in probabilities is the mean of \(D_i=I_{Bi}-I_{Ai}\). A pair that attains under both scenarios contributes zero to the difference, as does a pair that attains under neither. Discordant pairs contribute \(+1\) or \(-1\). The standard error is calculated from these paired differences, not by assuming independent estimated proportions.
A relative contrast may be more naturally expressed on a log scale. If both scenario outputs are strictly positive, define \(D_i=\log(Y_{Bi})-\log(Y_{Ai})=\log(Y_{Bi}/Y_{Ai})\). The mean log contrast describes an average log ratio; exponentiating it gives the geometric mean ratio under the simulated pairing. The scale and interpretation should be selected before examining the results, because a difference and a ratio answer different questions.
Worked Example: Paired Threshold-Attainment Probabilities
In six hypothetical matched replicates, suppose the binary outcomes for a target exposure threshold are as follows. Values of 1 indicate attainment.
| Replicate | Attains under A | Attains under B | \(D_i=B-A\) |
|---|---|---|---|
| 1 | 1 | 1 | 0 |
| 2 | 0 | 1 | 1 |
| 3 | 1 | 1 | 0 |
| 4 | 0 | 0 | 0 |
| 5 | 1 | 0 | -1 |
| 6 | 0 | 0 | 0 |
Each scenario attains the threshold in three of six replicates, so each estimated probability is \(3/6=0.500\), and the estimated difference is \(0.500-0.500=0\). The paired differences are \(0,1,0,0,-1,0\), with mean zero. Their squared deviations sum to 2, giving \(s_D^2=2/(6-1)=0.400\) and \(SE(\bar{D})=\sqrt{0.400/6}=0.258\). The simulation gives no estimated difference in attainment probability, but this small example has appreciable Monte Carlo error. Equal marginal estimates do not imply that the paired outcomes were identical: two replicates are discordant in opposite directions.
Checking Whether Pairing Is Working
After a run, inspect the structure and the effect of pairing rather than assuming CRN improved precision. Useful checks include:
- Verify alignment. Each scenario output should carry the same replicate ID for the shared inputs that generated it. Sorting, filtering, or joining outputs incorrectly can silently pair unrelated replicates.
- Check the marginal targets. Confirm that each scenario uses its intended parameter, covariate, and variability distributions. A smaller standard error is not evidence that those distributions were implemented correctly.
- Estimate within-pair dependence. Compare the scenario outputs within pairs, or compare the variance of paired differences with the independent-simulation variance estimate. Weak or negative covariance means little gain or a possible loss in precision.
- Inspect the contrast distribution. Look for a few unusually large differences or patterns that suggest scenario-specific failures. A mean and standard error can conceal influential replicates.
- Check implementation under changes in execution. Branches that consume different numbers of random draws can break alignment if the code relies on sequential stream position. Explicitly indexed draws or pre-generated shared inputs can make the intended pairing clearer.
A useful empirical summary is the variance ratio: the estimated variance of paired differences divided by the sum of the two scenario sample variances. Values below 1 indicate reduced variance relative to an independent-scenario difference with the same number of replicates, while values above 1 indicate increased variance. This is a diagnostic for the chosen coupling and scenario pair, not a guarantee that CRN will help other contrasts.
Comparison, Limitations, and Common Mistakes
| Approach | How the contrast is formed | Main implication |
|---|---|---|
| Paired CRN simulation | Difference within each shared-input replicate, then average differences | Uses within-pair covariance; often reduces Monte Carlo error when outputs move together |
| Independent scenario simulations | Difference between two separately simulated means | Simple to implement, but has no induced within-pair covariance to reduce contrast variance |
| Scenario-specific coupling | Shares only selected inputs that are genuinely common across scenarios | Can better represent differing mechanisms, but may provide less variance reduction |
The primary limitation is that variance reduction depends on the covariance created by the coupling. Some scenario outputs respond in opposite directions to shared inputs, producing negative covariance and noisier paired differences. A second limitation is interpretive: the Monte Carlo standard error quantifies finite simulation error for a specified target, not uncertainty about whether the model structure, parameter estimates, or chosen target reflect reality. Those sources of uncertainty require their own representation in the analysis.
Common mistake: “I used the same random seed for both scenarios, so the results are paired and must be more precise.” A shared seed does not by itself establish that corresponding outputs used corresponding draws; scenario branches may consume random numbers differently, and even valid pairs can have negative covariance. The correction is to align and verify the underlying replicate-level draws, calculate \(D_i\) within those pairs, and assess the observed covariance or difference variance. Pairing can reduce Monte Carlo noise, but it does not guarantee a reduction.
For a clean comparison, keep the target contrast, number of replicate pairs, and definitions of shared inputs explicit. If a scenario changes a mechanism that should not be shared, preserve its correct marginal distribution rather than forcing identical randomness. If paired and independent estimates differ substantially, investigate alignment and implementation before attributing the difference to simulation noise alone.
Worked Example: Predicting the Effect of Shared Draws
Consider a simple hypothetical model with a random input \(\theta\) having mean 10 and variance 4, and residual variability \(\varepsilon\) having mean 0 and variance 9. Suppose the scenario outputs are \(Y_A=\theta+\varepsilon\) and \(Y_B=1.1\theta+\varepsilon\). A CRN design uses the same \(\theta\) and \(\varepsilon\) in both outputs. These equations are a simplified illustration of coupling, not a pharmacokinetic model.
The expected outputs are \(E(Y_A)=10\) and \(E(Y_B)=11\), so the target contrast \(E(Y_B-Y_A)\) is 1 output unit. The marginal variances are \(\operatorname{Var}(Y_A)=4+9=13\) and \(\operatorname{Var}(Y_B)=1.1^2(4)+9=13.84\). Under shared draws, the covariance is \(1.1(4)+9=13.4\), because both the parameter and residual terms are shared. Thus,
Equivalently, \(Y_B-Y_A=0.1\theta\), whose variance is \(0.1^2(4)=0.040\). If the scenario outputs instead came from independent draws, their covariance would be zero and the difference variance would be \(13.84+13=26.84\). With \(N=1{,}000\), the corresponding standard errors are \(\sqrt{0.040/1{,}000}=0.00632\) output units for paired CRN and \(\sqrt{26.84/1{,}000}=0.164\) output units for independent simulation. The large reduction occurs because the shared terms nearly cancel in the contrast. In a real application, the degree of cancellation depends on the model and must be checked rather than assumed.
Check Your Understanding
Answer each question, then compare with the solutions.
- For \(N\) paired replicates, what is the replicate-level quantity used to estimate the difference in scenario means, and how is its Monte Carlo standard error calculated?
- Suppose two scenarios have standard deviations 5 and 4, and their within-pair correlation is 0.8. What is the variance of the paired difference?
- Why can the same random seed fail to create correctly aligned pairs?
- Two paired binary outcomes have equal estimated attainment probabilities. Does that imply all outcomes agree within pairs? Explain.
- What does a reduced Monte Carlo standard error establish, and what does it not establish?
Solutions. (1) Calculate \(D_i=Y_{Bi}-Y_{Ai}\) in each aligned pair; then calculate \(SE(\bar{D})=s_D/\sqrt{N}\). (2) The covariance is \(0.8(5)(4)=16\), so \(\operatorname{Var}(Y_B-Y_A)=4^2+5^2-2(16)=9\) output-units squared. (3) Different scenario code paths can consume random numbers differently, so the same seed may not assign corresponding draws to the same replicate. Alignment should be checked through replicate IDs or explicitly shared inputs. (4) No. Equal marginal proportions can coexist with discordant pairs; opposite-direction discordances can cancel in the estimated difference. (5) It establishes greater simulation precision for that estimated target under the simulated design. It does not establish that the model or target is scientifically correct, nor does it remove parameter or structural uncertainty.