Skip to main content

qualia_core_db/solvers/statistics/hypothesis/
chi_square.rs

1//! χ² tests — goodness-of-fit and independence (contingency table), with a real
2//! upper-tail p-value from [`chi_squared`](super::super::distributions::chi_squared).
3
4use super::super::distributions::chi_squared;
5
6/// χ² test result.
7#[derive(Debug, Clone, Copy, PartialEq)]
8pub struct ChiSquareResult {
9    pub statistic: f64,
10    pub p_value: f64,
11    pub dof: f64,
12}
13
14/// Pearson goodness-of-fit: `Σ (Oᵢ − Eᵢ)² / Eᵢ` with `dof = k − 1`. `None` if the
15/// vectors differ in length, are shorter than 2, or any expected count is ≤ 0.
16pub fn chi_square_gof(observed: &[f64], expected: &[f64]) -> Option<ChiSquareResult> {
17    if observed.len() != expected.len() || observed.len() < 2 {
18        return None;
19    }
20    if expected.iter().any(|&e| e <= 0.0) {
21        return None;
22    }
23    let stat: f64 = observed
24        .iter()
25        .zip(expected.iter())
26        .map(|(&o, &e)| (o - e).powi(2) / e)
27        .sum();
28    let dof = (observed.len() - 1) as f64;
29    Some(ChiSquareResult {
30        statistic: stat,
31        p_value: chi_squared::upper_p(stat, dof),
32        dof,
33    })
34}
35
36/// χ² test of independence on an `R×C` contingency table of counts. Expected
37/// `Eᵢⱼ = rowᵢ·colⱼ / total`, `dof = (R−1)(C−1)`. `None` if the table is not at
38/// least 2×2, is ragged, or the grand total is 0.
39pub fn chi_square_independence(table: &[&[f64]]) -> Option<ChiSquareResult> {
40    let rows = table.len();
41    if rows < 2 {
42        return None;
43    }
44    let cols = table[0].len();
45    if cols < 2 || table.iter().any(|r| r.len() != cols) {
46        return None;
47    }
48    let row_sums: Vec<f64> = table.iter().map(|r| r.iter().sum()).collect();
49    let mut col_sums = vec![0.0; cols];
50    for r in table {
51        for (j, &v) in r.iter().enumerate() {
52            col_sums[j] += v;
53        }
54    }
55    let total: f64 = row_sums.iter().sum();
56    if total <= 0.0 {
57        return None;
58    }
59    let mut stat = 0.0;
60    for (i, r) in table.iter().enumerate() {
61        for (j, &o) in r.iter().enumerate() {
62            let e = row_sums[i] * col_sums[j] / total;
63            if e > 0.0 {
64                stat += (o - e).powi(2) / e;
65            }
66        }
67    }
68    let dof = ((rows - 1) * (cols - 1)) as f64;
69    Some(ChiSquareResult {
70        statistic: stat,
71        p_value: chi_squared::upper_p(stat, dof),
72        dof,
73    })
74}
75
76#[cfg(test)]
77mod tests {
78    use super::*;
79
80    #[test]
81    fn gof_fair_die_is_not_rejected() {
82        // 60 rolls of a fair die, observed close to expected 10 each.
83        let observed = [9.0, 11.0, 10.0, 12.0, 8.0, 10.0];
84        let expected = [10.0; 6];
85        let r = chi_square_gof(&observed, &expected).unwrap();
86        assert_eq!(r.dof, 5.0);
87        assert!(
88            r.p_value > 0.5,
89            "fair die should not be rejected: p={}",
90            r.p_value
91        );
92    }
93
94    #[test]
95    fn gof_loaded_die_is_rejected() {
96        let observed = [5.0, 5.0, 5.0, 5.0, 5.0, 35.0]; // sixes way over
97        let expected = [10.0; 6];
98        let r = chi_square_gof(&observed, &expected).unwrap();
99        assert!(r.statistic > 11.07, "statistic {}", r.statistic); // > χ²_{0.95,5}
100        assert!(r.p_value < 0.001);
101    }
102
103    #[test]
104    fn independence_known_example() {
105        // 2×2 table with a clear association.
106        let table: [&[f64]; 2] = [&[90.0, 60.0], &[30.0, 120.0]];
107        let r = chi_square_independence(&table).unwrap();
108        assert_eq!(r.dof, 1.0);
109        // Strong association → large statistic, tiny p.
110        assert!(r.statistic > 30.0, "statistic {}", r.statistic);
111        assert!(r.p_value < 1e-6);
112    }
113
114    #[test]
115    fn independence_of_independent_table() {
116        // Rows proportional → no association → small statistic, large p.
117        let table: [&[f64]; 2] = [&[10.0, 20.0], &[20.0, 40.0]];
118        let r = chi_square_independence(&table).unwrap();
119        assert!(r.statistic < 1e-9);
120        assert!((r.p_value - 1.0).abs() < 1e-6);
121    }
122
123    #[test]
124    fn guards_bad_shapes() {
125        assert!(chi_square_gof(&[1.0], &[1.0]).is_none());
126        assert!(chi_square_gof(&[1.0, 2.0], &[1.0, 0.0]).is_none()); // zero expected
127        let ragged: [&[f64]; 2] = [&[1.0, 2.0], &[3.0]];
128        assert!(chi_square_independence(&ragged).is_none());
129    }
130}