Statistics¶
Agentic AI · Fundamentals
The standard error of a success rate, observed in simulation and derived from the Bernoulli distribution; the normal-based 95 % interval; the Wald test for the difference between two success rates.
Overview¶
A success rate measured on a finite test set is itself a random quantity: a second test set of the same size, drawn from the same population of tasks, gives a different rate. The standard error states how far the measured rate typically lies from the true one.
| Section | The question it answers |
|---|---|
| Bernoulli Trial | What is the outcome of one task, and what are its expected value and variance? |
| Sample Proportion | What is the measured success rate $\hat p$, and what is its expected value? |
| Sampling Distribution | How does $\hat p$ vary between test sets of the same size? |
| Variance of the Sample Proportion | Why does the variance of $\hat p$ equal $p(1-p)/n$ for independent tasks? |
| Clustered Tasks | By which factor does the variance grow when several tasks share one source, such as one text passage? |
| Repeated Runs | How much does running each task $K$ times reduce the variance of $\hat p$, and which standard error does such a test set carry? |
| Standard Error | How far does $\hat p$ typically lie from $p$, and how is that distance estimated from the test set alone? |
| Central Limit Theorem & the 95 % Interval | Which interval $\hat p \pm 1.96\,\widehat{\text{se}}$ contains $p$ in 95 % of test sets, and under which condition does that hold? |
| Wald Test for Two Success Rates | Do two systems, measured on separate test sets, differ in their success rates, and what does the p-value of the test state? |
| Paired Comparison | How much smaller is the standard error of the difference when both systems solve the same tasks? |
| Minimum Detectable Effect | Which difference can a test set of $n$ tasks detect, and how many tasks does a given difference require? |
| Standard Error of Other Metrics | Which standard error does accuracy, precision, recall or a mean score carry? |
| Bootstrap | How is the standard error of a metric such as $F_1$ estimated when no formula is at hand? |
Method: derivation first, simulation second. Each section states and derives a result; a simulation then checks it. The simulated agent has a true success rate $p$ that is set in the code, so that the true value is known, which it never is in a real evaluation. The code uses only NumPy and Matplotlib, and the fixed random seed makes every run print the same numbers.
What it builds on: the Bernoulli distribution, the rules for expectation and variance, and the central limit theorem. The notation for random variables and distributions is that of the recap Probability.
Bernoulli Trial¶
One task of a test set is a Bernoulli trial: the agent solves it or does not. The outcome is a random variable $X$ with
$$X = \begin{cases} 1 & \text{with probability } p \quad\text{(success)}\\ 0 & \text{with probability } 1-p \quad\text{(failure)}\end{cases}$$
where $p$ is the true success rate, the probability that the agent solves a task drawn from the population of tasks.
Expected value — the probability-weighted mean of the values: $$\mathbb{E}[X] = 1\cdot p + 0\cdot(1-p) = p.$$
Variance — the expected squared distance from the mean, $\operatorname{Var}(X)=\mathbb{E}[(X-\mathbb{E}[X])^2]$:
$$ \begin{aligned} \operatorname{Var}(X) &= \mathbb{E}[(X-\mathbb{E}[X])^2] \\ &= \mathbb{E}[X^2 - 2X \mathbb{E}[X] + \mathbb{E}[X]^2] \\ &= \mathbb{E}[X^2] - \mathbb{E}[2X \mathbb{E}[X]] + \mathbb{E}\left[\mathbb{E}[X]^2\right] \\ &= \mathbb{E}[X^2] - 2\mathbb{E}[X]\mathbb{E}[X] + \mathbb{E}[X]^2 \\ &= \mathbb{E}[X^2] - \mathbb{E}[X]^2 \end{aligned} $$
with $\mathbb{E}[X^2] = 1^2 \cdot p + 0^2 \cdot (1-p) = p$
Substituting these back into the variance formula:
$$ \begin{aligned} \operatorname{Var}(X) &= p - p^2 \\ &= p(1 - p) \end{aligned} $$
Simulation. A uniform random number in $[0, 1)$ is smaller than $p$ with probability $p$, so the
comparison rng.random() < p is one Bernoulli trial. Over many trials the share of ones approaches
$p$, and the empirical variance approaches $p(1-p)$.
import numpy as np
from statistics import NormalDist
from ipynb.fs.defs.plot_statistics import ( # the figures; their code is in plot_statistics.ipynb
plot_bernoulli_variance, plot_sampling_distribution, plot_clustered, plot_repeated_runs,
plot_standard_error,
plot_standard_normal, plot_clt, plot_intervals, plot_coverage, plot_wald_null, plot_power,
plot_paired_se, plot_bootstrap)
rng = np.random.default_rng(0) # fixed seed: every run gives the same numbers
p = 0.7 # true success rate of the simulated agent
x = (rng.random(20) < p).astype(int) # 20 tasks: 1 = solved, 0 = failed
print("outcomes of 20 tasks:", x)
x = (rng.random(100_000) < p).astype(int) # 100 000 tasks
print(f"mean of X: {x.mean():.4f} E[X] = p = {p:.4f}")
print(f"variance of X: {x.var():.4f} Var(X) = p(1-p) = {p * (1 - p):.4f}")
rng = np.random.default_rng(42)
p = 0.7
N = 100_000
x = (rng.random(N) < p).astype(int)
sample_sizes = np.logspace(1, 5, num=500, dtype=int)
running_vars = [x[:n].var() for n in sample_sizes] # variance of the first n outcomes
plot_bernoulli_variance(p, sample_sizes, running_vars)
Sample Proportion¶
The measured success rate is the mean of $n$ Bernoulli outcomes. A test set of $n$ tasks gives the outcomes $X_1,\dots,X_n$, assumed independent and identically distributed (iid): each $X_i\sim\text{Bernoulli}(p)$, and no outcome changes the distribution of another. The sample proportion is
$$\hat p = \frac{X_1 + X_2 + \dots + X_n}{n} = \frac{1}{n}\sum_{i=1}^{n} X_i .$$
$\hat p$ is a random variable: a different test set gives different $X_i$ and a different $\hat p$.
The expectation is linear. For two random variables $X, Y$ and constants $a, b$, summing the joint probability $P(X=x, Y=y)$ over $y$ gives the marginal $P(X=x)$, and over $x$ the marginal $P(Y=y)$; hence
$$ \begin{aligned} \mathbb{E}[aX + bY] &= \sum_{x,y} (ax + by)\,P(X=x, Y=y) \\ &= a\sum_x x \sum_y P(X=x, Y=y) + b\sum_y y \sum_x P(X=x, Y=y) \\ &= a\sum_x x\,P(X=x) + b\sum_y y\,P(Y=y) \\ &= a\,\mathbb{E}[X] + b\,\mathbb{E}[Y] . \end{aligned} $$
Applied to the $n$ terms of $\hat p$, with $\mathbb{E}[X_i] = p$ from Bernoulli Trial:
$$\mathbb{E}[\hat p] = \frac{1}{n}\sum_{i=1}^{n}\mathbb{E}[X_i] = \frac{1}{n}\cdot n\cdot p = p .$$
The expected value of the measured rate is the true rate: $\hat p$ is an unbiased estimator of $p$ (Wasserman, 2004, Ex. 6.8).
The number of solved tasks follows the binomial distribution. The count $K = X_1 + \dots + X_n = n\hat p$ takes the values $k = 0, 1, \dots, n$ with the probabilities
$$P(K = k) = \binom{n}{k}\,p^k\,(1-p)^{n-k} ,$$
written $K \sim \text{Binomial}(n, p)$ (Wasserman, 2004, §2.3). The factor $p^k(1-p)^{n-k}$ is the probability of one particular sequence of $k$ successes and $n-k$ failures, because the tasks are independent; $\binom{n}{k}$ is the number of such sequences. The distribution of $\hat p$ follows from it: $P(\hat p = k/n) = P(K = k)$.
Simulation. Five test sets of $n = 100$ tasks, drawn from the same agent with $p = 0.7$, give five different values of $\hat p$; none of them needs to equal $p$.
p, n = 0.7, 100
for k in range(1, 6):
x = (rng.random(n) < p).astype(int)
print(f"test set {k}: {x.sum():3d} of {n} tasks solved, p_hat = {x.mean():.2f}")
Sampling Distribution¶
The sampling distribution of $\hat p$ is the distribution of the measured success rate over all test sets of $n$ tasks that could have been drawn (Wasserman, 2004, §6.3). One evaluation observes one value from this distribution. The variance, the standard error and the interval of the following sections are properties of the sampling distribution, and each of them is computed from a single test set.
The sampling distribution of $\hat p$ is known exactly. By Sample Proportion, $P(\hat p = k/n) = P(K = k)$ with $K \sim \text{Binomial}(n, p)$.
Simulation. For an agent with $p = 0.7$, 10 000 test sets are drawn task by task at each of the
sizes $n = 25$, $100$ and $400$; the true binomial probabilities are drawn as a solid line, and the
histogram of the simulated $\hat p$ as a dashed line. A histogram bar has the width $1/n$, so the
solid line is drawn at the height $n\,P(K = k)$, the probability per unit of $\hat p$. The two
lines agree, so the later
simulations draw $K$ directly from the binomial distribution (rng.binomial) instead of drawing
every task. The three histograms show two properties that the following sections derive:
- every histogram is centred at the true rate $p$ — $\hat p$ is unbiased (see Sample Proportion);
- the width shrinks as $n$ grows, but more slowly than $n$ grows: at $n = 400$ the histogram is about half as wide as at $n = 100$.
from math import comb
p, reps = 0.7, 10_000
p_hats, pmfs = {}, {}
for n in (25, 100, 400):
K = (rng.random((reps, n)) < p).sum(axis=1) # solved tasks, drawn task by task
p_hats[n] = K / n
pmfs[n] = np.array([comb(n, j) * p**j * (1 - p)**(n - j) for j in range(n + 1)]) # P(K = k) of Binomial(n, p)
plot_sampling_distribution(p, p_hats, pmfs)
def sample_rates(p, n, reps):
# Success rates p_hat = K/n of `reps` simulated test sets of n tasks each. K is drawn directly
# from Binomial(n, p), which the plot shows to be the distribution of the task-by-task count.
return rng.binomial(n, p, size=reps) / n
Variance of the Sample Proportion¶
The variance of $\hat p$ falls in proportion to $1/n$. The derivation uses two rules of the variance, both obtained from its definition $\operatorname{Var}(X) = \mathbb{E}[(X-\mathbb{E}[X])^2]$ and the linearity of the expectation (see Sample Proportion).
Rule (1): a constant factor is squared. With $\mathbb{E}[aX] = a\,\mathbb{E}[X]$:
$$ \begin{aligned} \operatorname{Var}(aX) &= \mathbb{E}[(aX - a\,\mathbb{E}[X])^2] \\ &= \mathbb{E}[a^2 (X - \mathbb{E}[X])^2] \\ &= a^2\,\operatorname{Var}(X) . \end{aligned} $$
Rule (2): the variance of a sum contains the covariance. The covariance of $X$ and $Y$ is $\operatorname{Cov}(X,Y) = \mathbb{E}[(X-\mathbb{E}[X])(Y-\mathbb{E}[Y])]$. Expanding the square $(u+v)^2 = u^2 + v^2 + 2uv$ with $u = X-\mathbb{E}[X]$ and $v = Y-\mathbb{E}[Y]$:
$$ \begin{aligned} \operatorname{Var}(X+Y) &= \mathbb{E}[((X-\mathbb{E}[X]) + (Y-\mathbb{E}[Y]))^2] \\ &= \mathbb{E}[(X-\mathbb{E}[X])^2] + \mathbb{E}[(Y-\mathbb{E}[Y])^2] + 2\,\mathbb{E}[(X-\mathbb{E}[X])(Y-\mathbb{E}[Y])] \\ &= \operatorname{Var}(X) + \operatorname{Var}(Y) + 2\operatorname{Cov}(X,Y) . \end{aligned} $$
Applied repeatedly to $n$ terms: $\operatorname{Var}\!\left(\sum_i X_i\right) = \sum_i \operatorname{Var}(X_i) + 2\sum_{i<j}\operatorname{Cov}(X_i,X_j)$.
Independence removes the covariance. Expanding the product in the covariance and applying the linearity of the expectation term by term, with the constants $\mathbb{E}[X]$ and $\mathbb{E}[Y]$ taken out of the expectation, gives the same shortcut as for $\operatorname{Var}(X)$ in Bernoulli Trial:
$$ \begin{aligned} \operatorname{Cov}(X,Y) &= \mathbb{E}[(X-\mathbb{E}[X])(Y-\mathbb{E}[Y])] \\ &= \mathbb{E}\big[XY - X\,\mathbb{E}[Y] - \mathbb{E}[X]\,Y + \mathbb{E}[X]\,\mathbb{E}[Y]\big] \\ &= \mathbb{E}[XY] - \mathbb{E}[X]\,\mathbb{E}[Y] - \mathbb{E}[X]\,\mathbb{E}[Y] + \mathbb{E}[X]\,\mathbb{E}[Y] \\ &= \mathbb{E}[XY] - \mathbb{E}[X]\,\mathbb{E}[Y] . \end{aligned} $$
For independent $X$ and $Y$ the joint probability factorises, $P(X=x, Y=y) = P(X=x)\,P(Y=y)$, and therefore
$$ \mathbb{E}[XY] = \sum_{x,y} x\,y\,P(X=x)\,P(Y=y) = \Big(\sum_x x\,P(X=x)\Big)\Big(\sum_y y\,P(Y=y)\Big) = \mathbb{E}[X]\,\mathbb{E}[Y] , $$
so $\operatorname{Cov}(X,Y) = 0$, and the variance of a sum of independent random variables is the sum of their variances (the two rules: Wasserman, 2004, Thm. 3.15, Thm. 3.20).
Applied to $\hat p = \frac{1}{n}\sum_i X_i$:
$$\operatorname{Var}(\hat p) = \operatorname{Var}\!\left(\frac{1}{n}\sum_{i=1}^{n} X_i\right) \overset{(1)}{=} \frac{1}{n^2}\operatorname{Var}\!\left(\sum_{i=1}^{n} X_i\right) \overset{(2)}{=} \frac{1}{n^2}\sum_{i=1}^{n}\operatorname{Var}(X_i) = \frac{1}{n^2}\cdot n\cdot p(1-p) = \frac{p(1-p)}{n} .$$
The last step uses $\operatorname{Var}(X_i) = p(1-p)$ for every $i$ (see Bernoulli Trial). Step (2) is where the independence of the tasks enters. For dependent outcomes the terms $2\operatorname{Cov}(X_i,X_j)$ of rule (2) remain; positively correlated tasks, for example several questions on one text, therefore make $\operatorname{Var}(\hat p)$ larger than $p(1-p)/n$ (see Clustered Tasks).
Simulation. The table compares the variance of 10 000 simulated values of $\hat p$ with $p(1-p)/n$ for test sets of $n = 10$ to $n = 800$ independent tasks.
p, reps = 0.7, 10_000
ns = np.array([10, 25, 50, 100, 200, 400, 800])
var_sim = np.array([sample_rates(p, n, reps).var() for n in ns])
print(" n Var(p_hat), simulated p(1-p)/n")
for n, v in zip(ns, var_sim):
print(f"{n:4d} {v:.6f} {p * (1 - p) / n:.6f}")
Clustered Tasks¶
Tasks are clustered when several of them share a source that makes them easy or hard together, for example several questions on one text passage. Their outcomes are then positively correlated, and the variance of $\hat p$ exceeds $p(1-p)/n$, although the expected value of $\hat p$ is still $p$.
A model in two stages. A test set has $J$ passages with $m$ questions each, $n = Jm$ tasks.
- Each passage $j$ receives its own success probability $\Theta_j$, drawn independently from a Beta Distribution $\text{Beta}(\alpha, \beta)$: a distribution on $[0, 1]$ with the density $f(\theta) \propto \theta^{\alpha-1}(1-\theta)^{\beta-1}$, the mean $\mathbb{E}[\Theta] = \alpha/(\alpha+\beta)$ and the variance $\operatorname{Var}(\Theta) = \alpha\beta/\big((\alpha+\beta)^2(\alpha+\beta+1)\big)$ (Wasserman, 2004, §2.4, §3.4.)
- Given $\Theta_j = \theta$, the $m$ questions of passage $j$ are independent $\text{Bernoulli}(\theta)$ trials.
The two stages correspond to a thumbtack machine: a machine with the parameters $\alpha, \beta$ produces thumbtacks with different probabilities $\theta$ of landing on the head, and each thumbtack is tossed $m$ times. A passage is a thumbtack, and a question is a toss. The two-coin example of Conditional Independence in the recap Probability has the same structure, with a coin drawn from two kinds instead of a probability drawn from a Beta distribution.
Notation. $X_{jk}$ is the outcome of question $k$ of passage $j$, with $j = 1, \dots, J$ and $k = 1, \dots, m$. The outcomes $X_1, \dots, X_n$ of the earlier sections are these $n = Jm$ outcomes numbered consecutively, so $\hat p = \frac1n \sum_{j=1}^{J} \sum_{k=1}^{m} X_{jk}$. The $\Theta_1, \dots, \Theta_J$ are independent and take different values, but all follow $\text{Beta}(\alpha, \beta)$; $\mathbb{E}[\Theta_j]$ and $\operatorname{Var}(\Theta_j)$ are therefore the same numbers for every passage and are written $\mathbb{E}[\Theta]$ and $\operatorname{Var}(\Theta)$.
Each question alone is still a $\text{Bernoulli}(p)$ trial. Averaging $P(X_{jk} = 1 \mid \Theta_j) = \Theta_j$ over the passage probability (the rule of iterated expectations, Wasserman, 2004, Thm. 3.24) gives
$$P(X_{jk} = 1) = \mathbb{E}[\Theta_j] = \frac{\alpha}{\alpha+\beta} =: p .$$
Two questions of one passage are correlated, although they are independent given $\Theta_j$. Two questions $k \ne l$ of passage $j$ are conditionally independent given $\Theta_j$, so $P(X_{jk} = 1, X_{jl} = 1 \mid \Theta_j) = \Theta_j^2$. Without conditioning on $\Theta_j$ they are dependent: a solved question $k$ makes a high value of $\Theta_j$ more probable, and with it the success of question $l$. Averaging over $\Theta_j$ gives the covariance:
$$ \begin{aligned} \operatorname{Cov}(X_{jk}, X_{jl}) &= \mathbb{E}[X_{jk} X_{jl}] - \mathbb{E}[X_{jk}]\,\mathbb{E}[X_{jl}] = P(X_{jk}{=}1, X_{jl}{=}1) - p^2 \\ &= \mathbb{E}[\Theta_j^2] - \mathbb{E}[\Theta_j]^2 = \operatorname{Var}(\Theta) . \end{aligned} $$
Questions of different passages remain independent: $\operatorname{Cov}(X_{jk}, X_{j'l}) = 0$ for $j \ne j'$. The covariance $\operatorname{Var}(\Theta)$ within a passage is positive unless $\Theta$ is a constant, that is, unless every passage has the same success probability. For a Beta distribution with finite $\alpha, \beta$ this does not occur: $\operatorname{Var}(\Theta) > 0$.
The variance of $\hat p$. Rule (2) of Variance of the Sample Proportion keeps the covariance terms. Each passage $j$ contributes $m(m-1)$ ordered pairs $(k, l)$, $k \ne l$, with $\operatorname{Cov}(X_{jk}, X_{jl}) = \operatorname{Var}(\Theta)$, the $J$ passages $J\,m(m-1)$ pairs:
$$ \begin{aligned} \operatorname{Var}(\hat p) &= \frac{1}{n^2}\Big(\sum_{j=1}^{J}\sum_{k=1}^{m}\operatorname{Var}(X_{jk}) + \sum_{j=1}^{J}\sum_{k \ne l}\operatorname{Cov}(X_{jk}, X_{jl})\Big) = \frac{1}{n^2}\Big(n\,p(1-p) + J\,m(m-1)\operatorname{Var}(\Theta)\Big) \\ &= \frac{p(1-p)}{n} + \frac{(m-1)\operatorname{Var}(\Theta)}{n} . \end{aligned} $$
For the Beta distribution, $\alpha\beta/(\alpha+\beta)^2 = p(1-p)$, so $\operatorname{Var}(\Theta) = p(1-p)/(\alpha+\beta+1)$ and
$$\operatorname{Var}(\hat p) = \frac{p(1-p)}{n}\left(1 + \frac{m-1}{\alpha+\beta+1}\right).$$
The factor in parentheses is the design effect $\operatorname{deff}$ of the test set: the variance of $\hat p$ divided by the variance $p(1-p)/n$ of $n$ independent tasks, that is, the factor by which the clustering inflates the variance. It grows with the number $m$ of questions per passage and shrinks as $\alpha + \beta$ grows: a large $\alpha + \beta$ concentrates the passage probabilities around $p$, the passages become alike, and $\operatorname{deff}$ approaches $1$, the case of independent tasks.
Example. $\alpha = 2.1$, $\beta = 0.9$ give $p = 2.1/3 = 0.7$; with $m = 5$ questions per passage the design effect is $\operatorname{deff} = 1 + 4/4 = 2$. A test set of $J = 20$ passages of $m = 5$ questions has $n = 100$ tasks, one task per question. Its variance is $\operatorname{Var}(\hat p) = 2 \cdot p(1-p)/100 = p(1-p)/50$, the variance of a test set of only 50 independent tasks: the clustering halves the number of independent tasks to which the test set is equivalent. The formula $p(1-p)/n$ with $n = 100$ understates the variance by the factor $2$.
Simulation. The left panel shows the Beta density of the passage probabilities together with the histogram of the drawn values. The right panel compares the histogram of $\hat p$ for 10 000 test sets of 20 passages of 5 questions with that for test sets of 100 independent tasks; the code prints both variances and their ratio next to the factor derived above.
from math import gamma
p, reps, J, m = 0.7, 10_000, 20, 5 # J passages of m questions, n = J * m tasks
a, b = 2.1, 0.9 # Beta parameters: mean a / (a + b) = 0.7
n = J * m
theta = rng.beta(a, b, size=(reps, J)) # success probability of each passage
x = rng.random((reps, J, m)) < theta[:, :, None] # each question solved with its passage's theta
p_hat_clustered = x.mean(axis=(1, 2))
p_hat_independent = sample_rates(p, n, reps)
print(f"mean of p_hat, clustered: {p_hat_clustered.mean():.3f} (p = {p})")
print(f"Var(p_hat), clustered: {p_hat_clustered.var():.6f}")
print(f"Var(p_hat), independent: {p_hat_independent.var():.6f} (p(1-p)/n = {p * (1 - p) / n:.6f})")
print(f"ratio, simulated: {p_hat_clustered.var() / p_hat_independent.var():.2f}")
print(f"factor 1 + (m-1)/(a+b+1): {1 + (m - 1) / (a + b + 1):.2f}")
grid = np.linspace(0.001, 0.999, 500)
density = gamma(a + b) / (gamma(a) * gamma(b)) * grid**(a - 1) * (1 - grid)**(b - 1) # Beta(a, b) density
plot_clustered(theta, grid, density, a, b, p, p_hat_independent, p_hat_clustered, J, m)
The design effect depends on the passage probabilities only through the correlation of two questions of one passage. This correlation is the intra-class correlation $\rho$ of the test set, for two questions $k \ne l$ of passage $j$. The correlation of two random variables is their covariance divided by the product of their standard deviations (Wasserman, 2004, Def. 3.18); with the variance and the covariance derived above:
$$\rho = \frac{\operatorname{Cov}(X_{jk}, X_{jl})}{\sqrt{\operatorname{Var}(X_{jk})\,\operatorname{Var}(X_{jl})}} = \frac{\operatorname{Var}(\Theta)}{p(1-p)} = \frac{1}{\alpha+\beta+1} .$$
$\rho$ is therefore also the share of the variance $p(1-p)$ of one question that comes from the differences between the passage probabilities. It is $0$ when every passage has the same success probability, the case of independent tasks, and $1$ when every $\Theta_j$ is $0$ or $1$, so that all questions of a passage have the same outcome. With $\operatorname{Var}(\Theta) = \rho\,p(1-p)$, the design effect is $\operatorname{deff} = 1 + (m-1)\rho$.
The design effect states how many independent tasks a clustered test set is worth. The effective sample size $n_{\text{eff}} = n/\operatorname{deff}$ is the number of independent tasks with the same variance:
$$\operatorname{Var}(\hat p) = \frac{p(1-p)}{n}\,\operatorname{deff} = \frac{p(1-p)}{n_{\text{eff}}} .$$
A clustered test set therefore reaches the standard error of $n$ independent tasks only with $n \cdot \operatorname{deff}$ tasks. In the example, $\operatorname{deff} = 2$: the standard error of 100 independent tasks requires 200 tasks, for example 40 passages of 5 questions.
Only further passages reduce the variance without limit. Substituting $n = Jm$:
$$ \begin{aligned} \operatorname{Var}(\hat p) &= \frac{p(1-p)}{Jm}\big(1 + (m-1)\rho\big) && n = Jm \\ &= \frac{p(1-p)}{J}\left(\rho + \frac{1-\rho}{m}\right) && \frac{1 + (m-1)\rho}{m} = \rho + \frac{1-\rho}{m} \\ &= \frac{\operatorname{Var}(\Theta)}{J} + \frac{p(1-p) - \operatorname{Var}(\Theta)}{Jm} && \rho\,p(1-p) = \operatorname{Var}(\Theta) . \end{aligned} $$
More questions per passage ($m$ grows, $J$ fixed) reduce only the second term, and the variance does not fall below $\operatorname{Var}(\Theta)/J$, the variance of the mean of the $J$ passage probabilities. More passages ($J$ grows) reduce both terms in proportion to $1/J$. With $\rho = 1/4$ and $J = 20$ passages, as in the example, no number of questions per passage reduces the variance below $p(1-p)/80$, the variance of 80 independent tasks.
Evidence. On DROP (Dua et al., 2019), a reading-comprehension benchmark in which several questions refer to one text paragraph, the standard error of the mean question score of an Anthropic model is 3.05 times larger when it is computed with the paragraph as the unit than when the questions are treated as independent (Miller, 2024, §2.2, Tab. 4). A ratio of 3.05 between the standard errors corresponds to an estimated design effect of $3.05^2 \approx 9.3$.
Consequence for an evaluation. A clustered test set is enlarged by passages, not by questions per passage; its standard error is computed with the passage as the unit, and the number of passages is reported beside the number of questions (Miller, 2024, §2.2).
Repeated Runs¶
The output of an agent is sampled from a model, so the same task, run $K$ times, need not give the same outcome. The mean of the $K$ runs of each task reduces the variance of $\hat p$, but only the part that comes from the sampling of runs; the part that comes from the choice of the tasks remains (Miller, 2024, §3, §3.1). Repeated runs therefore lower the standard error towards a limit that only more tasks can lower further. The model is that of Clustered Tasks, with a task in the place of the passage and a run in the place of the question.
A model in two stages. A test set has $n$ tasks, and each task is run $K$ times.
- Each task $i$ receives its own success probability $\Theta_i$, the probability that one run solves it. $\Theta_1, \dots, \Theta_n$ are independent, follow one distribution on $[0, 1]$, and $\mathbb{E}[\Theta] = p$.
- Given $\Theta_i = \theta$, the runs $X_{i1}, \dots, X_{iK}$ of task $i$ are independent $\text{Bernoulli}(\theta)$ trials.
The task score $\bar X_i = \frac1K \sum_{k=1}^{K} X_{ik}$ is the share of the $K$ runs that solve task $i$, and $\hat p = \frac1n \sum_{i=1}^{n} \bar X_i$. Because every task has the same $K$ runs, $\hat p$ is also the share of all $nK$ runs that succeed.
One run and two runs of a task. The derivation in Clustered Tasks of $P(X_{jk} = 1) = \mathbb{E}[\Theta]$ and $\operatorname{Cov}(X_{jk}, X_{jl}) = \operatorname{Var}(\Theta)$ uses only the two stages, not the Beta form of the distribution of $\Theta$. For two runs $k \ne l$ of one task it therefore gives $\operatorname{Var}(X_{ik}) = p(1-p)$ and $\operatorname{Cov}(X_{ik}, X_{il}) = \operatorname{Var}(\Theta)$.
The variance of one run splits into two parts. Averaged over the tasks, a run is a $\text{Bernoulli}(p)$ trial, since $P(X_{ik} = 1) = \mathbb{E}[\Theta] = p$; its variance is $\operatorname{Var}(X_{ik}) = p(1-p)$. Given $\Theta_i$, the same run is a $\text{Bernoulli}(\Theta_i)$ trial with the conditional variance $\Theta_i(1-\Theta_i)$ (see Bernoulli Trial). The mean of this conditional variance over the tasks is
$$ \begin{aligned} \mathbb{E}[\Theta(1-\Theta)] &= \mathbb{E}[\Theta] - \mathbb{E}[\Theta^2] && \text{linearity of } \mathbb{E} \\ &= p - \big(\operatorname{Var}(\Theta) + p^2\big) && \mathbb{E}[\Theta^2] = \operatorname{Var}(\Theta) + \mathbb{E}[\Theta]^2 \\ &= p(1-p) - \operatorname{Var}(\Theta) , \end{aligned} $$
so that
$$\operatorname{Var}(X_{ik}) = p(1-p) = \underbrace{\operatorname{Var}(\Theta)}_{\text{between tasks}} + \underbrace{\mathbb{E}[\Theta(1-\Theta)]}_{\text{within a task}} .$$
The first part is the variance of the success probability from task to task, the second the mean variance of the runs of one task. With $\mathbb{E}[X_{ik} \mid \Theta_i] = \Theta_i$ and $\operatorname{Var}(X_{ik} \mid \Theta_i) = \Theta_i(1-\Theta_i)$, the two parts are the two terms of the law of total variance, $\operatorname{Var}(X_{ik}) = \operatorname{Var}\big(\mathbb{E}[X_{ik} \mid \Theta_i]\big) + \mathbb{E}\big[\operatorname{Var}(X_{ik} \mid \Theta_i)\big]$, which Miller (2024, §3) states for any score.
The variance of $\hat p$. Rules (1) and (2) of Variance of the Sample Proportion, applied to the $K$ runs of one task with their $K$ variances and $K(K-1)$ covariances, and then to the $n$ independent task scores:
$$ \begin{aligned} \operatorname{Var}(\bar X_i) &= \operatorname{Var}\Big(\frac1K \sum_{k=1}^{K} X_{ik}\Big) && \text{definition of } \bar X_i \\ &= \frac{1}{K^2}\operatorname{Var}\Big(\sum_{k=1}^{K} X_{ik}\Big) && \text{rule (1) with } a = 1/K \\ &= \frac{1}{K^2}\Big(\sum_{k=1}^{K}\operatorname{Var}(X_{ik}) + \sum_{k \ne l}\operatorname{Cov}(X_{ik}, X_{il})\Big) && \text{rule (2), over the } K(K-1) \text{ ordered pairs } (k, l) \\ &= \frac{1}{K^2}\Big(K\,p(1-p) + K(K-1)\operatorname{Var}(\Theta)\Big) && \operatorname{Var}(X_{ik}) = p(1-p);\ \operatorname{Cov}(X_{ik}, X_{il}) = \operatorname{Var}(\Theta) \text{ for } k \ne l \\ &= \frac{p(1-p)}{K} + \frac{K-1}{K}\operatorname{Var}(\Theta) \\ &= \operatorname{Var}(\Theta) + \frac{p(1-p) - \operatorname{Var}(\Theta)}{K} && \tfrac{K-1}{K} = 1 - \tfrac{1}{K} \\ &= \operatorname{Var}(\Theta) + \frac{\mathbb{E}[\Theta(1-\Theta)]}{K} && \text{the decomposition} \\[4pt] \operatorname{Var}(\hat p) &= \frac{\operatorname{Var}(\bar X_1)}{n} = \frac{\operatorname{Var}(\Theta)}{n} + \frac{\mathbb{E}[\Theta(1-\Theta)]}{nK} && \text{independent, identically distributed task scores, rules (1) and (2)} . \end{aligned} $$
The first term depends only on the number $n$ of tasks; the second falls with the number $K$ of runs. As $K$ grows, $\operatorname{Var}(\hat p)$ approaches $\operatorname{Var}(\Theta)/n$ and not $0$: repeated runs do not replace tasks. A system that gives the same outcome on every run of a task has $\Theta_i \in \{0, 1\}$, so $\Theta_i(1-\Theta_i) = 0$, and repeated runs do not reduce its variance at all.
The variance of $\hat p$ depends on the number of runs through the correlation of two runs of one task. The intra-class correlation $\rho$ is the correlation of two runs $k \ne l$ of task $i$. With $\operatorname{Var}(X_{ik}) = p(1-p)$ and $\operatorname{Cov}(X_{ik}, X_{il}) = \operatorname{Var}(\Theta)$ of one and two runs of a task:
$$\rho = \frac{\operatorname{Cov}(X_{ik}, X_{il})}{\sqrt{\operatorname{Var}(X_{ik})\,\operatorname{Var}(X_{il})}} = \frac{\operatorname{Var}(\Theta)}{p(1-p)} ,$$
the share of the variance $p(1-p)$ of one run that lies between tasks. It is the intra-class correlation of Clustered Tasks, with the task in place of the passage and the run in place of the question. With $\operatorname{Var}(\Theta) = \rho\,p(1-p)$ and, by the decomposition, $\mathbb{E}[\Theta(1-\Theta)] = (1-\rho)\,p(1-p)$:
$$\operatorname{Var}(\hat p) = \frac{\rho\,p(1-p)}{n} + \frac{(1-\rho)\,p(1-p)}{nK} = \frac{p(1-p)}{n}\left(\rho + \frac{1-\rho}{K}\right) .$$
With $K = 1$ this is $p(1-p)/n$, the variance of $n$ tasks with one run each; as $K$ grows it approaches $\rho\,p(1-p)/n = \operatorname{Var}(\Theta)/n$. For task probabilities $\Theta_i \sim \text{Beta}(\alpha, \beta)$, $\operatorname{Var}(\Theta) = p(1-p)/(\alpha+\beta+1)$ as for the passage probabilities, so $\rho = 1/(\alpha+\beta+1)$; the parameters $\alpha = 2.1$, $\beta = 0.9$ of the passage example give $\rho = 1/4$:
| Runs per task $K$ | $1$ | $2$ | $4$ | $K \to \infty$ |
|---|---|---|---|---|
| $\operatorname{Var}(\hat p)$ in units of $p(1-p)/n$ | $1$ | $5/8$ | $7/16$ | $1/4$ |
| share of the variance at $K = 1$ removed | $0$ | $3/8$ | $9/16$ | $3/4$ |
For task probabilities uniform on $[0, 1]$, $\text{Beta}(1, 1)$ with $\rho = 1/3$, $K = 2$ removes a third of the variance, $K = 4$ half, and no $K$ more than two thirds (Miller, 2024, §3.1).
The $nK$ runs are not $nK$ independent tasks. As $nK$ independent tasks, the runs would give $\hat p$ the variance $p(1-p)/(nK)$. The variance derived above is larger by the factor
$$\frac{p(1-p)}{n}\left(\rho + \frac{1-\rho}{K}\right) \Big/ \frac{p(1-p)}{nK} = K\rho + 1 - \rho = 1 + (K-1)\rho ,$$
which is $1.75$ for $K = 4$ and $\rho = 1/4$: the $4n$ runs have the variance of $4n/1.75 \approx 2.3\,n$ independent tasks.
Consequence for an evaluation. The runs of a task are averaged into one task score, and the task, not the run, is the unit of the evaluation. Further runs reduce the variance of $\hat p$ only as long as $\mathbb{E}[\Theta(1-\Theta)]/K$ is not small against $\operatorname{Var}(\Theta)$; a variance below $\operatorname{Var}(\Theta)/n$ requires more tasks.
Simulation. Test sets of $n = 100$ tasks with $\Theta_i \sim \text{Beta}(2.1, 0.9)$ are run $K = 1$ to $K = 32$ times per task; given $\Theta_i$, the number of solved runs of task $i$ is $\text{Binomial}(K, \Theta_i)$. For each $K$, 10 000 test sets give the variance of $\hat p$, which the code prints next to the formula. The figure compares the simulated variance with the formula, with its limit $\operatorname{Var}(\Theta)/n$ and with the variance $p(1-p)/(nK)$ of $nK$ independent tasks. The last two printed columns, the mean estimated standard error of the pooled runs and of the task scores, are explained in Standard Error of Other Metrics.
p, n, reps = 0.7, 100, 10_000
a, b = 2.1, 0.9 # Beta parameters of the task probabilities, mean 0.7
rho = 1 / (a + b + 1) # share of p(1-p) between tasks: Var(Theta) = rho * p(1-p)
Ks = np.array([1, 2, 4, 8, 16, 32]) # runs per task
var_runs, se_pooled, se_tasks = [], [], []
for K in Ks:
theta = rng.beta(a, b, size=(reps, n)) # success probability of each task
scores = rng.binomial(K, theta) / K # task score: share of the K runs that solve the task
p_hat = scores.mean(axis=1)
var_runs.append(p_hat.var())
se_pooled.append(np.sqrt(p_hat * (1 - p_hat) / (n * K)).mean()) # all n*K runs treated as independent
se_tasks.append((scores.std(axis=1, ddof=1) / np.sqrt(n)).mean()) # S / sqrt(n) over the n task scores
var_runs = np.array(var_runs)
var_formula = p * (1 - p) / n * (rho + (1 - rho) / Ks)
print(" K Var(p_hat), simulated formula sd(p_hat) mean se, pooled runs mean se, task scores")
for K, v, f, sp, st in zip(Ks, var_runs, var_formula, se_pooled, se_tasks):
print(f"{K:2d} {v:.6f} {f:.6f} {np.sqrt(v):.4f} {sp:.4f} {st:.4f}")
print(f"limit Var(Theta)/n = {rho * p * (1 - p) / n:.6f}")
plot_repeated_runs(Ks, var_runs, var_formula, p, n, rho)
Standard Error¶
The standard error of an estimator is the standard deviation of its sampling distribution, the distribution of the estimator over repeated samples (Wasserman, 2004, §6.3, eq. 6.5). For the sample proportion it is the square root of the variance:
$$\text{se}(\hat p) = \sqrt{\operatorname{Var}(\hat p)} = \sqrt{\frac{p(1-p)}{n}} .$$
The true rate $p$ is unknown, so the measured rate takes its place; the result is the estimated standard error (Wasserman, 2004, Ex. 6.8):
$$\widehat{\text{se}}(\hat p) = \sqrt{\frac{\hat p(1-\hat p)}{n}} .$$
The plug-in estimate $\hat p(1-\hat p)$ underestimates the variance by the factor $(n-1)/n$. For 0/1 outcomes $X_i^2 = X_i$, so $\hat p(1-\hat p)$ is the mean squared deviation with the divisor $n$:
$$ \begin{aligned} \frac1n\sum_{i=1}^{n} (X_i - \hat p)^2 &= \frac1n\sum_{i=1}^{n} X_i^2 - \hat p^2 && \text{the shortcut of } \textit{Bernoulli Trial} \\ &= \hat p - \hat p^2 = \hat p(1-\hat p) && X_i^2 = X_i . \end{aligned} $$
Its expected value follows from $\mathbb{E}[\hat p] = p$ and $\operatorname{Var}(\hat p) = p(1-p)/n$:
$$ \begin{aligned} \mathbb{E}[\hat p(1-\hat p)] &= \mathbb{E}[\hat p] - \mathbb{E}[\hat p^2] && \text{linearity of } \mathbb{E} \\ &= p - \big(\operatorname{Var}(\hat p) + p^2\big) && \mathbb{E}[\hat p^2] = \operatorname{Var}(\hat p) + \mathbb{E}[\hat p]^2 \\ &= p(1-p) - \frac{p(1-p)}{n} = \frac{n-1}{n}\,p(1-p) . \end{aligned} $$
The sample variance $S^2 = \frac{1}{n-1}\sum_i (X_i - \hat p)^2 = \frac{n}{n-1}\,\hat p(1-\hat p)$ removes this bias, $\mathbb{E}[S^2] = p(1-p)$ (Wasserman, 2004, Thm. 3.17); it is the estimate that Standard Error of Other Metrics uses for a mean score.
Neither divisor gives an unbiased standard error. The variance of $S$ is not negative, so $\operatorname{Var}(S) = \mathbb{E}[S^2] - \mathbb{E}[S]^2 \ge 0$ gives $\mathbb{E}[S] \le \sqrt{\mathbb{E}[S^2]}$, with equality only when $S$ is constant: $S/\sqrt n$, too, underestimates the standard error on average. The two estimates of the standard error differ by the factor $\sqrt{n/(n-1)}$, which is $1.005$ for $n = 100$ and approaches $1$ as $n$ grows: at $n = 100$ the plug-in estimate is $0.995$ times the estimate with $S$. The recap uses $\hat p(1-\hat p)/n$ for a success rate, as Wasserman (2004, Ex. 6.8) and Miller (2024), §2.1, Eq. 2, do.
Clustered tasks enlarge the standard error by the square root of the factor. The standard error is the square root of the variance, so the factor $1 + (m-1)/(\alpha+\beta+1)$ of Clustered Tasks enters it under the root. In the example there, the factor is $2$: for 20 passages of 5 questions, $\sqrt{p(1-p)/n}$ with $n = 100$ understates the standard error by the factor $\sqrt 2 \approx 1.41$.
The standard error of a success rate is the standard error of the mean of 0/1 outcomes. $\hat p$ is the sample mean of the outcomes $X_i$, so $\text{se}(\hat p)$ is the standard error of the mean (SEM) $\sigma/\sqrt n$ for Practical Reporting of Measurement Uncertainty in Estimators and Uncertainties , with $\sigma^2 = \operatorname{Var}(X_i) = p(1-p)$. One difference follows from the Bernoulli distribution: the variance is a function of the mean, so it is not estimated separately from the spread of the data; $\hat p$ is inserted into $p(1-p)$.
A reported uncertainty names its multiplier. A result is reported as best estimate ± uncertainty, and the same $\widehat{\text{se}}$ gives two common reports: $\hat p \pm \widehat{\text{se}}$ contains $p$ with probability about $0.68$, and $\hat p \pm 1.96\,\widehat{\text{se}}$ with probability about $0.95$ (derived in Central Limit Theorem & the 95 % Interval). This recap reports the 95 % interval $\hat p \pm 1.96\,\widehat{\text{se}}$. A report that states only ± 0.05 is ambiguous: the 0.05 may be one standard error or 1.96 standard errors.
Example. $\hat p = 0.6$ on $n = 100$ tasks: $\widehat{\text{se}} = \sqrt{0.6\cdot 0.4/100} \approx 0.049$. On $n = 400$ tasks: $\widehat{\text{se}} \approx 0.0245$. Because $n$ stands under the square root, the standard error halves only when $n$ quadruples.
Simulation. The left panel plots the standard deviation of the simulated $\hat p$ (from the variance table in Variance of the Sample Proportion) against $n$ on logarithmic axes: the points lie on the line $\sqrt{p(1-p)/n}$, whose slope is $-1/2$. The right panel computes $\widehat{\text{se}}$ on each of 10 000 test sets of $n = 100$ tasks. $\widehat{\text{se}}$ is a function of $\hat p$, so it varies from test set to test set because $\hat p$ does, and it varies around the true $\text{se}$.
p, n, reps = 0.7, 100, 10_000
p_hat = sample_rates(p, n, reps)
se_hat = np.sqrt(p_hat * (1 - p_hat) / n)
se = np.sqrt(p * (1 - p) / n)
plot_standard_error(ns, var_sim, p, se_hat, se, n)
Central Limit Theorem & the 95 % Interval¶
For large $n$, $\hat p$ is approximately normally distributed around $p$, and this turns the standard error into an interval. The central limit theorem states that the mean of $n$ iid random variables with mean $\mu$ and variance $\sigma^2$, standardised to $\sqrt{n}\,(\bar X_n-\mu)/\sigma$, converges in distribution, written $\xrightarrow{d}$, to the standard normal distribution $\mathcal{N}(0,1)$ (Wasserman, 2004, Thm. 5.8). The theorem is not derived here; the simulation shows it for the sample proportion.
Standardising $\hat p$. For the sample proportion $\bar X_n = \hat p$, $\mu = p$ and $\sigma = \sqrt{p(1-p)}$, so the standardised mean is
$$\frac{\sqrt{n}\,(\hat p - p)}{\sqrt{p(1-p)}} = \frac{\hat p - p}{\sqrt{p(1-p)/n}} = \frac{\hat p - p}{\text{se}} \;\xrightarrow{\;d\;}\; \mathcal{N}(0,1) \qquad (n \to \infty) .$$
For a fixed, large $n$ the limit is used as an approximation, written $\overset{\text{approx.}}{\sim}$:
$$\frac{\hat p - p}{\text{se}} \;\overset{\text{approx.}}{\sim}\; \mathcal{N}(0,1), \qquad\text{equivalently}\qquad \hat p \;\overset{\text{approx.}}{\sim}\; \mathcal{N}\!\left(p,\ \text{se}^2\right) .$$
Simulation. The code standardises 10 000 simulated values of $\hat p$, at $p = 0.7$, to $(\hat p - p)/\text{se}$, and the figure overlays the density of $\mathcal{N}(0,1)$: at $n = 10$ the histogram is coarse and asymmetric, at $n = 100$ it follows the normal density.
p, reps = 0.7, 10_000
z = {}
for n in (10, 100):
se = np.sqrt(p * (1 - p) / n)
z[n] = (sample_rates(p, n, reps) - p) / se # standardised p_hat
plot_clt(p, z)
The constant 1.96. Let $Z \sim \mathcal{N}(0,1)$ with density $\varphi$ and cumulative distribution function (CDF) $\Phi$:
$$\varphi(t) = \frac{1}{\sqrt{2\pi}}\,e^{-t^2/2}, \qquad \Phi(z) = P(Z \le z) = \int_{-\infty}^{z} \varphi(t)\,dt .$$
The density is symmetric around $0$, $\varphi(-t) = \varphi(t)$, so the substitution $t = -u$ gives
$$ \begin{aligned} \Phi(-z) &= \int_{-\infty}^{-z} \varphi(t)\,dt \\ &= \int_{z}^{\infty} \varphi(-u)\,du && t = -u,\ dt = -du \\ &= \int_{z}^{\infty} \varphi(u)\,du && \varphi \text{ symmetric} \\ &= 1 - \Phi(z) && \textstyle\int_{-\infty}^{\infty} \varphi(u)\,du = 1 \end{aligned} $$
and
$$P(-z \le Z \le z) = \Phi(z) - \Phi(-z) = 2\,\Phi(z) - 1 .$$
Setting this probability to $0.95$ gives $\Phi(z) = 0.975$, and $z = \Phi^{-1}(0.975) \approx 1.96$. For $z = 1$ the same formula gives $2\,\Phi(1) - 1 \approx 0.68$, the probability for the report $\hat p \pm \widehat{\text{se}}$.
Simulation. The code computes $\Phi^{-1}(0.975)$. The left panel of the figure shows the probability $P(-z \le Z \le z)$ as the area under $\varphi$: $0.95$ for $z = 1.96$ (shaded), $0.68$ for $z = 1$ (hatched) and $0.025$ in each tail beyond $\pm 1.96$. The right panel reads $\Phi(1.96) = 0.975$ and $\Phi(-1.96) = 0.025$ off the CDF.
z_star = NormalDist().inv_cdf(0.975)
print(f"Phi^-1(0.975) = {z_star:.4f}")
plot_standard_normal(z_star)
The interval. The event $-1.96 \le (\hat p - p)/\text{se} \le 1.96$ is rearranged into an interval around $\hat p$:
$$ \begin{aligned} 0.95 &\approx P\!\left(-1.96 \le \frac{\hat p - p}{\text{se}} \le 1.96\right) \\ &= P\!\left(\hat p - 1.96\,\text{se} \le p \le \hat p + 1.96\,\text{se}\right) . \end{aligned} $$
The random quantity in the second line is the interval, not $p$: over repeated test sets, the interval contains the true rate $p$ in about 95 % of the cases. The unknown $\text{se}$ is replaced by $\widehat{\text{se}}$, which gives the normal-based 95 % confidence interval (Wasserman, 2004, Thm. 6.16)
$$\hat p \pm 1.96\,\widehat{\text{se}}(\hat p) .$$ Example. $\hat p = 0.6$ on $n=100$: $0.6 \pm 1.96\cdot 0.049 = 0.6 \pm 0.096$, i.e. $60 \pm 9.6$ percentage points; on $n = 400$, $60 \pm 4.8$ percentage points.
Simulation. The figure draws the 95 % intervals of 100 test sets of $n = 100$ tasks at $p = 0.7$: an interval that misses $p$ is red, and about 5 of the 100 miss it. Below the figure, the code prints the coverage, the share of 100 000 test sets whose interval contains $p$.
def wald_interval(p_hat, n):
# Normal-based 95 % interval p_hat +- 1.96 * se_hat (see "Central Limit Theorem & the 95 % Interval").
se_hat = np.sqrt(p_hat * (1 - p_hat) / n)
return p_hat - 1.96 * se_hat, p_hat + 1.96 * se_hat
p, n, sets = 0.7, 100, 100
p_hat = sample_rates(p, n, sets)
lo, hi = wald_interval(p_hat, n)
plot_intervals(p_hat, lo, hi, p, n)
lo, hi = wald_interval(sample_rates(p, n, 100_000), n)
print(f"coverage over 100 000 test sets: {np.mean((lo <= p) & (p <= hi)):.3f}")
The interval is valid only when both successes and failures are numerous. The normal approximation describes the count of successes well only when neither $n\hat p$ nor $n(1-\hat p)$ is small; the interval is a large-sample approximation (Wasserman, 2004, §6.3). At $\hat p = 0.02$ on $n = 100$ tasks, $n\hat p = 2$, and the interval $0.02 \pm 0.027$ extends below zero, which no success rate can.
Simulation. The coverage — the share of 100 000 test sets of $n = 100$ tasks whose interval contains $p$ — is computed for every true rate from $p = 0.01$ to $p = 0.99$. It is close to 0.95 in the middle of the range and falls far below 0.95 near $0$ and $1$. At $\hat p = 0$ the formula gives $\widehat{\text{se}} = 0$ and an interval of width zero, which cannot contain $p$.
n, reps = 100, 100_000
ps = np.linspace(0.01, 0.99, 99)
coverage = []
for p in ps:
lo, hi = wald_interval(sample_rates(p, n, reps), n)
coverage.append(np.mean((lo <= p) & (p <= hi)))
plot_coverage(ps, coverage, n)
Wald Test for Two Success Rates¶
Two systems differ measurably only when the difference of their success rates exceeds about two standard errors of that difference. The Wald test of the null hypothesis $\delta = p_1 - p_2 = 0$ uses the statistic
$$W = \frac{\hat\delta}{\widehat{\text{se}}(\hat\delta)}, \qquad \hat\delta = \hat p_1 - \hat p_2 ,$$
and rejects the null hypothesis at level $0.05$ when $|W| > 1.96$ (Wasserman, 2004, Def. 10.3, Ex. 10.7).
The distribution of $W$. On separate test sets of sizes $m$ and $n$, $\hat p_1$ and $\hat p_2$ are independent, and each is approximately normal (see Central Limit Theorem & the 95 % Interval). By the linearity of the expectation $\mathbb{E}[\hat\delta] = p_1 - p_2 = \delta$; by rule (1) with $a = -1$ and rule (2) with zero covariance (see Variance of the Sample Proportion),
$$\operatorname{Var}(\hat p_1 - \hat p_2) = \operatorname{Var}(\hat p_1) + (-1)^2\operatorname{Var}(\hat p_2) = \frac{p_1(1-p_1)}{m} + \frac{p_2(1-p_2)}{n} .$$
A difference of two independent normal random variables is normal, so $\hat\delta \overset{\text{approx.}}{\sim} \mathcal{N}(\delta,\ \text{se}(\hat\delta)^2)$. When $\delta = 0$, $W$ is $\hat\delta$ standardised, and $W \overset{\text{approx.}}{\sim} \mathcal{N}(0,1)$.
On separate test sets, $\hat p_1$ and $\hat p_2$ in place of $p_1$ and $p_2$ give the estimate $\widehat{\text{se}}(\hat\delta) = \sqrt{\hat p_1(1-\hat p_1)/m + \hat p_2(1-\hat p_2)/n}$. When both systems solve the same $n$ tasks, $\hat p_1$ and $\hat p_2$ are not independent, the covariance term of rule (2) does not vanish, and the standard error is the one derived in Paired Comparison.
Example. Two systems measured on separate test sets of $m = n = 100$ tasks reach $\hat p_1 = 0.02$ and $\hat p_2 = 0.04$, so $\hat\delta = \hat p_1 - \hat p_2 = -0.02$ and
$$\widehat{\text{se}}(\hat\delta) = \sqrt{\frac{\hat p_1(1-\hat p_1)}{m} + \frac{\hat p_2(1-\hat p_2)}{n}} = \sqrt{\frac{0.02\cdot0.98}{100} + \frac{0.04\cdot0.96}{100}} \approx 0.024 ,$$
$|W| = 0.02/0.024 \approx 0.83 < 1.96$ — no evidence of a difference. With $n\hat p_1 = 2$ and $n\hat p_2 = 4$ successes, these rates violate the condition of the normal approximation (see Central Limit Theorem & the 95 % Interval), so the value of $|W|$ is only approximate.
Size. The size of a test is the probability that it rejects the null hypothesis when the null hypothesis is true; such a rejection is a type I error, a difference reported where none exists. When $\delta = 0$, $W \overset{\text{approx.}}{\sim} \mathcal{N}(0,1)$, and by the derivation of the constant 1.96
$$P(|W| > 1.96) = 1 - P(-1.96 \le W \le 1.96) \approx 1 - (2\,\Phi(1.96) - 1) = 0.05 .$$
The Wald test therefore has size about $0.05$, and the approximation becomes exact as $n \to \infty$ (Wasserman, 2004, Thm. 10.4): of 100 comparisons of two equally good systems, about 5 report a difference.
p-value. A test that reports only reject or retain discards how far $|W|$ lies beyond the threshold. The p-value is the smallest level $\alpha$ at which the test rejects the null hypothesis (Wasserman, 2004, Def. 10.11). For the Wald test with the observed value $w$ of $W$, it is the probability, when $\delta = 0$, of a statistic at least as far from $0$ as $w$ (Wasserman, 2004, Thm. 10.13):
$$ \begin{aligned} p\text{-value} &= P_{\delta = 0}(|W| \ge |w|) \approx P(|Z| \ge |w|) && W \overset{\text{approx.}}{\sim} \mathcal{N}(0,1) \text{ when } \delta = 0;\ Z \sim \mathcal{N}(0,1) \\ &= P(Z \ge |w|) + P(Z \le -|w|) && \text{disjoint events} \\ &= 2\,\Phi(-|w|) && P(Z \ge a) = P(Z \le -a) \text{ by the symmetry of } \mathcal{N}(0,1) . \end{aligned} $$
$\Phi$ is increasing, so $|w| > 1.96$ holds exactly when $2\,\Phi(-|w|) < 2\,\Phi(-1.96) = 0.05$: the test rejects at level $0.05$ exactly when the p-value is below $0.05$. In the example of two systems with $\hat p_1 = 0.02$ and $\hat p_2 = 0.04$, $|w| \approx 0.83$ gives the p-value $2\,\Phi(-0.83) \approx 0.41$.
What a p-value does not state (Wasserman, 2004, §10.1, §10.2):
| Reading | Why it is wrong | |---|---| | the probability that the two systems are equally good | the p-value is computed under the assumption $\delta = 0$; it is a probability of the data given $\delta = 0$, not of $\delta = 0$ given the data | | a large p-value shows that the systems are equally good | a large p-value also arises when $\delta \ne 0$ and the power of the test is low, as for small test sets | | a small p-value shows a large difference | $|W| \approx |\hat\delta|/\text{se}$ and $\text{se}$ falls with $1/\sqrt n$: for a fixed $\delta \neq 0$, $|W|$ grows with $\sqrt n$, and on enough tasks a difference of 0.5 percentage points gives a small p-value |
The size of $0.05$ also holds for each of several comparisons: of 20 comparisons of equally good systems, on average one has a p-value below $0.05$.
Power. The power of a test is the probability that it rejects the null hypothesis when the null hypothesis is false; a missed rejection is a type II error, a real difference that the test does not report. When $\delta \neq 0$, replacing $\widehat{\text{se}}$ by $\text{se}$ gives $W \approx \hat\delta/\text{se}$, and
$$\frac{\hat\delta}{\text{se}} = \frac{\hat\delta - \delta}{\text{se}} + \frac{\delta}{\text{se}}$$
is a standard normal variable shifted by $\delta/\text{se}$: $W \overset{\text{approx.}}{\sim} \mathcal{N}(\delta/\text{se},\ 1)$. Let $Z = W - \delta/\text{se}$; then $Z \overset{\text{approx.}}{\sim} \mathcal{N}(0,1)$, and
$$ \begin{aligned} P(|W| > 1.96) &= P(W > 1.96) + P(W < -1.96) && \text{disjoint events} \\ &= P\!\left(Z > 1.96 - \tfrac{\delta}{\text{se}}\right) + P\!\left(Z < -1.96 - \tfrac{\delta}{\text{se}}\right) && Z = W - \tfrac{\delta}{\text{se}} \overset{\text{approx.}}{\sim} \mathcal{N}(0,1) \\ &\approx 1 - \Phi\!\left(1.96 - \frac{\delta}{\text{se}}\right) + \Phi\!\left(-1.96 - \frac{\delta}{\text{se}}\right) && P(Z > a) = 1 - \Phi(a),\ P(Z < b) = \Phi(b) . \end{aligned} $$
The power grows with $|\delta|/\text{se}$: with the true difference $|\delta| = |p_1 - p_2|$ and with the size of the test sets. For $m = n$, $\text{se} = \sqrt{\big(p_1(1-p_1) + p_2(1-p_2)\big)/n}$ is proportional to $1/\sqrt n$, so $|\delta|/\text{se}$ grows with $\sqrt n$ (Wasserman, 2004, Thm. 10.6).
Example. Two systems with true rates $p_1 = 0.70$ and $p_2 = 0.65$, so $\delta = 0.05$, are measured on separate test sets of $m = n = 100$ tasks:
$$\text{se}(\hat\delta) = \sqrt{\frac{p_1(1-p_1)}{m} + \frac{p_2(1-p_2)}{n}} = \sqrt{\frac{0.70\cdot0.30}{100} + \frac{0.65\cdot0.35}{100}} \approx 0.066 ,$$
$\delta/\text{se} \approx 0.76$, and the power is about $0.12$: of 100 such comparisons, about 12 detect the difference and about 88 miss it. At $m = n = 1600$, $\text{se} \approx 0.0165$, $\delta/\text{se} \approx 3.02$, and the power is about $0.86$: of 100 such comparisons, about 86 detect the difference and about 14 miss it.
Simulation. The figure compares two systems with the same true rate $p = 0.7$ on separate test sets of $n = 100$ tasks, 10 000 times: $W$ follows $\mathcal{N}(0,1)$ approximately, and about 5 % of the comparisons exceed $|W| = 1.96$.
p, n, reps = 0.7, 100, 10_000
p1_hat, p2_hat = sample_rates(p, n, reps), sample_rates(p, n, reps)
se_hat = np.sqrt(p1_hat * (1 - p1_hat) / n + p2_hat * (1 - p2_hat) / n)
W = (p1_hat - p2_hat) / se_hat
plot_wald_null(W)
Simulation. The figure gives the simulated power for a true difference of 5 percentage points, 70 % against 65 %, as a function of the size $n$ of each test set, together with the approximation derived above.
p1, p2, reps = 0.70, 0.65, 10_000
ns = np.array([50, 100, 200, 400, 800, 1600, 3200])
power = []
for n in ns:
a, b = sample_rates(p1, n, reps), sample_rates(p2, n, reps)
se_hat = np.sqrt(a * (1 - a) / n + b * (1 - b) / n)
power.append(np.mean(np.abs((a - b) / se_hat) > 1.96))
se = np.sqrt(p1 * (1 - p1) / ns + p2 * (1 - p2) / ns)
d = (p1 - p2) / se # delta / se
approx = [1 - NormalDist().cdf(1.96 - x) + NormalDist().cdf(-1.96 - x) for x in d]
plot_power(ns, power, approx, p1, p2)
Paired Comparison¶
A paired comparison evaluates two systems on the same tasks and tests the per-task differences; it removes the variation in task difficulty that both systems share (Wasserman, 2004, Ex. 10.7; Miller, 2024, §4.2, Eq. 7).
Both systems solve the same $n$ tasks. For task $i$, $X_i = 1$ if the first system solves it and $X_i = 0$ otherwise, and $Y_i = 1$ if the second system solves it and $Y_i = 0$ otherwise (Wasserman, 2004, Ex. 10.7). The per-task difference is $D_i = X_i - Y_i \in \{-1, 0, 1\}$, and the null hypothesis $\delta = 0$ concerns the difference of the two success probabilities, $\delta = \mathbb{E}[D_i] = P(X_i = 1) - P(Y_i = 1)$. The mean of the differences is the difference of the two success rates:
$$\bar D = \frac{1}{n}\sum_{i=1}^{n} (X_i - Y_i) = \frac{1}{n}\sum_{i=1}^{n} X_i - \frac{1}{n}\sum_{i=1}^{n} Y_i = \hat p_1 - \hat p_2 = \hat\delta .$$
The standard error for separate test sets does not apply. The standard error of Wald Test for Two Success Rates adds the variances of $\hat p_1$ and $\hat p_2$, because rule (2) loses its covariance term when the two rates are independent. On the same tasks, $X_i$ and $Y_i$ are dependent: a hard task tends to be failed by both, so $\operatorname{Cov}(X_i, Y_i) > 0$. This dependence enters the variance of $\hat\delta$ only through the covariance term, and the derivation below therefore keeps it.
The variance of a difference contains the covariance. Rules (1) and (2) of Variance of the Sample Proportion give
$$ \begin{aligned} \operatorname{Var}(D_i) &= \operatorname{Var}\big(X_i + (-Y_i)\big) && D_i = X_i + (-Y_i) \\ &= \operatorname{Var}(X_i) + \operatorname{Var}(-Y_i) + 2\operatorname{Cov}(X_i, -Y_i) && \text{rule (2)} \\ &= \operatorname{Var}(X_i) + \operatorname{Var}(Y_i) - 2\operatorname{Cov}(X_i, Y_i) && \text{rule (1) with } a = -1;\ \operatorname{Cov}(X_i, -Y_i) = -\operatorname{Cov}(X_i, Y_i) . \end{aligned} $$
With $\operatorname{Var}(X_i) = p_1(1-p_1)$ and $\operatorname{Var}(Y_i) = p_2(1-p_2)$ (see Bernoulli Trial), the first two terms are those of the separate-set formula; the term $-2\operatorname{Cov}(X_i, Y_i)$ is what the shared tasks change.
The paired standard error. The pairs of different tasks are independent. The tasks are drawn independently, and the difficulty of a task enters only the two outcomes $X_i, Y_i$ of that task; the $D_i$ are therefore independent and identically distributed, while $X_i$ and $Y_i$ of one task are dependent. Tasks that share a passage share a difficulty as well, and their $D_i$ are then dependent, as in Clustered Tasks. For independent tasks, the derivation of $\operatorname{Var}(\hat p)$ applies unchanged to the mean of the $D_i$:
$$ \begin{aligned} \operatorname{Var}(\bar D) &= \operatorname{Var}\Big(\frac{1}{n}\sum_{i=1}^{n} D_i\Big) && \text{definition of } \bar D \\ &= \frac{1}{n^2}\operatorname{Var}\Big(\sum_{i=1}^{n} D_i\Big) && \text{rule (1) with } a = 1/n \\ &= \frac{1}{n^2}\sum_{i=1}^{n}\operatorname{Var}(D_i) && \text{rule (2)},\ \operatorname{Cov}(D_i, D_j) = 0 \text{ for } i \neq j \\ &= \frac{\operatorname{Var}(D_1)}{n} && D_i \text{ identically distributed} \\ &= \frac{p_1(1-p_1) + p_2(1-p_2)}{n} - \frac{2\operatorname{Cov}(X_1, Y_1)}{n} && \text{the variance of } D_i . \end{aligned} $$
The first fraction is the variance of the separate-set formula for $m = n$; the shared tasks reduce it by $2\operatorname{Cov}(X_1, Y_1)/n$ (Miller, 2024, §4.2).
The paired test. $\operatorname{Var}(D_1)$ is unknown and, unlike $p(1-p)$, is not a function of $\delta$: since $D_i^2 = 1$ exactly when the two systems disagree on task $i$, $\operatorname{Var}(D_1) = \mathbb{E}[D_1^2] - \delta^2 = P(X_1 \ne Y_1) - \delta^2$. With $p_{10} = P(X_1 = 1, Y_1 = 0)$ and $p_{01} = P(X_1 = 0, Y_1 = 1)$, the difference is $\delta = p_{10} - p_{01}$, while $P(X_1 \ne Y_1) = p_{10} + p_{01}$: one value of $\delta$ admits every probability of disagreement between $|\delta|$ and $1$, and with it different variances. $\operatorname{Var}(D_1)$ is estimated by the variance of the observed differences, the plug-in estimate
$$S^2 = \frac{1}{n}\sum_{i=1}^{n} (D_i - \bar D)^2 = \frac{1}{n}\sum_{i=1}^{n} D_i^2 - \bar D^2 ,$$
the share of tasks on which the two systems disagree minus $\bar D^2$ (Wasserman, 2004, Ex. 10.7). By the weak law of large numbers, $\frac1n\sum_i D_i^2$ approaches $\mathbb{E}[D_1^2]$ and $\bar D$ approaches $\delta$ as $n$ grows (Wasserman, 2004, Thm. 5.6), so $S^2$ approaches $\operatorname{Var}(D_1)$; $S/\sqrt n$ in place of $\sqrt{\operatorname{Var}(D_1)/n}$ therefore leaves the approximate normal distribution of the standardised mean unchanged (Wasserman, 2004, Thm. 5.10). The test statistic is
$$W = \frac{\bar D}{S/\sqrt n} .$$
$S^2$ has the divisor $n$ of Wasserman (2004, Ex. 10.7). Miller (2024), Eq. 7, divides by $n-1$; the two standard errors differ by the factor $\sqrt{n/(n-1)}$, which is negligible for an evaluation of a few hundred tasks (see Standard Error).
The distribution of $W$. The $D_i$ are independent and identically distributed with finite variance, so $\bar D$ is approximately normal by the central limit theorem (see Central Limit Theorem & the 95 % Interval), with $\mathbb{E}[\bar D] = \mathbb{E}[D_i] = \delta$ by the linearity of the expectation. When $\delta = 0$, $W$ is $\bar D$ standardised, and $W \overset{\text{approx.}}{\sim} \mathcal{N}(0,1)$. The paired test rejects the null hypothesis $\delta = 0$ in favour of $\delta \ne 0$ when $|W| > 1.96$ and, like the test for separate test sets, has size about $0.05$ (see Wald Test for Two Success Rates).
With $\operatorname{Cov}(X_i,Y_i) > 0$ the differences $D_i$ vary less than the two outcomes separately, the paired standard error is smaller, and the same test set separates smaller differences. The formula for separate test sets omits the covariance term; applied to two systems on the same tasks, it overstates the standard error.
Example. Two systems are run on the same $n = 400$ tasks. The contingency table counts the tasks by the joint outcome $(X_i, Y_i)$ of both systems:
| second solves ($Y_i = 1$) | second fails ($Y_i = 0$) | total | |
|---|---|---|---|
| first solves ($X_i = 1$) | 242 | 38 | 280 |
| first fails ($X_i = 0$) | 18 | 102 | 120 |
| total | 260 | 140 | 400 |
The row total of the first row gives $\hat p_1 = 280/400 = 0.70$, the column total of the first column $\hat p_2 = 260/400 = 0.65$. The two off-diagonal cells are the discordant tasks, on which the systems disagree: $D_i = 1$ on the 38 that only the first system solves, $D_i = -1$ on the 18 that only the second solves, and $D_i = 0$ on the 344 on the diagonal. Hence $\bar D = (38 - 18)/400 = 0.05$ and $\frac1n\sum_i D_i^2 = 56/400 = 0.14$. The shortcut $\operatorname{Var}(X) = \mathbb{E}[X^2] - \mathbb{E}[X]^2$ of Bernoulli Trial, applied to the 400 values $D_i$, gives
$$S^2 = \frac{1}{n}\sum_{i=1}^{n} D_i^2 - \bar D^2 = 0.14 - 0.05^2 = 0.1375, \qquad \frac{S}{\sqrt n} = \sqrt{\frac{0.1375}{400}} \approx 0.0185, \qquad |W| = \frac{0.05}{0.0185} \approx 2.70 > 1.96 .$$
The paired test detects the difference. The formula for separate test sets, applied to the same outcomes, gives $\sqrt{0.70\cdot0.30/400 + 0.65\cdot0.35/400} \approx 0.0331$ and $|W| \approx 1.51 < 1.96$: no evidence of a difference. Only the 56 discordant tasks contribute to $S^2$; the 344 tasks on the diagonal have $D_i = 0$.
The standard error of $\hat\delta$ depends on how the two rates were measured:
| Measurement | $\widehat{\text{se}}(\hat\delta)$ |
|---|---|
| separate test sets of sizes $m$ and $n$ | $\sqrt{\hat p_1(1-\hat p_1)/m + \hat p_2(1-\hat p_2)/n}$ — the variances of independent estimates add, derived in Wald Test for Two Success Rates |
| the same $n$ tasks (paired comparison) | $S/\sqrt{n}$, where $D_i = X_i - Y_i$ is the difference of the two outcomes on task $i$, $\bar D$ their mean and $S^2 = \frac{1}{n}\sum_i (D_i - \bar D)^2$ — derived above |
Simulation. Each task receives its own difficulty, drawn independently of the other tasks and shared by both systems on that task, and each system adds its own independent variation. The true rates are 70 % and 65 %, each test set has $n = 400$ tasks, and 10 000 test sets are drawn. On the same outcomes, the simulation computes both standard errors and the share of test sets in which each test detects the difference of 5 percentage points. It also counts the test sets in which the paired test detects the difference although the two 95 % intervals overlap, so that the overlap rule would report no difference.
p1, p2, n, reps = 0.70, 0.65, 400, 10_000
shared = 0.6 # share of the variation that both systems share
t1, t2 = NormalDist().inv_cdf(p1), NormalDist().inv_cdf(p2) # thresholds for success rates p1, p2
difficulty = rng.standard_normal((reps, n)) # difficulty of task i: drawn independently for each task, shared by both systems on that task
own1, own2 = rng.standard_normal((reps, n)), rng.standard_normal((reps, n))
X = (np.sqrt(shared) * difficulty + np.sqrt(1 - shared) * own1 < t1).astype(int) # system 1 solves task i
Y = (np.sqrt(shared) * difficulty + np.sqrt(1 - shared) * own2 < t2).astype(int) # system 2 solves task i
p1_hat, p2_hat = X.mean(axis=1), Y.mean(axis=1)
D = X - Y
se_separate = np.sqrt(p1_hat * (1 - p1_hat) / n + p2_hat * (1 - p2_hat) / n)
se_paired = D.std(axis=1) / np.sqrt(n) # S / sqrt(n), S^2 = mean of (D_i - D_bar)^2
print(f"correlation of X_i and Y_i: {np.corrcoef(X.ravel(), Y.ravel())[0, 1]:.2f}")
print(f"mean se, separate-sets formula: {se_separate.mean():.4f}")
print(f"mean se, paired formula: {se_paired.mean():.4f}")
print(f"power, separate-sets formula: {np.mean(np.abs((p1_hat - p2_hat) / se_separate) > 1.96):.2f}")
print(f"power, paired formula: {np.mean(np.abs(D.mean(axis=1) / se_paired) > 1.96):.2f}")
half1 = 1.96 * np.sqrt(p1_hat * (1 - p1_hat) / n) # arms of the two 95 % intervals
half2 = 1.96 * np.sqrt(p2_hat * (1 - p2_hat) / n)
overlap = np.abs(p1_hat - p2_hat) < half1 + half2
detected = np.abs(D.mean(axis=1) / se_paired) > 1.96
print(f"detected by the paired test although the 95 % intervals overlap: "
f"{np.mean(overlap & detected):.2f} of all test sets")
plot_paired_se(se_separate, se_paired, n)
Minimum Detectable Effect¶
The minimum detectable effect of a test set is the smallest true difference $\delta$ of two success rates that a test at level $0.05$ detects with a given power; the sample-size formula inverts it and gives the number of tasks that a given $\delta$ requires (Miller, 2024, §5, Eq. 9, Eq. 10). Both are computed before an evaluation is run: whether a test set of $n$ tasks can separate the differences that matter, and how many tasks a new test set needs. The power used below is $0.8$, the value of Miller's (2024, §5) example.
Derivation. The power of the Wald test (see Wald Test for Two Success Rates) is $1 - \Phi(1.96 - \delta/\text{se}) + \Phi(-1.96 - \delta/\text{se})$ for $\delta > 0$. At a power of $0.8$, $\delta/\text{se} \approx 2.8$, and the second term is $\Phi(-4.76) < 10^{-5}$; it is dropped:
$$ \begin{aligned} 1 - \Phi\!\left(1.96 - \frac{\delta}{\text{se}}\right) &= 0.8 && \text{power } 0.8 \\ \Phi\!\left(1.96 - \frac{\delta}{\text{se}}\right) &= 0.2 \\ 1.96 - \frac{\delta}{\text{se}} &= \Phi^{-1}(0.2) \approx -0.84 && \Phi \text{ is increasing and invertible} \\ \delta &= (1.96 + 0.84)\,\text{se} = 2.80\,\text{se} . \end{aligned} $$
Both standard errors of the difference have the form $\text{se} = \sigma/\sqrt n$, with $\sigma^2$ the variance per task of the measured difference. Substituting it gives the minimum detectable effect and, solved for $n$, the sample size:
$$\delta_{\min} = 2.80\,\frac{\sigma}{\sqrt n}, \qquad n = \frac{2.80^2\,\sigma^2}{\delta^2} \approx \frac{7.84\,\sigma^2}{\delta^2} .$$
| Measurement | $\sigma^2$ | derived in |
|---|---|---|
| separate test sets of $n$ tasks each | $p_1(1-p_1) + p_2(1-p_2)$ | Wald Test for Two Success Rates |
| the same $n$ tasks | $\operatorname{Var}(D_1) = P(X_1 \ne Y_1) - \delta^2$ | Paired Comparison |
For one run per task, Eq. 9 of Miller (2024) is this formula for the paired case; it adds terms for repeated runs, which divide the variance within a task by the number $K$ of runs (see Repeated Runs). On a clustered test set, $\sigma^2$ is multiplied by the design effect, and the number of tasks by the same factor (see Clustered Tasks).
Example. The rates $p_1 = 0.70$ and $p_2 = 0.65$ of the two preceding sections differ by $\delta = 0.05$. On separate test sets, $\sigma^2 = 0.21 + 0.2275 = 0.4375$ and $n \approx 7.84 \cdot 0.4375 / 0.05^2 \approx 1370$ tasks per test set; at $n = 1600$ the power of the Wald test is therefore above $0.8$, about $0.86$. On the same tasks, with the share $P(X_1 \ne Y_1) = 0.14$ of discordant tasks of the contingency table, $\sigma^2 = 0.14 - 0.05^2 = 0.1375$ and $n \approx 7.84 \cdot 0.1375 / 0.05^2 \approx 430$ tasks; the code, with the unrounded quantiles $\Phi^{-1}(0.975)$ and $\Phi^{-1}(0.8)$, gives 1374 and 432.
With the variances of the example held fixed, a test set of $n = 500$ tasks has
| Measurement | $\delta_{\min}$ at $n = 500$ |
|---|---|
| separate test sets | $2.80\sqrt{0.4375/500} \approx 0.083$, i.e. 8.3 percentage points |
| the same tasks | $2.80\sqrt{0.1375/500} \approx 0.046$, i.e. 4.6 percentage points |
A difference of 1 or 2 percentage points between two systems on 500 tasks lies far below both values: the test detects it with a power well below $0.8$, and a p-value above $0.05$ is then no evidence that the two systems are equally good.
Evidence. Miller (2024, §5) computes the sample size of a paired comparison for $\delta = 0.03$, power $0.8$, level $0.05$ and a variance $1/9$ of the per-question difference, and obtains $n = 7.84 \cdot (1/9)/0.03^2 \approx 969$ questions; the parameters are his own assumptions, not measurements.
Consequence for an evaluation. Before two systems are compared, $\delta_{\min}$ is computed for the available number of tasks. A test set whose $\delta_{\min}$ exceeds the differences of interest is enlarged, or the comparison is reported with $\delta_{\min}$ beside the p-value.
Simulation. The code computes $\delta_{\min}$ for several $n$ and the sample size for several $\delta$, both for separate test sets and for the same tasks. It then simulates 10 000 comparisons on separate test sets of the computed size for $\delta = 0.05$ and prints the share that the Wald test detects, which is close to $0.8$.
z = NormalDist().inv_cdf(0.975) + NormalDist().inv_cdf(0.80) # 1.96 + 0.84 = 2.80
p1, p2 = 0.70, 0.65
var_separate = p1 * (1 - p1) + p2 * (1 - p2) # sigma^2 per task, separate test sets
var_paired = 0.14 - (p1 - p2) ** 2 # P(X != Y) - delta^2, the contingency table
print(f"{'n':>6} {'delta_min separate':>18} {'delta_min paired':>16}")
for n in (100, 200, 500, 1000, 2000):
print(f"{n:6d} {z * np.sqrt(var_separate / n):18.3f} {z * np.sqrt(var_paired / n):16.3f}")
print(f"\n{'delta':>6} {'n separate':>10} {'n paired':>8}")
for delta in (0.10, 0.05, 0.03, 0.01):
print(f"{delta:6.2f} {z**2 * var_separate / delta**2:10.0f} {z**2 * var_paired / delta**2:8.0f}")
n = int(np.ceil(z**2 * var_separate / (p1 - p2) ** 2))
rng_mde = np.random.default_rng(1) # own generator: the later cells keep their numbers
a, b = rng_mde.binomial(n, p1, 10_000) / n, rng_mde.binomial(n, p2, 10_000) / n
se_hat = np.sqrt(a * (1 - a) / n + b * (1 - b) / n)
print(f"\nseparate test sets of n = {n}: simulated power {np.mean(np.abs((a - b) / se_hat) > 1.96):.2f}")
Standard Error of Other Metrics¶
Every metric computed on a finite test set is an estimate and is reported with its standard error, in the same way as a success rate. A second test set drawn from the same population of tasks gives a different accuracy, a different precision and a different mean score. The formula for the standard error depends on what the metric averages:
| Metric | Examples | Estimated standard error |
|---|---|---|
| proportion over all $n$ items | success rate, accuracy, pass rate, share of tool calls that fit the schema | $\sqrt{\hat p(1-\hat p)/n}$, as in Standard Error |
| proportion over a subset of the items | recall, precision, the recall of one class | the same formula, with $n$ the size of the subset |
| mean of a score per item | progress rate, mean reward, mean grade of a judge, mean cost per task | $S/\sqrt n$, with $S^2$ the sample variance of the scores |
| function of several counts or means | $F_1$, macroaverage over classes, cost per completed task | the Bootstrap |
A mean score has the standard error of the sample proportion, with $\sigma^2$ in place of $p(1-p)$. For independent scores $X_1, \dots, X_n$ with variance $\sigma^2$, rules (1) and (2) of Variance of the Sample Proportion give $\operatorname{Var}(\bar X) = \sigma^2/n$ for the mean $\bar X = \frac1n\sum_i X_i$; only the last step of that derivation inserts $\operatorname{Var}(X_i) = p(1-p)$. For a score that is not 0/1, $\sigma^2$ is not a function of the mean and is estimated by the sample variance $S^2 = \frac{1}{n-1}\sum_i (X_i - \bar X)^2$, for which $\mathbb{E}[S^2] = \sigma^2$ (Wasserman, 2004, Thm. 3.17). The estimated standard error is $S/\sqrt n$, the standard error of the mean of Estimators and Uncertainties, and the standard error that Miller (2024), §2.1, Eq. 1, recommends for the mean score of an evaluation.
A proportion over a subset uses the size of the subset. The four counts of a binary decision are the true positives $tp$ (actual positives that the system flags), the false negatives $fn$ (actual positives that it misses), the false positives $fp$ (actual negatives that it flags) and the true negatives. Recall $R = tp/(tp+fn)$ is the success rate of the system on the $n_{\text{pos}} = tp + fn$ actual positives; the other items do not enter it. Given $n_{\text{pos}}$, these items are independent Bernoulli trials, and the formula of the sample proportion applies with $n_{\text{pos}}$ in place of $n$. Precision $P = tp/(tp+fp)$ is the success rate on the $n_{\text{pred}} = tp + fp$ items that the system flags. Both subset sizes are themselves random — $n_{\text{pred}}$ depends on the system — and treating them as fixed is an approximation.
Example. Of $n = 1000$ actions, $40$ are harmful, and a guardrail blocks $32$ of them: $\hat R = 0.8$ and $\widehat{\text{se}} = \sqrt{0.8 \cdot 0.2/40} \approx 0.063$, a 95 % interval of $0.8 \pm 0.12$. With $n = 1000$ in the formula, $\widehat{\text{se}} \approx 0.013$ would be smaller by the factor $\sqrt{1000/40} = 5$.
With repeated runs, the estimated standard error is computed over the task scores. The task score $\bar X_i$ of Repeated Runs, the share of the $K$ runs that solve task $i$, is a score per task. The $n$ task scores are independent and identically distributed, $\hat p$ is their mean, and its estimated standard error is $S/\sqrt n$, with $S^2$ the sample variance of the $n$ task scores (Miller, 2024, §3.1). The formula $\sqrt{\hat p(1-\hat p)/(nK)}$ of Standard Error, applied to all $nK$ runs, treats them as independent: it understates the standard error by the square root of the factor $1 + (K-1)\rho$ derived in Repeated Runs, by $\sqrt{1.75} \approx 1.32$ for $K = 4$ and $\rho = 1/4$. The simulation of Repeated Runs prints both estimates for $K = 1$ to $K = 32$.
Bootstrap¶
The bootstrap estimates the standard error of any metric by recomputing the metric on test sets that are resampled from the one test set at hand (Wasserman, 2004, §8.2). The standard error is the standard deviation of a metric over repeated test sets (see Standard Error). A real evaluation has only one test set, and the bootstrap draws the repeated test sets from it:
- Draw $n$ items from the test set with replacement: an item may be drawn several times or not at all.
- Compute the metric $T^*$ on the drawn items.
- Repeat steps 1 and 2 $B$ times, which gives $T^*_1, \dots, T^*_B$.
- The bootstrap standard error $\widehat{\text{se}}_{\text{boot}}$ is the standard deviation of $T^*_1, \dots, T^*_B$.
Why the resamples stand in for new test sets. Drawing one item at random from the test set is a draw from the empirical distribution, which puts the probability $1/n$ on each item of the test set. The empirical distribution takes the place of the unknown population of tasks, so the variation of $T^*$ over resamples approximates the variation of the metric $T$ over test sets (Wasserman, 2004, §8.2). The 95 % interval is $T \pm 1.96\,\widehat{\text{se}}_{\text{boot}}$ when the distribution of $T$ is close to normal (Wasserman, 2004, eq. 8.2); otherwise the percentile interval, between the 2.5 % and the 97.5 % quantile of the $T^*_b$, is used (Wasserman, 2004, §8.3).
The resampled unit is the unit that is independent. Items drawn one by one presuppose independent items; where items are dependent, the unit that is drawn is the group within which they depend:
| Test set | Resampled unit | Reason |
|---|---|---|
| several questions per passage | the passage, with all its questions | passages are independent, questions of one passage are not (see Clustered Tasks) |
| several runs of each task | the task, with all its runs | the runs of one task share its difficulty (see Repeated Runs) |
| two systems on the same tasks | the task; both metrics and their difference are computed on the same resample | the shared difficulty enters both metrics, as in Paired Comparison |
For a proportion or a mean the formulas of Standard Error of Other Metrics give the standard error without resampling; Miller (2024), §2.1, calls the bootstrap "unnecessary unless a complicated sampling scheme or estimator is being used".
Simulation. A guardrail is evaluated on test sets of $n = 1000$ actions, of which 4 % are harmful; it flags a harmful action with probability 0.8 and a harmless one with probability 0.03. The standard deviation of recall and $F_1$ over 5 000 simulated test sets is the true standard error. On one of these test sets, the code computes the binomial standard error of recall with $n_{\text{pos}}$ and with $n$, and the bootstrap standard error of recall and $F_1$ from 2 000 resamples. The figure compares the $F_1$ values over test sets with the $F_1$ values over resamples of one test set: the second histogram is centred at the $F_1$ of that test set, not at the true value, and only its spread is used.
n, prevalence, tpr, fpr = 1000, 0.04, 0.8, 0.03 # actions per test set; share of harmful ones; P(flag | harmful); P(flag | harmless)
reps, B = 5_000, 2_000 # simulated test sets; bootstrap resamples
def draw_test_sets(reps):
y = rng.random((reps, n)) < prevalence # action is harmful
flag = np.where(y, rng.random((reps, n)) < tpr, rng.random((reps, n)) < fpr) # guardrail flags it
return y, flag
def recall_f1(y, flag):
tp = (y & flag).sum(axis=-1)
fn = (y & ~flag).sum(axis=-1)
fp = (~y & flag).sum(axis=-1)
R, P = tp / (tp + fn), tp / (tp + fp)
return R, 2 * P * R / (P + R)
R_sets, F1_sets = recall_f1(*draw_test_sets(reps)) # the metrics over many test sets
y, flag = draw_test_sets(1) # the one test set of a real evaluation
y, flag = y[0], flag[0]
R_hat, F1_hat = recall_f1(y, flag)
n_pos = y.sum()
idx = rng.integers(0, n, size=(B, n)) # B resamples: n item indices drawn with replacement
R_boot, F1_boot = recall_f1(y[idx], flag[idx])
print(f"one test set: {n_pos} harmful actions of n = {n}, recall {R_hat:.3f}, F1 {F1_hat:.3f}")
print(f"{'':32s}{'recall':>8s}{'F1':>8s}")
print(f"{'sd over test sets (true se)':32s}{R_sets.std():8.4f}{F1_sets.std():8.4f}")
print(f"{'bootstrap se':32s}{R_boot.std():8.4f}{F1_boot.std():8.4f}")
print(f"{'binomial se, n_pos':32s}{np.sqrt(R_hat * (1 - R_hat) / n_pos):8.4f}")
print(f"{'binomial se, n (wrong subset)':32s}{np.sqrt(R_hat * (1 - R_hat) / n):8.4f}")
plot_bootstrap(F1_sets, F1_boot, F1_hat)
References¶
Peer-reviewed papers
Dua, Dheeru, et al. "DROP: A Reading Comprehension Benchmark Requiring Discrete Reasoning Over Paragraphs." Proceedings of the 2019 Conference of the North American Chapter of the Association for Computational Linguistics: Human Language Technologies, Volume 1 (Long and Short Papers). 2019, pp. 2368–2378. arXiv:1903.00161. — the reading-comprehension benchmark, several questions per paragraph, on which Miller (2024) measures the clustered standard error.
Preprints & other non-peer-reviewed papers
Miller, Evan. "Adding Error Bars to Evals: A Statistical Approach to Language Model Evaluations." arXiv preprint arXiv:2411.00640 (2024). — the standard error of a mean evaluation score (§2.1, Eq. 1) and of a 0/1 score (Eq. 2), the bootstrap as unnecessary for simple estimators (§2.1), the paired standard error of two models on the same questions (§4.2, Eq. 7), the variance of a mean score split into a part between and a part within questions, and its reduction by resampling (§3, §3.1), the clustered standard error and its ratio to the naive one on DROP (§2.2, Tab. 4), the minimum detectable effect and the sample-size formula (§5, Eq. 9, Eq. 10).
Textbooks & monographs
Wasserman, Larry. All of Statistics: A Concise Course in Statistical Inference. Springer, 2004. doi:10.1007/978-0-387-21736-9. — the binomial distribution (§2.3), the Beta distribution (§2.4) and its mean and variance (§3.4), expectation, variance and covariance (Ex. 3.2, Def. 3.14, Thm. 3.15, Ex. 3.16, Thm. 3.20), correlation (Def. 3.18), the sample mean and the sample variance (Thm. 3.17), the rule of iterated expectations (Thm. 3.24), the weak law of large numbers (Thm. 5.6), the central limit theorem (Thm. 5.8) and its form with the estimated standard deviation (Thm. 5.10), the standard error of a proportion and the normal-based confidence interval (§6.3, eq. 6.5, Ex. 6.8, Thm. 6.16), the Wald test and the comparison of two algorithms on separate and on the same test set (Def. 10.3, Ex. 10.7), the p-value, the p-value of the Wald test and its misreadings (§10.1, §10.2, Def. 10.11, Thm. 10.13), bootstrap variance estimation and bootstrap intervals (§8.2, eq. 8.2, §8.3).