Vectorized fast paths: trailing stack, moving-reference stretching, filter-once ensembles - #33
Closed
mdenolle wants to merge 2 commits into
Closed
Vectorized fast paths: trailing stack, moving-reference stretching, filter-once ensembles#33mdenolle wants to merge 2 commits into
mdenolle wants to merge 2 commits into
Conversation
…-once ensembles Three performance fast paths, each verified to reproduce the per-day loop it replaces to ~1e-15 in dv/v (regression tests at rtol=0, atol=1e-12): - _trailing_stack: difference of float64 cumulative sums, O(ndays*nlag) independent of the stack length instead of O(ndays*k*nlag) (~2.3x at k=45). - measure_stretching_trailing: vectorized stretching against a trailing reference. The stretched sample positions t/(1+eps) are data-independent, so the linear-interpolation gather indices/weights are computed once per epsilon and applied to all days at once; trailing references come from a cumulative sum and the band-pass runs once over the whole matrix. deviations._moving_reference dispatches to it for "stretching (TS)" (measured 4.7x on the 3-year volcano synthetic, 12.2 -> 2.6 s), keeping the generic per-day loop for other estimators. collect_cc behavior from the return_cc branch is preserved. - run_pipeline(..., prefiltered=True): callers evaluating several stack/reference variants at the same band can band-pass the raw CCF matrix once; the estimators skip their internal band-pass. Exact because the band-pass is linear and commutes with linear stacking; restricted to the estimators whose band usage is that one linear filter (stretching, WCC, DTW, MWCS), ValueError otherwise. Combined, a 5-member same-band ensemble drops ~4x per pair-band. Co-Authored-By: Claude Fable 5 <noreply@anthropic.com>
There was a problem hiding this comment.
Pull request overview
This PR adds several vectorized “fast paths” to speed up dv/v processing while preserving bit-level agreement (to floating-point rounding) with the existing per-day loop implementations. It also introduces an optional prefiltered=True mode so ensembles at the same band can apply the band-pass once and skip redundant internal filtering.
Changes:
- Reworked
_trailing_stackto use float64 cumulative sums for an O(ndays × nlag) trailing mean. - Added a vectorized moving-reference stretching implementation (
measure_stretching_trailing) and dispatched to it fromdeviations._moving_referencefor"stretching (TS)". - Added
run_pipeline(..., prefiltered=True)plus estimator-level support to skip internal band-pass where valid, and updated tests/changelog accordingly.
Reviewed changes
Copilot reviewed 4 out of 4 changed files in this pull request and generated 1 comment.
| File | Description |
|---|---|
tests/test_deviations.py |
Adds regression tests for the new fast paths and prefiltered=True behavior. |
src/codameter/synthetic_demo.py |
Implements vectorized trailing stretching and adds prefiltered plumbing to supported estimators. |
src/codameter/deviations.py |
Routes moving-reference stretching to the vectorized path; adds prefiltered support + validation in run_pipeline. |
CHANGELOG.md |
Documents the new fast paths and the prefiltered option. |
💡 Add Copilot custom instructions for smarter, more guided reviews. Learn how to get started.
| ) | ||
| np.testing.assert_allclose(dvv_b, dvv_a, rtol=0, atol=1e-12) | ||
| np.testing.assert_array_equal(val_b, val_a) | ||
| np.testing.assert_allclose(cc_b, cc_a, rtol=0, atol=1e-12) |
Copilot flagged assert_allclose(cc_b, cc_a, ...) as failing on the all-NaN cc arrays (inversion/mwcs cases, where no CC is collected) because equal_nan defaults to False. That's true of plain np.allclose, but numpy.testing.assert_allclose already defaults equal_nan=True -- verified all 4 parametrized cases (fixed/moving/inversion/mwcs) were already passing. Not a real bug, but making the default explicit (with a comment on why) so the next reader doesn't hit the same false alarm.
Contributor
Author
|
Superseded by #34 — this PR's base branch ( |
This file contains hidden or bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment
Add this suggestion to a batch that can be applied as a single commit.This suggestion is invalid because no changes were made to the code.Suggestions cannot be applied while the pull request is closed.Suggestions cannot be applied while viewing a subset of changes.Only one suggestion per line can be applied in a batch.Add this suggestion to a batch that can be applied as a single commit.Applying suggestions on deleted lines is not supported.You must change the existing code in this line in order to create a valid suggestion.Outdated suggestions cannot be applied.This suggestion has been applied or marked resolved.Suggestions cannot be applied from pending reviews.Suggestions cannot be applied on multi-line comments.Suggestions cannot be applied while the pull request is queued to merge.Suggestion cannot be applied right now. Please check back later.
Implements the three performance optimizations verified bit-exact by prototypes during the 2026-08-02 audit (max dv/v difference ~1e-15 vs the current code on the 1095-day volcano synthetic).
Stacked on #32 (
feat/run-pipeline-return-cc), because the moving-reference fast path must preserve thecollect_ccbehavior added there. When #32 merges, this PR retargets tomasterautomatically.Changes
_trailing_stackvia cumulative sums (synthetic_demo.py). Difference of float64 cumsums: O(ndays x nlag) independent of stack length, replacing the O(ndays x k x nlag) per-day loop. Measured 2.2x at k=45 (roughly at parity below k~15, where the loop was already cheap).Vectorized moving-reference stretching (
measure_stretching_trailing, dispatched fromdeviations._moving_referencefor"stretching (TS)"). The stretched sample positionst/(1+eps)are data-independent, so the linear-interpolation gather indices/weights are computed once per epsilon and applied to every day's band-passed trailing reference at once; trailing references are built by cumsum and the band-pass runs once over the whole matrix (2 FFTs total instead of 2 per day). Measured 4.7x on the 3-year volcano synthetic (12.2 -> 2.6 s), max |dv/v| difference 1.4e-15. Other estimators keep the generic per-day loop;collect_ccfrom Add return_cc option to run_pipeline #32 is preserved (regression-tested).Filter-once ensembles:
run_pipeline(..., prefiltered=True). Band-passing is linear, so it commutes exactly with linear stacking; callers evaluating several stack/reference variants at the same band can band-pass the raw CCF matrix once and the estimators skip their internal band-pass. Only valid at an identical band and only for estimators whose band usage is that one linear filter (stretching, WCC, DTW, MWCS); the wavelet estimators raiseValueError. Measured 1.7x on a 5-member same-band ensemble on top of the other fast paths.Combined effect on a 5-member ensemble: ~4x per pair-band.
Tests
TestFastPathRegressionsintests/test_deviations.py: each fast path vs the replaced per-day loop (replicated inline), andprefiltered=Truevs the internal band-pass for fixed / moving / inversion references and MWCS — all asserted atrtol=0, atol=1e-12(observed ~1e-15), plus theValueErrorfor wavelet estimators.tests/test_deviations.py,tests/test_golden.py, plustest_synthetic_demo,test_use_cases,test_uq_bayes,test_frugalmind_export,test_private_golden,test_bench: 93 passed (Python 3.13 env; the repo pixi env has a pre-3.10 interpreter that cannot import the package).python_version = "3.10"); a current mypy run on the changed files shows no new errors vs the branch baseline.CHANGELOG updated under Unreleased.
🤖 Generated with Claude Code