From 2d0148f1992a2f11968fdc4b24ef809462d5c395 Mon Sep 17 00:00:00 2001 From: Mohammad Abdul Sahil <127765312+abdulsaheel@users.noreply.github.com> Date: Thu, 20 Aug 2026 07:41:20 +0530 Subject: [PATCH 1/3] cardio stager: one beat-gather, not three copies of it MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit _windowRmssd and _windowRemFeatures each hand-copied _cleanBeatsInWindow, which is exactly what that function's docstring says can't happen. folded both back into it; the rem one passes halfWinMs: 90s because its window is deliberately different (spanSec, the 240-point ls grid and the <16 abstain gate are all specified against ±90s, not ±2.5min) — noted that in the doc so nobody "simplifies" it later. also fixed the docstring rationale: it claimed the beats had to match so the z-scores stayed comparable, but rmssd isn't z-scored any more. real-night gate unchanged: wake=52 rem=300 nrem=716 deep=137, same confidence to the last digit. --- lib/src/onehz/sleep/cardio_stager.dart | 110 +++++++++---------------- 1 file changed, 41 insertions(+), 69 deletions(-) diff --git a/lib/src/onehz/sleep/cardio_stager.dart b/lib/src/onehz/sleep/cardio_stager.dart index 7a7957f..e6639aa 100644 --- a/lib/src/onehz/sleep/cardio_stager.dart +++ b/lib/src/onehz/sleep/cardio_stager.dart @@ -950,30 +950,7 @@ CardioStagerResult _abstain(int epochSec) => CardioStagerResult( /// window centred on epoch [s,t). Returns NaN when too few clean beats. double _windowRmssd(List rrMs, List rrTsMs, List accel, int s, int t, int epochSec) { - if (rrMs.isEmpty || rrTsMs.length != rrMs.length) return double.nan; - final mid = (s + t) ~/ 2; - if (mid >= accel.length) return double.nan; - final centreMs = accel[mid].tsMs; - const halfWinMs = 150 * 1000; // ±2.5 min for a stable RMSSD on sparse RR - final lo = centreMs - halfWinMs, hi = centreMs + halfWinMs; - // Gather clean beats in window (300–2000 ms, drop big successive jumps). - final beats = []; - double? prev; - for (var i = 0; i < rrMs.length; i++) { - final ts = rrTsMs[i]; - if (ts < lo || ts > hi) continue; - final v = rrMs[i]; - if (v < _rrMin || v > _rrMax) { - prev = null; - continue; - } - if (prev != null && (v - prev).abs() > _rrMaxStep) { - prev = v; - continue; - } - beats.add(v); - prev = v; - } + final beats = _cleanBeatsInWindow(rrMs, rrTsMs, accel, s, t).beats; if (beats.length < 5) return double.nan; var ss = 0.0; for (var i = 1; i < beats.length; i++) { @@ -995,7 +972,7 @@ double _windowRmssd(List rrMs, List rrTsMs, /// flat RMSSD of 67.3 / 70.7 / 79.7. double _windowSdnn(List rrMs, List rrTsMs, List accel, int s, int t, int epochSec) { - final beats = _cleanBeatsInWindow(rrMs, rrTsMs, accel, s, t); + final beats = _cleanBeatsInWindow(rrMs, rrTsMs, accel, s, t).beats; if (beats.length < 5) return double.nan; final m = mean(beats)!; var ss = 0.0; @@ -1005,18 +982,33 @@ double _windowSdnn(List rrMs, List rrTsMs, return math.sqrt(ss / (beats.length - 1)); } -/// Clean RR beats (ms) inside the ±2.5-min window centred on epoch [s,t). -/// Shared by [_windowRmssd] and [_windowSdnn] so both see exactly the same -/// beats — a divergence there would make their z-scores incomparable. -List _cleanBeatsInWindow(List rrMs, List rrTsMs, - List accel, int s, int t) { - if (rrMs.isEmpty || rrTsMs.length != rrMs.length) return const []; +/// Clean RR beats (ms) inside a ±[halfWinMs] window centred on epoch [s,t), +/// with each beat's time rebased to the window start (seconds). +/// +/// The physiologic gate (300–2000 ms, drop big successive jumps) lives here and +/// only here, so every window feature sees exactly the same beats for a given +/// window — [_windowRmssd] and [_windowSdnn] are computed over one window and +/// must agree on its contents. +/// +/// [_windowRemFeatures] passes a DIFFERENT [halfWinMs] on purpose (±90 s, not +/// ±2.5 min): its LF/HF grid and its `beats.length` abstain gate are specified +/// against that shorter window. Do not hand it the default list. +({List beats, List tsSec}) _cleanBeatsInWindow( + List rrMs, + List rrTsMs, + List accel, + int s, + int t, { + int halfWinMs = 150 * 1000, +}) { + const empty = (beats: [], tsSec: []); + if (rrMs.isEmpty || rrTsMs.length != rrMs.length) return empty; final mid = (s + t) ~/ 2; - if (mid >= accel.length) return const []; + if (mid >= accel.length) return empty; final centreMs = accel[mid].tsMs; - const halfWinMs = 150 * 1000; final lo = centreMs - halfWinMs, hi = centreMs + halfWinMs; final beats = []; + final tsSec = []; double? prev; for (var i = 0; i < rrMs.length; i++) { final ts = rrTsMs[i]; @@ -1031,9 +1023,18 @@ List _cleanBeatsInWindow(List rrMs, List rrTsMs, continue; } beats.add(v); + // Rebase to the window start (lo), NOT absolute epoch ms. Lomb–Scargle is + // time-shift invariant (the τ phase reference cancels any offset), so the + // LF/HF output is unchanged — but this keeps beat times in [0, 180] s + // instead of ~1.75e9 s. Absolute epoch seconds force every sin/cos in the + // periodogram onto libm's __kernel_rem_pio2 multi-precision slow path + // (args ~9e9 rad), which — run per 30-s epoch over a full night on the + // main isolate — caused main-thread ANRs on Android (Crashlytics: libm.so + // __kernel_rem_pio2 / sin / cos, 0.9.13). + tsSec.add((ts - lo) / 1000.0); prev = v; } - return beats; + return (beats: beats, tsSec: tsSec); } /// Webster sleep-continuity rescore: brief wake bouts flanked by enough sleep @@ -1135,41 +1136,12 @@ void _websterRescore(List sm, int epochSec) { /// nulls when too few clean beats for a stable estimate. ({double? lfhf, double? rk}) _windowRemFeatures(List rrMs, List rrTsMs, List accel, int s, int t, int epochSec) { - if (rrMs.isEmpty || rrTsMs.length != rrMs.length) { - return (lfhf: null, rk: null); - } - final mid = (s + t) ~/ 2; - if (mid >= accel.length) return (lfhf: null, rk: null); - final centreMs = accel[mid].tsMs; - const halfWinMs = 90 * 1000; // ±90 s per the REM feature spec - final lo = centreMs - halfWinMs, hi = centreMs + halfWinMs; - final beats = []; // clean RR (ms) - final beatTsSec = []; // matching beat times (s) - double? prev; - for (var i = 0; i < rrMs.length; i++) { - final ts = rrTsMs[i]; - if (ts < lo || ts > hi) continue; - final v = rrMs[i]; - if (v < _rrMin || v > _rrMax) { - prev = null; - continue; - } - if (prev != null && (v - prev).abs() > _rrMaxStep) { - prev = v; - continue; - } - beats.add(v); - // Rebase to the window start (lo), NOT absolute epoch ms. Lomb–Scargle is - // time-shift invariant (the τ phase reference cancels any offset), so the - // LF/HF output is unchanged — but this keeps beat times in [0, 180] s - // instead of ~1.75e9 s. Absolute epoch seconds force every sin/cos in the - // periodogram onto libm's __kernel_rem_pio2 multi-precision slow path - // (args ~9e9 rad), which — run per 30-s epoch over a full night on the - // main isolate — caused main-thread ANRs on Android (Crashlytics: libm.so - // __kernel_rem_pio2 / sin / cos, 0.9.13). - beatTsSec.add((ts - lo) / 1000.0); - prev = v; - } + // ±90 s per the REM feature spec — a DIFFERENT window from the RMSSD/SDNN + // one; see [_cleanBeatsInWindow]. + final win = _cleanBeatsInWindow(rrMs, rrTsMs, accel, s, t, + halfWinMs: 90 * 1000); + final beats = win.beats; // clean RR (ms) + final beatTsSec = win.tsSec; // matching beat times (s), rebased to window if (beats.length < 16) return (lfhf: null, rk: null); // spectral stability gate // R(k): mean absolute successive difference of instantaneous HR (bpm). From de52888f5c090cc90bc685aa1e2f1d6ee1138d81 Mon Sep 17 00:00:00 2001 From: Mohammad Abdul Sahil <127765312+abdulsaheel@users.noreply.github.com> Date: Thu, 20 Aug 2026 07:42:51 +0530 Subject: [PATCH 2/3] util: percentileSorted, so the local hr window sorts once not twice MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit the rolling hr baseline took median(win) and percentile(win, 25) of the same 361-epoch window, and each of those sorted a fresh copy — 25 lines above the comment bragging about saving sorts. sort win in place once (it's built right there) and take both off it. also hoisted the median(sleepHr) that ran twice in the night-observation record. percentileSorted keeps the double? / null-on-empty contract. the two private _percentileSorted copies in advanced_stager and load_trimp return 0 on empty — not promoting that, a 0.0 bpm hr floor recorded as a measurement is a bug we've already had. deduping those two is a separate job. real-night gate unchanged: wake=52 rem=300 nrem=716 deep=137. --- lib/src/onehz/sleep/cardio_stager.dart | 14 +++++++++----- lib/src/onehz/util.dart | 16 ++++++++++++---- test/onehz/util_test.dart | 9 +++++++++ 3 files changed, 30 insertions(+), 9 deletions(-) diff --git a/lib/src/onehz/sleep/cardio_stager.dart b/lib/src/onehz/sleep/cardio_stager.dart index e6639aa..36e72ba 100644 --- a/lib/src/onehz/sleep/cardio_stager.dart +++ b/lib/src/onehz/sleep/cardio_stager.dart @@ -674,11 +674,14 @@ CardioStagerResult classifyCardioEpochs( for (var k = lo; k < hi; k++) if (still(k) && !hr[k].isNaN) hr[k] ]; - final m = median(win); + // One sort for both percentiles of this 361-epoch window (stddev doesn't + // care about order). `win` is freshly built here, so sorting it is local. + win.sort(); + final m = percentileSorted(win, 50); if (m != null) { hrMedLocal[e] = m; hrArousalLocal[e] = m + math.max(6.0, (stddev(win) ?? 6)); - hrP25Local[e] = percentile(win, 25) ?? m; + hrP25Local[e] = percentileSorted(win, 25) ?? m; } } // Blend the per-epoch LOCAL HR gates toward the sleeper's rolling profile. @@ -872,13 +875,14 @@ CardioStagerResult classifyCardioEpochs( // Only when the edge armed recording AND the staged span is a real sleep // (≥60 min) — naps must not pollute the sleeper's night profile. if (cardioRecordObservations && nEpoch >= 120 && sleepHr.isNotEmpty) { + final hrSleepMed = median(sleepHr); _cardioObservations.add(SleepNightObservation( epochs: nEpoch, hrFloorP5: percentile(sleepHr, 5), hrFloorP25: percentile(sleepHr, 25), - hrSleepMedian: median(sleepHr), - hrArousal: (median(sleepHr) ?? hrMedGlobal) + - math.max(6.0, stddev(sleepHr) ?? 6.0), + hrSleepMedian: hrSleepMed, + hrArousal: + (hrSleepMed ?? hrMedGlobal) + math.max(6.0, stddev(sleepHr) ?? 6.0), rmssdMed: rmssdMed, rmssdMad: mad(sleepRmssd), enmoStillCut: motMed + 1.5 * motMad, diff --git a/lib/src/onehz/util.dart b/lib/src/onehz/util.dart index ff584fc..5b73896 100644 --- a/lib/src/onehz/util.dart +++ b/lib/src/onehz/util.dart @@ -40,10 +40,14 @@ double? stddevPop(List xs) { return math.sqrt(s / xs.length); } -/// Linear-interpolated percentile (p in [0,100]); null if empty. -double? percentile(List values, double p) { - if (values.isEmpty) return null; - final sorted = [...values]..sort(); +/// Linear-interpolated percentile (p in [0,100]) of an ALREADY-ASCENDING list; +/// null if empty. For when you hold a sorted list already, or want several +/// percentiles of one window without re-sorting it each time. +/// +/// Null on empty, never 0 — a 0.0 is a measurement, and an empty window is the +/// absence of one. +double? percentileSorted(List sorted, double p) { + if (sorted.isEmpty) return null; if (sorted.length == 1) return sorted[0]; final rank = (p / 100) * (sorted.length - 1); final lo = rank.floor(); @@ -53,6 +57,10 @@ double? percentile(List values, double p) { return sorted[lo] + (sorted[hi] - sorted[lo]) * frac; } +/// Linear-interpolated percentile (p in [0,100]); null if empty. +double? percentile(List values, double p) => + percentileSorted([...values]..sort(), p); + /// Median (linear-interpolated), or null if empty. double? median(List xs) => percentile(xs, 50); diff --git a/test/onehz/util_test.dart b/test/onehz/util_test.dart index 177c50a..49d6e55 100644 --- a/test/onehz/util_test.dart +++ b/test/onehz/util_test.dart @@ -16,6 +16,15 @@ void main() { expect(stddev([5]), isNull); }); + test('percentileSorted matches percentile and stays null on empty', () { + final sorted = [1.0, 2.0, 3.0, 4.0]; + expect(percentileSorted(sorted, 50), percentile(sorted, 50)); + expect(percentileSorted(sorted, 25), percentile(sorted, 25)); + expect(percentileSorted([7.0], 25), 7.0); + // NOT 0 — an empty window is no measurement, not a 0 bpm floor. + expect(percentileSorted([], 25), isNull); + }); + test('MAD and robust z; MAD=0 guard on quantized data', () { // [1,1,1,1,1] -> median 1, MAD 0 => robustZ null (not div-by-zero). expect(mad([1, 1, 1, 1, 1], scaled: false), 0); From 8bf3be886a1fbeb0b4df30d4d86071e865a04a87 Mon Sep 17 00:00:00 2001 From: Mohammad Abdul Sahil <127765312+abdulsaheel@users.noreply.github.com> Date: Thu, 20 Aug 2026 07:45:56 +0530 Subject: [PATCH 3/3] cardio stager: binary-search the rr window instead of scanning the night MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit the per-epoch window gather walked all ~40k beats of the night to find the few hundred inside ±2.5 min, once per epoch per feature. lower-bound in, break at the far edge. the bound has to be a TRUE lower bound. rr_ts_ms is rec_ts * 1000 so beats tie on the second boundary, and "first index strictly greater than lo" drops the tied beats sitting on the edge — i tried it: deep goes 137 -> 132 epochs on the real night, silently. test pins ties at both edges and fails on that exact mutation. no cursor threaded through the epoch loop on purpose — centreMs comes from accel, not rrTsMs, so a cross-epoch cursor would assume an ordering the caller was never asked for. cardioStager on the real night: ~950ms -> ~805ms, so about 15%. not more — lomb-scargle's trig is the rest of it and this doesn't touch it. counts unchanged: wake=52 rem=300 nrem=716 deep=137. --- lib/src/onehz/sleep/cardio_stager.dart | 38 +++++++++++++++++-- test/onehz/real_night_cardio_stager_test.dart | 36 ++++++++++++++++++ 2 files changed, 70 insertions(+), 4 deletions(-) diff --git a/lib/src/onehz/sleep/cardio_stager.dart b/lib/src/onehz/sleep/cardio_stager.dart index 36e72ba..a6e97d1 100644 --- a/lib/src/onehz/sleep/cardio_stager.dart +++ b/lib/src/onehz/sleep/cardio_stager.dart @@ -366,8 +366,9 @@ void resetCardioObservations() => _cardioObservations.clear(); /// [hr1hz] per-second HR (bpm; 0 = off-skin) over the in-bed window. /// [accel] per-second gravity vectors, SAME length/time base as [hr1hz]. /// [rrMs] / [rrTsMs] beat-to-beat RR (ms) and their ABSOLUTE times (ms), same -/// clock as `accel[i].tsMs`. Sparse/empty is fine — REM/deep just lean more on -/// HR and confidence drops. +/// clock as `accel[i].tsMs`, ASCENDING in [rrTsMs] (the per-epoch window +/// gather lower-bounds into it). Sparse/empty is fine — REM/deep just lean +/// more on HR and confidence drops. CardioStagerResult cardioStager( List hr1hz, List accel, { @@ -1011,12 +1012,28 @@ double _windowSdnn(List rrMs, List rrTsMs, if (mid >= accel.length) return empty; final centreMs = accel[mid].tsMs; final lo = centreMs - halfWinMs, hi = centreMs + halfWinMs; + // [rrTsMs] is ascending, so lower-bound into it and stop at the far edge + // instead of walking every beat of the night once per epoch per feature. + // + // TRUE lower bound — the first index with ts >= lo, NOT the first index + // strictly greater than lo. `rr_ts_ms` is `rec_ts * 1000`, so beats TIE on + // the second boundary; dropping the tied beats at a window edge changes + // RMSSD / SDNN / LF-HF and with them the REM and deep epoch counts, silently. + var a = 0, b = rrTsMs.length; + while (a < b) { + final m = (a + b) >> 1; + if (rrTsMs[m] < lo) { + a = m + 1; + } else { + b = m; + } + } final beats = []; final tsSec = []; double? prev; - for (var i = 0; i < rrMs.length; i++) { + for (var i = a; i < rrMs.length; i++) { final ts = rrTsMs[i]; - if (ts < lo || ts > hi) continue; + if (ts > hi) break; final v = rrMs[i]; if (v < _rrMin || v > _rrMax) { prev = null; @@ -1041,6 +1058,19 @@ double _windowSdnn(List rrMs, List rrTsMs, return (beats: beats, tsSec: tsSec); } +/// [_cleanBeatsInWindow]'s beats, exposed so the window-edge regression test +/// can drive the gather directly. Not a staging entry point — call +/// [cardioStager]. +List cleanBeatsInWindowForTest( + List rrMs, + List rrTsMs, + List accel, + int s, + int t, { + int halfWinMs = 150 * 1000, +}) => + _cleanBeatsInWindow(rrMs, rrTsMs, accel, s, t, halfWinMs: halfWinMs).beats; + /// Webster sleep-continuity rescore: brief wake bouts flanked by enough sleep /// are re-labelled sleep (NREM). This is the published actigraphy step that /// prevents normal in-sleep repositioning from inflating WASO. diff --git a/test/onehz/real_night_cardio_stager_test.dart b/test/onehz/real_night_cardio_stager_test.dart index c1d470c..7ac34cf 100644 --- a/test/onehz/real_night_cardio_stager_test.dart +++ b/test/onehz/real_night_cardio_stager_test.dart @@ -45,6 +45,8 @@ import 'package:test/test.dart'; import 'package:openstrap_analytics/onehz.dart'; void main() { + _windowEdgeTies(); + test('cardioStager on the real 2026-07 overnight capture matches Apple ' 'Watch ground truth within a wide band (regression: was wake=294min ' 'rem=41min pre-fix)', () { @@ -135,3 +137,37 @@ void main() { expect(tot, closeTo(stages.length * epochMin, 0.01)); }); } + +// The per-epoch RR window is found by binary search. `rr_ts_ms` is +// `rec_ts * 1000`, so beat timestamps TIE — several beats share one second, +// and a second boundary can land exactly on a window edge. The search must be +// a TRUE lower bound (first index >= lo) and the far edge must be inclusive; +// "first index strictly greater than lo" drops the tied beats at the edge, +// which moves RMSSD / SDNN / LF-HF and with them the REM and deep counts. +void _windowEdgeTies() { + test('window gather keeps tied beats at BOTH edges', () { + final accel = [ + for (var i = 0; i < 30; i++) AccelSample(i * 1000.0, 0, 0, 1) + ]; + // epoch [0,30) -> mid 15 -> centre 15000 ms; ±5 s -> window [10000, 20000]. + final ts = [ + 9000, // just before -> out + 10000, 10000, 10000, // tie ON the low edge -> in + 15000, + 20000, 20000, 20000, // tie ON the high edge -> in + 21000, // just after -> out + ]; + final rr = [for (var i = 0; i < ts.length; i++) 900.0]; + + final beats = + cleanBeatsInWindowForTest(rr, ts, accel, 0, 30, halfWinMs: 5000); + // 7, not 4 (low ties dropped) and not 5 (both edges exclusive). + expect(beats.length, 7); + + // and the beats really are the in-window ones, in order. + final marked = [900, 901, 902, 903, 904, 905, 906, 907, 908]; + final got = + cleanBeatsInWindowForTest(marked, ts, accel, 0, 30, halfWinMs: 5000); + expect(got, [901, 902, 903, 904, 905, 906, 907]); + }); +}