Skip to main content

regression_diagnostics/residuals/
white_test.rs

1use super::{lm_heteroskedasticity_test, non_intercept_columns};
2use crate::OlsFit;
3
4/// Result of White's heteroskedasticity test (LM form).
5#[derive(Debug, Clone, Copy, PartialEq)]
6pub struct WhiteTest {
7    /// Lagrange-multiplier statistic `n·R²_aux`.
8    pub statistic: f64,
9    /// Degrees of freedom (number of auxiliary regressors).
10    pub df: usize,
11    /// Upper-tail p-value under χ²(df).
12    pub p_value: f64,
13}
14
15/// White's test for heteroskedasticity.
16///
17/// Like [`breusch_pagan`](super::breusch_pagan) it regresses the squared
18/// residuals on auxiliary regressors, but the auxiliary set is richer: the
19/// original regressors, their **squares**, and all **pairwise cross-products**.
20/// This makes White's test sensitive to more general (including nonlinear) forms
21/// of heteroskedasticity and, because it includes cross-products, it doubles as a
22/// specification check.
23///
24/// # Tradeoff versus Breusch-Pagan
25///
26/// The extra terms cost degrees of freedom. In **small samples** that spreads the
27/// test's power thin, so White's can fail to flag heteroskedasticity that the
28/// tighter Breusch-Pagan catches. White's is the more general test; Breusch-Pagan
29/// has more power when you already suspect a linear variance relationship. The
30/// crate offers both rather than picking one for you.
31///
32/// With `k` non-intercept regressors the auxiliary regression has `k` linear
33/// terms, `k` squares, and `k(k−1)/2` cross-products; duplicate columns that can
34/// arise (e.g. the square of a 0/1 dummy equals the dummy) are dropped to keep
35/// the auxiliary design full rank.
36///
37/// Because that regressor count grows quadratically in `k`, small samples can
38/// leave the auxiliary regression without enough residual degrees of freedom to
39/// estimate. When there are too few observations the statistic and p-value are
40/// returned as `NaN` rather than a fabricated number.
41pub fn white_test(fit: &OlsFit) -> WhiteTest {
42    let n = fit.n_observations();
43    let cols = non_intercept_columns(fit);
44    let x = fit.design_matrix();
45
46    let base: Vec<Vec<f64>> = cols.iter().map(|&j| x.column(j).to_vec()).collect();
47    let k = base.len();
48
49    let mut regressors: Vec<Vec<f64>> = Vec::new();
50    // Linear terms.
51    for col in &base {
52        regressors.push(col.clone());
53    }
54    // Squares.
55    for col in &base {
56        regressors.push(col.iter().map(|v| v * v).collect());
57    }
58    // Pairwise cross-products.
59    for a in 0..k {
60        for b in (a + 1)..k {
61            let cross: Vec<f64> = base[a]
62                .iter()
63                .zip(base[b].iter())
64                .map(|(u, v)| u * v)
65                .collect();
66            regressors.push(cross);
67        }
68    }
69
70    // Drop near-duplicate / near-constant auxiliary columns to protect the
71    // auxiliary design's rank (e.g. the square of a 0/1 dummy).
72    let regressors = dedup_regressors(n, regressors);
73
74    let resid_sq: Vec<f64> = fit.residuals().iter().map(|e| e * e).collect();
75
76    let (statistic, df, p_value) = lm_heteroskedasticity_test(n, &regressors, &resid_sq);
77    WhiteTest {
78        statistic,
79        df,
80        p_value,
81    }
82}
83
84/// Remove constant columns and columns that duplicate an earlier one (up to a
85/// tight relative tolerance), which would otherwise make the auxiliary design
86/// rank-deficient.
87fn dedup_regressors(n: usize, cols: Vec<Vec<f64>>) -> Vec<Vec<f64>> {
88    let mut kept: Vec<Vec<f64>> = Vec::new();
89    for col in cols {
90        let first = col[0];
91        let is_constant = col
92            .iter()
93            .all(|&v| (v - first).abs() <= 1e-12 * first.abs().max(1.0));
94        if is_constant {
95            continue;
96        }
97        let dup = kept
98            .iter()
99            .any(|k| (0..n).all(|i| (k[i] - col[i]).abs() <= 1e-12 * col[i].abs().max(1.0)));
100        if !dup {
101            kept.push(col);
102        }
103    }
104    kept
105}