gam_math/
quantile.rs

1/// Linear-interpolation quantile matching numpy.quantile default (method='linear').
2pub fn quantile_from_sorted(sorted: &[f64], q: f64) -> f64 {
3    let n = sorted.len();
4    if n == 0 {
5        return f64::NAN;
6    }
7    if n == 1 {
8        return sorted[0];
9    }
10    let pos = q.clamp(0.0, 1.0) * (n - 1) as f64;
11    let lo = pos.floor() as usize;
12    let hi = (lo + 1).min(n - 1);
13    let frac = pos - lo as f64;
14    sorted[lo] * (1.0 - frac) + sorted[hi] * frac
15}
16
17/// Exact 1-based order statistic from an already sorted slice.
18///
19/// Returns `NaN` when the slice is empty, `rank` is zero, or `rank` exceeds
20/// the number of observations. This is intentionally not an interpolating
21/// quantile: split-conformal calibration needs the observed `k`-th value to
22/// preserve the finite-sample coverage proof.
23pub fn order_statistic_from_sorted(sorted: &[f64], rank: usize) -> f64 {
24    if sorted.is_empty() || rank == 0 || rank > sorted.len() {
25        return f64::NAN;
26    }
27    sorted[rank - 1]
28}
29
30/// Exact 1-based order statistic from an unsorted slice.
31///
32/// This centralizes the sort+select path for code that needs an observed
33/// sample value rather than `quantile_from_sorted`'s linear interpolation.
34pub fn order_statistic(values: &[f64], rank: usize) -> f64 {
35    let mut sorted = values.to_vec();
36    sorted.sort_by(f64::total_cmp);
37    order_statistic_from_sorted(&sorted, rank)
38}
39
40#[cfg(test)]
41mod tests {
42    use super::*;
43
44    // ── quantile_from_sorted ───────────────────────────────────────────────
45
46    #[test]
47    fn quantile_from_sorted_empty_returns_nan() {
48        assert!(quantile_from_sorted(&[], 0.5).is_nan());
49    }
50
51    #[test]
52    fn quantile_from_sorted_single_element_returns_it_for_any_q() {
53        assert_eq!(quantile_from_sorted(&[7.0], 0.0), 7.0);
54        assert_eq!(quantile_from_sorted(&[7.0], 0.5), 7.0);
55        assert_eq!(quantile_from_sorted(&[7.0], 1.0), 7.0);
56    }
57
58    #[test]
59    fn quantile_from_sorted_q0_returns_min_q1_returns_max() {
60        let v = [1.0, 2.0, 3.0, 4.0, 5.0];
61        assert_eq!(quantile_from_sorted(&v, 0.0), 1.0);
62        assert_eq!(quantile_from_sorted(&v, 1.0), 5.0);
63    }
64
65    #[test]
66    fn quantile_from_sorted_median_five_element() {
67        // numpy.quantile([1,2,3,4,5], 0.5, method='linear') == 3.0
68        let v = [1.0, 2.0, 3.0, 4.0, 5.0];
69        assert_eq!(quantile_from_sorted(&v, 0.5), 3.0);
70    }
71
72    #[test]
73    fn quantile_from_sorted_interpolates_between_adjacent_elements() {
74        // [1.0, 3.0]: q=0.5 → pos=0.5, lo=0, hi=1, frac=0.5 → 1.0*0.5 + 3.0*0.5 = 2.0
75        let v = [1.0, 3.0];
76        assert_eq!(quantile_from_sorted(&v, 0.5), 2.0);
77        // q=0.25 → pos=0.25, lo=0, frac=0.25 → 1.0*0.75 + 3.0*0.25 = 1.5
78        assert!((quantile_from_sorted(&v, 0.25) - 1.5).abs() < 1e-14);
79        // q=0.75 → pos=0.75, lo=0, frac=0.75 → 1.0*0.25 + 3.0*0.75 = 2.5
80        assert!((quantile_from_sorted(&v, 0.75) - 2.5).abs() < 1e-14);
81    }
82
83    #[test]
84    fn quantile_from_sorted_q_clamped_below_zero() {
85        // q < 0 should clamp to q=0 → return minimum
86        let v = [2.0, 4.0, 6.0];
87        assert_eq!(quantile_from_sorted(&v, -1.0), 2.0);
88        assert_eq!(quantile_from_sorted(&v, -0.0001), 2.0);
89    }
90
91    #[test]
92    fn quantile_from_sorted_q_clamped_above_one() {
93        // q > 1 should clamp to q=1 → return maximum
94        let v = [2.0, 4.0, 6.0];
95        assert_eq!(quantile_from_sorted(&v, 2.0), 6.0);
96        assert_eq!(quantile_from_sorted(&v, 1.0001), 6.0);
97    }
98
99    #[test]
100    fn quantile_from_sorted_two_element_boundary_ranks() {
101        // q=0 → first, q=1 → last
102        let v = [10.0, 20.0];
103        assert_eq!(quantile_from_sorted(&v, 0.0), 10.0);
104        assert_eq!(quantile_from_sorted(&v, 1.0), 20.0);
105    }
106
107    #[test]
108    fn quantile_from_sorted_matches_numpy_linear_reference() {
109        // Reference values from numpy.quantile([0,1,2,3,4,5,6,7,8,9], q, method='linear')
110        let v: Vec<f64> = (0..10).map(|i| i as f64).collect();
111        // q=0.1 → pos=0.9, lo=0, hi=1, frac=0.9 → 0*0.1 + 1*0.9 = 0.9
112        assert!((quantile_from_sorted(&v, 0.1) - 0.9).abs() < 1e-14);
113        // q=0.9 → pos=8.1, lo=8, hi=9, frac=0.1 → 8*0.9 + 9*0.1 = 8.1
114        assert!((quantile_from_sorted(&v, 0.9) - 8.1).abs() < 1e-14);
115        // q=0.5 → pos=4.5, lo=4, hi=5, frac=0.5 → 4*0.5 + 5*0.5 = 4.5
116        assert!((quantile_from_sorted(&v, 0.5) - 4.5).abs() < 1e-14);
117    }
118
119    // ── order_statistic / order_statistic_from_sorted ─────────────────────
120
121    #[test]
122    fn order_statistic_returns_nan_on_empty() {
123        assert!(order_statistic(&[], 1).is_nan());
124        assert!(order_statistic_from_sorted(&[], 1).is_nan());
125    }
126
127    #[test]
128    fn order_statistic_returns_nan_on_zero_rank() {
129        let v = [3.0, 1.0, 2.0];
130        assert!(order_statistic(&v, 0).is_nan());
131        let sorted = [1.0, 2.0, 3.0];
132        assert!(order_statistic_from_sorted(&sorted, 0).is_nan());
133    }
134
135    #[test]
136    fn order_statistic_returns_nan_when_rank_exceeds_len() {
137        let v = [3.0, 1.0, 2.0];
138        assert!(order_statistic(&v, 4).is_nan());
139        let sorted = [1.0, 2.0, 3.0];
140        assert!(order_statistic_from_sorted(&sorted, 4).is_nan());
141    }
142
143    #[test]
144    fn order_statistic_hits_mid_rank() {
145        // Unsorted input; the 1-based 3rd-smallest of {1,2,3,4,5} is 3.
146        let v = [5.0, 1.0, 4.0, 2.0, 3.0];
147        assert_eq!(order_statistic(&v, 3), 3.0);
148        let sorted = [1.0, 2.0, 3.0, 4.0, 5.0];
149        assert_eq!(order_statistic_from_sorted(&sorted, 3), 3.0);
150        // Boundary ranks: first and last.
151        assert_eq!(order_statistic(&v, 1), 1.0);
152        assert_eq!(order_statistic(&v, 5), 5.0);
153    }
154}
gam_math/quantile.rs

gam_math/
quantile.rs
