mirror of
https://github.com/wshobson/agents
synced 2026-06-21 14:13:58 +00:00
feat(plugin-eval): add statistics module (bootstrap, Wilson, Clopper-Pearson, kappa)
This commit is contained in:
Seth Hobson
@@ -0,0 +1,181 @@
|
||||
"""Statistical methods for PluginEval confidence intervals and reliability metrics."""
|
||||
|
||||
from __future__ import annotations
|
||||
|
||||
import math
|
||||
import random
|
||||
from collections import Counter
|
||||
|
||||
|
||||
def wilson_score_ci(
|
||||
successes: int,
|
||||
trials: int,
|
||||
confidence: float = 0.95,
|
||||
) -> tuple[float, float]:
|
||||
"""Wilson score interval for binomial proportion."""
|
||||
if trials == 0:
|
||||
raise ValueError("trials must be > 0")
|
||||
if successes > trials:
|
||||
raise ValueError(f"successes ({successes}) cannot exceed trials ({trials})")
|
||||
|
||||
z = _z_score(confidence)
|
||||
n = trials
|
||||
p_hat = successes / n
|
||||
|
||||
denominator = 1 + z**2 / n
|
||||
center = (p_hat + z**2 / (2 * n)) / denominator
|
||||
margin = (z / denominator) * math.sqrt(p_hat * (1 - p_hat) / n + z**2 / (4 * n**2))
|
||||
|
||||
lower = max(0.0, center - margin)
|
||||
upper = min(1.0, center + margin)
|
||||
return lower, upper
|
||||
|
||||
|
||||
def bootstrap_ci(
|
||||
data: list[float],
|
||||
confidence: float = 0.95,
|
||||
n_resamples: int = 1000,
|
||||
seed: int | None = None,
|
||||
) -> tuple[float, float]:
|
||||
"""Bootstrap confidence interval via percentile method."""
|
||||
if len(data) == 0:
|
||||
raise ValueError("data must not be empty")
|
||||
if len(data) == 1:
|
||||
return data[0], data[0]
|
||||
|
||||
rng = random.Random(seed)
|
||||
n = len(data)
|
||||
means = []
|
||||
for _ in range(n_resamples):
|
||||
sample = [rng.choice(data) for _ in range(n)]
|
||||
means.append(sum(sample) / n)
|
||||
|
||||
means.sort()
|
||||
alpha = 1 - confidence
|
||||
lower_idx = int(math.floor(alpha / 2 * n_resamples))
|
||||
upper_idx = int(math.ceil((1 - alpha / 2) * n_resamples)) - 1
|
||||
|
||||
lower_idx = max(0, min(lower_idx, n_resamples - 1))
|
||||
upper_idx = max(0, min(upper_idx, n_resamples - 1))
|
||||
|
||||
return means[lower_idx], means[upper_idx]
|
||||
|
||||
|
||||
def clopper_pearson_ci(
|
||||
failures: int,
|
||||
trials: int,
|
||||
confidence: float = 0.95,
|
||||
) -> tuple[float, float]:
|
||||
"""Clopper-Pearson exact confidence interval for failure rate."""
|
||||
if trials == 0:
|
||||
raise ValueError("trials must be > 0")
|
||||
if failures > trials:
|
||||
raise ValueError(f"failures ({failures}) cannot exceed trials ({trials})")
|
||||
|
||||
alpha = 1 - confidence
|
||||
|
||||
lower = 0.0 if failures == 0 else _beta_ppf(alpha / 2, failures, trials - failures + 1)
|
||||
upper = 1.0 if failures == trials else _beta_ppf(1 - alpha / 2, failures + 1, trials - failures)
|
||||
|
||||
return lower, upper
|
||||
|
||||
|
||||
def coefficient_of_variation(data: list[float]) -> float:
|
||||
"""Coefficient of variation (std / mean). Lower = more consistent."""
|
||||
if len(data) == 0:
|
||||
raise ValueError("data must not be empty")
|
||||
mean = sum(data) / len(data)
|
||||
if mean == 0:
|
||||
return 0.0
|
||||
variance = sum((x - mean) ** 2 for x in data) / len(data)
|
||||
return math.sqrt(variance) / abs(mean)
|
||||
|
||||
|
||||
def cohens_kappa(rater1: list[int], rater2: list[int]) -> float:
|
||||
"""Cohen's kappa for inter-rater agreement."""
|
||||
if len(rater1) != len(rater2):
|
||||
raise ValueError("Rater lists must have the same length")
|
||||
|
||||
n = len(rater1)
|
||||
categories = sorted(set(rater1) | set(rater2))
|
||||
|
||||
observed_agreement = sum(a == b for a, b in zip(rater1, rater2, strict=True)) / n
|
||||
|
||||
counts1 = Counter(rater1)
|
||||
counts2 = Counter(rater2)
|
||||
expected_agreement = sum((counts1.get(c, 0) / n) * (counts2.get(c, 0) / n) for c in categories)
|
||||
|
||||
if expected_agreement == 1.0:
|
||||
return 1.0
|
||||
|
||||
return (observed_agreement - expected_agreement) / (1 - expected_agreement)
|
||||
|
||||
|
||||
def _z_score(confidence: float) -> float:
|
||||
"""Approximate z-score for common confidence levels."""
|
||||
z_table = {0.90: 1.645, 0.95: 1.960, 0.99: 2.576}
|
||||
if confidence in z_table:
|
||||
return z_table[confidence]
|
||||
p = (1 + confidence) / 2
|
||||
t = math.sqrt(-2 * math.log(1 - p))
|
||||
c0, c1, c2 = 2.515517, 0.802853, 0.010328
|
||||
d1, d2, d3 = 1.432788, 0.189269, 0.001308
|
||||
return t - (c0 + c1 * t + c2 * t**2) / (1 + d1 * t + d2 * t**2 + d3 * t**3)
|
||||
|
||||
|
||||
def _beta_ppf(p: float, a: float, b: float) -> float:
|
||||
"""Approximate beta distribution percent point function via Newton's method."""
|
||||
if p <= 0:
|
||||
return 0.0
|
||||
if p >= 1:
|
||||
return 1.0
|
||||
|
||||
mu = a / (a + b)
|
||||
sigma = math.sqrt(a * b / ((a + b) ** 2 * (a + b + 1)))
|
||||
z = _z_score(2 * p - 1) if p > 0.5 else -_z_score(2 * (1 - p) - 1)
|
||||
x = max(0.001, min(0.999, mu + sigma * z))
|
||||
|
||||
for _ in range(50):
|
||||
cdf = _beta_cdf(x, a, b)
|
||||
pdf = _beta_pdf(x, a, b)
|
||||
if pdf < 1e-12:
|
||||
break
|
||||
x_new = x - (cdf - p) / pdf
|
||||
x_new = max(0.001, min(0.999, x_new))
|
||||
if abs(x_new - x) < 1e-10:
|
||||
break
|
||||
x = x_new
|
||||
|
||||
return x
|
||||
|
||||
|
||||
def _beta_pdf(x: float, a: float, b: float) -> float:
|
||||
"""Beta distribution probability density function."""
|
||||
if x <= 0 or x >= 1:
|
||||
return 0.0
|
||||
log_pdf = (a - 1) * math.log(x) + (b - 1) * math.log(1 - x) - _log_beta(a, b)
|
||||
return math.exp(log_pdf)
|
||||
|
||||
|
||||
def _beta_cdf(x: float, a: float, b: float, steps: int = 200) -> float:
|
||||
"""Beta CDF via numerical integration (Simpson's rule)."""
|
||||
if x <= 0:
|
||||
return 0.0
|
||||
if x >= 1:
|
||||
return 1.0
|
||||
|
||||
h = x / steps
|
||||
total = _beta_pdf(0.001, a, b) + _beta_pdf(x, a, b)
|
||||
for i in range(1, steps):
|
||||
xi = i * h
|
||||
if xi <= 0 or xi >= 1:
|
||||
continue
|
||||
weight = 4 if i % 2 == 1 else 2
|
||||
total += weight * _beta_pdf(xi, a, b)
|
||||
|
||||
return total * h / 3
|
||||
|
||||
|
||||
def _log_beta(a: float, b: float) -> float:
|
||||
"""Log of the beta function."""
|
||||
return math.lgamma(a) + math.lgamma(b) - math.lgamma(a + b)
|
||||
@@ -0,0 +1,100 @@
|
||||
import pytest
|
||||
|
||||
from plugin_eval.stats import (
|
||||
bootstrap_ci,
|
||||
clopper_pearson_ci,
|
||||
cohens_kappa,
|
||||
coefficient_of_variation,
|
||||
wilson_score_ci,
|
||||
)
|
||||
|
||||
|
||||
class TestWilsonScore:
|
||||
def test_perfect_activation(self):
|
||||
lower, upper = wilson_score_ci(successes=50, trials=50, confidence=0.95)
|
||||
assert lower > 0.90
|
||||
assert upper == pytest.approx(1.0, abs=0.01)
|
||||
|
||||
def test_half_activation(self):
|
||||
lower, upper = wilson_score_ci(successes=25, trials=50, confidence=0.95)
|
||||
assert lower < 0.50
|
||||
assert upper > 0.50
|
||||
assert lower > 0.35
|
||||
assert upper < 0.65
|
||||
|
||||
def test_zero_trials_raises(self):
|
||||
with pytest.raises(ValueError):
|
||||
wilson_score_ci(successes=0, trials=0)
|
||||
|
||||
def test_successes_exceed_trials_raises(self):
|
||||
with pytest.raises(ValueError):
|
||||
wilson_score_ci(successes=10, trials=5)
|
||||
|
||||
|
||||
class TestBootstrapCI:
|
||||
def test_tight_data(self):
|
||||
data = [0.80, 0.82, 0.81, 0.83, 0.79, 0.80, 0.82, 0.81]
|
||||
lower, upper = bootstrap_ci(data, confidence=0.95, n_resamples=1000, seed=42)
|
||||
assert lower > 0.78
|
||||
assert upper < 0.84
|
||||
assert lower < upper
|
||||
|
||||
def test_single_value(self):
|
||||
lower, upper = bootstrap_ci([0.5], confidence=0.95, n_resamples=100, seed=42)
|
||||
assert lower == pytest.approx(0.5)
|
||||
assert upper == pytest.approx(0.5)
|
||||
|
||||
def test_empty_raises(self):
|
||||
with pytest.raises(ValueError):
|
||||
bootstrap_ci([], confidence=0.95)
|
||||
|
||||
|
||||
class TestClopperPearson:
|
||||
def test_zero_failures(self):
|
||||
lower, upper = clopper_pearson_ci(failures=0, trials=50, confidence=0.95)
|
||||
assert lower == 0.0
|
||||
assert upper < 0.10
|
||||
|
||||
def test_some_failures(self):
|
||||
lower, upper = clopper_pearson_ci(failures=2, trials=50, confidence=0.95)
|
||||
assert lower < 0.04
|
||||
assert upper > 0.04
|
||||
assert upper < 0.15
|
||||
|
||||
def test_zero_trials_raises(self):
|
||||
with pytest.raises(ValueError):
|
||||
clopper_pearson_ci(failures=0, trials=0)
|
||||
|
||||
|
||||
class TestCoefficientOfVariation:
|
||||
def test_low_variation(self):
|
||||
data = [0.80, 0.82, 0.81, 0.83, 0.79]
|
||||
cv = coefficient_of_variation(data)
|
||||
assert cv < 0.05
|
||||
|
||||
def test_high_variation(self):
|
||||
data = [0.20, 0.90, 0.10, 0.95, 0.50]
|
||||
cv = coefficient_of_variation(data)
|
||||
assert cv > 0.40
|
||||
|
||||
def test_empty_raises(self):
|
||||
with pytest.raises(ValueError):
|
||||
coefficient_of_variation([])
|
||||
|
||||
|
||||
class TestCohensKappa:
|
||||
def test_perfect_agreement(self):
|
||||
rater1 = [1, 2, 3, 4, 5]
|
||||
rater2 = [1, 2, 3, 4, 5]
|
||||
k = cohens_kappa(rater1, rater2)
|
||||
assert k == pytest.approx(1.0)
|
||||
|
||||
def test_no_agreement(self):
|
||||
rater1 = [1, 2, 3, 4, 5]
|
||||
rater2 = [5, 4, 3, 2, 1]
|
||||
k = cohens_kappa(rater1, rater2)
|
||||
assert k < 0.0
|
||||
|
||||
def test_mismatched_length_raises(self):
|
||||
with pytest.raises(ValueError):
|
||||
cohens_kappa([1, 2], [1, 2, 3])
|
||||
Reference in New Issue
Block a user