Skip to main content

r2rs_stats/tests/
linear_estimator.rs

1// "Whatever you do, work at it with all your heart, as working for the Lord,
2// not for human masters, since you know that you will receive an inheritance
3// from the Lord as a reward. It is the Lord Christ you are serving."
4// (Col 3:23-24)
5
6use std::{
7    fmt::{Debug, Formatter},
8    sync::Arc,
9};
10
11use nalgebra::DMatrix;
12use strafe_trait::{Assumption, Manipulator, StatisticalEstimate, StatisticalEstimator};
13use strafe_type::{Alpha64, ModelMatrix, Rational64};
14
15use crate::funcs::{confidence_interval, predict, prediction_interval};
16
17#[derive(Copy, Clone, Debug)]
18pub struct LinearEstimator {
19    alpha: Alpha64,
20    scale: Option<f64>,
21    df: Option<Rational64>,
22}
23
24impl StatisticalEstimator for LinearEstimator {
25    type Input = (ModelMatrix, DMatrix<f64>, DMatrix<f64>, DMatrix<f64>);
26    type Output = Result<LinearEstimate, Box<dyn std::error::Error>>;
27
28    fn assumptions() -> Vec<Box<dyn Assumption>> {
29        todo!()
30    }
31
32    fn estimate(&self, (x, y, b, w): &Self::Input) -> Self::Output {
33        let b_temp = b.clone();
34        let x = x.matrix();
35
36        Ok(LinearEstimate {
37            predictor_estimate: Manipulator::new(Arc::new(move |x: &ModelMatrix| {
38                ModelMatrix::from(predict(&x.matrix().insert_column(0, 1.0), &b_temp))
39            })),
40            predictor_confidence_interval: {
41                let (lower, upper) =
42                    prediction_interval(&x, &y, &b, &w, self.alpha, self.scale, self.df);
43                let mut pi: ModelMatrix = DMatrix::<f64>::from_iterator(
44                    x.nrows(),
45                    2,
46                    lower.into_iter().chain(upper.into_iter()).cloned(),
47                )
48                .into();
49                pi.set_name_index(0, "Lower Prediction Interval");
50                pi.set_name_index(1, "Upper Prediction Interval");
51
52                pi
53            },
54            predictions_estimate: ModelMatrix::from(predict(&x, b)),
55            predictions_confidence_interval: {
56                let (lower, upper) =
57                    confidence_interval(&x, &y, &b, &w, self.alpha, self.scale, self.df);
58                let mut pi: ModelMatrix = DMatrix::<f64>::from_iterator(
59                    x.nrows(),
60                    2,
61                    lower.into_iter().chain(upper.into_iter()).cloned(),
62                )
63                .into();
64                pi.set_name_index(0, "Lower Confidence Interval");
65                pi.set_name_index(1, "Upper Confidence Interval");
66
67                pi
68            },
69            alpha: self.alpha,
70        })
71    }
72}
73
74impl Default for LinearEstimator {
75    fn default() -> Self {
76        Self {
77            alpha: 0.05.into(),
78            scale: None,
79            df: None,
80        }
81    }
82}
83
84impl LinearEstimator {
85    pub fn new() -> Self {
86        Self::default()
87    }
88
89    pub fn with_alpha<A: Into<Alpha64>>(self, alpha: A) -> Self {
90        Self {
91            alpha: alpha.into(),
92            ..self
93        }
94    }
95
96    pub fn with_scale(self, scale: f64) -> Self {
97        Self {
98            scale: Some(scale),
99            ..self
100        }
101    }
102
103    pub fn with_df<R: Into<Rational64>>(self, df: R) -> Self {
104        Self {
105            df: Some(df.into()),
106            ..self
107        }
108    }
109}
110
111#[derive(Clone)]
112pub struct LinearEstimate {
113    pub predictor_estimate: Manipulator<ModelMatrix, ModelMatrix>,
114    pub predictor_confidence_interval: ModelMatrix,
115    pub predictions_estimate: ModelMatrix,
116    pub predictions_confidence_interval: ModelMatrix,
117    alpha: Alpha64,
118}
119
120impl Debug for LinearEstimate {
121    fn fmt(&self, f: &mut Formatter<'_>) -> std::fmt::Result {
122        f.debug_struct("RankLinearEstimate")
123            // .field("predictor_estimate", &self.predictor_estimate)
124            .field(
125                "predictor_confidence_interval",
126                &self.predictor_confidence_interval,
127            )
128            .field("predictions_estimate", &self.predictions_estimate)
129            .field(
130                "predictions_confidence_interval",
131                &self.predictions_confidence_interval,
132            )
133            .field("alpha", &self.alpha)
134            .finish()
135    }
136}
137
138impl StatisticalEstimate<Manipulator<ModelMatrix, ModelMatrix>, ModelMatrix> for LinearEstimate {
139    fn alpha(&self) -> Alpha64 {
140        self.alpha
141    }
142
143    fn estimate(&self) -> Manipulator<ModelMatrix, ModelMatrix> {
144        self.predictor_estimate.clone()
145    }
146
147    fn confidence_interval(&self) -> ModelMatrix {
148        self.predictor_confidence_interval.clone()
149    }
150}
151
152impl StatisticalEstimate<ModelMatrix, ModelMatrix> for LinearEstimate {
153    fn alpha(&self) -> Alpha64 {
154        self.alpha
155    }
156
157    fn estimate(&self) -> ModelMatrix {
158        self.predictions_estimate.clone()
159    }
160
161    fn confidence_interval(&self) -> ModelMatrix {
162        self.predictions_confidence_interval.clone()
163    }
164}