RQA2

The RQA2 module provides a modern, object-oriented interface for Recurrence Quantification Analysis (RQA). It supersedes the procedural functions in SMdRQA.RQA_functions and adds surrogate-data testing and chaotic-system simulation in a single coherent package.

RQA2

class RQA2(data=None, normalize=True, **kwargs)

Bases: object

Comprehensive Recurrence Quantification Analysis class that handles all RQA computations, visualizations, and batch processing in an object-oriented manner.

Fixed version with proper 0-based indexing throughout.

classmethod batch_process(input_path, output_path, group_level=False, group_level_estimates=None, **kwargs)

Analyse all .npy files in input_path and write results to output_path.

Two CSV files are written: rqa_results.csv (one row per file with all RQA measures) and error_report.csv (files that could not be processed). Recurrence plots are saved as .npy arrays alongside the CSVs.

Parameters:

• input_path (str or path-like) – Directory containing .npy time-series files.

• output_path (str or path-like) – Directory where outputs are written (created if absent).

• group_level (bool, default False) – If True, replace per-file parameters with group averages for the parameters listed in group_level_estimates.

• group_level_estimates (list of {‘tau’, ‘m’, ‘eps’}, optional) – Which parameters to estimate at the group level.

• **kwargs – Passed to RQA2.__init__ for each file.

Returns:

• results (list of dict) – Per-file dictionaries containing file name, RQA parameters, and measure values.

• error_files (list of dict) – Per-file dictionaries for failed files, containing 'file' and 'error' keys.

compute_embedded_signal()

Build the time-delay embedding tensor.

Constructs delay vectors of the form [x(t), x(t+τ), …, x(t+(m-1)τ)] for each valid time index t.

Returns:

ndarray, shape (N_embedded, m, D) – Delay-embedded signal where N_embedded = N - (m-1)*τ, m is the embedding dimension, and D is the number of original signal dimensions.

Raises:

ValueError – If no data have been loaded or if the data are too short for the chosen (τ, m) combination.

compute_embedding_dimension()

Estimate the optimal embedding dimension m via False Nearest Neighbours.

Uses the FNN criterion: the embedding dimension is the smallest m for which the fraction of false nearest neighbours drops below self.config['bound'].

Returns:

int – Embedding dimension m ≥ 1.

Raises:

ValueError – If no data have been loaded.

compute_neighborhood_radius(reqrr=None)

Find the neighbourhood radius ε that achieves a target recurrence rate.

The search scans self.config['epsdiv'] evenly-spaced candidate values between epsmin and epsmax and returns the first ε for which |RR − reqrr| < rr_delta.

Parameters:

reqrr (float, optional) – Target recurrence rate in (0, 1). Clamped to [0.01, 0.99]. Defaults to self.config['reqrr'].

Returns:

float – Neighbourhood radius ε > 0.

Raises:

ValueError – If no data have been loaded.

compute_recurrence_plot()

Build the binary recurrence plot matrix from the current parameters.

Uses the delay-embedded signal with the stored (or lazily computed) τ, m, and ε. Two delay vectors are considered recurrent when their Euclidean distance is less than ε.

Returns:

ndarray of int, shape (N_embedded, N_embedded) – Symmetric binary matrix; 1 indicates recurrence, 0 otherwise.

Raises:

ValueError – If no data have been loaded or if the embedding parameters yield an empty embedded signal.

compute_rqa_measures(lmin=None)

Compute the full set of RQA measures from the recurrence plot.

Parameters:

lmin (int, optional) – Minimum line length for DET and LAM computation. Defaults to self.config['lmin'] (typically 2).

Returns:

dict – Dictionary with the following keys:

recurrence_ratefloat

Fraction of recurrent points (RR).

determinismfloat

Fraction of recurrent points on diagonal lines ≥ lmin (DET).

laminarityfloat

Fraction of recurrent points on vertical lines ≥ lmin (LAM).

diagonal_entropyfloat

Shannon entropy of diagonal line-length distribution (L_entr).

vertical_entropyfloat

Shannon entropy of vertical line-length distribution (V_entr).

average_diagonal_lengthfloat

Mean length of diagonal lines ≥ lmin (L).

average_vertical_lengthfloat

Mean length of vertical lines ≥ lmin (TT – trapping time).

max_diagonal_lengthfloat

Maximum diagonal line length (L_max).

max_vertical_lengthfloat

Maximum vertical line length (V_max).

diagonal_modefloat

Mode of the diagonal line-length distribution.

vertical_modefloat

Mode of the vertical line-length distribution.

compute_time_delay(method=None, mi_method=None)

Estimate the optimal time delay τ from the mutual-information curve.

Parameters:

• method ({‘default’, ‘polynomial’}, optional) – Algorithm used to locate the minimum of the MI curve. 'default' selects the first local minimum; 'polynomial' fits a cross-validated polynomial and uses its first minimum. Defaults to self.config['tau_method'].

• mi_method ({‘histdd’, ‘avg’}, optional) – Mutual-information estimator. Defaults to self.config['mi_method'].

Returns:

int – Optimal time delay τ ≥ 1.

Raises:

ValueError – If no data have been loaded.

determinism(lmin=None)

Return the determinism (DET) measure.

Parameters:

lmin (int, optional) – Minimum diagonal line length threshold.

Returns:

float – DET in [0, 1].

property embedded_signal

Time-delayed embedded signal.

property eps

Neighborhood radius.

get_summary()

Return a nested dictionary summarising data info, parameters, and (if already computed) RQA measures.

Returns:

dict – Top-level keys: 'Data Info', 'Parameters', and optionally 'RQA Measures'.

laminarity(lmin=None)

Return the laminarity (LAM) measure.

Parameters:

lmin (int, optional) – Minimum vertical line length threshold.

Returns:

float – LAM in [0, 1].

load_data(data, normalize=True)

Load (and optionally normalise) time-series data, resetting all caches.

Parameters:

• data (array_like) – Time-series array of shape (N,) or (N, D).

• normalize (bool, default True) – Apply z-score normalisation column-wise.

Returns:

None

load_results(filepath)

Restore RQA state from a pickle file previously written by save_results().

Parameters:

filepath (str or path-like) – Path to a .pkl file created by save_results().

Returns:

None

property m

Embedding dimension.

plot_fnn_curve(max_m=None, figsize=(10, 6), save_path=None)

Plot the false nearest-neighbours ratio as a function of embedding dimension m.

A vertical dashed line marks the automatically selected m value.

Parameters:

• max_m (int, optional) – Maximum embedding dimension to evaluate.

• figsize (tuple of float, default (10, 6)) – Figure dimensions in inches.

• save_path (str or path-like, optional) – Destination path for saving the figure at 300 dpi.

Returns:

None

plot_recurrence_plot(figsize=(8, 8), title=None, save_path=None)

Display the recurrence plot.

Parameters:

• figsize (tuple of float, default (8, 8)) – Width and height of the figure in inches.

• title (str, optional) – Figure title. Defaults to 'Recurrence Plot (RR=<value>)'.

• save_path (str or path-like, optional) – If given, the figure is saved to this path at 300 dpi before being displayed.

Returns:

None

plot_rqa_measures_summary(figsize=(12, 8), save_path=None)

Display a 2×3 panel summarising the main RQA measures and parameters.

Panels: main measures (RR, DET, LAM), entropy measures, average line lengths, maximum line lengths, mode line lengths, and RQA parameters (τ, m, ε).

Parameters:

• figsize (tuple of float, default (12, 8)) – Figure dimensions in inches.

• save_path (str or path-like, optional) – Destination path for saving the figure at 300 dpi.

Returns:

None

plot_tau_mi_curve(max_tau=None, figsize=(10, 6), save_path=None)

Plot the mutual information as a function of time delay τ.

A vertical dashed line marks the automatically selected τ value.

Parameters:

• max_tau (int, optional) – Maximum τ to evaluate. Defaults to min(100, N/4).

• figsize (tuple of float, default (10, 6)) – Figure dimensions in inches.

• save_path (str or path-like, optional) – Destination path for saving the figure at 300 dpi.

Returns:

None

plot_time_series(figsize=(12, 6), save_path=None)

Plot the original (unnormalised) time-series data.

For univariate data a single line plot is drawn; for multivariate data each of the first five dimensions is shown in a stacked subplot.

Parameters:

• figsize (tuple of float, default (12, 6)) – Figure dimensions in inches.

• save_path (str or path-like, optional) – Destination path for saving the figure at 300 dpi.

Returns:

None

Raises:

ValueError – If no data have been loaded.

property recurrence_plot

Recurrence plot matrix.

property recurrence_rate

Recurrence rate of the RP.

save_results(filepath)

Serialise all computed RQA state to a pickle file.

Saved fields: data, tau, m, eps, recurrence_plot, rqa_measures, config.

Parameters:

filepath (str or path-like) – Destination file path (e.g. 'analysis.pkl').

Returns:

None

property tau

Time delay parameter.

trapping_time(lmin=None)

Return the trapping time (TT) – mean vertical line length.

Parameters:

lmin (int, optional) – Minimum vertical line length threshold.

Returns:

float – Average vertical line length TT ≥ 0.

RQA2_simulators

class RQA2_simulators(seed: int | None = None)

Bases: object

Generators for chaotic dynamical systems used to test RQA2 and surrogate methods.

All ODE systems are integrated with scipy.integrate.solve_ivp using the RK45 solver at tight tolerances (rtol=1e-9, atol=1e-12). Discrete maps are iterated directly.

Parameters:

seed (int or None, optional) – Seed for the internal numpy.random.Generator. Use an integer for reproducible results.

Examples

>>> sim = RQA2_simulators(seed=42)
>>> x, y, z = sim.rossler(n=1000)
>>> battery = sim.generate_test_battery()

chua(tmax: int = 10000, n: int = 2000, Xi: Tuple[float, float, float] = (0.1, 0.1, 0.1), alpha: float = 15.6, beta: float = 28.0, m0: float = -1.143, m1: float = -0.714, dt: float = 0.01) → Tuple[ndarray, ndarray, ndarray]

Simulate Chua’s circuit attractor.

Uses the piecewise-linear Chua diode characteristic h(x) and the three-dimensional ODE:

dx/dt = α·(y - x - h(x))
dy/dt = x - y + z
dz/dt = -β·y

Parameters:

• tmax (int, default 10000) – Number of integration steps.

• n (int, default 2000) – Number of output points.

• Xi (tuple of float, default (0.1, 0.1, 0.1)) – Initial conditions (x₀, y₀, z₀).

• alpha, beta (float) – Circuit parameters.

• m0, m1 (float) – Slopes of the piecewise-linear Chua diode.

• dt (float, default 0.01) – Integration step size.

Returns:

x, y, z (ndarray, shape (n,)) – State variables at n equidistant times.

generate_test_battery() → dict

Generate a standard battery of chaotic and periodic test signals.

Returns:

dict – Keys: 'rossler_chaotic', 'rossler_sync', 'lorenz', 'henon', 'chua'. Each value is a dict with keys 'x', 'y', (optionally 'z'), and 'regime'.

henon(n: int = 2000, Xi: Tuple[float, float] = (0.1, 0.1), a: float = 1.4, b: float = 0.3) → Tuple[ndarray, ndarray]

Simulate the Hénon map.

Defined by the 2-D discrete-time recurrence:

x_{n+1} = 1 - a·x_n² + y_n
y_{n+1} = b·x_n

The classic chaotic attractor uses a=1.4, b=0.3.

Parameters:

• n (int, default 2000) – Number of iterations.

• Xi (tuple of float, default (0.1, 0.1)) – Initial conditions (x₀, y₀).

• a, b (float) – Map parameters.

Returns:

X, Y (ndarray, shape (n,)) – Trajectory of the x and y coordinates.

lorenz(tmax: int = 10000, n: int = 2000, Xi: Tuple[float, float, float] = (1.0, 1.0, 1.0), sigma: float = 10.0, rho: float = 28.0, beta: float = 2.6666666666666665, dt: float = 0.01) → Tuple[ndarray, ndarray, ndarray]

Simulate the Lorenz attractor.

The system is defined by:

dx/dt = σ·(y - x)
dy/dt = ρ·x - y - x·z
dz/dt = x·y - β·z

The classic chaotic butterfly attractor uses σ=10, ρ=28, β=8/3.

Parameters:

• tmax (int, default 10000) – Number of integration steps.

• n (int, default 2000) – Number of output points.

• Xi (tuple of float, default (1.0, 1.0, 1.0)) – Initial conditions (x₀, y₀, z₀).

• sigma, rho, beta (float) – System parameters.

• dt (float, default 0.01) – Integration step size.

Returns:

x, y, z (ndarray, shape (n,)) – State variables at n equidistant times.

rossler(tmax: int = 10000, n: int = 2000, Xi: Tuple[float, float, float] = (1.0, 1.0, 1.0), a: float = 0.2, b: float = 0.2, c: float = 5.7, dt: float = 0.01) → Tuple[ndarray, ndarray, ndarray]

Simulate the Rössler attractor.

The system is defined by:

dx/dt = -y - z
dy/dt =  x + a·y
dz/dt =  b + z·(x - c)

With b=0.2, c=5.7: chaotic for a ≲ 0.2, periodic for a ≳ 0.2.

Parameters:

• tmax (int, default 10000) – Number of integration steps.

• n (int, default 2000) – Number of output points (sub-sampled from the solution).

• Xi (tuple of float, default (1.0, 1.0, 1.0)) – Initial conditions (x₀, y₀, z₀).

• a, b, c (float) – System parameters. Classic chaotic: a=0.2, b=0.2, c=5.7.

• dt (float, default 0.01) – Integration step size.

Returns:

x, y, z (ndarray, shape (n,)) – State variables sampled at n equidistant times.

seed: int | None = None

RQA2_tests

class RQA2_tests(signal: ndarray, seed: int | None = None, max_workers: int = 1)

Bases: object

Surrogate-data generation and statistical validation for nonlinear dynamics testing.

Implements six surrogate algorithms and a comprehensive validation framework that tests each algorithm against multiple nonlinear metrics (Lyapunov exponent, sample entropy, correlation dimension, etc.).

Parameters:

• signal (ndarray of float) – 1-D floating-point time series to be tested.

• seed (int or None, optional) – Seed for reproducible surrogate generation.

• max_workers (int, default 1) – Number of parallel workers for surrogate generation (parallelism kicks in when n_surrogates >= 50).

Raises:

TypeError – If signal is not a floating-point array.

Examples

>>> import numpy as np
>>> rng = np.random.default_rng(0)
>>> signal = rng.standard_normal(512).astype(float)
>>> tester = RQA2_tests(signal, seed=42)
>>> surrogates = tester.generate('IAAFT', n_surrogates=100)
>>> surrogates.shape
(100, 512)

comprehensive_validation(systems_data: Dict[str, Dict[str, ndarray]], n_surrogates: int = 200, save_path: str = 'surrogate_validation_results.png') → Dict[str, Dict[str, float]]

Run all six surrogate algorithms against multiple dynamical systems and return rank-based two-sided p-values for six nonlinear metrics.

A heatmap of the results is saved to save_path.

Parameters:

• systems_data (dict) – Mapping from system name (str) to a dict with at least an 'x' key containing a 1-D ndarray (the x-coordinate of the attractor). Typically produced by RQA2_simulators.generate_test_battery().

• n_surrogates (int, default 200) – Number of surrogates per algorithm per system.

• save_path (str, default ‘surrogate_validation_results.png’) – Path where the validation heatmap is saved.

Returns:

dict – Nested dict: results[system_name][surrogate_method][metric] → p-value (float, or nan if the metric could not be computed).

generate(kind: Literal['FT', 'AAFT', 'IAAFT', 'IDFS', 'WIAAFT', 'PPS'] = 'FT', *, n_surrogates: int = 200, **kwargs) → ndarray

Generate an ensemble of surrogate time series.

Each surrogate is produced with an independent random seed to ensure statistical independence. When n_surrogates >= 50 and max_workers > 1 the ensemble is generated in parallel.

Parameters:

• kind ({‘FT’, ‘AAFT’, ‘IAAFT’, ‘IDFS’, ‘WIAAFT’, ‘PPS’}, default ‘FT’) – Surrogate algorithm:

  • FT – Fourier-Transform phase randomisation.

  • AAFT – Amplitude-Adjusted Fourier Transform.

  • IAAFT – Iterative AAFT (n_iter=100 by default).

  • IDFS – Iterative Digitally-Filtered Shuffled.

  • WIAAFT – Wavelet-based IAAFT.

  • PPS – Pseudo-Periodic Surrogate.

• n_surrogates (int, default 200) – Number of surrogates to generate.

• **kwargs – Passed to the underlying surrogate algorithm (e.g. n_iter=200 for IAAFT, wavelet='db4' for WIAAFT).

Returns:

ndarray, shape (n_surrogates, N) – Array of surrogate time series, one per row.

Raises:

KeyError – If kind is not a recognised algorithm name.

max_workers: int = 1

seed: int | None = None

signal: ndarray