| 134 | } |
| 135 | |
| 136 | pub fn _get_betas( |
| 137 | _x: &[Vec<f64>], |
| 138 | _y: &[Vec<f64>], |
| 139 | ) -> StructuralBreakResult<(Vec<f64>, Vec<Vec<f64>>)> { |
| 140 | let x_matrix = to_matrix(_x)?; |
| 141 | let y_matrix = to_matrix(_y)?; |
| 142 | |
| 143 | let rows = x_matrix.nrows(); |
| 144 | let cols = x_matrix.ncols(); |
| 145 | let y_cols = y_matrix.ncols(); |
| 146 | |
| 147 | let xy = x_matrix.transpose() * &y_matrix; |
| 148 | let xx = x_matrix.transpose() * &x_matrix; |
| 149 | |
| 150 | let Some(xx_inv) = xx.try_inverse() else { |
| 151 | let b_mean = vec![f64::NAN; cols]; |
| 152 | let b_var = vec![vec![f64::NAN; cols]; cols]; |
| 153 | return Ok((b_mean, b_var)); |
| 154 | }; |
| 155 | |
| 156 | let b_mean = &xx_inv * xy; |
| 157 | let err = y_matrix - x_matrix * &b_mean; |
| 158 | let err_t_err = err.transpose() * err; |
| 159 | let denom = rows as f64 - cols as f64; |
| 160 | let scale = err_t_err / denom; |
| 161 | |
| 162 | let b_var_matrix = if y_cols == 1 { |
| 163 | let scalar = scale[(0, 0)]; |
| 164 | xx_inv * scalar |
| 165 | } else if scale.nrows() == cols && scale.ncols() == cols { |
| 166 | xx_inv.component_mul(&scale) |
| 167 | } else { |
| 168 | let scalar = scale[(0, 0)]; |
| 169 | xx_inv * scalar |
| 170 | }; |
| 171 | |
| 172 | let mut b_mean_vec = Vec::with_capacity(cols); |
| 173 | for i in 0..cols { |
| 174 | b_mean_vec.push(b_mean[(i, 0)]); |
| 175 | } |
| 176 | |
| 177 | Ok((b_mean_vec, matrix_to_vec(b_var_matrix))) |
| 178 | } |
| 179 | |
| 180 | fn get_y_x( |
| 181 | series: &[f64], |