cs.CLSep 20, 2026
SaveHow to Run Statistics over LLM Judges and Trust the Results: Calibrated Inference for Small-Sample AI Evaluation with evalstats
Organizations: Department of Computer Science and Operations Research (DIRO) Université de Montréal Montréal, Quebec, Canada
Abstract
Researchers across academia increasingly base significance claims on LLM judge scores and small-sample AI evaluations. Yet without well-calibrated confidence intervals (CIs), hypothesis tests, and judge-bias corrections, such claims are unreliable. We address these issues in several contributions. First, we find that running statistics over raw LLM judge scores leads to inflated false positives: counterintuitively, for many inter-rater agreement metrics, false positive risk peaks at "almost perfect" human-LLM agreement. To help researchers understand how to run statistics over LLM judges responsibly, we present guidance and tooling for the statistical analysis of mixed human-AI judge designs, and implement nine hypothesis tests via prediction-powered inference (PPI), including the first known PPI corrections for four rank-based tests (Wilcoxon signed-rank, Mann-Whitney U, and omnibus variants). To keep PPI++ stable with small human-labeled calibration sets, we introduce bootstrap-adaptive power tuning, which shrinks the estimated weight toward a target estimated from the labeled data, and accounts for that weight's own sampling variance. Second, through Monte Carlo simulations, we derive recommendations for what CI, p-value, and FWER correction methods to use for small-sample AI evaluations (N<100), and warn researchers against bootstrap CIs. We package these recommendations into evalstats, an open-source Python package that selects calibrated methods automatically, and demonstrate it in three scenarios, including one where a real LLM judge validated at "substantial agreement" would have led a researcher to publish a spurious finding. evalstats is publicly available at https://github.com/ianarawjo/evalstats.
Figures & tables
Figure 1. Running evalstats from the command line on a spreadsheet of AI evaluation results prints 95% confidence intervals, -values, gradient plots ( Correll and Gleicher, 2014 ) , pairwise tests, and an executive summary, with conditions grouped into rank bands ( Demšar, 2006 ) . Each method was chosen automatically for the sample size and data type, backed by our simulations (§ 5 - 6 ). Shown are results from BiGGen Bench ( Kim et al., 2025 ) , a benchmark that uses LLM judges as evaluators, sampling three models at items each with human labels as calibration. Every estimate on screen—means, CIs, tests, -values, and effect size—is corrected for LLM judge bias using prediction-powered inference ( Angelopoulos et al., 2023a ) . The PPI-Friedman and PPI-Wilcoxon tests shown are contributions of this paper.
| What users need to know | Common practice without evalstats | With evalstats |
| Is the best-scoring condition really best? | Read overlapping error bars by eye, which does not answer the question; or run pairwise tests and correct them via a familiar correction like Bonferroni when comparisons share items and are therefore strongly dependent. | compare(factors=...) : pairwise CIs already widened for simultaneous coverage, plus rank bands marking statistical ties. |
| Does the result hold if I run it again? | Averaging or pass@k ( Hariri et al., 2026 ) across runs, which hides per-item disagreement that caused the instability. Even if trying to accommodate runs, struggle to pick an estimator that respects the nesting. | Multi-run data is detected; per-input variance is reported alongside the means, or in isolation via stability() . |
| Which condition wins once cost or latency counts? | Rank on the primary metric and inspect the secondary one informally, eyeballing means. | Simply passing secondary_metric= on compare() to get an uncertainty-aware Pareto front labeling each condition frontier , dominated , or ambiguous . |
| Which combination of two factors is best? | Test every factor level against every other—e.g., grid of model-prompt pairs, for 780 comparisons—then choose a basic FWER correction that is overly conservative under strong correlation between conditions. | The same call with two factors: Šidák simultaneous CIs and step-down FWER correction (Shaffer or Romano–Wolf) over the full grid. |
| My metric is an LLM judge. Can I trust the comparison? | Compute pooled inter-rater reliability once against a hand-labeled sample collected by convenience. | label from the CLI to randomly sample items for labeling, then run judge_alignment() . |
| Only part of my data is human-labeled, and I need to run statistics over LLM judge scores. | Analyze the labeled subset alone and discard the judge scores, losing power—or keep the judge scores, inherit their bias, and argue with reviewers or managers that high IRR with human raters warrants running familiar statistics over LLM judge scores, potentially leading to spurious findings. | alignment= on compare() propagates the labeled sample into every estimate, CI and -value via PPI correction. PPI tests recover power from the unlabeled items while remaining calibrated, worth the human labels at (§ 6.3 ). |
Table 1. How evalstats ’ features can help developers and researchers running statistics over AI evaluation results.
Figure 3. Evidence for our recommendations in Figure 2 , showing coverage rates, interval score ( Gneiting and Raftery, 2007 ) , and worst-case coverage for CI methods; and error, power, and worst-case error for -values and FWER corrections. Blue marks the recommended method, with orange its strongest runner-up; gray is every other method we tested there (9-14 per panel). Shading marks differences indistinguishable from Monte Carlo error. Per-method values in Tables 4 - 6 and 10 , and full reasoning in § A.2 . Replicates per panel: 2,000 (single-sample CIs), 300 (paired CIs, pairwise -values, simultaneous CIs), 500 (multi-arm).
Figure 4. Type I error remains calibrated across all nine PPI-corrected tests under an aggressive factorial sweep compounding judge bias magnitude and judge noise with sample size , labeled-item count , and data type (2,046 null scenarios/test, 200 MC replicates each, 3.68M trials). PPI calibration holds even when multiple judge biases compound simultaneously. These are null scenarios under MCAR labeling; MNAR appears in Fig. 10 . Classical tests on uncorrected judge scores instead yield false positives ( gray dots). The shaded band is a perfectly-calibrated test’s scatter at this replicate count. The -axis is square-root scaled to resolve the region near while keeping the uncorrected spread visible. § B.3 has full details.
Figure 5. PPI-corrected tests ( orange ) recover power under biased LLM judges, approaching a fully-human oracle ( green ) and improving on the small human-only subset ( gray dotted) across all data types. Pink : running tests over uncorrected scores show false-positive rates at zero effect; purple : swapping in some human gold labels still miscalibrates. Bottom: power across human-only size at Cohen’s . Continuous/Likert averages five PPI-corrected tests and their classical counterparts (independent, Welch’s, paired -tests, Wilcoxon, Mann-Whitney U); binary averages the two -tests. For more details, see § B.6 . 200 Monte Carlo replicates per condition (Appendix B ).
Figure 6. False-positive rate against human-LLM judge alignment on MCAR data, with five IRR metrics across three data types. Uncorrected, the false-positive rate peaks where agreement is high on Likert (c-e) and continuous data (b). PPI (blue) and a no-bias control (dashed) hold nominal throughout. Shading marks alignment our sweep did not reach. A wider sweep across 11 IRR metrics shows similar patterns (full-page plots in Appendix F ). 200 MC replicates per scenario (§ B ).
Figure 7. How many human ratings can LLM judges save, once judge bias is corrected for? Top : the label-efficiency multiplier ( ) plotted across eval types, the factor by which a PPI-corrected test with human labels out of judge-scored items behaves like a larger human-only sample , for . Judge quality is , the squared judge-human correlation (Pearson for parametric, Spearman for rank-based tests). Each point pools labeling budgets, effect sizes, two-sample tests, and judge-error shapes. Dotted lines show the predicted efficiency . Below : a table form of the same data, expressed in terms of the multiplier. means PPI doubled the effective sample size, acting like double the human raters. Each cell is the median over four effect sizes with a bootstrapped 95% CI. Likert and Binary exceed their predictions slightly at high through a power-curve inversion artifact. See Appendix D for simulation details.
Figure 8. End-to-end validation of evalstats ’ compare() with every method choice left to "auto" , across 7,176 cells with 100 MC repetitions. Top row : CI method coverage at increasing multiplicity: per-arm CIs, pairwise CIs at , and FWER-corrected simultaneous CIs at . Bottom row : family-wise CI coverage and -value Type-I error split by non-PPI and PPI scenarios (dashed). Bottom-right shows the power that PPI tests recover after FWER-correction for groups, over analyzing the human-labeled subset alone. Dashed black lines mark nominal ( coverage, , and zero gain). Coverage is nominal at every level, Type-I is controlled with and without PPI, and the PPI gains hold. For sweep details behind this plot, see § E .
Appendix figures & tables18 assets
Supplementary material from the paper’s appendix.
Appendix
| Model | Benchmark | Model | Benchmark |
|---|---|---|---|
| falcon-7b-instruct | mmlu-pro | gpt-4.1-mini | hi-tom |
| gpt-4o | culturalbench | gemma-3-27b-it | ifeval c |
| o4-mini | opentom | falcon-40b-instruct | imdb |
| llama-65b | bbq | qwen-2.5-72b-instruct | omni-math |
| vicuna-13b-v1.3 | cnndm c | kimi-k2 | salad-bench |
| DeepSeek-V3-0324 | do-not-answer c | llama-2-70b | xsum c |
Table 3. The two real evals corpora. OpenEval pairs one benchmark with each model; Inspect AI crosses every model with every benchmark. Unmarked benchmarks are binary, c marks continuous.
| Overall | Type-I error by | Real evals data | ||||||||||||
|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
| Method | Type-I | Max | Power | Type-I | Max | Power | ||||||||
| Binary | ||||||||||||||
| bootstrap | .029 | .110 | .370 | .019 | .025 | .026 | .040 | .040 | .047 | .127 | .368 | .037 | .045 | .049 |
| bayes_bootstrap | .043 | .143 | .408 | .047 | .041 | .043 | .046 | .046 | .064 | .173 | .402 | .074 | .059 | .056 |
| smooth_bootstrap | .020 | .070 | .342 | .013 | .017 | .022 | .029 | .030 | .033 | .093 | .341 | .028 | .035 | .040 |
| bootstrap_t | .062 | .190 | .374 | .050 | .033 | .076 | .069 | .065 | .059 | .213 | .368 | .052 | .051 | .050 |
Table 6. Pairwise p-value methods (nominal , 300 MC reps per cell), on synthetic and real evals data.
| Correction | FWER | Best-arm power | Time (ms) | Type | n=15 | n=30 | n=50 | n=100 | n=200 | n=500 | n=1000 | k=3 | k=5 | k=10 | k=20 |
|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
| holm (bin) | 0.023 | 0.219 | 0.0 | bin | 0.002 | 0.010 | 0.018 | 0.026 | 0.031 | 0.036 | 0.038 | 0.029 | 0.023 | 0.022 | 0.018 |
| bonferroni (bin) | 0.023 | 0.210 | 0.0 | bin | 0.002 | 0.010 | 0.018 | 0.026 | 0.031 | 0.036 | 0.038 | 0.029 | 0.023 | 0.022 | 0.018 |
| fdr_bh (bin) | 0.025 | 0.260 | 0.0 | bin | 0.003 | 0.011 | 0.020 | 0.028 | 0.033 | 0.040 | 0.041 | 0.031 | 0.026 | 0.024 | 0.021 |
| hochberg (bin) | 0.023 | 0.220 | 0.0 | bin | 0.002 | 0.010 | 0.018 | 0.026 | 0.031 | 0.037 | 0.038 | 0.029 | 0.024 | 0.022 | 0.018 |
| shaffer (bin) | 0.023 | 0.228 | 0.1 | bin | 0.002 | 0.010 | 0.018 | 0.026 | 0.031 | 0.036 | 0.038 | 0.029 | 0.023 | 0.022 | 0.018 |
| friedman_nemenyi (bin) | 0.000 | 0.078 | 12.7 | bin | 0.000 | 0.000 | 0.000 | 0.000 | 0.000 | 0.000 | 0.000 | 0.000 | 0.000 | 0.000 | 0.000 |
Table 9. p-value FWER correction methods’ Type I error rate and true effect discovery power across our synthetic data suite (500 MC reps per cell). Corrections are applied to the base -value evalstats reports for the data type—mid- on binary, Wilcoxon otherwise—except max-T, Romano–Wolf and Westfall–Young, which build their own joint statistic.
Figure 9. Family-wise error rate on real benchmark data, one gray line per set. Shading is the Monte-Carlo interval for a perfectly calibrated method at the reps behind this plot. Romano–Wolf is unstable below ; Shaffer is conservative throughout. 500 replicates per cell.
| Eval type | Population SD |
|---|---|
| Binary | 0.5000 |
| Continuous | 0.1206 |
| Likert (1–5) | 1.1447 |
Table 11. Population standard deviation of each eval type’s representative truth distribution, used as the standardization unit for all judge-bias/noise/effect-size magnitudes in Table 12 and § 6.1 . This ensures that we can compare power across data types.
Figure 10. PPI-corrected Type-I error for hypothesis tests under missing-not-at-random (MNAR) labeling, where items for human labeling were not randomly selected, violating PPI’s assumptions (e.g., labeling only the “easy” items). Computed from our aggressive factorial synthetic data sweep, which compounds multiple bias mechanisms per scenario. Spans 4,092 null scenarios (2,046 each at mild/strong MNAR strength) across all nine tests and 200 Monte Carlo repetitions each, for a total of 7.3M simulated trials. Rejection rates for paired t-test, Wilcoxon, and RM-ANOVA remain near nominal even as judge biases compound, while others exhibit instability. Friedman, MWU, and K-W rank-based tests appear to be the least reliable if MCAR labeling is violated.
| Synthetic, MCAR | Synthetic, MNAR | Real judge data | Power at | ||||||||
|---|---|---|---|---|---|---|---|---|---|---|---|
| Test | PPI | Unc. | Max | PPI | Unc. | Max | PPI | Unc. | Max | PPI | Subset |
| -test (indep) | .047 | .153 | .100 | .063 | .153 | .210 | .059 | .571 | .120 | .438 | .165 |
| Welch’s | .048 | .153 | .110 | .063 | .153 | .210 | .058 | .571 | .115 | .363 | .148 |
| paired | .053 | .154 | .120 | .054 | .154 | .220 | .053 | .788 | .110 | .477 | .232 |
| Mann-Whitney | .049 | .154 | .120 | .051 | .154 | .255 | .052 | .637 | .110 | .383 | .165 |
| Wilcoxon | .044 | .156 | .125 | .047 | .157 | .180 | .043 | .855 | .080 | .328 | .145 |
Table 13. Type I error of the nine PPI-corrected tests at nominal (200 MC reps per cell), the numbers behind Figures 4 (MCAR) and 10 (MNAR). PPI and Unc. are the mean rejection rate across null cells with and without correction. Max is the worst single null cell, the analogue of MinCov in the CI tables. Synthetic columns pool 2,046 MCAR and 4,092 MNAR null cells per test over continuous and Likert data; real columns pool 264 cells per test over the four multi-judge corpora of Table 14 . Power is at Cohen’s against a human-only subset of the same labeling budget. Per-factor breakdowns are in Supplementary Table S1.
Figure 11. Type I error control of PPI-corrected tests on real LLM judge data (counterpart to Figure 4 ) under an MCAR labeling regime. Uncorrected rejection rates (gray) reach 1.0 for every test; PPI-corrected rates (colored) cluster tightly around nominal . 200 Monte Carlo replicates per test and corpus.
Figure 12. PPI correction of the four CI methods evalstats returns for method="auto" , on the synthetic suite (red circles) and real judge data (blue squares). Small marks are individual scenarios, large marks their median; the dashed lines mark zero bias and the nominal 95% target. Bonett-Price’s real-data column runs conservative (99.7%) because our real data simulation gives both arms the same human labels, leaving its rectifier no disagreements to learn from. 2,000 replicates per synthetic scenario, 200 per real-data scenario.
Figure 13. Omnibus-test (3+ condition) analogue of Figure 5 , pooling ANOVA-independent, ANOVA-repeated, Friedman, and Kruskal-Wallis rejection rates across the same five-estimator comparison (all-human oracle, PPI-corrected, LLM-plus-label-overwrite, raw uncorrected LLM-only judgment, and a human-only subset matched in size to the labeled budget) for continuous and Likert data types. Top row: rejection rate as a function of true effect size at a fixed 20% labeled-item budget. Bottom row: rejection rate as a function of labeled-item count at effect size Cohen’s . PPI-corrected power tracks close behind the oracle and clearly ahead of both the human-only subset and the uncorrected LLM baselines in every panel. 200 Monte Carlo replicates per condition.
Figure 14. The variance multiplier is a property of judge quality only for parametric tests, and degrades with effect size for rank-based ones. Each panel sweeps the true effect at fixed judge quality ( , , continuous, near-centred truth marginal), plotting the implied by the measured multiplier (solid, MC standard error over 20,000 replicates). Parametric tests stay flat, as the algebra requires ( Angelopoulos et al., 2023a ) , while rank-based tests fall ( Mann-Whitney, Wilcoxon). Dashed lines show that evalstats ’ influence-function , output from its judge_alignment method, stays within of the measured value.
| Planned test | Use | Computed between LLM judge & humans on… |
| Estimation and two-condition tests | ||
| Mean, proportion, CI | raw scores, per measure | |
| Independent-samples -test | raw scores, within each condition | |
| Paired -test | difference scores ( ) | |
| Mann–Whitney | each item’s mid-rank placement against the other condition | |
| Wilcoxon signed-rank | difference scores, via their Hájek projection | |
Table 15. Which judge-human correlation to use in the formula, and what to compute it on. evalstats ’ judge_alignment computes the correct value for a given test from labeled data. For Mann–Whitney and Wilcoxon it correlates influence-function values, which reduce to Spearman’s at a null effect. (Note: For within-subjects/repeated-measures designs, randomly sample rows (subjects) with all conditions filled, and have human raters label all per selected subject.)
Figure 15. Label-efficiency multiplier against judge-human agreement , one line per data type, with bootstrap 95% CIs on the median. The shaded band marks savings under : statistically above , but perhaps not worth the cost. By all three data types have cleared it cleanly. Results are averaged over parametric and nonparametric tests using their respective , showing the threshold is invariant to test type, data type, effect size, and setup. A per-family breakdown is in Figure 16 . 500 replicates per cell (same run as Figure 7 ).
| CI coverage (target ) | PPI power gain over labels-only (points) | ||||||
| Data type | Per-arm | Pairwise | Family | Type-I | Peak | at | at |
| Binary | 0.957 | 0.952 | 0.969 | 0.031 | +18.9 | 100 | +0.5 |
| Continuous | 0.952 | 0.949 | 0.962 | 0.040 | +11.7 | 40 | +0.0 |
| Likert | 0.952 | 0.950 | 0.958 | 0.040 | +13.1 | 40 | +4.1 |
| All | 0.954 | 0.950 | 0.963 | 0.037 | +14.0 | 40 | +1.6 |
Table 16. End-to-end results of compare() by data type (100 MC reps per cell), a descriptive summary of Figure 8 .
Figure 17. 5-point Likert: uncorrected false-positive rate (red) against each IRR metric (rows), by judge bias (columns, in population SDs). Faint lines are the eleven judge variants; bold is their mean, with the dot marking its peak. PPI-corrected rates in blue hold at nominal ( dashed) in every cell.
Figure 18. 7-point Likert, same layout as Fig. 17 . The peaks sit higher (0.92–0.95) because half a scale point, the rounding threshold, is a smaller fraction of the scale’s SD.
Figure 19. Continuous data IRR metrics. False positive risk rises monotonically across IRR metrics.
Figure 20. Binary data IRR metrics as judge bias grows. False positive risk falls as agreement rises.
More in this topic
Statistical Hypothesis Testing
Statistical hypothesis testing develops methods for assessing evidence against a stated null while controlling false-positive and false-negative rates.
See all 101 papers