From 05e0a9c6274aa9ddac6bda4f42ea8993e2ddaaef Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Carlson=20B=C3=BCth?= Date: Thu, 30 Jul 2026 17:55:58 +0200 Subject: [PATCH 1/2] perf(ordinal-patterns): serial path for small data and bypass decorator overhead --- .../continuous_ordinal_patterns.py | 91 ++++++++++++++----- 1 file changed, 66 insertions(+), 25 deletions(-) diff --git a/delaynet/connectivities/continuous_ordinal_patterns.py b/delaynet/connectivities/continuous_ordinal_patterns.py index cc38479..877070c 100644 --- a/delaynet/connectivities/continuous_ordinal_patterns.py +++ b/delaynet/connectivities/continuous_ordinal_patterns.py @@ -4,7 +4,10 @@ from numba import njit, prange from ..decorators import connectivity -from .granger import gt_multi_lag +from ..utils.lag_steps import find_optimal_lag +from .granger import gt_multi_lag, gt_single_lag + +_SERIAL_THRESHOLD = 200 @connectivity @@ -53,7 +56,9 @@ def random_patterns( t_ts2 = pattern_transform(np.copy(ts2), rnd_patterns) for i in range(num_rnd_patterns): - p_v, pv_idx = gt_multi_lag(t_ts1[i, :], t_ts2[i, :], lag_steps=lag_steps) + p_v, pv_idx = find_optimal_lag( + gt_single_lag, t_ts1[i, :], t_ts2[i, :], lag_steps=lag_steps + ) if best_pv > p_v: best_pv = p_v best_lag = pv_idx @@ -99,7 +104,6 @@ def pattern_transform(ts: np.ndarray, patterns: np.ndarray) -> np.ndarray: return pattern_transform_2d(ts, patterns).squeeze() -@njit(nogil=True, parallel=True) def pattern_transform_2d( ts: np.ndarray, patterns: np.ndarray ) -> np.ndarray: # pragma: no cover @@ -128,34 +132,74 @@ def pattern_transform_2d( return transformed -@njit(nogil=True, parallel=True) -def norm_windows(ts: np.ndarray, window_size: int) -> np.ndarray: # pragma: no cover - """Normalise sliding windows of a time series to values between -1 and 1. +@njit(nogil=True, cache=True) +def _norm_windows_serial( + ts: np.ndarray, window_size: int +) -> np.ndarray: # pragma: no cover + windows = np.lib.stride_tricks.as_strided( + x=ts, + strides=(ts.strides[0], ts.strides[0]), + shape=(ts.shape[0] - window_size + 1, window_size), + ) + normed_windows = np.zeros_like(windows) + for i in range(windows.shape[0]): + normed_windows[i] = norm_window(windows[i]) + return normed_windows - :param ts: Time series. - :type ts: numpy.ndarray - :param window_size: Size of the window. - :type window_size: int - :return: Normalised windows. - :rtype: numpy.ndarray - """ - # Create a sliding window view of the input array - # windows = np.lib.stride_tricks.sliding_window_view( - # x=ts, window_shape=window_size, writeable=False - # ) + +@njit(nogil=True, parallel=True, cache=True) +def _norm_windows_parallel( + ts: np.ndarray, window_size: int +) -> np.ndarray: # pragma: no cover windows = np.lib.stride_tricks.as_strided( x=ts, strides=(ts.strides[0], ts.strides[0]), shape=(ts.shape[0] - window_size + 1, window_size), ) normed_windows = np.zeros_like(windows) - # Normalise each window to [-1, 1] for i in prange(windows.shape[0]): normed_windows[i] = norm_window(windows[i]) return normed_windows -@njit(nogil=True, parallel=True) +def norm_windows(ts: np.ndarray, window_size: int) -> np.ndarray: + """Normalise sliding windows of a time series to values between -1 and 1. + + :param ts: Time series. + :type ts: numpy.ndarray + :param window_size: Size of the window. + :type window_size: int + :return: Normalised windows. + :rtype: numpy.ndarray + """ + n_windows = ts.shape[0] - window_size + 1 + if n_windows < _SERIAL_THRESHOLD: + return _norm_windows_serial(ts, window_size) + return _norm_windows_parallel(ts, window_size) + + +@njit(nogil=True, cache=True) +def _pattern_distance_serial( + windows: np.ndarray, pattern: np.ndarray +) -> np.ndarray: # pragma: no cover + distances = np.zeros(windows.shape[0]) + for i in range(windows.shape[0]): + for j in range(pattern.shape[0]): + distances[i] += np.abs(windows[i, j] - pattern[j]) + return distances / pattern.shape[0] / 2.0 + + +@njit(nogil=True, parallel=True, cache=True) +def _pattern_distance_parallel( + windows: np.ndarray, pattern: np.ndarray +) -> np.ndarray: # pragma: no cover + distances = np.zeros(windows.shape[0]) + for i in prange(windows.shape[0]): + for j in prange(pattern.shape[0]): + distances[i] += np.abs(windows[i, j] - pattern[j]) + return distances / pattern.shape[0] / 2.0 + + def pattern_distance( windows: np.ndarray, pattern: np.ndarray ) -> np.ndarray: # pragma: no cover @@ -168,9 +212,6 @@ def pattern_distance( :return: Distance between the windows and the pattern. :rtype: numpy.ndarray """ - distances = np.zeros(windows.shape[0]) - for i in prange(windows.shape[0]): - for j in prange(pattern.shape[0]): - distances[i] += np.abs(windows[i, j] - pattern[j]) - return distances / pattern.shape[0] / 2.0 - # equiv. to np.sum(np.abs(windows - pattern), axis=1) / pattern.shape[0] / 2.0 + if windows.shape[0] < _SERIAL_THRESHOLD: + return _pattern_distance_serial(windows, pattern) + return _pattern_distance_parallel(windows, pattern) From 9dd0b4240db054c8c09bf47012e9a197b65a5697 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Carlson=20B=C3=BCth?= Date: Thu, 30 Jul 2026 18:03:05 +0200 Subject: [PATCH 2/2] test(ordinal-patterns): random-noise tests and serial-vs-parallel equivalence --- .../test_continuous_ordinal_patterns.py | 75 +++++++++++++++++++ 1 file changed, 75 insertions(+) diff --git a/tests/connectivities/test_continuous_ordinal_patterns.py b/tests/connectivities/test_continuous_ordinal_patterns.py index d7e2db2..e00a650 100644 --- a/tests/connectivities/test_continuous_ordinal_patterns.py +++ b/tests/connectivities/test_continuous_ordinal_patterns.py @@ -1,10 +1,17 @@ """Tests for the continuous ordinal patterns connectivity measure.""" +import random as _random +import numpy as np import pytest from numpy import array, allclose, linspace, random, roll +from numpy.lib.stride_tricks import as_strided from delaynet import connectivity from delaynet.connectivities.continuous_ordinal_patterns import ( + _norm_windows_parallel, + _norm_windows_serial, + _pattern_distance_parallel, + _pattern_distance_serial, norm_window, norm_windows, pattern_distance, @@ -12,6 +19,35 @@ random_patterns, ) +_SEED = _random.randint(0, 10000) # deterministic per-session + + +def _norm_windows_expected(ts, window_size): + """Pure-numpy reference for sliding window normalisation (no numba).""" + windows = as_strided( + ts, + strides=(ts.strides[0], ts.strides[0]), + shape=(ts.shape[0] - window_size + 1, window_size), + ) + result = np.zeros_like(windows) + for i in range(windows.shape[0]): + w = windows[i] + nw = w - w.min() + nw = nw / nw.max() + nw = (nw - 0.5) * 2.0 + nw[np.isnan(nw)] = 0.0 + result[i] = nw + return result + + +def _pattern_distance_expected(windows, pattern): + """Pure-numpy reference for pattern distance (no numba).""" + distances = np.zeros(windows.shape[0]) + for i in range(windows.shape[0]): + for j in range(pattern.shape[0]): + distances[i] += np.abs(windows[i, j] - pattern[j]) + return distances / pattern.shape[0] / 2.0 + @pytest.mark.parametrize( "ts, expected", @@ -58,6 +94,26 @@ def test_norm_windows(ts, window_size, expected): assert allclose(norm_windows(array(ts), window_size), array(expected)) +@pytest.mark.parametrize("window_size", [2, 5, 10]) +def test_norm_windows_random_noise(window_size): + """norm_windows on non-monotonic random data matches independent reference.""" + rng = np.random.default_rng(_SEED) + ts = rng.normal(0, 1, size=200) + assert allclose( + norm_windows(ts, window_size), _norm_windows_expected(ts, window_size) + ) + + +def test_norm_windows_serial_parallel_equivalence(): + """Serial and parallel paths produce identical output on large data.""" + rng = np.random.default_rng(_SEED + 1) + ts = rng.normal(0, 1, size=500) + window_size = 2 # 499 windows > _SERIAL_THRESHOLD = 200 + serial = _norm_windows_serial(ts, window_size) + parallel = _norm_windows_parallel(ts, window_size) + assert allclose(serial, parallel) + + @pytest.mark.parametrize( "windows, pattern, expected", [ @@ -73,6 +129,25 @@ def test_pattern_distance(windows, pattern, expected): assert allclose(pattern_distance(array(windows), array(pattern)), array(expected)) +def test_pattern_distance_random_noise(): + """pattern_distance on random data matches independent reference.""" + rng = np.random.default_rng(_SEED + 2) + windows = rng.uniform(-1, 1, (50, 5)) + pattern = rng.uniform(-1, 1, (5,)) + expected = _pattern_distance_expected(windows, pattern) + assert allclose(pattern_distance(windows, pattern), expected) + + +def test_pattern_distance_serial_parallel_equivalence(): + """Serial and parallel paths produce identical output on large data.""" + rng = np.random.default_rng(_SEED + 3) + windows = rng.uniform(-1, 1, (500, 3)) + pattern = rng.uniform(-1, 1, (3,)) + serial = _pattern_distance_serial(windows, pattern) + parallel = _pattern_distance_parallel(windows, pattern) + assert allclose(serial, parallel) + + @pytest.mark.parametrize( "ts, patterns, expected", [