reach.uncertainty

Statistical uncertainty estimation, Wilson score confidence intervals, and power analysis.

Calculate statistical confidence intervals and power for rate metrics.

Interval

Bases: BaseModel

Represent a statistical confidence interval with lower and upper bounds.

Source code in src/reach/uncertainty.py

class Interval(BaseModel):
    """Represent a statistical confidence interval with lower and upper bounds."""

    model_config = ConfigDict(frozen=True, extra="forbid")

    low: Annotated[float, Field(ge=-1.0, le=1.0)]
    high: Annotated[float, Field(ge=-1.0, le=1.0)]
    confidence: Annotated[float, Field(gt=0.0, lt=1.0)] = DEFAULT_CONFIDENCE

    @classmethod
    def zero(cls, confidence: float = DEFAULT_CONFIDENCE) -> Self:
        """Return a zero-width interval centered at zero."""
        return cls(low=0.0, high=0.0, confidence=confidence)

    @classmethod
    def unit(cls, confidence: float = DEFAULT_CONFIDENCE) -> Self:
        """Return a unit interval [0.0, 1.0]."""
        return cls(low=0.0, high=1.0, confidence=confidence)

    @classmethod
    def from_tuple(
        cls,
        bounds: Interval | Sequence[float],
        confidence: float = DEFAULT_CONFIDENCE,
    ) -> Self:
        """Construct an Interval from an existing Interval or a 2-element sequence."""
        if isinstance(bounds, Interval):
            if bounds.confidence == confidence and type(bounds) is cls:
                return cast(Self, bounds)
            return cls(low=bounds.low, high=bounds.high, confidence=confidence)
        if len(bounds) != _INTERVAL_BOUNDS_LEN:
            msg = f"expected 2 elements for interval bounds, got {len(bounds)}"
            raise ValueError(msg)
        return cls(low=float(bounds[0]), high=float(bounds[1]), confidence=confidence)

    @model_validator(mode="after")
    def _bounds_are_ordered(self) -> Self:
        """Validate that lower bound does not exceed upper bound."""
        if self.low > self.high:
            msg = f"interval bounds are inverted: [{self.low}, {self.high}]"
            raise ValueError(msg)
        return self

    @property
    def width(self) -> float:
        """Return the span (high - low) of the confidence interval."""
        return self.high - self.low

    def excludes(self, rate: float) -> bool:
        """Return True if rate falls strictly outside the interval bounds."""
        return not self.low <= rate <= self.high

    def overlaps(self, other: Interval) -> bool:
        """Return True if this interval intersects with another interval."""
        return self.low <= other.high and other.low <= self.high

    def format_percent(self, digits: int = 1, separator: str = " - ") -> str:
        """Format the interval as a percentage range string '[low% - high%]'."""
        return f"[{self.low * 100:.{digits}f}%{separator}{self.high * 100:.{digits}f}%]"

width property

width: float

Return the span (high - low) of the confidence interval.

excludes

excludes(rate: float) -> bool

Return True if rate falls strictly outside the interval bounds.

Source code in src/reach/uncertainty.py

def excludes(self, rate: float) -> bool:
    """Return True if rate falls strictly outside the interval bounds."""
    return not self.low <= rate <= self.high

format_percent

format_percent(
    digits: int = 1, separator: str = " - "
) -> str

Format the interval as a percentage range string '[low% - high%]'.

Source code in src/reach/uncertainty.py

def format_percent(self, digits: int = 1, separator: str = " - ") -> str:
    """Format the interval as a percentage range string '[low% - high%]'."""
    return f"[{self.low * 100:.{digits}f}%{separator}{self.high * 100:.{digits}f}%]"

from_tuple classmethod

from_tuple(
    bounds: Interval | Sequence[float],
    confidence: float = DEFAULT_CONFIDENCE,
) -> Self

Construct an Interval from an existing Interval or a 2-element sequence.

Source code in src/reach/uncertainty.py

@classmethod
def from_tuple(
    cls,
    bounds: Interval | Sequence[float],
    confidence: float = DEFAULT_CONFIDENCE,
) -> Self:
    """Construct an Interval from an existing Interval or a 2-element sequence."""
    if isinstance(bounds, Interval):
        if bounds.confidence == confidence and type(bounds) is cls:
            return cast(Self, bounds)
        return cls(low=bounds.low, high=bounds.high, confidence=confidence)
    if len(bounds) != _INTERVAL_BOUNDS_LEN:
        msg = f"expected 2 elements for interval bounds, got {len(bounds)}"
        raise ValueError(msg)
    return cls(low=float(bounds[0]), high=float(bounds[1]), confidence=confidence)

overlaps

overlaps(other: Interval) -> bool

Return True if this interval intersects with another interval.

Source code in src/reach/uncertainty.py

def overlaps(self, other: Interval) -> bool:
    """Return True if this interval intersects with another interval."""
    return self.low <= other.high and other.low <= self.high

unit classmethod

unit(confidence: float = DEFAULT_CONFIDENCE) -> Self

Return a unit interval [0.0, 1.0].

Source code in src/reach/uncertainty.py

@classmethod
def unit(cls, confidence: float = DEFAULT_CONFIDENCE) -> Self:
    """Return a unit interval [0.0, 1.0]."""
    return cls(low=0.0, high=1.0, confidence=confidence)

zero classmethod

zero(confidence: float = DEFAULT_CONFIDENCE) -> Self

Return a zero-width interval centered at zero.

Source code in src/reach/uncertainty.py

@classmethod
def zero(cls, confidence: float = DEFAULT_CONFIDENCE) -> Self:
    """Return a zero-width interval centered at zero."""
    return cls(low=0.0, high=0.0, confidence=confidence)

cluster_wilson_interval

cluster_wilson_interval(
    hits: int,
    probes: int,
    attempts: int = 1,
    confidence: float = DEFAULT_CONFIDENCE,
    intra_cluster_correlation: float = 0.6,
    *,
    queries_per_skill: float = 1.0,
    skill_icc: float = 0.0,
) -> Interval | None

Calculate Wilson score confidence interval adjusted for cluster design effect.

Source code in src/reach/uncertainty.py

def cluster_wilson_interval(
    hits: int,
    probes: int,
    attempts: int = 1,
    confidence: float = DEFAULT_CONFIDENCE,
    intra_cluster_correlation: float = 0.6,
    *,
    queries_per_skill: float = 1.0,
    skill_icc: float = 0.0,
) -> Interval | None:
    """Calculate Wilson score confidence interval adjusted for cluster design effect."""
    deff = _design_effect(attempts, intra_cluster_correlation, queries_per_skill, skill_icc)
    if deff <= 1.0:
        return wilson_interval(hits, probes, confidence=confidence)
    return _wilson_from_counts(hits, probes, max(1.0, probes / deff), confidence)

detectable_delta

detectable_delta(
    probes: int,
    confidence: float = DEFAULT_CONFIDENCE,
    power: float = DEFAULT_POWER,
) -> float | None

Calculate the minimum detectable rate delta for a given sample size per arm.

Source code in src/reach/uncertainty.py

def detectable_delta(
    probes: int,
    confidence: float = DEFAULT_CONFIDENCE,
    power: float = DEFAULT_POWER,
) -> float | None:
    """Calculate the minimum detectable rate delta for a given sample size per arm."""
    if probes < 0:
        msg = f"probes cannot be negative, got {probes}"
        raise ValueError(msg)
    if not probes:
        return None
    alpha = 1.0 - _checked_confidence(confidence)
    if not 0.0 < power < 1.0:
        msg = f"power must lie strictly between 0 and 1, got {power}"
        raise ValueError(msg)
    return min(1.0, math.sqrt(_two_proportion_constant(alpha, power) / probes))

effective_sample_size

effective_sample_size(
    sample_size: int,
    attempts: int = 1,
    intra_cluster_correlation: float = 0.6,
    *,
    queries_per_skill: float = 1.0,
    skill_icc: float = 0.0,
) -> int

Calculate survey-style effective sample size adjusting for cluster correlation.

Source code in src/reach/uncertainty.py

def effective_sample_size(
    sample_size: int,
    attempts: int = 1,
    intra_cluster_correlation: float = 0.6,
    *,
    queries_per_skill: float = 1.0,
    skill_icc: float = 0.0,
) -> int:
    """Calculate survey-style effective sample size adjusting for cluster correlation."""
    if sample_size <= 0:
        return 0
    deff = _design_effect(attempts, intra_cluster_correlation, queries_per_skill, skill_icc)
    return max(1, round(sample_size / deff))

estimate_skill_icc

estimate_skill_icc(
    outcomes_by_skill: Mapping[str, Sequence[float]],
) -> float | None

Estimate intra-skill correlation via one-way random-effects ANOVA across skills.

Parameters:

Name	Type	Description	Default
outcomes_by_skill	Mapping[str, Sequence[float]]	Mapping of skill names to sequences of per-query pass rates.	required

Returns:

Type	Description
float | None	Estimated intra-class correlation coefficient in [0.0, 1.0], or None when
float | None	fewer than two skills or zero within-skill degrees of freedom are available.

Source code in src/reach/uncertainty.py

def estimate_skill_icc(
    outcomes_by_skill: Mapping[str, Sequence[float]],
) -> float | None:
    """Estimate intra-skill correlation via one-way random-effects ANOVA across skills.

    Args:
        outcomes_by_skill: Mapping of skill names to sequences of per-query pass rates.

    Returns:
        Estimated intra-class correlation coefficient in [0.0, 1.0], or None when
        fewer than two skills or zero within-skill degrees of freedom are available.
    """
    groups = [[float(v) for v in vals] for vals in outcomes_by_skill.values() if len(vals) >= 1]
    s_skills = len(groups)
    if s_skills < _MIN_ANOVA_GROUPS:
        return None
    total_q = sum(len(g) for g in groups)
    if total_q <= s_skills:
        return None

    grand_mean = sum(sum(g) for g in groups) / total_q
    ssb = 0.0
    ssw = 0.0
    for g in groups:
        q_s = len(g)
        mean_s = sum(g) / q_s
        ssb += q_s * ((mean_s - grand_mean) ** 2)
        ssw += sum((y - mean_s) ** 2 for y in g)

    msb = ssb / (s_skills - 1)
    msw = ssw / (total_q - s_skills)
    if msb == 0.0 and msw == 0.0:
        return 0.0

    q_bar = (total_q - sum(len(g) ** 2 for g in groups) / total_q) / (s_skills - 1)
    denom = msb + (q_bar - 1.0) * msw
    if denom <= 0.0:
        return 0.0
    rho = (msb - msw) / denom
    return round(max(0.0, min(1.0, rho)), 4)

required_probes

required_probes(
    delta: float,
    confidence: float = DEFAULT_CONFIDENCE,
    power: float = DEFAULT_POWER,
) -> int

Calculate sample size per arm required to detect a rate delta with power.

Source code in src/reach/uncertainty.py

def required_probes(
    delta: float,
    confidence: float = DEFAULT_CONFIDENCE,
    power: float = DEFAULT_POWER,
) -> int:
    """Calculate sample size per arm required to detect a rate delta with power."""
    if not 0.0 < delta <= 1.0:
        msg = f"delta must lie in (0, 1], got {delta}"
        raise ValueError(msg)
    alpha = 1.0 - _checked_confidence(confidence)
    if not 0.0 < power < 1.0:
        msg = f"power must lie strictly between 0 and 1, got {power}"
        raise ValueError(msg)
    return math.ceil(_two_proportion_constant(alpha, power) / (delta * delta))

wilson_interval

wilson_interval(
    hits: int,
    probes: int,
    confidence: float = DEFAULT_CONFIDENCE,
) -> Interval | None

Calculate the Wilson score confidence interval for a binomial proportion.

Source code in src/reach/uncertainty.py

def wilson_interval(
    hits: int,
    probes: int,
    confidence: float = DEFAULT_CONFIDENCE,
) -> Interval | None:
    """Calculate the Wilson score confidence interval for a binomial proportion."""
    return _wilson_from_counts(hits, probes, float(probes), confidence)