Skip to main content

qualia_core_db/solvers/learning/regression/
pcr.rs

1//! Principal Components Regression (ISL ch 6.3.1) — regress the response on the
2//! first `m` principal components of the predictors. Reuses
3//! `dimensionality::pca` for the projection and `regression::linear` for the OLS
4//! (no duplicated math): it is literally PCA followed by least squares on the
5//! component scores, which tames collinearity by discarding low-variance directions.
6
7use crate::solvers::learning::dimensionality::pca::{self, Pca};
8use crate::solvers::learning::regression::linear::{self, LinearModel};
9use crate::solvers::learning::LearningError;
10
11/// A fitted PCR model: the PCA projection plus an OLS fit on the component scores.
12#[derive(Debug, Clone)]
13pub struct PcrModel {
14    pca: Pca,
15    ols: LinearModel,
16    n_components: usize,
17    p: usize,
18}
19
20impl PcrModel {
21    /// Fit PCR with `n_components` principal components (clamped to `p`). Fails
22    /// closed via the PCA / OLS solvers.
23    pub fn fit(
24        x: &[f64],
25        y: &[f64],
26        n: usize,
27        p: usize,
28        n_components: usize,
29    ) -> Result<Self, LearningError> {
30        if n == 0 || p == 0 || x.len() != n * p || y.len() != n {
31            return Err(LearningError::InvalidDimension);
32        }
33        let m = n_components.clamp(1, p);
34        let pca = pca::fit(x, n, p)?;
35        let scores = pca
36            .transform(x, n, m)
37            .ok_or(LearningError::InvalidDimension)?;
38        let ols = linear::fit(&scores, y, n, m, true)?;
39        Ok(Self {
40            pca,
41            ols,
42            n_components: m,
43            p,
44        })
45    }
46
47    /// Predict for one predictor row (length `p`).
48    pub fn predict_row(&self, x_row: &[f64]) -> f64 {
49        // Project the (single) row onto the components, then apply the OLS fit.
50        let scores = self
51            .pca
52            .transform(x_row, 1, self.n_components)
53            .unwrap_or_default();
54        self.ols.predict_row(&scores)
55    }
56
57    pub fn predict(&self, x: &[f64], m: usize) -> Vec<f64> {
58        (0..m)
59            .map(|i| self.predict_row(&x[i * self.p..(i + 1) * self.p]))
60            .collect()
61    }
62
63    pub fn n_components(&self) -> usize {
64        self.n_components
65    }
66}
67
68#[cfg(test)]
69mod tests {
70    use super::*;
71    use crate::solvers::learning::metrics::regression::r2_score;
72
73    #[test]
74    fn full_components_matches_ols_fit_quality() {
75        // With all components retained, PCR explains the data as well as OLS.
76        let x = [1.0, 2.0, 2.0, 1.0, 3.0, 0.0, 4.0, 5.0, 5.0, 4.0, 6.0, 1.0];
77        let y = [3.0, 5.0, 4.0, 9.0, 13.0, 8.0];
78        let pcr = PcrModel::fit(&x, &y, 6, 2, 2).unwrap();
79        let preds = pcr.predict(&x, 6);
80        assert!(r2_score(&y, &preds).unwrap() > 0.5);
81    }
82
83    #[test]
84    fn one_component_captures_dominant_direction() {
85        // y depends mostly on the high-variance direction; 1 component suffices.
86        let n = 20;
87        let mut x = vec![0.0; n * 2];
88        let mut y = vec![0.0; n];
89        for i in 0..n {
90            let t = i as f64;
91            x[i * 2] = t; // high variance
92            x[i * 2 + 1] = 0.01 * ((i % 3) as f64); // tiny variance
93            y[i] = 2.0 * t + 1.0;
94        }
95        let pcr = PcrModel::fit(&x, &y, n, 2, 1).unwrap();
96        assert_eq!(pcr.n_components(), 1);
97        let preds = pcr.predict(&x, n);
98        assert!(r2_score(&y, &preds).unwrap() > 0.99);
99    }
100
101    #[test]
102    fn guards() {
103        assert_eq!(
104            PcrModel::fit(&[1.0, 2.0, 3.0], &[1.0, 2.0], 2, 2, 1).unwrap_err(),
105            LearningError::InvalidDimension
106        );
107    }
108}