pyhctsa.operations.stationarity.slow_feature_analysis

pyhctsa.operations.stationarity.slow_feature_analysis(y, min_window_length=20)

Multi-scale slow feature analysis of windowed statistics: a calibrated test for slow drift.

At each scale, the time series is split into W non-overlapping windows and five statistics are computed in every window: mean, variance, skewness, lag-1 autocorrelation (AC1), and lag-1 trev, a normalized third-order moment of the increments measuring time-reversal asymmetry, mean(dz^3) / mean(dz^2)^(3/2) with dz = diff(window) (AC1 is symmetric under time reversal, so trev lets a slow drift in the irreversibility of the dynamics show up). Each statistic is replaced by the normal scores of its ranks across the W windows, and Slow Feature Analysis (SFA) [1] finds the linear combination of the five that varies most slowly across the sequence of windows: after PCA-whitening, the eigenvector of the covariance of the first differences with the smallest eigenvalue, eta1 (about 2 for windows in random order, near 0 for a smooth trend). Non-overlapping windows are used because overlap would make adjacent window statistics correlated and eta1 spuriously small.

The scales are W = 10, 20, 40, ..., 1280 windows, using every W with at least min_window_length samples per window, so the feature is defined from N = 10 * min_window_length (200 by default). Few, long windows detect slow, smooth drift; many, short windows are needed for a slow modulation with several cycles in the record (about eight windows per cycle).

Each scale’s eta1 is converted to a one-sided p-value against a Monte Carlo null for W windows with exchangeable (iid) rows of statistics (scripts/gen_sfa_null_table.py). Because the statistics enter as rank normal scores, this null does not depend on their marginal distributions, so it holds for heavy-tailed series as well as Gaussian ones. Serial correlation in the series correlates the statistics of adjacent windows, which would make eta1 smaller than under the null. To stay calibrated, the first difference of each statistic is divided by sqrt(1 - rho) before the SFA eigenvalue is computed, where rho is the correlation between that statistic in adjacent windows of a noise model fitted to the series. The noise model is an ARFIMA(p, d, 0) process, fitted by Whittle’s method to the pooled periodogram of ten windows, each demeaned and scaled to unit variance and without its lowest frequency, so that a slow drift is not mistaken for noise; white noise gives rho = 0. For the window mean, rho follows from the model’s autocorrelation function (including long memory, d > 0, which correlates the means of windows of every length); for the other four statistics it is the adjacent-window correlation of the window variance of the model’s short-memory part (Gaussian approximation), which captures narrowband (oscillatory) noise, whose amplitude varies slowly. Neither rho is allowed below the adjacent-window correlation of means of an AR(1) process with the series’ (bias-corrected) lag-1 autocorrelation.

The scales are combined by the Sidak-corrected minimum p-value. In simulations (2,000 series of each process and length, N = 200, 400, 800, 1600 and 4000), the false-positive rate stayed within sampling error of the nominal level or below it for white noise, AR(1) with coefficient 0.8, 0.9 or 0.95, AR(2) with coefficients (1.0, -0.5) or (1.9, -0.95), and Student-t noise with 2 or 3 degrees of freedom or Cauchy noise. It was at most 0.08 at the 0.05 level and 0.035 at the 0.01 level for ARFIMA(0, d, 0) noise with d = 0.1, 0.25 or 0.4, fractional Gaussian noise with Hurst exponent 0.9, and narrowband AR(2) noise with coefficients (1.6, -0.9), the excess growing with N. Long memory and slow drift cannot be fully told apart in a single record.

This operation is a compact variant of feature-based SFA (f-SFA) [2]_, which applies SFA to sliding-window catch22/catch24 features; here five fixed window statistics are used.

Cost: about 1 ms at N = 200 (one scale), 8 ms at N = 4000 (five scales) and 17 ms at N = 20000 (eight scales).

References

invariances.” Neural Computation 14(4), 715-770 (2002). .. [2] Owens, K.S., Tamaki, M. & Fulcher, B.D. “Parameter inference from a non-stationary unknown process using statistical feature-based slow feature analysis”, arXiv:2609.01651 (2026). The full feature-based method (f-SFA: sliding-window catch22/catch24 features + SFA) is implemented in the Python package fsfa: https://github.com/KieranOwens/fsfa

Parameters:
y : array-like

The input time series.

min_window_length : int, optional

The fewest samples a window may have for its scale to be used (at least 10). Window skewness, AC1 and trev are unreliable in very short windows. Default is 20.

Returns:

  • driftScore: -log10 of the Sidak-corrected minimum p-value over scales: 0 to about 1.3 for a stationary series, above 2 (p < 0.01) for clear slow variation of the window statistics.

  • driftWinLength: the window length (samples) at the reported scale: among the scales that are significant on their own (Sidak p < 0.05), the one where eta1 is smallest relative to its null median (the strongest slow structure); else the scale with the smallest p-value. Long windows for a smooth drift, about an eighth of the period or shorter for a modulation.

  • slowPeriod: the period (samples) of the slow component at that scale, from the peak of the periodogram of the slowest SFA component across windows (refined between frequency bins). A monotone drift gives the record length. NaN unless driftScore passes p < 0.05: the period of nothing significant cannot be estimated.

At the reported scale:

  • etaEnd: the largest (fastest) SFA eigenvalue.

  • etaStd: the standard deviation of the five SFA eigenvalues (spread of the slowness spectrum).

  • pc1VarFrac: the fraction of the variance of the five (rank-scored) statistics on the leading principal component.

  • slowPCA1corr: the absolute correlation between the slowest SFA component and the leading principal component: near 1 when the slow direction is just the dominant (highest-variance) direction, near 0 when SFA isolates a separate, low-variance slow mode.

  • slowCorrMean, slowCorrVar, slowCorrSkew, slowCorrAC1, slowCorrTrev: the absolute correlation between each (rank-scored) window statistic and the slowest SFA component, indicating which statistic drives the slow mode.

Returns NaN if the time series is shorter than 10 * min_window_length, or if no scale has defined, non-degenerate window statistics (e.g. a constant series).

Return type:

dict