Numerical Integration Of Conditional Error

Why Numerical Integration Is Needed

A conditional error function assigns an error probability to the interim data already observed. It answers: under the null model, given this interim result, what is the probability that the planned design will eventually reject? The answer is needed at the observed value, not merely on average across all possible interim results. It can also be used when the future statistic has a distribution for which a convenient closed-form probability is unavailable.

This tutorial develops the computation for a specified joint distribution and a fixed future rejection rule. It focuses on evaluating the probability integral reliably. It does not determine a new continuation boundary. In particular, extreme interim results deserve special attention: the conditional error can be so small that a numerical routine reports zero, or so close to one that rounding hides a meaningful difference.

Model and Conditional Error Integral

Definition: Let \(Z_I\) be the interim statistic and \(Z_F\) the final statistic. Under the null hypothesis, suppose \((Z_I,Z_F)\) has a specified joint distribution. If the planned final rejection region is \(R\), the conditional error at observed interim value \(z\) is \(CE(z)=P(Z_F\in R\mid Z_I=z,H_0)\). It is a conditional probability under the null model and the planned design.

For the examples, assume the null joint distribution is bivariate standard normal: \(Z_I\sim N(0,1)\), \(Z_F\sim N(0,1)\), and their correlation is \(\rho\), with \(-1<\rho<1\). The conditional distribution is \(Z_F\mid Z_I=z,H_0\sim N(\rho z,1-\rho^2)\). Thus, for a one-sided final rule that rejects when \(Z_F\ge c\),

CE(z)=\int_c^\infty \frac{1}{\sqrt{1-\rho^2}}\, \phi\!\left(\frac{x-\rho z}{\sqrt{1-\rho^2}}\right)\,dx =\int_{(c-\rho z)/\sqrt{1-\rho^2}}^\infty \phi(u)\,du,

where \(\phi(u)=(2\pi)^{-1/2}\exp(-u^2/2)\) is the standard normal density. The change of variable is \(u=(x-\rho z)/\sqrt{1-\rho^2}\). For a two-sided rule, the rejection region has two tails, and the conditional error is the sum of the integrals over those tails. These formulas require the stated joint model and a rejection region fixed by the planned design. If either is misspecified, accurate integration still gives the wrong conditional error.

The conditional normal result follows by writing \(Z_F=\rho Z_I+\sqrt{1-\rho^2}\,E\), where \(E\sim N(0,1)\) is independent of \(Z_I\). Given \(Z_I=z\), the first term is fixed and the second is normal with mean zero and variance \(1-\rho^2\). This establishes the integrand and its scale; it also shows why simply treating the final statistic as independent of the interim statistic is generally incorrect.

Composite Simpson Quadrature

Composite Simpson quadrature approximates an integral on \([a,b]\) by dividing the interval into an even number \(n\) of equal subintervals, with \(h=(b-a)/n\). If \(f_j=f(a+jh)\), then

\int_a^b f(u)\,du\ \approx\ \frac{h}{3}\left[ f_0+f_n+ 4\sum_{\substack{j=1\\j\ {\rm odd}}}^{n-1}f_j+ 2\sum_{\substack{j=2\\j\ {\rm even}}}^{n-2}f_j \right].

The rule is accurate for smooth integrands when the grid is sufficiently fine, but it does not automatically account for an infinite integration limit. For a normal tail, one practical approach is to integrate from the lower limit to a finite upper limit \(U\), then bound or separately calculate the remaining tail. A second approach is a change of variable mapping the infinite interval to a finite one; the transformed integrand must include the Jacobian and should be checked near both endpoints.

Worked Example: A One-Sided Conditional Error

Suppose the null interim and final statistics are standard normal with correlation \(\rho=0.50\). The planned final rule rejects for \(Z_F\ge c\), where \(c=1.645\). At the observed interim value \(z=0\), compute the conditional error by quadrature.

The model gives conditional mean \(\rho z=0\) and conditional standard deviation \(\sqrt{1-\rho^2}=\sqrt{0.75}=0.8660\). The standardized lower integration limit is \(a=1.645/\sqrt{0.75}\approx1.89948\), approximately 1.90. For the Simpson calculation below, round the lower limit to 1.90; the exact conditional error is \(\int_{1.89948}^\infty\phi(u)\,du\).

Use Simpson quadrature from 1.90 to 4.00 with \(h=0.35\), giving six subintervals. Rounded density values at \(1.90,2.25,2.60,2.95,3.30,3.65,4.00\) are \(0.06562,0.03174,0.01358,0.005143,0.001723,0.0005105,0.000134\). The Simpson calculation is

\frac{0.35}{3}\left[ 0.06562+0.000134+ 4(0.03174+0.005143+0.0005105)+ 2(0.01358+0.001723) \right] \approx 0.02869.

The remaining standard normal tail above 4 is approximately \(0.0000317\), so adding it to the truncated-grid Simpson estimate gives about \(0.02872\). Evaluating the normal tail directly at the exact lower limit \(1.645/\sqrt{0.75}\) gives approximately \(0.0287505\), or \(0.02875\). Thus the conditional error at \(Z_I=0\) is approximately 0.02875, or 2.875%. It is a null rejection probability conditional on this interim statistic, not the probability that the null hypothesis is true.

As a check, evaluating the normal tail directly gives a value close to the quadrature result. Agreement is reassuring, but it does not replace checking the model, the integration limits, and the quadrature grid.

Extreme Interim Values and Numerical Stability

The standardized lower limit \(a=(c-\rho z)/\sqrt{1-\rho^2}\) moves with the interim observation. Large positive \(z\) can move it far below zero, making the conditional error close to one; large negative \(z\) can move it far into the upper tail, making the conditional error tiny. The same absolute quadrature tolerance is not equally informative in both cases. For example, an absolute error of \(10^{-6}\) is negligible for a probability near 0.5 but is substantial relative to a probability near \(10^{-8}\).

Worked Example: Checking Both Tails of the Interim Range

Keep \(\rho=0.50\) and \(c=1.645\), but evaluate the planned one-sided conditional error at \(z=-6\) and \(z=6\). These are deliberately extreme, hypothetical interim values used to test numerical behavior.

At \(z=-6\), the conditional mean is \(-3\), and the standard deviation remains 0.8660. Hence \(a=(1.645-(-3))/0.8660=5.364\). Integrating the standard normal density from 5.363 to infinity gives approximately \(4.1\times10^{-8}\). A routine that returns exactly zero has not established that the conditional error is zero; it may have underflowed or stopped with a tolerance too large for this tail.

At \(z=6\), the conditional mean is 3, giving \(a=(1.645-3)/0.8660=-1.565\). The integral from \(-1.564\) to infinity is approximately 0.941. A stable calculation can evaluate the complementary lower tail, \(1-P(U<-1.564)\), rather than numerically accumulating a nearly full area and losing precision in its small complement.

The two answers make sense under the model: an interim statistic far below the null mean makes a future positive crossing unlikely, while a large positive interim statistic makes it likely. Check both extremes in an implementation, not just values near zero.

For very small probabilities, useful safeguards include relative as well as absolute error tolerances, evaluating log densities, and comparing the computed tail with a trusted normal-tail calculation when the model is normal. For probabilities near one, compute the smaller complementary probability directly and subtract it from one. Refining the grid should lead to stable digits; if the answer changes materially as the grid is refined, the earlier grid was inadequate.

Normalization Is Not Optional

Some numerical procedures integrate an unnormalized kernel and then divide by its total integral. The raw integral is not itself a probability. A simple grid calculation makes this distinction explicit and provides an arithmetic check for Simpson weights.

Worked Example: A Normalized Conditional Distribution

Let \(v\), conditional on an interim result, have a uniform distribution on \([0,12]\). Use the unnormalized constant kernel \(k(v)=0.09725\). With Simpson quadrature using \(h=0.50\), there are 24 subintervals: 12 odd-indexed interior points and 11 even-indexed interior points. The normalizing integral is

D\approx \frac{0.50}{3} \left[2(0.09725)+4\{12(0.09725)\}+2\{11(0.09725)\}\right] =1.167.

The value 1.167 is an unnormalized integral, not an error probability. The normalized density is \(k(v)/D=0.09725/1.167\), approximately \(1/12\). If continuation results in rejection exactly when \(v\ge6\), the numerator integral is over half the interval:

N\approx 6(0.09725)=0.5835,\qquad CE=\frac{N}{D}=\frac{0.5835}{1.167}=0.500.

The uniform distribution assigns conditional probability 0.500 to \(v\ge6\), as expected. This example shows why an integral greater than one can be valid as an intermediate calculation only when the kernel is explicitly unnormalized and the final probability is normalized correctly.

Diagnostics, Alternatives, and Limitations

A useful numerical audit separates model checks from integration checks. Verify that the conditional distribution uses the specified correlation, variance, and observed interim value; verify that the integration region matches the rejection rule; then assess numerical convergence. In particular:

Refine the grid or tighten the adaptive tolerance and compare results at the displayed precision.
Check that a computed probability lies between 0 and 1, and that disjoint rejection tails are added without overlap.
For a truncated infinite integral, bound the omitted tail and ensure it is small relative to the required accuracy.
Evaluate representative interim values near the center and in both extremes; check for underflow, cancellation, and accidental reporting of a rounded zero or one.
Confirm that any density kernel is normalized, or explicitly use a numerator-to-denominator ratio.

For the bivariate normal model, direct normal-tail evaluation is usually simpler and more accurate than numerical quadrature. Quadrature is valuable when the conditional density or rejection region is more complicated, or when an integral is part of a larger computation. Monte Carlo simulation is another alternative: it estimates the conditional probability by the fraction of simulated future statistics in the rejection region. Its sampling error can be large for very small conditional errors, so a simulation that sees no rejections does not demonstrate a zero probability. Adaptive quadrature allocates more evaluations where the integrand changes rapidly; composite Simpson quadrature is easier to inspect and reproduce for a smooth, one-dimensional integral.

All these approaches share a limitation: they compute the probability implied by the supplied model. They cannot repair an incorrect joint distribution, an incorrectly coded rejection region, or an omitted design feature such as a planned stopping rule. Numerical accuracy and statistical validity are distinct requirements.

Common mistake: Treating the integral of an unnormalized kernel as a conditional error, or trusting a quadrature output of zero in an extreme tail. The correction is to normalize the kernel when necessary, check the integration region and omitted tail, and evaluate small probabilities with tolerances and methods appropriate to their scale.

Check Your Understanding

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

Under the bivariate normal model, what is the conditional variance of \(Z_F\) given \(Z_I=z\)?
For \(\rho=0.40\), \(z=1\), and one-sided boundary \(c=1.645\), write the standardized lower limit for the normal-tail integral.
Why does a Simpson approximation to an integral over an infinite interval require an additional truncation or transformation step?
A numerical routine reports conditional error zero at an extreme interim value. Give two checks before interpreting this as zero.
An unnormalized kernel has integral 1.167 and its rejection-region integral is 0.5835. What is the normalized conditional error?

Solutions. (1) \(1-\rho^2\), which does not depend on \(z\). (2) The conditional standard deviation is \(\sqrt{1-0.40^2}=\sqrt{0.84}\), so the lower limit is \((1.645-0.40)/\sqrt{0.84}\approx1.358\). (3) Simpson's rule operates on a finite interval; the omitted tail must be bounded or the infinite range mapped to a finite one with its Jacobian included. (4) Refine the calculation and use an appropriately small absolute and relative tolerance; also evaluate the tail using a stable method, such as log-density or a trusted tail calculation. (5) Divide numerator by denominator: \(0.5835/1.167=0.500\).