fdars_core/detrend/
linear.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::TrendResult;
8
9/// Remove linear trend from functional data using least squares.
10///
11/// # Arguments
12/// * `data` - Matrix (n x m): n samples, m evaluation points
13/// * `argvals` - Time/argument values of length m
14///
15/// # Returns
16/// TrendResult with trend, detrended data, and coefficients (intercept, slope)
17#[must_use = "returns the detrending result which should not be discarded"]
18pub fn detrend_linear(data: &FdMatrix, argvals: &[f64]) -> TrendResult {
19    let (n, m) = data.shape();
20    if n == 0 || m < 2 || argvals.len() != m {
21        return TrendResult::empty(data, n, m, Cow::Borrowed("linear"), 2);
22    }
23
24    let mean_t: f64 = argvals.iter().sum::<f64>() / m as f64;
25    let ss_t: f64 = argvals.iter().map(|&t| (t - mean_t).powi(2)).sum();
26
27    // Compute only scalar coefficients in parallel (no Vec allocations per sample)
28    let scalars: Vec<(f64, f64, f64)> = iter_maybe_parallel!(0..n)
29        .map(|i| {
30            let mut sum_y = 0.0;
31            for j in 0..m {
32                sum_y += data[(i, j)];
33            }
34            let mean_y = sum_y / m as f64;
35            let mut sp = 0.0;
36            for j in 0..m {
37                sp += (argvals[j] - mean_t) * (data[(i, j)] - mean_y);
38            }
39            let slope = if ss_t.abs() > 1e-15 { sp / ss_t } else { 0.0 };
40            let intercept = mean_y - slope * mean_t;
41            let mut rss = 0.0;
42            for j in 0..m {
43                let residual = data[(i, j)] - (intercept + slope * argvals[j]);
44                rss += residual * residual;
45            }
46            (intercept, slope, rss)
47        })
48        .collect();
49
50    // Write directly into pre-allocated output matrices
51    let mut trend = FdMatrix::zeros(n, m);
52    let mut detrended = FdMatrix::zeros(n, m);
53    let mut coefficients = FdMatrix::zeros(n, 2);
54    let mut rss = vec![0.0; n];
55
56    for (i, &(intercept, slope, r)) in scalars.iter().enumerate() {
57        coefficients[(i, 0)] = intercept;
58        coefficients[(i, 1)] = slope;
59        rss[i] = r;
60        for j in 0..m {
61            let trend_val = intercept + slope * argvals[j];
62            trend[(i, j)] = trend_val;
63            detrended[(i, j)] = data[(i, j)] - trend_val;
64        }
65    }
66
67    TrendResult {
68        trend,
69        detrended,
70        method: Cow::Borrowed("linear"),
71        coefficients: Some(coefficients),
72        rss,
73        n_params: 2,
74    }
75}
fdars_core/detrend/linear.rs

fdars_core/detrend/
linear.rs
