From e8f1e3b964847240a24a47fd74f4b1b9767cd23e Mon Sep 17 00:00:00 2001 From: Mohammad Abdul Sahil <127765312+abdulsaheel@users.noreply.github.com> Date: Sun, 4 Oct 2026 19:45:28 +0530 Subject: [PATCH 1/4] rmssd: when the jitter gate refuses, accept a breathing line that stays put in Hz while hr drifts --- lib/src/onehz/clinical/hrv_time.dart | 282 +++++++++++++++++++++++---- test/onehz/clinical_test.dart | 54 +++++ 2 files changed, 293 insertions(+), 43 deletions(-) diff --git a/lib/src/onehz/clinical/hrv_time.dart b/lib/src/onehz/clinical/hrv_time.dart index 5b7607a..c51d537 100644 --- a/lib/src/onehz/clinical/hrv_time.dart +++ b/lib/src/onehz/clinical/hrv_time.dart @@ -256,6 +256,157 @@ List _clearedWindows(double? nightAcf1, List rmssds, ]; } +final _hann128 = [ + for (var i = 0; i < 128; i++) 0.5 - 0.5 * math.cos(2 * math.pi * i / 127) +]; + +/// Welch power (128-beat Hann, 50 % overlap) of one window's RR, rebuilt from +/// its difference runs, at cycles-per-beat frequencies [f]. Null when no run +/// holds a full segment. +List? _windowPsd(List> diffRuns, List f) { + const n = 128; + final p = List.filled(f.length, 0.0); + var segs = 0; + for (final r in diffRuns) { + final x = List.filled(r.length + 1, 0.0); + for (var i = 0; i < r.length; i++) { + x[i + 1] = x[i] + r[i]; + } + for (var s = 0; s + n <= x.length; s += n ~/ 2) { + var m = 0.0; + for (var i = 0; i < n; i++) { + m += x[s + i]; + } + m /= n; + for (var j = 0; j < f.length; j++) { + final c = math.cos(2 * math.pi * f[j]), sn = math.sin(2 * math.pi * f[j]); + var cr = 1.0, ci = 0.0, re = 0.0, im = 0.0; + for (var i = 0; i < n; i++) { + final v = (x[s + i] - m) * _hann128[i]; + re += v * cr; + im += v * ci; + final t = cr * c - ci * sn; + ci = cr * sn + ci * c; + cr = t; + } + p[j] += re * re + im * im; + } + segs++; + } + } + return segs == 0 ? null : p; +} + +/// Last resort for a night the jitter gate refused: the windows whose +/// breathing line holds still in Hz while the heart rate moves. +/// +/// RSA follows breathing (Hirsch & Bishop 1981), a rate in breaths per +/// minute, so its line sits at a fixed frequency in Hz; a heart rate that +/// drifts across the night slides it in cycles per beat. Beat-timing jitter, +/// alternation and grid artifacts are tied to the beat, not the clock. So each +/// 5-min window's spectrum is pooled twice, once on a cycles-per-beat axis and +/// once rescaled by its own mean RR onto a Hz axis, and the night passes only +/// when the Hz pooling shows a line in 8–30 br/min that is sharper than the +/// beat pooling. If the heart rate barely moves the two poolings coincide and +/// the night stays refused: stability in Hz says nothing there. +/// +/// Ceiling: needs >= 12 windows and a heart rate that wanders; breathing at +/// half the heart rate is still an alternation and still refused. +/// Returns the indices of the windows whose own peak sits on the line. +List? _steadyBreathingWindows( + List>> winRuns, List meanRrMs) { + var ssd = 0.0, nd = 0; + final all = [for (final w in winRuns) ...w]; + for (final r in all) { + for (final d in r) { + ssd += d * d; + nd++; + } + } + if (nd == 0 || _onCoarseLattice(all, ssd / nd)) return null; + final ref = median([for (final m in meanRrMs) if (m > 0) m]); + if (ref == null) return null; + // 0.15–0.5 cycles/beat at the night's median beat. + final fc = [for (var j = 38; j <= 128; j++) j / 256]; + final byBeat = List.filled(fc.length, 0.0); + final byHz = List.filled(fc.length, 0.0); + final cover = List.filled(fc.length, 0); + final peakHz = {}; // window -> its own peak, cycles per ms + for (var w = 0; w < winRuns.length; w++) { + if (meanRrMs[w] <= 0) continue; + final pb = _windowPsd(winRuns[w], fc); + if (pb == null) continue; + final norm = median(pb)!; + if (norm <= 0) continue; + final k = meanRrMs[w] / ref; + final fh = [for (final f in fc) f * k]; + final ph = _windowPsd(winRuns[w], fh)!; + var best = -1; + for (var j = 0; j < fc.length; j++) { + byBeat[j] += pb[j] / norm; + if (fh[j] > 0.5) continue; + byHz[j] += ph[j] / norm; + cover[j]++; + if (best < 0 || ph[j] > ph[best]) best = j; + } + if (best >= 0) peakHz[w] = fc[best] / ref; + } + final nw = peakHz.length; + if (nw < 12) return null; + final bins = [for (var j = 0; j < fc.length; j++) if (cover[j] == nw) j]; + if (bins.length < 20) return null; + final hz = [for (final j in bins) byHz[j]]; + final bt = [for (final j in bins) byBeat[j]]; + var pk = 0; + for (var i = 1; i < hz.length; i++) { + if (hz[i] > hz[pk]) pk = i; + } + if (pk == 0 || pk == hz.length - 1) return null; + final f0 = fc[bins[pk]] / ref; + final brpm = f0 * 60000; + if (brpm < 8 || brpm > 30) return null; + // Same bins, same normalisation: the line must stand taller aligned in Hz + // than the beat pooling stands anywhere within ±12 bins of it. + var near = 0.0; + for (var i = math.max(0, pk - 12); i <= math.min(bt.length - 1, pk + 12); i++) { + near = math.max(near, bt[i]); + } + // The pooled noise ripple shrinks as 1/√windows, so the bar does too. + final promHz = hz[pk] / median(hz)!; + if (promHz < math.max(2.5, 1.5 + 7 / math.sqrt(nw))) return null; + if (promHz < 1.1 * near / median(bt)!) return null; + // Alternation lives at Nyquist on the beat axis, outside these bins. + if (byBeat.last / median(bt)! >= 0.5 * promHz) return null; + final keep = [ + for (final e in peakHz.entries) + if (math.log(e.value / f0).abs() <= 0.15) e.key + ]; + return keep.length >= 6 ? keep : null; +} + +/// The window RMSSDs a windowed headline publishes, or null for none: the +/// usual gate first, and only when it leaves nothing, the windows on a +/// breathing line that holds still in Hz ([_steadyBreathingWindows]), which +/// publish at the confidence floor. +({List rmssds, bool byBreathing, String note})? _keptWindows( + double? acf1, + bool refused, + List rmssds, + List>> winRuns, + List meanRrMs) { + final kept = refused ? const [] : _clearedWindows(acf1, rmssds, winRuns); + if (kept.isNotEmpty) return (rmssds: kept, byBreathing: false, note: ''); + if (acf1 == null || acf1 >= kNnDiffAcf1Floor) return null; + final idx = _steadyBreathingWindows(winRuns, meanRrMs); + if (idx == null) return null; + return ( + rmssds: [for (final i in idx) rmssds[i]], + byBreathing: true, + note: ' Jitter gate failed; kept the ${idx.length} windows on a breathing ' + 'line steady in Hz while heart rate drifted, at floor confidence.', + ); +} + /// Confidence multiplier for a measured [acf1]: 1.0 on a smooth tachogram, /// falling linearly to 0 at [kNnDiffAcf1Floor] so confidence bottoms out /// exactly where RMSSD is refused. 1.0 when ACF1 could not be measured. @@ -362,8 +513,12 @@ Metric hrvTime( // always advised. final acf1 = nnDiffAcf1(runs); final jittery = _jitterRefused(acf1, runs); - final rmssd = (pairs > 0 && !jittery) ? math.sqrt(ssd / pairs) : null; + var rmssd = (pairs > 0 && !jittery) ? math.sqrt(ssd / pairs) : null; final pnn50 = (pairs > 0 && !jittery) ? 100.0 * nn50 / pairs : null; + // Refused: a long record may still show a breathing line steady in Hz + // ([_steadyBreathingWindows]); then RMSSD alone comes from those 5-min windows. + if (jittery && gapAware) rmssd = _breathingRmssd(nnMs, nnTimesMs); + final byBreathing = jittery && rmssd != null; final sdnn = stddev(nnMs); double? sdann, sdnnIndex; @@ -384,11 +539,13 @@ Metric hrvTime( // were ~pure noise. The beat-count term is capped BEFORE the quality terms // multiply it; multiplying first let an all-night beat count (n/250 ≈ 100) // swallow any penalty and re-clamp to 0.95 regardless. - final conf = ((nnMs.length / 250.0).clamp(0.0, 1.0) // ~250 beats ≈ 5 min - * - _acf1Quality(acf1) * - (1 - artifactFraction)) - .clamp(0.3, 0.95); + final conf = byBreathing + ? 0.3 + : ((nnMs.length / 250.0).clamp(0.0, 1.0) // ~250 beats ≈ 5 min + * + _acf1Quality(acf1) * + (1 - artifactFraction)) + .clamp(0.3, 0.95); return Metric( value: HrvTime( rmssd: rmssd, @@ -404,12 +561,47 @@ Metric hrvTime( inputs_used: inputs, note: jittery ? '${_jitterNote(acf1!)}. SDNN/SDANN survive it and are the lead here. ' - 'PRV not ECG-HRV.' + '${byBreathing ? 'RMSSD only from 5-min windows on a breathing ' + 'line steady in Hz, at floor confidence. ' : ''}PRV not ECG-HRV.' : 'PRV not ECG-HRV; RMSSD/pNN50 are quantization-sensitive at 1 Hz ' '— lead with SDNN/SDANN', ); } +/// RMSSD pooled over the 5-min windows of [nn] that sit on a breathing line +/// steady in Hz, or null. Same seam rule as [hrvTime]. +double? _breathingRmssd(List nn, List times) { + final wins = >>{}; + final rrSum = {}, rrN = {}; + var prevWin = -1; + for (var i = 0; i < nn.length; i++) { + final w = ((times[i] - times.first) / 300000.0).floor(); + rrSum[w] = (rrSum[w] ?? 0) + nn[i]; + rrN[w] = (rrN[w] ?? 0) + 1; + final runs = wins[w] ??= >[]; + if (i > 0 && w == prevWin && times[i] - times[i - 1] <= nn[i] + 0.5) { + runs.last.add(nn[i] - nn[i - 1]); + } else { + runs.add([]); + } + prevWin = w; + } + final keys = wins.keys.toList()..sort(); + final idx = _steadyBreathingWindows( + [for (final k in keys) wins[k]!], [for (final k in keys) rrSum[k]! / rrN[k]!]); + if (idx == null) return null; + var ss = 0.0, n = 0; + for (final i in idx) { + for (final r in wins[keys[i]]!) { + for (final d in r) { + ss += d * d; + n++; + } + } + } + return n == 0 ? null : math.sqrt(ss / n); +} + /// Robust NOCTURNAL RMSSD (ms). /// /// A single whole-night RMSSD is dominated by the few high-Δ segments produced @@ -465,6 +657,7 @@ Metric nocturnalRmssd( // calmest-looking windows while the night pooled to −0.43/−0.51. final runs = >[]; final perWindow = >>[]; + final meanRr = []; final indices = buckets.keys.toList()..sort(); for (final idx in indices) { if (stageMaskPerSec != null) { @@ -504,47 +697,46 @@ Metric nocturnalRmssd( if (nd < minBeatsPerWindow) continue; runs.addAll(winRuns); perWindow.add(winRuns); + meanRr.add(mean([for (final i in seg) nnMs[i]])!); rmssds.add(math.sqrt(ssd / nd)); } final acf1 = nnDiffAcf1(runs); - if (_jitterRefused(acf1, runs)) { - return Metric.absent( - tier: Tier.high, - inputs_used: inputs, - note: _jitterNote(acf1!), - ); - } - if (rmssds.isEmpty) { + final refused = _jitterRefused(acf1, runs); + if (rmssds.isEmpty && !refused) { return const Metric.absent( tier: Tier.high, inputs_used: inputs, note: 'no usable 5-min windows for nocturnal RMSSD', ); } - final kept = _clearedWindows(acf1, rmssds, perWindow); - if (kept.isEmpty) { + final kept = _keptWindows(acf1, refused, rmssds, perWindow, meanRr); + if (kept == null) { return Metric.absent( tier: Tier.high, inputs_used: inputs, - note: '${_jitterNote(acf1!)}; no 5-min window clears it on its own', + note: refused + ? _jitterNote(acf1!) + : '${_jitterNote(acf1!)}; no 5-min window clears it on its own', ); } - final robust = median(kept)!; + final robust = median(kept.rmssds)!; // Confidence scales with how many windows we could median over, and with the // measured jitter level (see [kNnDiffAcf1Floor]). - final conf = - ((kept.length / 12.0).clamp(0.0, 1.0) * _acf1Quality(acf1)).clamp( - // 12 ≈ 1 h - 0.3, - 0.95); + final conf = kept.byBreathing + ? 0.3 + : ((kept.rmssds.length / 12.0).clamp(0.0, 1.0) * _acf1Quality(acf1)) + .clamp( + // 12 ≈ 1 h + 0.3, + 0.95); return Metric( value: robust, confidence: conf, tier: Tier.high, inputs_used: inputs, - note: 'robust nocturnal RMSSD = MEDIAN of ${kept.length} consecutive ' + note: 'robust nocturnal RMSSD = MEDIAN of ${kept.rmssds.length} consecutive ' '5-min-window RMSSDs (REM/arousal-robust). PRV not ECG-HRV; ' - 'RMSSD is quantization-sensitive at 1 Hz.', + 'RMSSD is quantization-sensitive at 1 Hz.${kept.note}', ); } @@ -607,11 +799,15 @@ Metric sleepSessionWindowedRmssd( final rmssds = []; final runs = >[]; // pooled jitter floor — see [nocturnalRmssd] final perWindow = >>[]; + final meanRr = []; final indices = buckets.keys.toList()..sort(); for (final idx in indices) { - final diffRuns = [ + final rrRuns = [ for (final r in _cleanWindowRuns(buckets[idx]!, bucketsTs[idx]!)) - if (r.length >= 2) [for (var i = 1; i < r.length; i++) r[i] - r[i - 1]] + if (r.length >= 2) r + ]; + final diffRuns = [ + for (final r in rrRuns) [for (var i = 1; i < r.length; i++) r[i] - r[i - 1]] ]; var ssd = 0.0; var nd = 0; @@ -626,6 +822,7 @@ Metric sleepSessionWindowedRmssd( if (nd < 5) continue; runs.addAll(diffRuns); perWindow.add(diffRuns); + meanRr.add(mean([for (final r in rrRuns) ...r])!); rmssds.add(math.sqrt(ssd / nd)); } @@ -633,38 +830,37 @@ Metric sleepSessionWindowedRmssd( // are noise, the honest output is no headline, not a plausible one — the // readiness composite already treats a null HRV driver as absent. final acf1 = nnDiffAcf1(runs); - if (_jitterRefused(acf1, runs)) { - return Metric.absent( - tier: Tier.high, - inputs_used: inputs, - note: _jitterNote(acf1!), - ); - } - if (rmssds.isEmpty) { + final refused = _jitterRefused(acf1, runs); + if (rmssds.isEmpty && !refused) { return const Metric.absent( tier: Tier.high, inputs_used: inputs, note: 'no valid 5-min windows for sleep-session RMSSD', ); } - final kept = _clearedWindows(acf1, rmssds, perWindow); - if (kept.isEmpty) { + final kept = _keptWindows(acf1, refused, rmssds, perWindow, meanRr); + if (kept == null) { return Metric.absent( tier: Tier.high, inputs_used: inputs, - note: '${_jitterNote(acf1!)}; no 5-min window clears it on its own', + note: refused + ? _jitterNote(acf1!) + : '${_jitterNote(acf1!)}; no 5-min window clears it on its own', ); } - final meanRmssd = mean(kept)!; - final conf = ((kept.length / 12.0).clamp(0.0, 1.0) * _acf1Quality(acf1)) - .clamp(0.3, 0.95); + final meanRmssd = mean(kept.rmssds)!; + final conf = kept.byBreathing + ? 0.3 + : ((kept.rmssds.length / 12.0).clamp(0.0, 1.0) * _acf1Quality(acf1)) + .clamp(0.3, 0.95); return Metric( value: meanRmssd, confidence: conf, tier: Tier.high, inputs_used: inputs, - note: 'sleep-session HRV: mean RMSSD over cleaned 5-min windows.', + note: 'sleep-session HRV: mean RMSSD over cleaned 5-min windows.' + '${kept.note}', ); } diff --git a/test/onehz/clinical_test.dart b/test/onehz/clinical_test.dart index 852ed61..1ba0573 100644 --- a/test/onehz/clinical_test.dart +++ b/test/onehz/clinical_test.dart @@ -311,6 +311,60 @@ void main() { expect(m.value!.rmssd, closeTo(30, 2), reason: '√2·30·sin(π/4)'); }); + // 6 h at 18 br/min with beat-time jitter loud enough that the pooled + // spectral exemption refuses; HR 50 ± [drift] bpm over a 3 h cycle. + (List, List) slowNight(int seed, double drift, double rsa) { + final rnd = math.Random(seed); + final rr = [], ts = []; + var t = 0.0, e0 = 0.0; + while (t < 6 * 3600e3) { + final base = 60000 / (50 + drift * math.sin(2 * math.pi * t / 10800e3)); + final e1 = (rnd.nextDouble() - 0.5) * 60; + final v = base + rsa * math.sin(2 * math.pi * 0.3 * t / 1000) + e1 - e0; + e0 = e1; + t += v; + rr.add(v); + ts.add(t); + } + return (rr, ts); + } + + test('HRV-02: slow-heart RSA steady in Hz while HR drifts publishes', () { + for (var seed = 0; seed < 2; seed++) { + final (rr, ts) = slowNight(seed, 8, 10); + final ss = sleepSessionWindowedRmssd(rr, ts, + startSec: 1, endSec: (ts.last / 1000).floor()); + expect(ss.present, isTrue, reason: 'seed $seed'); + expect(ss.note, contains('steady in Hz')); + expect(ss.confidence, 0.3); + expect(nocturnalRmssd(rr, ts).present, isTrue, reason: 'seed $seed'); + expect(hrvTime(rr, nnTimesMs: ts).value!.rmssd, isNotNull); + // Same night with the heart rate held flat: Hz and beats coincide, so + // stability in Hz proves nothing and the night stays refused. + final (fr, ft) = slowNight(seed, 0, 10); + expect( + sleepSessionWindowedRmssd(fr, ft, + startSec: 1, endSec: (ft.last / 1000).floor()) + .present, + isFalse, + reason: 'seed $seed'); + } + }); + + test('HRV-02: beat-time jitter with a drifting HR stays refused', () { + for (var seed = 0; seed < 4; seed++) { + final (rr, ts) = slowNight(seed, 8, 0); + expect(hrvTime(rr, nnTimesMs: ts).value!.rmssd, isNull); + expect(nocturnalRmssd(rr, ts).present, isFalse, reason: 'seed $seed'); + expect( + sleepSessionWindowedRmssd(rr, ts, + startSec: 1, endSec: (ts.last / 1000).floor()) + .present, + isFalse, + reason: 'seed $seed'); + } + }); + test('HRV-02: confidence carries jitter and artifact, not beat count alone', () { // It used to be clamp(n/250, .3, .95), which published 0.95 on all 13 From 9870ef4c8568aa2eaefc9bb312e70e6c52c988f1 Mon Sep 17 00:00:00 2001 From: Mohammad Abdul Sahil <127765312+abdulsaheel@users.noreply.github.com> Date: Sun, 4 Oct 2026 20:34:51 +0530 Subject: [PATCH 2/4] rmssd breathing line: refuse when beat-time jitter is most of the number --- lib/src/onehz/clinical/hrv_time.dart | 39 ++++++++++++++++++++++++---- test/onehz/clinical_test.dart | 17 +++++++++--- 2 files changed, 48 insertions(+), 8 deletions(-) diff --git a/lib/src/onehz/clinical/hrv_time.dart b/lib/src/onehz/clinical/hrv_time.dart index c51d537..11f68cf 100644 --- a/lib/src/onehz/clinical/hrv_time.dart +++ b/lib/src/onehz/clinical/hrv_time.dart @@ -259,10 +259,11 @@ List _clearedWindows(double? nightAcf1, List rmssds, final _hann128 = [ for (var i = 0; i < 128; i++) 0.5 - 0.5 * math.cos(2 * math.pi * i / 127) ]; +final _hann128Sq = _hann128.fold(0.0, (a, w) => a + w * w); -/// Welch power (128-beat Hann, 50 % overlap) of one window's RR, rebuilt from -/// its difference runs, at cycles-per-beat frequencies [f]. Null when no run -/// holds a full segment. +/// Welch power (128-beat Hann, 50 % overlap, averaged over segments) of one +/// window's RR, rebuilt from its difference runs, at cycles-per-beat +/// frequencies [f]. Null when no run holds a full segment. List? _windowPsd(List> diffRuns, List f) { const n = 128; final p = List.filled(f.length, 0.0); @@ -294,7 +295,7 @@ List? _windowPsd(List> diffRuns, List f) { segs++; } } - return segs == 0 ? null : p; + return segs == 0 ? null : [for (final v in p) v / segs]; } /// Last resort for a night the jitter gate refused: the windows whose @@ -332,6 +333,7 @@ List? _steadyBreathingWindows( final byHz = List.filled(fc.length, 0.0); final cover = List.filled(fc.length, 0); final peakHz = {}; // window -> its own peak, cycles per ms + final jit = >{}; // window -> Hz-axis PSD over jitter shape for (var w = 0; w < winRuns.length; w++) { if (meanRrMs[w] <= 0) continue; final pb = _windowPsd(winRuns[w], fc); @@ -350,6 +352,10 @@ List? _steadyBreathingWindows( if (best < 0 || ph[j] > ph[best]) best = j; } if (best >= 0) peakHz[w] = fc[best] / ref; + jit[w] = [ + for (var j = 0; j < fc.length; j++) + ph[j] / (2 - 2 * math.cos(2 * math.pi * math.min(fh[j], 0.5))) + ]; } final nw = peakHz.length; if (nw < 12) return null; @@ -381,7 +387,30 @@ List? _steadyBreathingWindows( for (final e in peakHz.entries) if (math.log(e.value / f0).abs() <= 0.15) e.key ]; - return keep.length >= 6 ? keep : null; + if (keep.length < 6) return null; + // The line proves breathing is there, not that it carries RMSSD: beat-time + // jitter loud enough to fail the gate would still be most of the number. + // Same ceiling as [nnDiffNoiseShare], over the kept windows. The floor is + // read as beat-time jitter σ², whose RR spectrum is σ²·(2 − 2cos ω) and + // whose differences carry 6σ², off the Hz pooling where the line is one + // narrow peak the median steps over. White RR noise reads louder here than + // it is, which can only refuse more. + var sq = 0.0; + final pooled = List.filled(bins.length, 0.0); + for (final w in keep) { + var n = 0; + for (final r in winRuns[w]) { + for (final d in r) { + sq += d * d; + n++; + } + } + for (var i = 0; i < bins.length; i++) { + pooled[i] += jit[w]![bins[i]] * n; + } + } + final noise = 6 * median(pooled)! / _hann128Sq; + return noise < kNnDiffNoiseShareCeiling * sq ? keep : null; } /// The window RMSSDs a windowed headline publishes, or null for none: the diff --git a/test/onehz/clinical_test.dart b/test/onehz/clinical_test.dart index 1ba0573..434077b 100644 --- a/test/onehz/clinical_test.dart +++ b/test/onehz/clinical_test.dart @@ -313,13 +313,14 @@ void main() { // 6 h at 18 br/min with beat-time jitter loud enough that the pooled // spectral exemption refuses; HR 50 ± [drift] bpm over a 3 h cycle. - (List, List) slowNight(int seed, double drift, double rsa) { + (List, List) slowNight(int seed, double drift, double rsa, + [double jitter = 60]) { final rnd = math.Random(seed); final rr = [], ts = []; var t = 0.0, e0 = 0.0; while (t < 6 * 3600e3) { final base = 60000 / (50 + drift * math.sin(2 * math.pi * t / 10800e3)); - final e1 = (rnd.nextDouble() - 0.5) * 60; + final e1 = (rnd.nextDouble() - 0.5) * jitter; final v = base + rsa * math.sin(2 * math.pi * 0.3 * t / 1000) + e1 - e0; e0 = e1; t += v; @@ -331,7 +332,7 @@ void main() { test('HRV-02: slow-heart RSA steady in Hz while HR drifts publishes', () { for (var seed = 0; seed < 2; seed++) { - final (rr, ts) = slowNight(seed, 8, 10); + final (rr, ts) = slowNight(seed, 8, 10, 20); final ss = sleepSessionWindowedRmssd(rr, ts, startSec: 1, endSec: (ts.last / 1000).floor()); expect(ss.present, isTrue, reason: 'seed $seed'); @@ -348,6 +349,16 @@ void main() { .present, isFalse, reason: 'seed $seed'); + // A real line under jitter that is most of the RMSSD (~42 of ~44 ms + // here, RSA alone ~13): the line is there, the number is not. + final (lr, lt) = slowNight(seed, 8, 10); + expect( + sleepSessionWindowedRmssd(lr, lt, + startSec: 1, endSec: (lt.last / 1000).floor()) + .present, + isFalse, + reason: 'seed $seed'); + expect(nocturnalRmssd(lr, lt).present, isFalse, reason: 'seed $seed'); } }); From 202aee3e9713ecdfadbad70cd5ebd5c4712c1866 Mon Sep 17 00:00:00 2001 From: Mohammad Abdul Sahil <127765312+abdulsaheel@users.noreply.github.com> Date: Sun, 4 Oct 2026 22:49:05 +0530 Subject: [PATCH 3/4] rmssd breathing line: read jitter off the floor outside the night's breathing band, not a median a wandering line spreads into --- lib/src/onehz/clinical/hrv_time.dart | 50 +++++++++++++++++++-------- test/onehz/clinical_test.dart | 51 ++++++++++++++++++++++++++++ 2 files changed, 86 insertions(+), 15 deletions(-) diff --git a/lib/src/onehz/clinical/hrv_time.dart b/lib/src/onehz/clinical/hrv_time.dart index 11f68cf..175572e 100644 --- a/lib/src/onehz/clinical/hrv_time.dart +++ b/lib/src/onehz/clinical/hrv_time.dart @@ -333,7 +333,7 @@ List? _steadyBreathingWindows( final byHz = List.filled(fc.length, 0.0); final cover = List.filled(fc.length, 0); final peakHz = {}; // window -> its own peak, cycles per ms - final jit = >{}; // window -> Hz-axis PSD over jitter shape + final psd = >{}; // window -> its cycles-per-beat PSD for (var w = 0; w < winRuns.length; w++) { if (meanRrMs[w] <= 0) continue; final pb = _windowPsd(winRuns[w], fc); @@ -352,10 +352,7 @@ List? _steadyBreathingWindows( if (best < 0 || ph[j] > ph[best]) best = j; } if (best >= 0) peakHz[w] = fc[best] / ref; - jit[w] = [ - for (var j = 0; j < fc.length; j++) - ph[j] / (2 - 2 * math.cos(2 * math.pi * math.min(fh[j], 0.5))) - ]; + psd[w] = pb; } final nw = peakHz.length; if (nw < 12) return null; @@ -390,13 +387,21 @@ List? _steadyBreathingWindows( if (keep.length < 6) return null; // The line proves breathing is there, not that it carries RMSSD: beat-time // jitter loud enough to fail the gate would still be most of the number. - // Same ceiling as [nnDiffNoiseShare], over the kept windows. The floor is - // read as beat-time jitter σ², whose RR spectrum is σ²·(2 − 2cos ω) and - // whose differences carry 6σ², off the Hz pooling where the line is one - // narrow peak the median steps over. White RR noise reads louder here than - // it is, which can only refuse more. - var sq = 0.0; - final pooled = List.filled(bins.length, 0.0); + // Same ceiling as [nnDiffNoiseShare], over the kept windows. Jitter σ² has + // an RR spectrum σ²·(2 − 2cos ω), flat once divided by that shape, and its + // differences carry 6σ²; physiology only adds power on top. Breathing + // wanders across the night and its line spreads with it, so each window + // drops the bins of every rate the kept windows peaked at (plus a Hann + // main lobe), mapped onto its own beats. σ² is the floor of what is left, + // pooled over windows and smoothed over two resolution cells. + var lo = double.infinity, hi = 0.0; + for (final w in keep) { + lo = math.min(lo, peakHz[w]!); + hi = math.max(hi, peakHz[w]!); + } + var sq = 0.0, nAll = 0; + final acc = List.filled(fc.length, 0.0); + final wt = List.filled(fc.length, 0.0); for (final w in keep) { var n = 0; for (final r in winRuns[w]) { @@ -405,11 +410,26 @@ List? _steadyBreathingWindows( n++; } } - for (var i = 0; i < bins.length; i++) { - pooled[i] += jit[w]![bins[i]] * n; + nAll += n; + final a = lo * meanRrMs[w] - 3 / 128, b = hi * meanRrMs[w] + 3 / 128; + for (var j = 0; j < fc.length; j++) { + if (fc[j] > a && fc[j] < b) continue; + acc[j] += n * psd[w]![j] / (2 - 2 * math.cos(2 * math.pi * fc[j])); + wt[j] += n; + } + } + var floor = double.infinity; + for (var j = 4; j < fc.length - 4; j++) { + var s = 0.0, m = 0; + for (var i = j - 4; i <= j + 4; i++) { + if (wt[i] < nAll / 2) break; + s += acc[i] / wt[i]; + m++; } + if (m == 9) floor = math.min(floor, s / 9); } - final noise = 6 * median(pooled)! / _hann128Sq; + if (floor == double.infinity) return null; + final noise = 6 * floor / _hann128Sq * nAll; return noise < kNnDiffNoiseShareCeiling * sq ? keep : null; } diff --git a/test/onehz/clinical_test.dart b/test/onehz/clinical_test.dart index 434077b..6d2a5a7 100644 --- a/test/onehz/clinical_test.dart +++ b/test/onehz/clinical_test.dart @@ -362,6 +362,57 @@ void main() { } }); + test('HRV-02: slow-heart RSA with wandering breathing publishes near truth', + () { + // Real breathing is not a fixed tone: rate wanders ±2–3 br/min around + // 16–19, depth swings ±30–50 %, HR drifts 42–58, plus small beat-time + // jitter. Truth = RMSSD of the jitter-free RR. + for (var seed = 0; seed < 8; seed++) { + final rnd = math.Random(seed); + final br0 = 16 + 3 * rnd.nextDouble(), brA = 2 + rnd.nextDouble(); + final dA = 0.3 + 0.2 * rnd.nextDouble(), ph = 6.28 * rnd.nextDouble(); + final rr = [], ts = [], clean = []; + var t = 0.0, phi = 0.0, e0 = 0.0, br = br0, depth = 40.0; + while (t < 6 * 3600e3) { + final hr = 50 + 8 * math.sin(2 * math.pi * t / 10800e3 + ph); + final c = 60000 / hr + depth * math.sin(phi); + final e1 = (rnd.nextDouble() - 0.5) * 40; + final v = c + e1 - e0; + e0 = e1; + phi += 2 * math.pi * br / 60 * c / 1000; + if (phi > 2 * math.pi) { + // Each breath its own rate and depth. + phi -= 2 * math.pi; + final g = math.sqrt(-2 * math.log(1 - rnd.nextDouble())) * + math.cos(2 * math.pi * rnd.nextDouble()); + br = br0 + brA * g.clamp(-2.0, 1.5); + depth = 40 * (1 + dA * (2 * rnd.nextDouble() - 1)); + } + t += v; + clean.add(c); + rr.add(v); + ts.add(t); + } + var ss = 0.0; + for (var i = 1; i < clean.length; i++) { + ss += math.pow(clean[i] - clean[i - 1], 2); + } + final truth = math.sqrt(ss / (clean.length - 1)); + final got = [ + hrvTime(rr, nnTimesMs: ts).value!.rmssd, + nocturnalRmssd(rr, ts).value, + sleepSessionWindowedRmssd(rr, ts, + startSec: 1, endSec: (ts.last / 1000).floor()) + .value, + ]; + for (final g in got) { + expect(g, isNotNull, reason: 'seed $seed'); + expect(g! / truth, inInclusiveRange(1 / 1.3, 1.3), + reason: 'seed $seed: $g vs $truth'); + } + } + }); + test('HRV-02: beat-time jitter with a drifting HR stays refused', () { for (var seed = 0; seed < 4; seed++) { final (rr, ts) = slowNight(seed, 8, 0); From 76834c8a9cf761b85931d08aed93eeecfcb6c274 Mon Sep 17 00:00:00 2001 From: Mohammad Abdul Sahil <127765312+abdulsaheel@users.noreply.github.com> Date: Sun, 4 Oct 2026 23:25:15 +0530 Subject: [PATCH 4/4] rmssd breathing line: smooth the jitter floor wider and undo the min's low bias --- lib/src/onehz/clinical/hrv_time.dart | 12 +++++++----- 1 file changed, 7 insertions(+), 5 deletions(-) diff --git a/lib/src/onehz/clinical/hrv_time.dart b/lib/src/onehz/clinical/hrv_time.dart index 175572e..b5c7438 100644 --- a/lib/src/onehz/clinical/hrv_time.dart +++ b/lib/src/onehz/clinical/hrv_time.dart @@ -393,7 +393,7 @@ List? _steadyBreathingWindows( // wanders across the night and its line spreads with it, so each window // drops the bins of every rate the kept windows peaked at (plus a Hann // main lobe), mapped onto its own beats. σ² is the floor of what is left, - // pooled over windows and smoothed over two resolution cells. + // pooled over windows and smoothed over four resolution cells. var lo = double.infinity, hi = 0.0; for (final w in keep) { lo = math.min(lo, peakHz[w]!); @@ -419,17 +419,19 @@ List? _steadyBreathingWindows( } } var floor = double.infinity; - for (var j = 4; j < fc.length - 4; j++) { + for (var j = 8; j < fc.length - 8; j++) { var s = 0.0, m = 0; - for (var i = j - 4; i <= j + 4; i++) { + for (var i = j - 8; i <= j + 8; i++) { if (wt[i] < nAll / 2) break; s += acc[i] / wt[i]; m++; } - if (m == 9) floor = math.min(floor, s / 9); + if (m == 17) floor = math.min(floor, s / 17); } if (floor == double.infinity) return null; - final noise = 6 * floor / _hann128Sq * nAll; + // The lowest stretch of a noisy curve reads below its mean; 1.1 puts the + // floor back at these pooled sizes, so jitter alone is not undercounted. + final noise = 1.1 * 6 * floor / _hann128Sq * nAll; return noise < kNnDiffNoiseShareCeiling * sq ? keep : null; }