Why Check a Canonical Model by Simulation?
A canonical joint distribution gives precise predictions about the behavior of cumulative score statistics across analysis looks. Simulation provides a practical way to check whether code intended to represent that distribution reproduces its key features. The check is especially useful when several information levels are involved: an implementation can produce plausible values at each look while getting the dependence between looks wrong.
The aim here is to check a specified canonical model, not to establish that a real study satisfies its assumptions. A simulation can reveal coding errors or discrepancies from the model’s mathematical predictions. It cannot, by itself, show that the model is appropriate for data from a particular trial.
The information levels need not be integer-valued. Their spacing determines how much new information each increment contributes. Write \(t_0=0\), so the increment between looks \(k-1\) and \(k\) has information length \(t_k-t_{k-1}\).
Generating a Canonical Score Path
A convenient simulation method starts with independent standard normal random variables \(\varepsilon_1,\ldots,\varepsilon_K\), each with mean zero and variance one. For each look, the cumulative score numerator includes the drift accumulated to that information level and the noise from all increments up to that look. Dividing by the square root of the information level gives the standardized statistic.
This construction is repeated independently for each simulated path. The same increment \(\varepsilon_j\) contributes to every later cumulative statistic, because later looks retain information accumulated at earlier looks. That reuse is essential: it creates the dependence between the cumulative statistics.
Why These Are the Canonical Predictions
The mean follows directly from the construction. Each \(\varepsilon_j\) has mean zero, so the only expected contribution to the numerator at look \(k\) is \(\theta t_k\). Dividing by \(\sqrt{t_k}\) gives \(E(Z_k)=\theta\sqrt{t_k}\).
For the variance and correlation, define the accumulated noise \(W_k=\sum_{j=1}^{k}\sqrt{t_j-t_{j-1}}\,\varepsilon_j\). Independent increments have zero covariance with one another. Therefore the variance of \(W_k\) is the sum of the increment variances:
Since \(Z_k\) is the sum of a fixed drift and \(W_k/\sqrt{t_k}\), its variance is \(t_k/t_k=1\). For \(i<j\), the noise \(W_i\) is part of \(W_j\); the additional increments from \(i+1\) through \(j\) are independent of \(W_i\). Thus \(\operatorname{Cov}(W_i,W_j)=\operatorname{Var}(W_i)=t_i\). Standardizing both statistics gives:
Because both variances are one, this covariance is also their correlation. The drift terms are fixed constants and do not affect covariance. This derivation also explains why the correlation is determined by information levels rather than by the value of \(\theta\).
A Practical Simulation Check
For a simulation with \(M\) independent paths, calculate the empirical mean and sample variance of the \(M\) values at each look. Calculate the sample correlation for each pair of looks. Compare those summaries with their theoretical targets, rather than expecting a finite simulation to reproduce them exactly.
Choose increasing information levels \(t_1,\ldots,t_K\) and a value of \(\theta\). State whether the check is under the null or an alternative value.
For each path, generate \(K\) independent standard normal values and use the canonical cumulative formula to calculate all \(K\) statistics.
At each look, compute the empirical mean and sample variance. For each pair of looks, compute the empirical correlation.
Check the empirical summaries against \(\theta\sqrt{t_k}\), one, and \(\sqrt{t_i/t_j}\), respectively. Interpret small discrepancies in light of Monte Carlo variability.
Worked Examples
Worked Example: Checking Three Looks Under the Null
Consider a hypothetical simulation with \(M=100{,}000\) paths, information levels \(t_1=0.25\), \(t_2=0.60\), and \(t_3=1.00\), and \(\theta=0\). The theoretical mean is zero and the theoretical variance is one at every look. The pairwise correlation targets are:
Suppose this clearly hypothetical run gives empirical means \((0.002,-0.004,0.001)\), sample variances \((0.998,1.006,0.995)\), and sample correlations \((0.644,0.500,0.776)\), listed in the order \((12,13,23)\). Each reported summary is close to its theoretical target.
For scale, under this null model the standard deviation of a simulated sample mean at any look is \(1/\sqrt{M}=0.003162\). The mean of \(-0.004\), for example, is about \(1.26\) such standard deviations below zero: \(-0.004/0.003162=-1.27\). That is a plausible Monte Carlo fluctuation, not evidence of a systematic mean shift. The sample variance estimates also fluctuate; for normal observations, their approximate standard error is \(\sqrt{2/(M-1)}\), or about \(0.004472\) here. These checks support the implementation for this run, but do not prove that it is error-free.
Worked Example: Checking Nonzero Means Under an Alternative
Use the same three information levels, but set \(\theta=0.80\). The predicted means are:
The variances remain one, and the correlations remain \(0.6455\), \(0.5000\), and \(0.7746\), because the drift changes the means but not the covariance structure. Suppose a hypothetical run with \(M=100{,}000\) paths reports means \((0.402,0.616,0.803)\), variances \((0.999,1.004,0.997)\), and correlations \((0.647,0.499,0.776)\). The results are close to all of the respective targets.
The expected standard deviation of a sample mean at look 2 is \(1/\sqrt{100{,}000}=0.003162\), since the variance of \(Z_2\) is one. The empirical difference from its target is \(0.616-0.6197=-0.0037\), roughly \(1.17\) standard errors. A nonzero mean does not justify checking against zero: doing so would incorrectly flag the intended alternative distribution as a failure.
Worked Example: Calculating One Path From Independent Increments
Suppose a single hypothetical path uses information levels \(t_1=0.40\), \(t_2=0.90\), and \(t_3=1.00\), with \(\theta=0.50\). Let its simulated independent increments be \(\varepsilon_1=0.50\), \(\varepsilon_2=-1.20\), and \(\varepsilon_3=0.30\). The information increments are \(0.40\), \(0.50\), and \(0.10\).
At the first look, the cumulative noise is \(\sqrt{0.40}(0.50)=0.3162\), and the drift is \(0.50(0.40)=0.20\). Therefore:
At the second look, the same first increment remains in the cumulative score, and the second increment contributes \(\sqrt{0.50}(-1.20)=-0.8485\). The drift is \(0.50(0.90)=0.45\), so:
At the third look, the additional noise is \(\sqrt{0.10}(0.30)=0.0949\), and the drift is \(0.50(1.00)=0.50\):
The statistics can rise or fall from one look to the next, even though information increases. The construction does not require cumulative \(Z\) statistics themselves to increase; it requires that their underlying score numerators retain earlier increments and add new ones.
Interpreting the Check and Avoiding a Common Error
Empirical means, variances, and correlations are useful summaries, but they are not exact acceptance criteria. Their values vary from run to run. Larger \(M\) generally reduces that Monte Carlo variation, while a mismatch that remains substantial across repeated runs merits investigation. Inspecting the information levels, increment variances, and reuse of earlier increments can help locate an implementation error.
A common mistake is to generate a separate independent standard normal value for each cumulative statistic. That procedure gives variance one at each look under the null, so a check of marginal variances alone may appear successful. But it gives zero correlation between distinct looks, rather than the canonical positive correlation \(\sqrt{t_i/t_j}\). Matching each marginal distribution is not sufficient to match the joint distribution.
There are two sound implementation approaches. Independent normal increments make the source of dependence explicit and are convenient for generating paths. Alternatively, one can generate a multivariate normal vector with the canonical mean vector and covariance matrix. That approach is mathematically equivalent when the covariance matrix is correctly specified, but may make it less obvious how earlier increments contribute to later looks. In either case, checking only means and variances misses dependence errors.
The check has limits. Agreement of these summaries is not a proof that every feature of the simulated joint distribution is correct, and a canonical-model simulation does not validate the assumptions that justify using that model for a real trial. The present check is best treated as an implementation verification step: it asks whether the code reproduces the stated model’s basic moment and dependence structure.
Check Your Understanding
Answer each question, then compare your reasoning with the solutions below.
- Under the null, what are the theoretical mean and variance of \(Z_k\) at an information level \(t_k\)?
- If \(t_1=0.36\) and \(t_2=0.81\), what is the canonical correlation between \(Z_1\) and \(Z_2\)?
- Under \(\theta=0.50\) and \(t=0.64\), what is the theoretical mean of the score statistic at that look?
- Why can independently simulated standard normal values at each look pass a variance check but fail a joint-distribution check?
Solutions. 1. Under the null, \(\theta=0\), so \(E(Z_k)=0\) and \(\operatorname{Var}(Z_k)=1\). 2. The correlation is \(\sqrt{0.36/0.81}=\sqrt{4/9}=2/3\), approximately \(0.6667\). 3. The mean is \(\theta\sqrt{t}=0.50\sqrt{0.64}=0.50(0.80)=0.4000\). 4. Each independently generated value has variance one, but distinct values have correlation zero. The canonical correlation is positive for increasing information levels, so the simulated dependence is wrong even though the marginal variances are correct.