Skip to main content

fdars_core/detrend/
diff.rs

1use crate::iter_maybe_parallel;
2use crate::matrix::FdMatrix;
3#[cfg(feature = "parallel")]
4use rayon::iter::ParallelIterator;
5use std::borrow::Cow;
6
7use super::{reassemble_polynomial_results, TrendResult};
8
9pub(super) fn diff_single_curve(
10    curve: &[f64],
11    m: usize,
12    order: usize,
13) -> (Vec<f64>, Vec<f64>, Vec<f64>, f64) {
14    let diff1: Vec<f64> = (0..m - 1).map(|j| curve[j + 1] - curve[j]).collect();
15    let detrended = if order == 2 {
16        (0..diff1.len() - 1)
17            .map(|j| diff1[j + 1] - diff1[j])
18            .collect()
19    } else {
20        diff1.clone()
21    };
22    let initial_values = if order == 2 {
23        vec![curve[0], curve[1]]
24    } else {
25        vec![curve[0]]
26    };
27    let rss: f64 = detrended.iter().map(|&x| x.powi(2)).sum();
28    let new_m = m - order;
29    let mut trend = vec![0.0; m];
30    trend[0] = curve[0];
31    if order == 1 {
32        for j in 1..m {
33            trend[j] = curve[j] - if j <= new_m { detrended[j - 1] } else { 0.0 };
34        }
35    } else {
36        trend = curve.to_vec();
37    }
38    let mut det_full = vec![0.0; m];
39    det_full[..new_m].copy_from_slice(&detrended[..new_m]);
40    (trend, det_full, initial_values, rss)
41}
42
43/// Remove trend by differencing.
44#[must_use = "returns the detrending result which should not be discarded"]
45pub fn detrend_diff(data: &FdMatrix, order: usize) -> TrendResult {
46    let (n, m) = data.shape();
47    if n == 0 || m <= order || order == 0 || order > 2 {
48        return TrendResult::empty(data, n, m, Cow::Owned(format!("diff{order}")), order);
49    }
50    let results: Vec<(Vec<f64>, Vec<f64>, Vec<f64>, f64)> = iter_maybe_parallel!(0..n)
51        .map(|i| {
52            let curve: Vec<f64> = (0..m).map(|j| data[(i, j)]).collect();
53            diff_single_curve(&curve, m, order)
54        })
55        .collect();
56    let (trend, detrended, coefficients, rss) = reassemble_polynomial_results(results, n, m, order);
57    TrendResult {
58        trend,
59        detrended,
60        method: Cow::Owned(format!("diff{order}")),
61        coefficients: Some(coefficients),
62        rss,
63        n_params: order,
64    }
65}