fdars_core/alignment/
srsf.rs

1//! SRSF (Square-Root Slope Function) transforms and warping utilities.
2
3use crate::fdata::deriv_1d;
4use crate::helpers::{cumulative_trapz, linear_interp};
5use crate::matrix::FdMatrix;
6
7// ─── SRSF Transform and Inverse ─────────────────────────────────────────────
8
9/// Compute the Square-Root Slope Function (SRSF) transform.
10///
11/// For each curve f, the SRSF is: `q(t) = sign(f'(t)) * sqrt(|f'(t)|)`
12///
13/// # Arguments
14/// * `data` — Functional data matrix (n × m)
15/// * `argvals` — Evaluation points (length m)
16///
17/// # Returns
18/// FdMatrix of SRSFs with the same shape as input.
19///
20/// # Examples
21///
22/// ```
23/// use fdars_core::matrix::FdMatrix;
24/// use fdars_core::alignment::srsf_transform;
25///
26/// let argvals: Vec<f64> = (0..20).map(|i| i as f64 / 19.0).collect();
27/// let data = FdMatrix::from_column_major(
28///     argvals.iter().map(|&t| (t * 6.0).sin()).collect(),
29///     1, 20,
30/// ).unwrap();
31/// let srsf = srsf_transform(&data, &argvals);
32/// assert_eq!(srsf.shape(), (1, 20));
33/// ```
34#[must_use = "expensive computation whose result should not be discarded"]
35pub fn srsf_transform(data: &FdMatrix, argvals: &[f64]) -> FdMatrix {
36    let (n, m) = data.shape();
37    if n == 0 || m == 0 || argvals.len() != m {
38        return FdMatrix::zeros(n, m);
39    }
40
41    let deriv = deriv_1d(data, argvals, 1);
42
43    let mut result = FdMatrix::zeros(n, m);
44    for i in 0..n {
45        for j in 0..m {
46            let d = deriv[(i, j)];
47            result[(i, j)] = d.signum() * d.abs().sqrt();
48        }
49    }
50    result
51}
52
53/// Reconstruct a curve from its SRSF representation.
54///
55/// Given SRSF q and initial value f0, reconstructs: `f(t) = f0 + ∫₀ᵗ q(s)|q(s)| ds`
56///
57/// # Arguments
58/// * `q` — SRSF values (length m)
59/// * `argvals` — Evaluation points (length m)
60/// * `f0` — Initial value f(argvals\[0\])
61///
62/// # Returns
63/// Reconstructed curve values.
64pub fn srsf_inverse(q: &[f64], argvals: &[f64], f0: f64) -> Vec<f64> {
65    let m = q.len();
66    if m == 0 {
67        return Vec::new();
68    }
69
70    // Integrand: q(s) * |q(s)|
71    let integrand: Vec<f64> = q.iter().map(|&qi| qi * qi.abs()).collect();
72    let integral = cumulative_trapz(&integrand, argvals);
73
74    integral.iter().map(|&v| f0 + v).collect()
75}
76
77// ─── Reparameterization ─────────────────────────────────────────────────────
78
79/// Reparameterize a curve by a warping function.
80///
81/// Computes `f(γ(t))` via linear interpolation.
82///
83/// # Arguments
84/// * `f` — Curve values (length m)
85/// * `argvals` — Evaluation points (length m)
86/// * `gamma` — Warping function values (length m)
87pub fn reparameterize_curve(f: &[f64], argvals: &[f64], gamma: &[f64]) -> Vec<f64> {
88    gamma
89        .iter()
90        .map(|&g| linear_interp(argvals, f, g))
91        .collect()
92}
93
94/// Compose two warping functions: `(γ₁ ∘ γ₂)(t) = γ₁(γ₂(t))`.
95///
96/// # Arguments
97/// * `gamma1` — Outer warping function (length m)
98/// * `gamma2` — Inner warping function (length m)
99/// * `argvals` — Evaluation points (length m)
100pub fn compose_warps(gamma1: &[f64], gamma2: &[f64], argvals: &[f64]) -> Vec<f64> {
101    gamma2
102        .iter()
103        .map(|&g| linear_interp(argvals, gamma1, g))
104        .collect()
105}
106
107/// Compute a single SRSF from a slice (single-row convenience).
108pub(crate) fn srsf_single(f: &[f64], argvals: &[f64]) -> Vec<f64> {
109    let m = f.len();
110    let mat = FdMatrix::from_slice(f, 1, m).expect("dimension invariant: data.len() == n * m");
111    let q_mat = srsf_transform(&mat, argvals);
112    q_mat.row(0)
113}
fdars_core/alignment/srsf.rs

fdars_core/alignment/
srsf.rs
