14. Severe Testing, Multiplicity, Martingales (e-Processes), and Conformal Prediction

The Scientist Who Peeked Every Morning Until the Coin Confessed

Suppose a quantitative researcher wants to prove that a new trading strategy (or clinical drug) beats a fair 50/50 coin toss ($H_0: p = 0.50$) at the classical significance threshold $\alpha = 0.05$.

• She starts an experiment, planning to collect up to $n = 1\text{,}000$ trades.
• However, instead of waiting until trade $1\text{,}000$, she checks the classical $p$-value on her screen after every single trade ($t = 10, 11, 12, \dots$), and the very first morning the $p$-value dips to $p = 0.048 < 0.05$, she yells “Statistically significant!”, stops the experiment, and publishes her paper.

If the strategy is actually a 100% useless fair coin ($H_0$ is true), what is the probability that a researcher who peeks after every trade up to $n = 1\text{,}000$ will see $p < 0.05$ at least once? It is not $5%$. Because the standardized random walk $\sum X_i / \sqrt{t}$ fluctuates back and forth by the Law of the Iterated Logarithm ($\limsup_{t \to \infty} \frac{S_t}{\sqrt{2t \ln\ln t}} = 1$), peeking continuously over $1\text{,}000$ steps inflates the false-positive rate from $5%$ to nearly $40%$—and with unlimited peeking ($t \to \infty$), a classical $p$-value is guaranteed to cross $p < 0.05$ with $100%$ probability!

Now combine optional stopping with two other everyday practices:

1. Trying $K = 20$ different feature combinations in the notebook and reporting only the winner (Multiplicity / $p$-Hacking).
2. Treating $10\text{,}000$ hourly rows from the same $10$ hospitals as if they were $10\text{,}000$ independent hospitals (Pseudoreplication / Cluster Dependence).

Algorithmic refutation of optional stopping, multiplicity fishing, pseudoreplication, and production drift rests on five inferential pillars: Severe Testing, Anytime-Valid Betting Martingales (e-processes), Cluster Intra-Class Correlation (ICC), Simulation-Based Calibration (SBC), and Conformal Prediction under Distribution Shift.

14.1 Severe Testing, Multiplicity ($K$), and the Five Statistical Errors

Karl Popper and Deborah Mayo’s Severity Principle

Why is finding a high correlation $R^2$ or a small training error not enough to validate a scientific claim $C$? Because an overparameterized curve can fit almost anything! Karl Popper and philosopher of statistics Deborah Mayo formulated the Severity Requirement:

Severe Test Criterion: Data $D$ provide genuine evidence for a claim $C$ if and only if $C$ passes a test procedure $T$ that had a very high probability of falsifying (rejecting) $C$ if $C$ were false.

The Five Canonical Statistical Errors

In introductory textbooks, students learn only Type I Error ($\alpha$, false positive: rejecting a true null $H_0$) and Type II Error ($\beta$, false negative: failing to reject a false null $H_0$). In operational reality, three additional errors cause far more damage:

1. Type III Error (A. W. Kimball, 1957): Computing the mathematically exact answer to the wrong operational question—for example, predicting annual average hospital bed demand with $99%$ precision when the nursing supervisor has to staff 8-hour shifts (enforced by Boundary Invariant $\Phi_5$ in Chapter 16).
2. Type S Error (Sign Error, Gelman & Tuerlinckx, 2000): Claiming with “statistical significance” that an intervention’s effect is positive when the true causal effect is negative!
3. Type M Error (Magnitude / Exaggeration Ratio, Gelman & Carlin, 2014): When the effective sample size $n_{\mathrm{eff}}$ is small and noise $\sigma$ is large, the only estimates $\hat{\theta}$ that cross the significance filter $|\hat{\theta}| > 1.96,\mathrm{SE}$ are those lucky draws where random noise happened to be huge. Consequently, condicioning on “statistically significant” findings exaggerates the true effect size by $3\times$ to $10\times$ (the Winner’s Curse)!

Multiplicity Inflation Across $K$ Attempted Specifications

Suppose a team tests $K$ independent candidate hypotheses or feature specifications at nominal level $\alpha = 0.05$, when in reality all $K$ candidates are pure noise ($H_0$ is true for all $K$).

• The probability that a single test does not produce a false alarm is $1 - \alpha = 0.95$.
• The probability that none of the $K$ independent tests produces a false alarm is $(1 - \alpha)^K$.
• Therefore, the Family-Wise Error Rate (FWER)—the probability of getting at least one false discovery across the $K$ attempts—is:

$$\alpha_{\text{family}}(K) = 1 - (1 - \alpha)^K$$

For $\alpha = 0.05$ and $K = 20$ silent attempts, $\alpha_{\text{family}}(20) = 1 - 0.95^{20} = \mathbf{64.15%}$! To keep the family-wise error rate below $\alpha$, Šidák / Bonferroni correction requires dividing the per-test threshold by $K$ ($\alpha_k = \alpha / K$) or multiplying the required e-value threshold by $K$ (which is why Chapter 15’s PROHIBITION_4 requires declaring specifications_tried_k and discarded_findings in every YAML Card!).

14.2 Testing by Betting: Martingales, e-Processes ($E_t$), and Ville’s Inequality

How do we monitor a live model or sequential experiment continuously over time $t = 1, 2, 3, \dots$ without ever invalidating our false-positive guarantee when we peek or stop early?

Game-theoretic statistics (Jean Ville, 1939; Glenn Shafer; Aaditya Ramdas) replaces fixed-sample $p$-values with Testing by Betting:

Definition 14.1 (e-Value and Cumulative e-Process)

An e-variable (evidence variable) at step $i$ is a non-negative random variable $L_i \ge 0$ whose expected value under the null hypothesis $H_0$ is at most $1$:

$$\mathbb{E}_{H_0}!\left[L_i ;\middle|; \mathcal{F}_{i-1}\right] \le 1$$

Think of $L_i$ as the one-step wealth multiplier of a fair ticket in a casino where the house ($H_0$) gives you zero edge: if you enter round $i$ with $$1$, your expected wealth after round $i$ under $H_0$ is at most $$1$.

If we reinvest our wealth round after round, starting from initial capital $E_0 = 1$, our cumulative wealth after $t$ observations is the product (cumulative e-process):

$$E_t = \prod_{i=1}^t L_i = E_{t-1} \cdot L_t$$

For example, if each step $i$ yields a bounded score $S_i \in [-1, 1]$ with $\mathbb{E}_{H_0}[S_i \mid \mathcal{F}_{i-1}] \le 0$ (and positive expectation under the alternative $H_1$), betting a fraction $\lambda_i \in [0, 1)$ of our current wealth produces the betting factor $L_i = 1 + \lambda_i S_i \ge 0$, or in likelihood-ratio form $L_i = \frac{p_{H_1}(Y_i)}{p_{H_0}(Y_i)}$.

Ville’s Inequality (Why You Can Peek at $E_t$ Whenever You Want)

Because $\mathbb{E}_{H_0}[E_t \mid \mathcal{F}_{t-1}] = E_{t-1},\mathbb{E}_{H_0}[L_t \mid \mathcal{F}_{t-1}] \le E_{t-1}$, the sequence ${E_t}_{t \ge 0}$ is a non-negative supermartingale starting at $E_0 = 1$.

In 1939, Jean Ville proved a sweeping generalization of Markov’s Inequality that holds simultaneously across the entire infinite trajectory $t = 1, 2, \dots, \infty$:

$$P_{H_0}!\left(\sup_{t \ge 1} E_t \ge \frac{1}{\alpha}\right) \le \alpha$$

Intuition of Ville’s Inequality: In a fair casino ($H_0$), if you walk in with $E_0 = $1$, the probability that your bankroll ever touches $$20$ ($1/\alpha$ for $\alpha = 0.05$) at any moment in your life—even if you peek after every single spin of the wheel and quit the instant you hit $$20$—is at most $1/20 = 5%$! Conversely, if the null hypothesis $H_0$ is false and the alternative $H_1$ is true ($\mathbb{E}_{H_1}[\ln L_i] > 0$), then $E_t = \exp(\sum_{i=1}^t \ln L_i)$ grows exponentially to $+\infty$, crossing $1/\alpha = 20$ rapidly!

This is why the Stochology Protocol (Chapters 15, 17, and 18) uses cumulative e-processes $E_t = \prod L_i$ as its official sequential evidence engine.

14.3 Cluster Dependence (Intra-Class Correlation $\rho$) and Simulation-Based Calibration (SBC)

Why Clustered Rows Are Not Independent: Variance Components and ICC ($\rho$)

Suppose we collect $n = 1\text{,}000$ observations grouped into $G = 10$ clusters of size $m = 100$ each (e.g., $100$ flights at each of $10$ airports, or $100$ patients in each of $10$ hospitals).

A two-level random-intercept model writes the outcome of unit $j$ inside cluster $g$ as:

$$Y_{gj} = \mu + u_g + e_{gj}$$

where $u_g \sim (0, \sigma_b^2)$ is the between-cluster shock shared by every row in cluster $g$ (e.g., a snowstorm at Airport $g$), and $e_{gj} \sim (0, \sigma_w^2)$ is the within-cluster individual noise. Look at the covariance between two different rows $j \ne k$ inside the same cluster $g$:

$$\mathrm{Cov}(Y_{gj}, Y_{gk}) = \mathrm{Cov}(u_g + e_{gj}, ; u_g + e_{gk}) = \mathrm{Var}(u_g) = \sigma_b^2 > 0$$

Dividing this shared covariance by the total variance $\mathrm{Var}(Y_{gj}) = \sigma_b^2 + \sigma_w^2$ gives the Intra-Class Correlation (ICC) $\rho \in [0, 1]$:

$$\rho = \frac{\sigma_b^2}{\sigma_b^2 + \sigma_w^2}$$

Now compute the variance of the grand average $\bar{Y} = \frac{1}{G}\sum_{g=1}^G \bar{Y}_g$ across all $n = Gm$ observations:

• Each cluster average $\bar{Y}_g = \mu + u_g + \frac{1}{m}\sum_{j=1}^m e_{gj}$ has variance $\sigma_b^2 + \frac{\sigma_w^2}{m}$.
• Averaging across $G$ independent clusters divides by $G$:

$$\mathrm{Var}(\bar{Y}) = \frac{1}{G}\left(\sigma_b^2 + \frac{\sigma_w^2}{m}\right) = \frac{\sigma_b^2 + \sigma_w^2}{G m}\left[1 + (m - 1)\rho\right] = \frac{\sigma_{\text{total}}^2}{n} \cdot \underbrace{\left[1 + (m - 1)\rho\right]}_{\text{Design Effect } \mathrm{DEFF}}$$

Look closely at the factor $\mathrm{DEFF} = 1 + (m - 1)\rho$:

• Even a tiny intra-cluster correlation $\rho = 0.05$ with cluster size $m = 100$ produces $\mathrm{DEFF} = 1 + 99(0.05) = 5.95$, shrinking your effective sample size from $n = 1\text{,}000$ down to $n_{\mathrm{eff}} = \frac{n}{1 + (m-1)\rho} = \frac{1000}{5.95} \approx \mathbf{168}$!
• In Chapter 16 ($\Phi_2$), we will extend this exact derivation to unequal cluster sizes $\tilde{m} = \bar{m}(1 + \mathrm{CV}^2)$ (Moulton-Campbell Correction) and tail threshold exceedances $\rho_I(p)$.

Simulation-Based Calibration (SBC, Cook, Gelman & Rubin / Talts et al.)

Before trusting any Bayesian or Monte Carlo inference engine on real data, how do we verify that the inference code itself is not buggy or biased?

1. Draw a synthetic ground-truth parameter from the prior: $\theta^{\text{true}} \sim \pi(\theta)$.
2. Simulate a synthetic dataset from the forward model: $y^{\text{sim}} \sim p(y \mid \theta^{\text{true}})$.
3. Run your inference algorithm on $y^{\text{sim}}$ to obtain $M$ posterior draws ${\theta^{(1)}, \dots, \theta^{(M)}}$, and count the rank $r \in {0, 1, \dots, M}$ of the true parameter among the posterior draws: $$r = \sum_{m=1}^M \mathbf{1}{\theta^{(m)} < \theta^{\text{true}}}$$

Theorem (Talts et al., 2018): If the inference algorithm is mathematically exact, the rank $r$ across repeated synthetic trials is strictly uniformly distributed on ${0, 1, \dots, M}$!

• If the rank histogram is $\cup$-shaped (piling up at the extremes $0$ and $M$), the posterior intervals are too narrow (overconfident).
• If the rank histogram is $\cap$-shaped (piling up in the middle), the posterior is underconfident.
• If it is tilted left or right, the estimator is systematically biased! Filter 3 of the Eliminatory Funnel (Chapter 18) enforces this uniformity via a $\chi^2$ goodness-of-fit test on the SBC rank histogram.

14.4 Distribution-Free Conformal Prediction and Shift Algebra

Split Conformal Prediction (Guaranteed Finite-Sample Coverage)

How can we wrap any point predictor $\hat{f}(x)$ in a prediction interval $\hat{C}(X_{n+1})$ that is guaranteed to contain the true future outcome $Y_{n+1}$ at least $1 - \alpha$ (say $90%$) of the time—without assuming Gaussian residuals or a correctly specified model?

Vladimir Vovk and colleagues discovered that exchangeability (every permutation of the calibration points $(X_1, Y_1), \dots, (X_n, Y_n)$ and the test point $(X_{n+1}, Y_{n+1})$ is equally likely) is all we need:

1. On a held-out calibration set of $n$ points, compute the non-conformity residuals $s_i = |Y_i - \hat{f}(X_i)|$ for $i = 1, \dots, n$.
2. Because the future test residual $s_{n+1} = |Y_{n+1} - \hat{f}(X_{n+1})|$ is exchangeable with $s_1, \dots, s_n$, its rank among all $n + 1$ residuals is uniformly distributed on ${1, 2, \dots, n + 1}$.
3. Therefore, setting $\hat{q}$ to the $\lceil (n + 1)(1 - \alpha) \rceil$-th smallest calibration residual $s_{(k)}$ guarantees that $s_{n+1} \le \hat{q}$ with probability at least $1 - \alpha$! Constructing the interval:

$$\hat{C}(x_{n+1}) = \left[\hat{f}(x_{n+1}) - \hat{q}, \quad \hat{f}(x_{n+1}) + \hat{q}\right] \implies P!\left(Y_{n+1} \in \hat{C}(X_{n+1})\right) \ge 1 - \alpha$$

What Happens When the Production Distribution Shifts? (Covariate vs. Label Shift)

In live deployment (Chapter 19), the target environment $Q(X, Y)$ often diverges from the training environment $P(X, Y)$ in one of two structured ways (corresponding to the two ways of factoring $P(X, Y)$):

1. Covariate Drift ($P(X, Y) = P(X)P(Y \mid X)$ where $Q(X) \ne P(X)$ but $Q(Y \mid X) = P(Y \mid X)$): The input demographics or weather distribution $Q(X)$ changes, while the physical conditional law $P(Y \mid X)$ stays invariant. Ryan Tibshirani et al. (2019) proved that weighting each calibration score $s_i$ by the density ratio $w(X_i) = \frac{dQ}{dP}(X_i)$ (Weighted Conformal Prediction) restores coverage, and Lihua Lei & Emmanuel Candès proved that under an estimated tilt weight $\hat{w}$, coverage degrades by at most the Total Variation distance $d_{\mathrm{TV}}$ (ExTRA-WCP-T in Chapter 19): $$P_Q!\left(Y_{n+1} \in \hat{C}(X_{n+1})\right) \ge 1 - \alpha - d_{\mathrm{TV}}$$
2. Label Shift ($P(X, Y) = P(Y)P(X \mid Y)$ where $Q(Y) \ne P(Y)$ but $Q(X \mid Y) = P(X \mid Y)$): In causal diagnostic problems ($Y \to X$, where disease $Y$ causes symptoms $X$), an epidemic can jump the base rate $Q(Y = 1)$ from $5%$ to $40%$ while the symptom profile given the disease $P(X \mid Y)$ stays identical! If we let $C_{jk} = P(\hat{Y} = j \mid Y = k)$ be the source confusion matrix and $\mathbf{q}_{\text{pred}} = Q(\hat{Y} = j)$ be the vector of predicted class frequencies observed in the new target environment, the Law of Total Probability gives the linear system: $$\mathbf{q}_{\text{pred}} = C , \boldsymbol{\pi}_Q \implies \hat{\boldsymbol{\pi}}_Q = C^{-1} \mathbf{q}_{\text{pred}}$$ Solving this equation (via BBSE, MLLS/SLD, DFM, or OTTER Optimal Transport in Chapter 19) recovers the new target prevalence $\boldsymbol{\pi}_Q$ without needing a single new ground-truth label in production!

Computational Audit in Rust

Every method introduced in this chapter is implemented in pure Rust across the workspace:

• Multiplicity & e-Processes (crates/stochology-core/src/card.rs, crates/stochology-siga/src/roadmap.rs): Audits specifications_tried_k, e_process_trajectory, and Ville’s threshold $E_t \ge K / \alpha$.
• Cluster ICC & Moulton-Campbell (crates/stochology-core/src/invariants.rs): Computes $\tilde{m}$, exceedance ICC $\rho_I(p)$, and $n_{\mathrm{eff}}(p)$.
• SBC Rank Uniformity (crates/stochology-siga/src/funnel.rs): Evaluates the $\chi^2$ rank statistic in Filter 3.
• Weighted Conformal (ExTRA-WCP-T) & Label Shift (BBSE, MLLS, DFM, OTTER) (crates/stochology-governance/src/drift.rs): Executes live drift adaptation and circuit breaking.

Run Labs 4, 5, and 6 to verify these mechanisms:

cargo run --release --bin labs -- --invariants
cargo run --release --bin labs -- --siga
cargo run --release --bin labs -- --governance

Exercises

1. Multiplicity Inflation Calculation: If a student runs $K = 14$ independent regressions at $\alpha = 0.05$ on pure white noise and reports only the best $p$-value, what is the probability $\alpha_{\text{family}} = 1 - (1 - 0.05)^{14}$ of getting at least one false positive? What Bonferroni threshold $\alpha / K$ must each test pass to keep the family error $\le 0.05$?
2. Ville’s Inequality and e-Processes: Suppose a sequential test produces likelihood ratios $L_1 = 2.0, ; L_2 = 1.8, ; L_3 = 1.5, ; L_4 = 4.0$. Compute the cumulative e-process trajectory $E_1, E_2, E_3, E_4$ ($E_t = \prod_{i=1}^t L_i$). At which step $t$ does $E_t$ cross the $\alpha = 0.05$ Ville threshold $1/\alpha = 20.0$?
3. Split Conformal Quantile by Hand: Suppose a held-out calibration set of $n = 9$ observations produces sorted absolute residuals $s_{(1..9)} = [0.2, 0.5, 0.8, 1.1, 1.4, 1.7, 2.1, 2.8, 4.5]$. For target coverage $1 - \alpha = 0.80$ ($\alpha = 0.20$), compute the index $k = \lceil (9 + 1)(0.80) \rceil = 8$, find $\hat{q} = s_{(8)}$, and write the $80%$ conformal prediction interval when $\hat{f}(x_{10}) = 10.0$.