pictor_eval/
bootstrap.rs

1//! Bootstrap confidence intervals.
2//!
3//! Given a sample of per-example metric values, [`bootstrap_ci`] produces a
4//! seed-deterministic confidence interval around the sample mean by
5//! resampling with replacement `n_resamples` times and reporting the
6//! percentile interval `[confidence_lo, confidence_hi]`.
7//!
8//! Uses a self-contained 64-bit xorshift PRNG (no external `rand` dependency)
9//! so results are bit-reproducible across platforms for a given seed.
10
11use serde::Serialize;
12
13use crate::error::EvalError;
14
15/// A bootstrap percentile confidence interval.
16#[derive(Debug, Clone, Serialize)]
17pub struct ConfidenceInterval {
18    /// Lower bound (confidence_lo percentile of resampled means).
19    pub lo: f32,
20    /// Upper bound (confidence_hi percentile of resampled means).
21    pub hi: f32,
22    /// Mean of the original sample.
23    pub mean: f32,
24    /// Confidence level used (e.g. 0.95 for 95% CI).
25    pub confidence: f32,
26    /// Number of resamples performed.
27    pub n_resamples: usize,
28}
29
30// ──────────────────────────────────────────────────────────────────────────────
31// Xorshift64*
32// ──────────────────────────────────────────────────────────────────────────────
33
34/// Marsaglia xorshift64* step. Avoids a zero seed by mixing in a constant.
35#[inline]
36fn xorshift64(state: u64) -> u64 {
37    // Ensure non-zero state — the transform hits a fixed point at 0.
38    let mut x = if state == 0 {
39        0x9E37_79B9_7F4A_7C15
40    } else {
41        state
42    };
43    x ^= x << 13;
44    x ^= x >> 7;
45    x ^= x << 17;
46    x
47}
48
49/// Draw one uniform `usize` in `[0, n)` from a 64-bit state (mutated).
50#[inline]
51fn next_index(state: &mut u64, n: usize) -> usize {
52    *state = xorshift64(*state);
53    if n == 0 {
54        0
55    } else {
56        (*state as usize) % n
57    }
58}
59
60// ──────────────────────────────────────────────────────────────────────────────
61// Public API
62// ──────────────────────────────────────────────────────────────────────────────
63
64/// Compute a bootstrap percentile confidence interval.
65///
66/// - `samples`: per-example metric values.
67/// - `n_resamples`: number of resamples (≥ 1; 0 yields a zero-width CI at the mean).
68/// - `confidence`: confidence level in `(0, 1)`, typically 0.95.
69/// - `seed`: PRNG seed for reproducibility.
70///
71/// Returns `Err(EvalError::DatasetEmpty)` when `samples` is empty, or
72/// `Err(EvalError::Numerical)` when `confidence` is out of range.
73pub fn bootstrap_ci(
74    samples: &[f32],
75    n_resamples: usize,
76    confidence: f32,
77    seed: u64,
78) -> Result<ConfidenceInterval, EvalError> {
79    if samples.is_empty() {
80        return Err(EvalError::DatasetEmpty);
81    }
82    if !(confidence > 0.0 && confidence < 1.0) {
83        return Err(EvalError::Numerical(format!(
84            "confidence must be in (0,1), got {}",
85            confidence
86        )));
87    }
88
89    let mean = mean_f32(samples);
90
91    if n_resamples == 0 {
92        return Ok(ConfidenceInterval {
93            lo: mean,
94            hi: mean,
95            mean,
96            confidence,
97            n_resamples: 0,
98        });
99    }
100
101    let mut state = if seed == 0 {
102        0xDEAD_BEEF_CAFE_F00D
103    } else {
104        seed
105    };
106    // Avoid a degenerate all-zero xorshift state on subsequent calls.
107    state = xorshift64(state);
108
109    let n = samples.len();
110    let mut resample_means: Vec<f32> = Vec::with_capacity(n_resamples);
111
112    for _ in 0..n_resamples {
113        let mut acc = 0.0f64;
114        for _ in 0..n {
115            let idx = next_index(&mut state, n);
116            acc += samples[idx] as f64;
117        }
118        resample_means.push((acc / n as f64) as f32);
119    }
120
121    resample_means.sort_by(|a, b| a.partial_cmp(b).unwrap_or(std::cmp::Ordering::Equal));
122
123    let alpha = 1.0 - confidence;
124    let lo_q = alpha / 2.0;
125    let hi_q = 1.0 - alpha / 2.0;
126    let lo = percentile(&resample_means, lo_q);
127    let hi = percentile(&resample_means, hi_q);
128
129    Ok(ConfidenceInterval {
130        lo,
131        hi,
132        mean,
133        confidence,
134        n_resamples,
135    })
136}
137
138fn mean_f32(xs: &[f32]) -> f32 {
139    if xs.is_empty() {
140        0.0
141    } else {
142        let s: f64 = xs.iter().map(|&v| v as f64).sum();
143        (s / xs.len() as f64) as f32
144    }
145}
146
147/// Percentile of a *sorted* slice via linear interpolation on rank.
148fn percentile(sorted: &[f32], q: f32) -> f32 {
149    if sorted.is_empty() {
150        return 0.0;
151    }
152    if sorted.len() == 1 {
153        return sorted[0];
154    }
155    let q = q.clamp(0.0, 1.0);
156    let rank = q * (sorted.len() - 1) as f32;
157    let lo = rank.floor() as usize;
158    let hi = rank.ceil() as usize;
159    if lo == hi {
160        sorted[lo]
161    } else {
162        let frac = rank - lo as f32;
163        sorted[lo] + (sorted[hi] - sorted[lo]) * frac
164    }
165}
pictor_eval/bootstrap.rs

pictor_eval/
bootstrap.rs
