Metrics and statistics

A metric turns thousands of graded outputs into one number. Statistics tells you how much to trust that number. Get the first wrong and you measure the wrong thing; get the second wrong and you celebrate noise. This module rebuilds the school statistics you need, with pictures first and formulas second.

Two questions behind every number

In Module 03 you built graders: each one turns a single output (or a whole agent trajectory) into a score. Now you have 500 scores. Two new questions appear:

1. What single number summarizes them? That is the metric. "Fraction of tasks solved" and "chance the agent solves a task every time a customer asks" sound similar but can differ by a factor of ten.
2. How much would that number move if I ran the eval again? That is statistics. Agents are random. If you rerun the same eval tomorrow, you get a different number. You need to know how different, or you cannot tell improvement from luck.

Success rate and accuracy

The simplest metric: each item is graded pass (1) or fail (0), and the metric is the average.

success rate = (number of passes) / (number of items)

People call it accuracy for question answering, resolve rate on SWE-bench (fraction of GitHub issues whose hidden tests pass), task success rate on OSWorld or WebArena, and pass@1 for code. All are the same arithmetic: the mean of 0/1 scores.

A few details that matter more than they look:

• What is the denominator? If the agent crashes, times out or the environment fails to start, is that a fail or excluded? Excluding crashes inflates the score. The honest default is "fail", and report the crash rate separately so infrastructure problems stay visible.
• Macro vs micro averaging. If your eval has 400 easy "billing" tasks and 50 hard "refund" tasks, the plain average (micro) is dominated by billing. The macro average first computes the rate per category, then averages the categories equally. Report both when categories are unbalanced; they answer different questions ("how often does a random ticket succeed?" vs "how good are we across the kinds of work?").
• Partial credit. Many agent tasks have several checks (the right row changed, the user was told, nothing else changed). You can score "all checks pass" (strict, usual for τ-bench style database checks) or "fraction of checks passed" (lenient). Strict success is what the user experiences; partial credit is a useful diagnostic for where it fails.

Precision, recall and F1

Some evals are really classification tests: the system says yes or no to something, and it can be wrong in two different directions. Examples: a guardrail that flags dangerous tool calls, a router that decides "escalate to a human", an LLM judge deciding "is this transcript a failure?" (Module 13 measures judges exactly this way).

Plain accuracy is misleading when the thing you are looking for is rare. Suppose 100 out of 1,000 agent transcripts contain a policy violation. A detector that never flags anything is 90% accurate and completely useless. So we split the errors apart.

• Precision (also "positive predictive value"): of the things you flagged, the fraction that were real. Low precision means alert fatigue: humans stop reading the alerts.
• Recall (also "sensitivity" or "true positive rate"): of the real things, the fraction you caught. Low recall means dangerous misses.
• F1: the harmonic mean 2PR/(P+R). The harmonic mean punishes imbalance: P = 1.0 and R = 0.1 gives F1 = 0.18, not 0.55. Use F_β if one side matters more (F2 weights recall higher).
• Specificity = TN / (TN + FP) = 860/900 = 95.6%. Together with recall it tells you how the detector behaves on each class separately, which is what you need when the class balance in production differs from your eval set.

Most detectors output a score, and you choose a threshold. Moving the threshold trades precision against recall. Plotting that trade-off gives a precision-recall curve; its area is a threshold-free summary. For a safety guardrail you usually fix a recall target first ("catch at least 95%") and then report the precision you get at that threshold.

Scores that are not pass or fail

Mean rubric score

Graders often output a number on a scale: a 1 to 5 helpfulness score from an LLM judge, or a 0 to 10 score from MT-Bench style single-answer grading (Zheng et al. 2023). The usual metric is the mean rubric score. It is still a mean, so it gets a standard error like anything else. Three cautions:

• Look at the distribution, not only the mean. A mean of 3.5 can be "everything is a 3 or 4" or "half 1s, half 5s". Those are very different products. Report the histogram or at least the share of scores at the bottom of the scale ("% rated 1 or 2").
• Scale points are not evenly spaced in meaning. Is the gap from 4 to 5 the same as from 1 to 2? Usually not. Treat means of rubric scores as a convenient summary, and prefer binary criteria ("did it cite the policy? yes/no") where you can. Binary checks are easier to grade consistently and easier to reason about statistically.
• Judges drift. If you change the judge model or prompt, the scale shifts. Never compare rubric means produced by different judge setups.

Win rate

Instead of an absolute score, you can ask a judge (human or model) "which of these two responses is better?". Pairwise judgments are often more reliable than absolute scores, because comparing is easier than calibrating. The metric is the win rate:

win rate of A over B = (wins + 0.5 × ties) / comparisons

A win rate of 50% means "indistinguishable". Pairwise judging has its own biases, notably position bias (the judge prefers whichever answer is shown first) and verbosity bias, documented for GPT-4 judges by Zheng et al. (2023). The standard mitigation is to judge each pair twice with the order swapped and count only consistent verdicts, or count inconsistent ones as ties.

Elo and Bradley-Terry: how LMArena ranks models

Win rates work for two systems. With fifty models you cannot run every pair equally often, and you want one number per model. That is what rating systems do, borrowed from chess.

Chatbot Arena (later LMArena, and simply "Arena" since January 2026; Chiang et al. 2024), shows a user two anonymous responses to their prompt and asks which is better (A, B, tie or both bad). From millions of such votes it computes a rating per model. It originally used online Elo and moved to a Bradley-Terry model fitted by maximum likelihood (in practice, logistic regression), with bootstrap confidence intervals.

The idea in one formula

Bradley-Terry assumes each model i has a hidden strength β_i, and the chance that i beats j depends only on the difference:

P(i beats j) = 1 / (1 + 10^(−(R_i − R_j) / 400))

where the rating R = 400 · log10(strength) + constant. The "400" is a chess convention: a 400-point gap means 10-to-1 odds. The only thing that matters is the gap. A model rated 1300 and one rated 1200 behave the same as 1100 vs 1000.

Elo vs Bradley-Terry. Online Elo updates ratings after every game: R_A ← R_A + K·(result − expected). It was designed for players whose skill changes over time, so recent games count more, and the final ratings depend on the order the votes arrived in. A model is a fixed artifact, so Bradley-Terry fits all votes at once and gives the same answer regardless of order, plus proper confidence intervals. Here is a tiny, dependency-free fit you can run:

import math
from collections import defaultdict

def fit_bradley_terry(battles, iters=200):
    """battles: list of (winner, loser). Returns Elo-style ratings.
    Uses the classic minorise-maximise update (Zermelo / Hunter 2004).
    Ties can be counted as half a win each way."""
    wins, games, players = defaultdict(float), defaultdict(float), set()
    for w, l in battles:
        wins[w] += 1
        games[frozenset((w, l))] += 1
        players |= {w, l}
    strength = {p: 1.0 for p in players}
    for _ in range(iters):
        new = {}
        for i in players:
            denom = sum(n / (strength[i] + strength[j])
                        for pair, n in games.items() if i in pair
                        for j in pair if j != i)
            new[i] = wins[i] / denom
        g = math.exp(sum(math.log(v) for v in new.values()) / len(new))   # fix the scale
        strength = {p: v / g for p, v in new.items()}
    # Put on the familiar Elo scale: 400 points = 10x odds, anchored at 1000.
    return {p: round(1000 + 400 * math.log10(s)) for p, s in strength.items()}

battles = ([("A", "B")] * 64 + [("B", "A")] * 36 +     # A beats B 64% of the time
           [("B", "C")] * 64 + [("C", "B")] * 36 +
           [("A", "C")] * 76 + [("C", "A")] * 24)
print(fit_bradley_terry(battles))

Output: {'A': 1100, 'B': 1000, 'C': 900}. A beat B 64% of the time, which is exactly the 100-point gap the curve above predicts.

pass@k: at least one success in k tries

Agents and code models are random: ask twice, get two different attempts. Sometimes you are allowed several attempts and only need one to work, for example when unit tests can automatically pick the passing solution, or a human reviews candidates. For that setting, Chen et al. (2021), the Codex / HumanEval paper, defined pass@k: the probability that at least one of k samples is correct.

The obvious estimator is biased

The obvious approach: run the task n times, observe c successes, estimate p̂ = c/n, and plug in 1 − (1 − p̂)^k. The paper points out that this is biased (for k > 1 it underestimates pass@k on average, because 1 − (1 − x)^k is curved and averaging over noisy p̂ drags it down). The alternative is very neat: instead of pretending you know p, simulate choosing k attempts out of the n you actually ran.

Why is this unbiased? Picking k of your n attempts uniformly at random is itself a fair way to get k independent-looking attempts, so the expected value of this fraction equals the true pass@k exactly. You need n ≥ k, and more n gives lower variance.

Numerical stability. C(200, 100) has 59 digits. Python integers cope, but floats and numpy overflow. The paper computes the ratio as a running product instead: C(n−c,k)/C(n,k) = ∏i=n−c+1..n (1 − k/i). That is what the code below does.

pass^k: success every time

pass@k answers "can it ever do this?". A customer-service agent needs the opposite question answered: "if 1,000 customers bring this same kind of request, will it get it right every time?". Yao et al. (2024), the τ-bench paper from Sierra, introduced pass^k (read "pass hat k"): the probability that all k independent trials of a task succeed. The unbiased estimator uses the same drawing trick, now counting all-pass draws:

pass^k (per task) = C(c, k) / C(n, k)

With the 10 attempts above but 7 passes instead of 3: C(7,3)/C(10,3) = 35/120 = 0.292. The task "passes 70% of the time" but passes three times in a row only 29% of the time.

At k = 1 the two metrics are identical (both equal the plain success rate). As k grows they separate in opposite directions:

from math import comb

def pass_at_k(n: int, c: int, k: int) -> float:
    """Unbiased pass@k for ONE task (Chen et al. 2021).
    n = attempts you ran, c = how many passed, k = budget you report.
    = 1 - C(n-c, k) / C(n, k): chance that k attempts drawn from your n
    (without replacement) contain at least one pass."""
    if n - c < k:          # fewer than k failures: every draw of k contains a pass
        return 1.0
    # Numerically stable product form from the paper, no giant factorials:
    # C(n-c,k)/C(n,k) = prod_{i=n-c+1}^{n} (1 - k/i)
    prob_all_fail = 1.0
    for i in range(n - c + 1, n + 1):
        prob_all_fail *= 1.0 - k / i
    return 1.0 - prob_all_fail

def pass_hat_k(n: int, c: int, k: int) -> float:
    """Unbiased pass^k for ONE task (Yao et al. 2024, tau-bench).
    = C(c, k) / C(n, k): chance that k attempts drawn from your n are ALL passes."""
    if c < k:
        return 0.0
    return comb(c, k) / comb(n, k)

# Benchmark score = average over tasks.
results = {"task_1": (10, 3), "task_2": (10, 10), "task_3": (10, 0), "task_4": (10, 7)}  # (n, c)
for k in (1, 3, 5):
    at  = sum(pass_at_k(n, c, k)  for n, c in results.values()) / len(results)
    hat = sum(pass_hat_k(n, c, k) for n, c in results.values()) / len(results)
    print(f"k={k}:  pass@k={at:.3f}   pass^k={hat:.3f}")

Output:

k=1:  pass@k=0.500   pass^k=0.500
k=3:  pass@k=0.675   pass^k=0.325
k=5:  pass@k=0.729   pass^k=0.271

Notice that task_2 (always passes) and task_3 (never passes) do not move with k. Only the flaky tasks (task_1 and task_4) spread the two metrics apart. The gap between pass@k and pass^k is therefore a direct measure of how much inconsistency your agent has.

Cost, latency, tokens and steps

An agent that solves 80% of tasks at $4 and 3 minutes each may be worse for your product than one that solves 75% at $0.20 and 20 seconds. These are not side notes; they are metrics, and they belong in the same table as accuracy.

Percentiles: p50 and p95

Latency distributions are right-skewed: most runs are quick, a few are very slow (a retry storm, a slow tool, a long reasoning chain). The mean is dragged up by the slow tail but still hides how bad the tail is. So we use percentiles: the p50 (median) is the time that half of runs beat; the p95 is the time that 95% of runs beat, so 1 in 20 users waits longer.

Percentile estimates are noisy too: the p95 of 100 runs is determined by roughly the 5 slowest runs. Use at least a few hundred runs before you trust a p95, and more for p99.

The Pareto frontier: accuracy vs cost

When you have two metrics you care about (higher accuracy, lower cost), one configuration dominates another if it is at least as good on both and better on one. The set of configurations that nobody dominates is the Pareto frontier. Every choice on the frontier is defensible; every choice off it is a mistake, because something else is cheaper and better.

Two practical rules. First, plot cost on a log scale: configurations span orders of magnitude. Second, compare at matched budgets: "our agent beats theirs" is meaningless if yours ran with 5 samples and a verifier. Research groups have argued that agent leaderboards should always show cost alongside accuracy for exactly this reason: simple baselines with retries can match elaborate agents at a fraction of the cost.

Cost-aware leaderboards: HAL

The Holistic Agent Leaderboard (HAL; Kapoor et al., Princeton; arXiv 2510.11977, accepted at ICLR 2026; hal.cs.princeton.edu) is the clearest example of treating cost as a first-class metric at leaderboard scale. Three of its design choices are worth copying in your own evals:

• One standardized harness. Every model × agent scaffold × benchmark cell runs under the same evaluation setup, so a difference between two cells is not secretly a difference in timeouts, sandboxes or graders. (Module 21 covers the harness side in depth.)
• Dollars next to accuracy. Every cell reports cost, and results are read as an accuracy versus cost frontier, exactly like the plot above. The HAL paper reports about 21,730 rollouts across 9 models and 9 benchmarks for about $40,000 in total (roughly $2 per rollout on average), and reports that higher reasoning effort reduced accuracy in the majority of runs. "More thinking" is a configuration to measure, not an upgrade to assume.
• Logs, not only scores. HAL keeps the full transcripts and inspects them with LLM help. The paper reports finding agents that searched online for the benchmark's answers instead of solving the task, and an agent that misused a credit card in a flight-booking task. Neither behavior moves the accuracy number in the right direction; both change what the number means.

Why scores wobble: sampling variance

Time to rebuild the statistics. Start with a thought experiment. Suppose your agent has a true success rate of exactly 70% on the kind of task in your eval: if you could run it on infinitely many such tasks, 70% would pass. You can't, so you run it on n tasks and get a measured rate. How far from 70% can the measured rate land?

For pass/fail items with true rate p, the number of passes follows a binomial distribution. You do not need the formula; you need the picture:

Look at the n = 50 panel and think about what it means for your work. If last week's version scored 64% and this week's scores 76% on a 50-task eval, both numbers are well inside the range that the same 70% agent produces by chance. Many "improvements" in prompt-engineering folklore are exactly this.

Two kinds of randomness

Miller (2024) makes a distinction that matters for agents. Your eval score is uncertain for two reasons:

1. Which questions you happened to include. Your 500 tasks are a sample from all the tasks you care about. A different sample of 500 would give a different score. This is present even if the agent were perfectly deterministic.
2. What the agent happened to do on each question. Run the same task twice and it may pass once and fail once. This is within-question randomness from sampling temperature, tool timing, simulated users and so on.

Running each question several times (several epochs) shrinks the second source but not the first. Adding more questions shrinks both. That is why "more questions" is usually the better investment, and why epochs still help when questions are expensive to write.

Standard error and confidence intervals

Standard deviation vs standard error

Two terms that everyone mixes up:

• The standard deviation (SD) describes how spread out the individual scores are. For 0/1 scores with rate p it is √(p(1−p)), which is about 0.46 at p = 0.7. It does not shrink with more data; individual tasks stay as variable as they are.
• The standard error (SE) describes how much the average would move between reruns of the whole eval. It is SD / √n. It shrinks as you add questions.

The reason averages are predictable is the central limit theorem (CLT): the average of many independent scores is approximately normally distributed (a bell curve) around the true value, with spread SE, no matter what the individual scores look like. The histograms in the previous figure are already nearly bell-shaped at n = 50.

The 95% confidence interval

A bell curve has about 95% of its mass within 1.96 standard errors of the center. So:

95% CI = mean ± 1.96 × SE

This is the Wald interval (also "normal approximation" or "CLT interval"). For a success rate: p̂ ± 1.96 √(p̂(1−p̂)/n).

What does "95% confidence" mean? Not "95% chance the true value is in this particular interval". It means: if you repeated the whole eval many times and built an interval each time, 95% of those intervals would contain the true rate. In practice, read it as "values inside the interval are all plausible; values far outside are not".

Worked example

Your agent passes 312 of 500 tasks.

• p̂ = 312/500 = 0.624
• SE = √(0.624 × 0.376 / 500) = √0.000469 = 0.0217, i.e. 2.17 points
• 95% CI = 62.4 ± 1.96 × 2.17 = 62.4 ± 4.25, i.e. [58.2%, 66.6%]

The Wilson interval (next) gives [58.1%, 66.5%]: nearly identical at this n.

The Wilson interval: use it for small n or extreme rates

The Wald interval misbehaves when n is small or p̂ is near 0% or 100%. At 10/10 passes it gives a width of zero ("100% ± 0"), which is absurd. At 9/10 it gives [71.4%, 108.6%], which goes above 100%. The Wilson score interval (Wilson, 1927) fixes both by solving for the values of p that are consistent with what you saw, instead of plugging in p̂:

centre = (p̂ + z²/2n) / (1 + z²/n)    half-width = z √(p̂(1−p̂)/n + z²/4n²) / (1 + z²/n)

with z = 1.96. At 9/10, Wilson gives [59.6%, 98.2%]. At 10/10, it gives [72.2%, 100%]: "we saw no failures, but with 10 tries the true rate could plausibly be as low as 72%". That is the honest answer, and a useful reality check for every "passed all 10 test cases!" claim.

How many questions do I need?

Because SE shrinks with √n, to halve your interval you need four times as many questions. Here are the actual Wilson 95% half-widths:

Keep this table in your head when you read benchmark tables. GPQA Diamond has 198 questions, the original τ-bench airline domain has 50 tasks, Terminal-Bench 2.0 has 89. A two-point gap between models on a 50 to 200 item benchmark is, on its own, not evidence of anything. (Paired comparisons, below, do better than this table suggests, which is one reason they matter so much.)

Adding error bars: Anthropic's five recommendations

In November 2024 Evan Miller at Anthropic published "Adding Error Bars to Evals: A Statistical Approach to Language Model Evaluations" (arXiv 2411.00640, with a companion blog post). It argues that the field reported eval results with almost no statistical hygiene, and gives five concrete recommendations. They are the backbone of how to report results in this course.

Clustered standard errors

The CLT assumes independent scores. Real evals often break that. A reading-comprehension benchmark asks 5 questions about each passage; an agent benchmark has 20 tasks built on the same website, or 30 SWE-bench tasks from one repository. If the agent misreads a passage, it tends to get all its questions wrong together. Those 5 questions carry less information than 5 unrelated ones.

The fix is simple to compute. Instead of treating each question's residual (score minus mean) independently, sum the residuals within each cluster first, then square:

SE_clustered = √( Σclusters c ( Σi in c (si − s̄) )² ) / n

If clusters are uncorrelated inside, this reduces to the ordinary SE. If they move together, it grows. A handy rule of thumb from survey statistics is the design effect 1 + (m − 1)ρ for clusters of size m with within-cluster correlation ρ: with 5 questions per passage and ρ = 0.5, the variance is 3× larger and the SE √3 ≈ 1.73× larger. In the simulation in the toolkit below (100 easy-or-hard passages of 5 questions), the naive SE is 0.022 and the clustered SE 0.043, nearly double.

For agents, cluster by whatever the tasks share: repository (SWE-bench), website (WebArena), app (AndroidWorld), task template, or source document. When you run several epochs, the epochs of one question are also a cluster; averaging per question first (recommendation 3) handles that automatically.

Reduce variance: resampling and epochs

If each question's result is random, run it K times and use the per-question average as that question's score. The between-question part of the variance stays, but the within-question part is divided by K. Tools support this directly: Inspect calls the repeats epochs and lets you choose how to combine them (mean, majority vote, or pass@k style reducers). METR runs each model about 8 times per task for its time-horizon measurements.

There is a floor, though. For 0/1 scores, even infinite epochs cannot shrink the SE below the between-question spread. With 500 questions, the 95% half-width near 62% is ±4.25 points with one epoch, and no less than about ±2.45 with three epochs (reached only if every question had the same success probability, so all variance is within-question). Beyond that, you need more questions.

The bootstrap: confidence intervals for anything

Formulas exist for means and proportions. What about the median latency, the F1 score, pass^k averaged over tasks, or a Bradley-Terry rating? The bootstrap (Efron, 1979) gives you a confidence interval for almost any statistic, with no formula at all.

The trick: you cannot rerun the world, but you can pretend your sample is the world. Draw n questions from your n questions with replacement (some appear twice, some not at all), recompute the metric, and repeat 10,000 times. The spread of those recomputed metrics approximates the spread you would see across real reruns.

Three rules that make bootstraps correct:

• Resample the independent unit. If you ran 3 epochs of each question, resample questions (carrying all their epochs along), not individual runs. If questions are clustered by passage or repo, resample clusters. This is the bootstrap version of clustered SEs.
• For comparisons, resample pairs. Draw the same question indices for both models, then compute the difference. This is the paired bootstrap.
• Use enough resamples (1,000 minimum; 10,000 is cheap) and fix the random seed so the report is reproducible.

LMArena uses exactly this idea for its rating intervals: resample the battles, refit Bradley-Terry, and read off the percentiles of each model's rating.

Comparing two models: the paired test

This is the most useful technique in the module. You ran Model A and Model B on the same 500 questions. A scored 68.0%, B scored 64.4%. Is A really better?

The wrong way: two separate error bars

Treat the two scores as independent. Each has SE ≈ 2.1 points, the difference 3.6 points has SE √(SE_A² + SE_B²) = 2.99 points, so the 95% CI for the difference is 3.6 ± 5.9, i.e. [−2.3, +9.5]. It includes zero. Conclusion: "no significant difference". That conclusion is wrong, because it throws away information.

The right way: look at each question

Most of the variation in scores comes from question difficulty: hard questions are hard for both models. Both models share that noise, so it cancels when you look at differences question by question. Only the questions where the models disagree carry information about which is better.

Step by step:

1. For every question i, compute d_i = score_A(i) − score_B(i). With 0/1 scores, d is +1, −1 or 0.
2. The mean of d is the difference in accuracy: 18/500 = 3.6 points.
3. The SE of that mean is SD(d)/√n. Here SD(d) = 0.351, so SE = 1.57 points.
4. 95% CI = 3.6 ± 3.07 = [+0.5, +6.7]. It excludes zero, so the difference is statistically significant at the 5% level (two-sided p ≈ 0.022). An exact McNemar test on the 40 vs 22 discordant pairs gives p ≈ 0.030, the same conclusion.
5. Report it honestly: "A beats B by 3.6 points (95% CI 0.5 to 6.7, paired, n = 500)". The interval says the true gap could be tiny, so do not oversell it.

Power analysis: how big an eval do I need?

Statistical power is the probability that your eval detects a real difference of a given size. The conventional target is 80% power at a 5% significance level. Power analysis runs the logic backwards: pick the smallest improvement you care about (the minimum detectable effect), and compute the number of questions required.

n ≈ ( (zα/2 + zpower) × SD(d) / δ )²   with zα/2 = 1.96, zpower = 0.84

Worked example: you want to detect a 3-point improvement from 65% to 68%.

• Unpaired (different questions for each model): SD(d) = √(0.65×0.35 + 0.68×0.32) = 0.667, so n ≈ (2.80 × 0.667 / 0.03)² ≈ 3,882 questions per model.
• Paired, assuming the models disagree on 12% of questions: SD(d) ≈ √(0.12 − 0.03²) = 0.345, so n ≈ 1,039 questions.

Two lessons. First, detecting small improvements needs far more questions than most teams have; a 100-item eval can only reliably detect differences of roughly 19 to 20 points or more unpaired (80% power). Second, the disagreement rate between models drives paired power, so you can estimate it from a pilot run and plan from there. If the math says you need 4,000 questions and you have 200, the honest move is to say "this eval can only detect large changes" rather than to read tea leaves in small ones.

Multiple comparisons and p-hacking

A 95% confidence level means that when there is truly no difference, you will still see a "significant" result 5% of the time. That is fine for one test. It is not fine when you run twenty.

Suppose you try 20 prompt variants, none actually better than the baseline. The chance that at least one looks significantly better at p < 0.05 is 1 − 0.95^20 = 64%. Pick the best of 20 and you have probably selected noise. This is the winner's curse, and it is exactly the leaderboard problem "The Leaderboard Illusion" raises about private testing of many variants.

Common forms in eval work, and what to do about each:

Where variance comes from in agents

Agent evals are noisier than single-turn evals because randomness enters at every step and compounds. Knowing the sources tells you which ones to control and which to average over.

• Model sampling. At temperature above zero, each token is sampled. Even at temperature 0, many hosted APIs are not bit-for-bit deterministic because of batching and floating-point effects, so do not assume two temperature-0 runs match. Evaluate with the sampling settings you will deploy with (temperature, or the effort level on newer models that fix temperature), and use epochs.
• Environment nondeterminism. Live websites change (a reason WebVoyager-style evals drift), clocks and dates change, random IDs differ. Prefer containerised, snapshotted environments; Module 05 covers this.
• Tool and infrastructure flakiness. A network timeout is not the agent's fault. Track infra errors separately, retry them (the run, not the agent's decision), and never let them silently count as agent failures or be silently dropped.
• Simulated users. τ-bench style evals use an LLM to play the customer. The research notes for this course record that its user simulator is itself noisy (sometimes hallucinating or ending early). That noise lands in your score.
• Judge noise. An LLM judge run twice on the same transcript can disagree with itself. Measure judge self-consistency and agreement with humans (Module 13), and consider averaging several judge calls.
• Compounding. In a 30-step trajectory, a small random difference at step 3 can send the agent down a completely different path. That is why agent success rates vary more across epochs than single-answer accuracy does, and why pass^k matters.

Reliability as its own science

pass^k told you whether an agent succeeds every time. That is one face of reliability, but only one. Engineers who build aircraft, bridges and power grids never describe reliability with a single number: they separately ask whether a system behaves the same way twice, whether it survives disturbances, whether its failures can be anticipated, and how bad its failures are. Agent evaluation is starting to borrow that habit.

The paper "Towards a Science of AI Agent Reliability" (arXiv 2602.16666, February 2026; presented at ICML 2026) comes from the Princeton group behind HAL, including Sayash Kapoor and Arvind Narayanan. It argues that a single accuracy score hides how agents fail, and proposes measuring reliability along four separate dimensions: consistency, robustness, predictability and safety. Its headline finding, as the authors and press coverage summarize it, is that agents have become much more capable while their reliability has improved far more slowly (Fortune's headline: "AI agents are getting more capable, but reliability is lagging"). HAL now hosts a reliability dashboard alongside its accuracy leaderboard.

What each dimension asks, and how to measure it yourself

The paper defines its own specific metrics under each dimension. You do not need its exact definitions to start; the table gives the question each dimension asks and a practical way to measure it in your own eval.

Here is a minimal "reliability report" over run records. It reuses the pass^k estimator from earlier, so you can see that consistency is literally the pass^k you already know, while the other three need extra information that a plain pass/fail log does not contain.

from math import comb

# One record per (task, run). Illustrative data for 4 tasks x 3 runs.
#   ok        : did the run succeed on the ORIGINAL task wording?
#   ok_perturb: did it succeed on a meaning-preserving paraphrase of the task?
#   conf      : the agent's own stated confidence (0..1) before finishing
#   severity  : 0 = no harm, 1 = minor (wrong but reversible), 3 = serious (destructive)
runs = [
    # task, ok, ok_perturb, conf, severity
    ("t1", 1, 1, 0.90, 0), ("t1", 1, 1, 0.85, 0), ("t1", 1, 0, 0.90, 0),
    ("t2", 1, 0, 0.80, 0), ("t2", 0, 0, 0.75, 1), ("t2", 1, 1, 0.80, 0),
    ("t3", 0, 0, 0.30, 0), ("t3", 0, 0, 0.40, 0), ("t3", 1, 0, 0.50, 0),
    ("t4", 0, 0, 0.90, 3), ("t4", 1, 1, 0.90, 0), ("t4", 0, 0, 0.85, 1),
]

def by_task(rows):
    out = {}
    for t, ok, *_ in rows:
        out.setdefault(t, []).append(ok)
    return out

def consistency(rows, k):
    """Mean unbiased pass^k over tasks: chance that k runs of a task ALL pass."""
    vals = [comb(sum(o), k) / comb(len(o), k) for o in by_task(rows).values()]
    return sum(vals) / len(vals)

def robustness_drop(rows):
    """Paired drop in success when the task is paraphrased (same tasks, same runs)."""
    return sum(r[2] - r[1] for r in rows) / len(rows)

def calibration_error(rows, bins=(0.0, 0.5, 0.8, 1.01)):
    """Expected calibration error: |stated confidence - actual success|, bin by bin."""
    ece = 0.0
    for lo, hi in zip(bins, bins[1:]):
        b = [r for r in rows if lo <= r[3] < hi]
        if b:
            gap = abs(sum(r[3] for r in b) / len(b) - sum(r[1] for r in b) / len(b))
            ece += len(b) / len(rows) * gap
    return ece

def harm_rate(rows):
    """Share of runs with any harmful failure, and severity-weighted harm per run."""
    return (sum(r[4] > 0 for r in rows) / len(rows), sum(r[4] for r in rows) / len(rows))

acc = sum(r[1] for r in runs) / len(runs)
print(f"accuracy (pass@1)       {acc:.0%}")
print(f"consistency  pass^3     {consistency(runs, 3):.0%}")
print(f"robustness   drop       {100 * robustness_drop(runs):+.0f} points")
print(f"predictability  ECE     {calibration_error(runs):.2f}")
share, weighted = harm_rate(runs)
print(f"safety  harmful runs    {share:.0%}  (severity-weighted {weighted:.2f} per run)")

Output:

accuracy (pass@1)       58%
consistency  pass^3     25%
robustness   drop       -25 points
predictability  ECE     0.15
safety  harmful runs    25%  (severity-weighted 0.42 per run)

A 58% agent whose tasks pass all three runs only 25% of the time, that loses 25 points when the task is reworded, and that is most confident on the task where it does the most damage (t4) is a very different product from what "58%" suggests. With 12 runs these numbers have enormous error bars; the point is the shape of the report, not the values.

The cost of evaluating: efficient benchmarking

Everything in this module pushes the eval bill up. Power analysis asks for more tasks. pass^k and consistency need several epochs. Robustness needs perturbed copies of every task. Comparisons need several configurations. And for agents, each "item" is not one model call but a whole multi-step rollout in a sandbox. The EvalEval Coalition's April 2026 essay "AI evals are becoming the new compute bottleneck" makes the case that evaluation compute is becoming a real constraint in its own right; HAL's roughly $40,000 for about 21,730 rollouts gives a sense of scale.

Option 1: subsample, and know what the ± means

The simplest saving is to run a random subset of n tasks out of the N in the benchmark. Your estimate of the full-benchmark score then has a standard error with a finite population correction (because you are sampling from a fixed, finite set of tasks):

SE = √( p(1−p)/n × (N−n)/(N−1) )

Worked numbers for a benchmark of N = 500 tasks and a score near 60%:

Three rules make subsampling work. Stratify: sample within categories or difficulty bands so each is represented, then weight back to the full mix. Pair: when comparing two versions, run both on the same subset and use the paired test from earlier, which is far tighter than two independent estimates. Keep the ladder in mind: a zero SE at n = 500 only means you have the full-benchmark score exactly; that score still has its own sampling uncertainty as an estimate of the agent's true skill, which is everything the first half of this module was about.

Option 2: anchor sets and Item Response Theory

Item Response Theory (IRT) comes from standardized testing, where a test must rank students fairly even if different students see different questions. It gives every task a difficulty b and a discrimination a (how sharply it separates weaker from stronger systems), and every system an ability θ on the same scale. The two-parameter logistic (2PL) model says:

P(system with ability θ solves task i) = 1 / (1 + e−ai(θ − bi))

When θ equals a task's difficulty, the chance is 50%. You fit a and b once, from many past models' full runs. Then a new system only needs to run a small anchor set of well-chosen tasks: you find the θ that best explains its anchor results, and predict its score on every task, including the ones it never ran. tinyBenchmarks (Polo et al., ICML 2024) applied this idea to popular LLM benchmarks and reports that about 100 curated examples per benchmark were enough to estimate full-benchmark accuracy within about 2% on average. Epoch AI's Capabilities Index uses IRT for a related job: putting results from many different benchmarks onto one scale (epoch.ai/eci).

Which tasks make good anchors? A task that every system solves, or that none solves, tells you nothing about the difference between them. Information is highest where a task's solve chance is near 50% for the systems you care about (the Fisher information of a 2PL item is a²·P·(1−P), which peaks at P = 0.5). The 2026 paper "Efficient Benchmarking of AI Agents" applies this intuition to agent benchmarks; its code release describes the approach as mid-range task selection for cost-efficient agent evaluation.

Worked example: predicting a full score from 10 anchors

An illustrative benchmark has 500 tasks in five difficulty tiers of 100 (b = −2, −1, 0, +1, +2), all with discrimination a = 1.5, fitted from past models. The anchor set is 10 tasks chosen mostly from the middle tiers, so its mix does not match the benchmark. A new agent passes 7 of the 10 anchors: 70% raw. Is its full-benchmark score 70%?

import math

def p_correct(theta, a, b):
    """2PL item response curve: chance an agent of ability theta solves a task
    with discrimination a and difficulty b."""
    return 1 / (1 + math.exp(-a * (theta - b)))

# Calibration (done once, from many past models' FULL runs):
# the 500-task benchmark has 5 difficulty tiers of 100 tasks, all with a = 1.5.
TIERS = [-2.0, -1.0, 0.0, 1.0, 2.0]
A = 1.5

# Anchor set: 10 tasks, chosen mostly near the middle (where tasks are most
# informative), so the mix does NOT match the benchmark. (difficulty b, result)
anchors = [(-2.0, 1),
           (-1.0, 1), (-1.0, 1),
           (0.0, 1), (0.0, 1), (0.0, 1), (0.0, 0),
           (1.0, 1), (1.0, 0),
           (2.0, 0)]

def log_lik(theta):
    return sum(math.log(p_correct(theta, A, b)) if y else
               math.log(1 - p_correct(theta, A, b)) for b, y in anchors)

# Maximum-likelihood ability: a simple grid search is enough for one number.
grid = [i / 1000 for i in range(-4000, 4001)]
theta_hat = max(grid, key=log_lik)

# Predicted score on the FULL benchmark = average solve chance over all 500 tasks.
def full_score(theta):
    return sum(p_correct(theta, A, b) for b in TIERS) / len(TIERS)

# Uncertainty: Fisher information of the 10 anchors -> SE of theta -> score range.
info = sum(A**2 * p_correct(theta_hat, A, b) * (1 - p_correct(theta_hat, A, b)) for b, _ in anchors)
se_theta = 1 / math.sqrt(info)
lo, hi = full_score(theta_hat - 1.96 * se_theta), full_score(theta_hat + 1.96 * se_theta)

print(f"raw anchor pass rate: {sum(y for _, y in anchors) / len(anchors):.0%}")
print(f"theta_hat = {theta_hat:.2f}  (SE {se_theta:.2f})")
print(f"predicted full-benchmark score: {full_score(theta_hat):.1%}  (95% range {lo:.1%} to {hi:.1%})")
for b in TIERS:
    print(f"  tier b={b:+.0f}: predicted solve rate {p_correct(theta_hat, A, b):.0%}")

Output:

raw anchor pass rate: 70%
theta_hat = 0.83  (SE 0.56)
predicted full-benchmark score: 65.8%  (95% range 45.1% to 84.1%)
  tier b=-2: predicted solve rate 99%
  tier b=-1: predicted solve rate 94%
  tier b=+0: predicted solve rate 78%
  tier b=+1: predicted solve rate 44%
  tier b=+2: predicted solve rate 15%

Two lessons. First, the IRT estimate corrects for an unrepresentative subset: the raw 70% overstates the full score because the anchors were concentrated in the middle, while the real benchmark has twice their share of hardest-tier tasks, where this agent mostly fails. Second, IRT does not create information from nothing. With 10 anchors the 95% range is 45% to 84%, far too wide to rank close systems. With 100 anchors in the same pattern, the same code gives a range of about 59% to 72%. Real anchor sets do better than this toy because they choose tasks with high discrimination, which is how tinyBenchmarks reaches small errors with about 100 examples.

How to report a result

Put everything together. A result that an expert reviewer will accept looks like this:

(The numbers are an illustrative example.) Why each part is there:

• ± and "95% CI": say what the ± means. "±3.4" could be one SE or a 95% half-width; they differ by about 2×.
• n and epochs: lets the reader sanity-check the interval. A reader who knows the table above can see that ±3.4 at n = 500 is plausible. It could not be ±2.1: with 500 pass/fail tasks the 95% half-width near 62% is ±4.25 with one epoch and cannot fall below about ±2.45 even with three.
• Clustering: tells the reader you accounted for correlated tasks and repeated epochs.
• Paired Δ with its own CI: the number decision-makers actually need ("is the new version better?"). Note that the Δ interval (±2.1) is narrower than the single-score interval (±3.4); that is pairing at work.
• pass^k, cost, latency, steps: the metrics that decide whether you can actually ship.
• Reliability profile (increasingly expected for agents): consistency across runs, the paired drop under meaning-preserving perturbations, calibration of the agent's confidence, and the harmful-failure rate (see "Reliability as its own science").
• How the number was obtained: full run or subsample / anchor estimate. An estimated score must say so and carry its own ±.
• Versions: evals are only comparable when system, grader and data are pinned (Module 02).

Your statistics toolkit (runnable code)

Everything in this module in one pure-Python file (standard library only, Python 3.8+). It reproduces every number in the worked examples above.

import math, random
from statistics import NormalDist

Z95 = NormalDist().inv_cdf(0.975)            # 1.95996...

def mean_se(scores):
    """Mean and CLT standard error of per-question scores (0/1 or fractional)."""
    n = len(scores)
    m = sum(scores) / n
    var = sum((s - m) ** 2 for s in scores) / (n - 1)   # sample variance
    return m, math.sqrt(var / n)

def wilson_ci(successes, n, z=Z95):
    """Wilson score interval for a binomial proportion. Safe at 0/n and n/n."""
    if n == 0:
        return (0.0, 1.0)
    p = successes / n
    denom = 1 + z * z / n
    centre = (p + z * z / (2 * n)) / denom
    half = z * math.sqrt(p * (1 - p) / n + z * z / (4 * n * n)) / denom
    return (max(0.0, centre - half), min(1.0, centre + half))

def bootstrap_ci(scores, n_boot=10_000, alpha=0.05, seed=0):
    """Percentile bootstrap CI for the mean: resample questions with replacement."""
    rng = random.Random(seed)
    n = len(scores)
    means = sorted(sum(rng.choice(scores) for _ in range(n)) / n for _ in range(n_boot))
    lo = means[int((alpha / 2) * n_boot)]
    hi = means[int((1 - alpha / 2) * n_boot) - 1]
    return lo, hi

def paired_diff(a, b, z=Z95):
    """Model A vs B on the SAME questions. a[i], b[i] = score on question i."""
    d = [x - y for x, y in zip(a, b)]
    m, se = mean_se(d)
    p_value = 2 * (1 - NormalDist().cdf(abs(m / se))) if se > 0 else float("nan")
    return m, se, (m - z * se, m + z * se), p_value

def clustered_se(scores, clusters):
    """CLT standard error when questions come in groups (Miller 2024, eq. for clustered SE).
    clusters[i] = id of the group question i belongs to (e.g. the passage)."""
    n = len(scores)
    m = sum(scores) / n
    by_c = {}
    for s, c in zip(scores, clusters):
        by_c[c] = by_c.get(c, 0.0) + (s - m)     # sum of residuals within each cluster
    return math.sqrt(sum(r * r for r in by_c.values())) / n

def n_for_power(sd_diff, delta, alpha=0.05, power=0.80):
    """Questions needed to detect a true mean difference `delta` (two-sided test)."""
    z_a = NormalDist().inv_cdf(1 - alpha / 2)
    z_b = NormalDist().inv_cdf(power)
    return math.ceil(((z_a + z_b) * sd_diff / delta) ** 2)

if __name__ == "__main__":
    print("Wilson 9/10:", [round(x, 3) for x in wilson_ci(9, 10)])
    print("Wilson 312/500:", [round(x, 3) for x in wilson_ci(312, 500)])

    # The worked example below: 500 questions, A right on 340, B right on 322.
    a = [1] * 300 + [1] * 40 + [0] * 22 + [0] * 138     # both, A only, B only, neither
    b = [1] * 300 + [0] * 40 + [1] * 22 + [0] * 138
    m, se, ci, p = paired_diff(a, b)
    print(f"paired: diff={m:.3f} se={se:.4f} CI=({ci[0]:.3f}, {ci[1]:.3f}) p={p:.3f}")
    print("bootstrap CI of A:", bootstrap_ci(a, n_boot=2000))

    # Clustering: 100 passages x 5 questions, questions in a passage tend to agree.
    rng = random.Random(1)
    scores, clus = [], []
    for c in range(100):
        p_passage = rng.choice([0.1, 0.9])            # easy or hard passage
        for _ in range(5):
            scores.append(1 if rng.random() < p_passage else 0); clus.append(c)
    print("naive SE:", round(mean_se(scores)[1], 4), " clustered SE:", round(clustered_se(scores, clus), 4))

    # Power: unpaired 65% vs 68% needs sd_diff = sqrt(p1q1 + p2q2); paired with 12% discordance.
    print("unpaired n:", n_for_power(math.sqrt(.65*.35 + .68*.32), 0.03))
    print("paired n:  ", n_for_power(math.sqrt(0.12 - 0.03**2), 0.03))

Output:

Wilson 9/10: [0.596, 0.982]
Wilson 312/500: [0.581, 0.665]
paired: diff=0.036 se=0.0157 CI=(0.005, 0.067) p=0.022
bootstrap CI of A: (0.64, 0.72)
naive SE: 0.0223  clustered SE: 0.0427
unpaired n: 3882
paired n:   1039

The numpy version (for real eval logs)

With a matrix of [questions × epochs] scores, numpy makes the bootstrap vectorised and fast:

import numpy as np

rng = np.random.default_rng(0)

# scores[q, e] = score of question q on epoch e (here: 500 questions x 3 epochs, 0/1)
scores = (rng.random((500, 3)) < rng.beta(1.3, 0.8, size=(500, 1))).astype(float)

per_q = scores.mean(axis=1)                        # average epochs FIRST (one number per question)
mean = per_q.mean()
se = per_q.std(ddof=1) / np.sqrt(len(per_q))       # CLT standard error over questions
print(f"mean={mean:.3f}  95% CI=({mean-1.96*se:.3f}, {mean+1.96*se:.3f})")

# Vectorised bootstrap: 10,000 resamples of question indices at once.
idx = rng.integers(0, len(per_q), size=(10_000, len(per_q)))
boot_means = per_q[idx].mean(axis=1)
lo, hi = np.percentile(boot_means, [2.5, 97.5])
print(f"bootstrap 95% CI=({lo:.3f}, {hi:.3f})")

# Paired bootstrap for a difference: resample the SAME indices for both models.
other = np.clip(per_q + rng.normal(0.03, 0.2, size=per_q.shape), 0, 1)   # a fake model B
diff = other - per_q
boot_diff = diff[idx].mean(axis=1)
print("paired diff CI:", np.percentile(boot_diff, [2.5, 97.5]).round(3))

If you have scipy, scipy.stats.binomtest(k, n).proportion_ci(method="wilson") gives the Wilson interval and scipy.stats.bootstrap does percentile or BCa bootstraps. For regression-style analyses with clustered errors, statsmodels supports cluster-robust covariance. The pure-Python versions above are there so you understand what those calls do.

Check yourself

Sources and further reading

• Chen et al. (OpenAI), 2021, Evaluating Large Language Models Trained on Code (HumanEval, unbiased pass@k estimator).
• Yao, Shinn, Razavi, Narasimhan (Sierra), 2024, τ-bench: A Benchmark for Tool-Agent-User Interaction in Real-World Domains (pass^k).
• Evan Miller (Anthropic), 2024, Adding Error Bars to Evals: A Statistical Approach to Language Model Evaluations, and the blog post A statistical approach to model evaluations.
• Chiang, Zheng et al., 2024, Chatbot Arena: An Open Platform for Evaluating LLMs by Human Preference (Elo, Bradley-Terry, bootstrap CIs).
• Zheng et al., 2023, Judging LLM-as-a-Judge with MT-Bench and Chatbot Arena (position and verbosity bias).
• Singh et al., 2025, The Leaderboard Illusion.
• Kwa, West et al. (METR), 2025, Measuring AI Ability to Complete Long Tasks (repeated runs per task).
• UK AI Security Institute, Inspect documentation (epochs and score reducers).
• Wikipedia, Binomial proportion confidence interval (Wald, Wilson and other intervals).
• Wikipedia, Bootstrapping (statistics), and Bradley-Terry model.
• Wikipedia, Multiple comparisons problem.
• Princeton researchers including Sayash Kapoor and Arvind Narayanan, 2026, Towards a Science of AI Agent Reliability (consistency, robustness, predictability, safety); authors' post New paper: Towards a science of AI agent reliability; Fortune, AI agents are getting more capable, but reliability is lagging (Mar 2026).
• Kapoor et al., Holistic Agent Leaderboard: The Missing Infrastructure for AI Agent Evaluation (ICLR 2026); HAL leaderboard and HAL reliability dashboard.
• EvalEval Coalition, AI evals are becoming the new compute bottleneck (Apr 29, 2026; also on the Hugging Face blog).
• Polo et al., 2024, tinyBenchmarks: evaluating LLMs with fewer examples (ICML 2024; code).
• Efficient Benchmarking of AI Agents (2026; code: mid-range task selection).
• Epoch AI, Epoch Capabilities Index (IRT across many benchmarks).
• Wikipedia, Item response theory.