Benchmarking Against Direct Multivariate Normal Probabilities

Why Use a Direct Benchmark?

A recursive calculation updates the probability density of the trial statistic among paths that have not stopped. A direct multivariate normal calculation instead asks for the probability that the complete vector of statistics lies inside the continuation region. For a small number of looks, the second calculation is an independent way to check the first: the methods express the same probability through different numerical routes.

The comparison is meaningful only when both routes use exactly the same null model, information times, stopping boundaries, and event definition. A numerical discrepancy can signal an integration error, but it can also arise because one calculation uses a different boundary convention or covariance matrix. The benchmark is therefore both a numerical check and a check that the probability being computed has been specified correctly.

Definition: Let $Z_1,\ldots,Z_K$ be the standardized test statistics at $K$ analysis looks under $H_0$, with information values $I_1<\cdots<I_K$. Under the canonical group-sequential model, $Z=(Z_1,\ldots,Z_K)$ has a multivariate normal distribution with mean zero, unit variances, and correlations $$ \operatorname{Corr}(Z_i,Z_j)=\sqrt{\frac{I_i}{I_j}},\qquad i<j. $$ For continuation intervals $(\ell_j,u_j)$, the probability of surviving through look $k$ is $S_k=P(\ell_j<Z_j<u_j\text{ for every }j\leq k)$. For a one-sided upper-boundary design, use $(\ell_j,u_j)=(-\infty,c_j)$. Set $S_0=1$. If crossing is the only stopping event, the probability of stopping at look $k$ is $S_{k-1}-S_k$.

This model assumes the canonical normal joint distribution is appropriate under the null and that the information values and boundaries used in the calculation are the planned ones. The correlation formula follows from independent information increments in the underlying score process: the statistic at an earlier information time shares its accumulated score with the later statistic. If the design has a different joint model, a multivariate normal benchmark must use that model rather than this canonical correlation structure.

What the Two Calculations Evaluate

The direct method evaluates a multivariate normal probability over a rectangle. For two looks, for example, the survival probability is the integral of the bivariate normal density over the product of the two continuation intervals. A recursive method obtains the same probability by integrating over the first look and then applying the conditional distribution of the next statistic. These are different representations of one probability, not different definitions of survival.

S_2=\int_{\ell_1}^{u_1}\phi(z_1) \left[ \Phi\!\left(\frac{u_2-\rho z_1}{\sqrt{1-\rho^2}}\right) - \Phi\!\left(\frac{\ell_2-\rho z_1}{\sqrt{1-\rho^2}}\right) \right]\,dz_1, \qquad \rho=\sqrt{\frac{I_1}{I_2}}.

Here $\phi$ and $\Phi$ are the standard normal density and distribution function. The formula uses $Z_2\mid Z_1=z_1\sim N(\rho z_1,1-\rho^2)$. The direct bivariate rectangle calculation integrates the joint density over the same two intervals.

The equivalence extends to any finite number of looks. The joint normal density factors into the first-look density and successive conditional densities. Integrating that factorization over the continuation region, one coordinate at a time, gives the recursive integration. Integrating the same joint density over the whole rectangle gives the direct calculation. Because both integrate the same nonnegative density over the same region, changing the order of integration does not change the exact probability. Any difference in computed answers is numerical or reflects mismatched inputs.

Benchmark principle: Compare probabilities for identical events under identical model inputs. Agreement to a prespecified numerical tolerance supports the recursive implementation for that case; it does not establish that every model assumption or every larger-look calculation is correct.

Two-Look Benchmarks

Worked Example: Two-Sided Boundaries at Two Looks

Suppose a hypothetical trial has two analyses at information values $I_1=1$ and $I_2=2$. Under $H_0$, the canonical statistics have correlation $\rho=\sqrt{1/2}\approx0.7071$. At both looks, stop for efficacy if $|Z_j|\geq2$. The continuation event through look 2 is therefore $-2<Z_1<2$ and $-2<Z_2<2$. The model assumptions are the stated canonical null model and the fixed symmetric boundaries.

At the first look, the survival probability is $$ S_1=P(-2<Z_1<2)=\Phi(2)-\Phi(-2)\approx0.9545. $$ For a direct bivariate calculation, integrate the bivariate normal density over the square $(-2,2)\times(-2,2)$. Equivalently, subtract the probability of leaving that square. Each marginal has two tails of probability $\Phi(-2)\approx0.02275$. The same-sign joint tail probability $P(Z_1\geq2,Z_2\geq2)$ is about $0.00751$; by symmetry the lower-lower tail has the same probability. The opposite-sign joint tails are very small. Thus, to the shown precision, $$ S_2\approx1-4\Phi(-2)+2(0.00751)\approx0.9240. $$

The recursive calculation integrates the conditional continuation probability in the displayed two-look formula, with both intervals equal to $(-2,2)$ and $\rho=\sqrt{1/2}$. It gives the same survival probability, approximately $0.9240$, when evaluated to comparable precision. Consequently, stopping at look 2 after surviving look 1 has probability $$ S_1-S_2\approx0.9545-0.9240=0.0305. $$ This is the total probability of crossing either boundary at look 2 after continuation at look 1. It is not the probability of crossing one particular boundary.

Worked Example: One-Sided Upper Boundaries

Consider another hypothetical two-look design with $I_1=1$, $I_2=2$, and the same null correlation $\rho=\sqrt{1/2}$. The upper efficacy boundaries are $c_1=1.8$ and $c_2=2.0$. The continuation event is $Z_1<1.8$ and $Z_2<2.0$. Under the canonical bivariate normal model, the first-look survival probability is $$ S_1=\Phi(1.8)\approx0.9641. $$

To compute the probability of surviving look 1 and crossing the upper boundary at look 2, start with $P(Z_2\geq2)=\Phi(-2)\approx0.02275$, then remove paths that had already crossed at look 1. The joint upper-tail probability $P(Z_1\geq1.8,Z_2\geq2)$, evaluated by bivariate normal integration, is approximately $0.00982$. Therefore, $$ P(Z_1<1.8,\ Z_2\geq2)\approx0.02275-0.00982=0.01293. $$

This is the stage-2 stopping probability. The survival probability through both looks is $S_2=S_1-0.01293\approx0.9511$, and the cumulative probability of crossing by look 2 is $1-S_2\approx0.0489$. A recursive calculation gets the stage-2 probability by integrating the conditional upper-tail probability over $Z_1<1.8$; a direct calculation gets it from the same bivariate distribution as a rectangle or tail difference. Their agreement checks the event definition as well as the arithmetic.

A Three-Look Check With a Closed-Form Benchmark

Worked Example: Zero Upper Boundaries at Three Looks

A deliberately simple hypothetical design has information values $I_1=1$, $I_2=2$, and $I_3=4$, and upper boundaries $c_1=c_2=c_3=0$. The null model is the canonical joint normal model. The continuation event through a look is that every statistic observed so far is below zero. Although a zero boundary is not offered here as a practical design recommendation, it creates an exact check on a three-dimensional probability calculation.

The correlations are $\rho_{12}=\sqrt{1/2}$, $\rho_{13}=\sqrt{1/4}=0.5$, and $\rho_{23}=\sqrt{2/4}=\sqrt{1/2}$. For a centered trivariate normal vector, the probability that all three components are negative is $$ S_3=\frac18+\frac{1}{4\pi} \left(\arcsin(\rho_{12})+\arcsin(\rho_{13})+\arcsin(\rho_{23})\right). $$ Here the sum of inverse sines is $\pi/4+\pi/6+\pi/4=2\pi/3$, so the direct benchmark is $S_3=1/8+1/6=7/24\approx0.2917$. For the first two looks, the bivariate zero-threshold formula gives $S_1=1/2$ and $S_2=1/4+\arcsin(\sqrt{1/2})/(2\pi)=3/8=0.3750$.

The recursive stagewise probabilities must therefore be stopping at look 1: $1-S_1=0.5000$; stopping at look 2: $S_1-S_2=0.1250$; stopping at look 3: $S_2-S_3=0.08333$; and continuing beyond look 3: $S_3=0.2917$. The stopping probabilities sum to $0.7083=1-S_3$. A recursive implementation for these same information times and boundaries should reproduce the three survival probabilities. This benchmark is especially useful because it checks the three-look dependence structure, not just a two-look marginal calculation.

How to Compare Results and Diagnose Disagreement

A useful comparison records each survival probability and, where relevant, each stagewise stopping probability. Checking only the final cumulative probability can conceal offsetting errors between stages. Keep sufficient internal precision during calculation; round only when reporting. Before attributing a discrepancy to numerical integration, confirm the following:

Same probability: Both methods must evaluate the same continuation rectangle or crossing event. For two-sided rules, include both tails at each look.
Same covariance: The correlation matrix must correspond to the same information times and null model. In the canonical model, $\rho_{ij}=\sqrt{I_i/I_j}$ for $i<j$.
Same boundary convention: Clarify whether equality at a boundary counts as crossing. Under a continuous normal model, equality has probability zero, but the convention still matters in implementation.
Probability accounting: Verify $S_{k-1}-S_k\geq0$, and verify that stagewise stopping probabilities plus final survival sum to one, up to numerical error.
Numerical stability: Repeat the calculation with tighter integration tolerances or a refined recursive grid. A discrepancy that shrinks with refinement is evidence of numerical error; a stable discrepancy suggests inconsistent inputs or an implementation defect.

The earlier material on truncating the Gaussian state space, grid refinement, and numerical error budgets is relevant when a recursive calculation uses a finite grid or a truncated range. The direct benchmark also has numerical error: it is not an error-free oracle simply because it evaluates a multivariate normal probability in one call. Record the method, tolerance, and reported precision for both routes, and compare beyond the rounding level rather than expecting identical last digits.

When Direct Integration Stops Being Practical

Direct multivariate normal integration is attractive for one, two, or a few looks because the event is easy to state as a low-dimensional rectangle and an independent implementation may expose errors in the recursion. Its computational burden generally increases with dimension. Generic deterministic quadrature needs many integration evaluations as dimensions are added, while specialized multivariate normal probability algorithms can still become costly or less stable for some covariance structures and tight tolerances. Closely spaced information times can also produce correlations near one, making the covariance matrix poorly conditioned for numerical work.

These limitations do not make direct probabilities useless. They make them most valuable as small-dimensional benchmarks, selected to cover representative boundary types and information patterns. Recursive integration is usually preferable when many looks are required because it exploits the sequential structure and reuses survivor-density calculations. Direct calculations are a useful independent comparison when the dimension is small; they are not a universal replacement for recursion or a proof of correctness for untested designs.

Common pitfall: It is incorrect to compare a recursive two-sided survival probability with a direct calculation that omits one tail, or to assume that the direct answer is correct because it is produced by a multivariate normal routine. Specify the rectangle, derive the correlations from the information times, and check the probability accounting before interpreting agreement or disagreement.

Check Your Understanding

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

For a two-look canonical design with $I_1=3$ and $I_2=12$, what is the null correlation between $Z_1$ and $Z_2$?
In the two-sided worked example, what is the approximate probability of stopping at look 2 after surviving look 1?
For the one-sided worked example, explain why the stage-2 crossing probability is not simply $\Phi(-2)$.
In the three-look zero-boundary example, what is the stopping probability at look 3, and what is the final continuation probability?
Name two reasons that a direct multivariate normal calculation may disagree numerically with a recursive calculation even when both are intended to represent the same design.

Solutions. 1. The correlation is $\sqrt{3/12}=\sqrt{1/4}=0.5000$. 2. It is $S_1-S_2\approx0.9545-0.9240=0.0305$, including crossing either boundary at look 2. 3. $\Phi(-2)$ includes paths that crossed the look-1 boundary already; subtracting $P(Z_1\geq1.8,Z_2\geq2)\approx0.00982$ gives the stage-2 probability $0.01293$. 4. The look-3 stopping probability is $S_2-S_3=0.3750-0.2917\approx0.08333$, and final continuation is $S_3=0.2917$. 5. Examples include using different boundaries or correlation matrices, different continuation-event definitions, inadequate grid refinement, truncation error, or insufficient integration tolerance. These should be checked before concluding that either method has been validated.