diff --git a/ALGORITHMS.md b/ALGORITHMS.md index 7529a12..fcd318b 100644 --- a/ALGORITHMS.md +++ b/ALGORITHMS.md @@ -71,7 +71,6 @@ Grouped by family (subdirectory under `lib/src/onehz/`). File paths are relative | `strainScoreFromSeries` | `clinical/load_trimp.dart` | **the headline 0-21 strain** from per-minute wake HR: `netTrimpAboveQuiet` log-mapped onto 0-21 | — | | `strainCurveFromSeries` | `clinical/load_trimp.dart` | cumulative 0-21 strain per minute, same arithmetic; last point == headline | — | | `strainScoreMetric` | `clinical/load_trimp.dart` | LUMP form (TRIMP minus wake-minutes × quiet cost) log-squashed onto 0-21; can debit exercise, kept only for rescaling a stored TRIMP with no per-minute series | — | -| `trimpStrain` | `clinical/load_trimp.dart` | TRIMP → 0-100 strain, honesty-wrapped (absent without real HRmax/RHR anchors) | — | | `ctlAtlTsb` | `clinical/load_trimp.dart` | Fitness-Fatigue-Form: EWMA CTL (42d) / ATL (7d) / TSB = CTL-ATL | Banister impulse-response model | | `baevskyStressIndex` | `clinical/stress_si.dart` | Baevsky Stress Index, 5-min sliding RR-histogram | Baevsky & Berseneva 2008 | | `irregularBeatScreen` | `clinical/irregular_rhythm.dart` | Poincaré SD1/SD2 + pNN70 irregular-rhythm screen | — (screen, not a diagnosis; no ECG) | @@ -133,7 +132,7 @@ Grouped by family (subdirectory under `lib/src/onehz/`). File paths are relative | `socialJetlag` / `chronotype` | `human/circadian_lifestyle.dart` | MSFsc-based social jetlag + chronotype | Wittmann & Roenneberg 2006; MCTQ | | `alcoholNightFlag` | `human/event_detection.dart` | dose-graded autonomic-stress signature — **reports state, never asserts a cause** | Pietilä et al. 2018 | | `roughNight` | `human/event_detection.dart` | neutral fallback descriptor when the signature is ambiguous | — | -| `percentileOfYou` / `personalRecord` | `human/percentile_of_you.dart` | percentile-vs-your-own-history, miss-tolerant personal-record streaks | — | +| `percentileOfYou` | `human/percentile_of_you.dart` | percentile-vs-your-own-history | — | | `glassBoxReadiness` | `human/readiness_glassbox.dart` | **deprecated** — kept only for its percentile-of-you breakdown + narrative and edge back-compat; `readinessComposite` is canonical | — | | `sleepNeed` / `strainTarget` / `recommendedBedtime` / `recommendedWake` / `sleepPerformance` | `human/coaching.dart` | deterministic coaching layer over the metrics above (naps come from `sleep/nap.dart`, not here) | — | | `journalCorrelations` | `human/coaching.dart` | per-tag mean-difference correlation vs. logged outcomes, on-device, personal | — | diff --git a/README.md b/README.md index 5edee80..efedba5 100644 --- a/README.md +++ b/README.md @@ -94,10 +94,10 @@ foundation layers: baselines, inverse-variance fusion. - **`clinical/`** (Tier-1) — HRV time/frequency domain (RMSSD/SDNN/pNNx, Lomb-Scargle LF/HF), PRSA (deceleration/acceleration capacity), nocturnal RHR/dip, an illness-risk - CUSUM state machine, Plews ln-RMSSD readiness, Baevsky stress index, Banister/Edwards - TRIMP + CTL/ATL/TSB training load (the 0-21 day strain, `strainScoreMetric`, is scored - above the user's own quiet waking heart-rate level rather than a fixed population - reference; the 0-100 `trimpStrain` effort score is raw TRIMP with no such subtraction), + CUSUM state machine, Plews ln-RMSSD readiness, Baevsky stress index, Banister + TRIMP + CTL/ATL/TSB training load (the 0-21 day strain is scored above the user's own + quiet waking heart-rate level rather than a fixed population reference, and it is the + only strain scale), a Poincaré irregular-rhythm screen, cosinor circadian fitting, and real-time cardiac coherence for guided breathing sessions. - **`sleep/`** — van Hees z-angle segmentation feeding a cardiac/motion stager (the single source of truth for sleep staging), AASM hypnogram metrics, cardiopulmonary @@ -114,8 +114,8 @@ foundation layers: - **`wellness/`** — the canonical composite readiness score, multivariate (Mahalanobis) anomaly detection, CUSUM changepoint detection, temperature-based illness flagging. - **`human/`** — sleep regularity index, social jetlag/chronotype, single-night event - detection (never names a specific cause — more on that below), percentile-of-you/ - personal records, and the deterministic coaching layer. + detection (never names a specific cause — more on that below), percentile-of-you, + and the deterministic coaching layer. ## The rule that matters most: never name a cause diff --git a/lib/src/onehz/clinical/cosinor.dart b/lib/src/onehz/clinical/cosinor.dart index aeda025..6d4ddfa 100644 --- a/lib/src/onehz/clinical/cosinor.dart +++ b/lib/src/onehz/clinical/cosinor.dart @@ -17,8 +17,15 @@ import '../util.dart'; class CosinorFit { final double mesor; // rhythm-adjusted mean final double amplitude; // half the peak-to-trough range - final double acrophaseRad; // phase of the peak (rad), in [-π, π] - final double acrophaseHours; // acrophase expressed as clock-hours of [period] + /// Phase of the peak (rad), in [-π, π]. Null when the zero-amplitude F-test + /// does not reject (Nelson 1979): with no resolvable rhythm the angle of + /// (β, γ) is set by noise or rounding, not by the data, and a caller would + /// otherwise present it as the user's own peak time. + final double? acrophaseRad; + final double? acrophaseHours; // acrophase expressed as clock-hours of [period] + + /// p-value of the zero-amplitude F-test, F(2, n−3). + final double amplitudeP; final double periodHours; final double r2; // raw goodness of fit (0..1) @@ -34,18 +41,34 @@ class CosinorFit { required this.periodHours, required this.r2, required this.r2Adj, + required this.amplitudeP, }); Map toJson() => { 'mesor': round6(mesor), 'amplitude': round6(amplitude), - 'acrophase_rad': round6(acrophaseRad), - 'acrophase_hours': round6(acrophaseHours), + 'acrophase_rad': acrophaseRad == null ? null : round6(acrophaseRad!), + 'acrophase_hours': + acrophaseHours == null ? null : round6(acrophaseHours!), 'period_hours': round6(periodHours), 'r2': round6(r2), 'r2_adj': round6(r2Adj), + 'amplitude_p': round6(amplitudeP), }; } +/// Significance level of the zero-amplitude test below which the acrophase is +/// reported. +const double cosinorAcrophaseAlpha = 0.05; + +/// Upper tail of F(2, d): P(F > f) = (1 + 2f/d)^(−d/2), exact for 2 numerator +/// degrees of freedom. +double _f2Survival(double f, int d) { + if (f.isNaN) return 1.0; + if (f.isInfinite) return 0.0; + if (f <= 0) return 1.0; + return math.pow(1 + 2 * f / d, -d / 2).toDouble(); +} + /// Minimum samples for a cosinor fit. The model has 3 free parameters /// (M, β, γ); with n = 4 there is a single residual degree of freedom, so R² /// is essentially an interpolation score (4 RANDOM points routinely fit at @@ -131,22 +154,41 @@ Metric cosinor( // score confidence 0.95 at tier HIGH. final double r2Adj = (1 - (1 - r2) * (n - 1) / (n - 3)).clamp(0, 1); + // Zero-amplitude F-test (Nelson 1979): the acrophase is only defined once + // the test rejects. A flat or rhythm-free series still yields some + // atan2(−γ, β) from residual noise or float rounding. + final dof = n - 3; + final double ampP; + if (ssTot <= 0) { + ampP = 1.0; + } else if (ssRes <= 0) { + ampP = 0.0; + } else { + final f = ((ssTot - ssRes) / 2) / (ssRes / dof); + ampP = _f2Survival(f, dof); + } + final phaseDefined = ampP < cosinorAcrophaseAlpha; + final conf = r2Adj.clamp(0.1, 0.95); return Metric( value: CosinorFit( mesor: mesor, amplitude: amplitude, - acrophaseRad: acro, - acrophaseHours: phaseHours, + acrophaseRad: phaseDefined ? acro : null, + acrophaseHours: phaseDefined ? phaseHours : null, periodHours: periodHours, r2: r2, r2Adj: r2Adj, + amplitudeP: ampP, ), confidence: conf, tier: Tier.high, inputs_used: inputs, - note: 'single-component cosinor; anchor circadian phase on the nocturnal ' - 'trough, not the daytime peak', + note: phaseDefined + ? 'single-component cosinor; anchor circadian phase on the nocturnal ' + 'trough, not the daytime peak' + : 'single-component cosinor; no resolvable rhythm ' + '(zero-amplitude test not rejected), acrophase undefined', ); } diff --git a/lib/src/onehz/clinical/hrv_time.dart b/lib/src/onehz/clinical/hrv_time.dart index 70cd001..d2b3510 100644 --- a/lib/src/onehz/clinical/hrv_time.dart +++ b/lib/src/onehz/clinical/hrv_time.dart @@ -785,8 +785,9 @@ double? _breathingRmssd(List nn, List times) { /// window into 5-min bins; required (returns absent without it). [windowMs] bin /// width (default 300 000 = 5 min). [minBeatsPerWindow] min NN diffs a window /// needs to contribute (default 5). [stageMaskPerSec] OPTIONAL per-second mask -/// (true = keep, e.g. NREM & immobile); a window is kept only when the mask is -/// true at the window's MIDPOINT second. +/// (true = keep, e.g. NREM & immobile) on the [nnTimesMs] clock — index s +/// covers nnTime s·1000 ms, NOT seconds since the first kept beat; a window is +/// kept only when the mask is true at the window's MIDPOINT second. /// /// Returns a Metric whose value is the median-of-windows RMSSD (ms). Keeps the /// PRV-not-ECG honesty note. Absent when there are too few usable windows, when @@ -840,7 +841,11 @@ Metric nocturnalRmssd( final indices = buckets.keys.toList()..sort(); for (final idx in indices) { if (stageMaskPerSec != null) { - final midSec = ((idx + 0.5) * windowMs / 1000.0).floor(); + // Windows are bucketed from the first KEPT beat, but `correctRr` keeps + // the elapsed time of a dropped leading artifact run, so t0 can sit + // minutes past the mask's origin. Read the mask on the beat clock or + // every window takes the stage label from t0 seconds earlier. + final midSec = ((t0 + (idx + 0.5) * windowMs) / 1000.0).floor(); final keep = midSec >= 0 && midSec < stageMaskPerSec.length && stageMaskPerSec[midSec]; diff --git a/lib/src/onehz/clinical/illness_cusum.dart b/lib/src/onehz/clinical/illness_cusum.dart index a8717e9..41fbe9a 100644 --- a/lib/src/onehz/clinical/illness_cusum.dart +++ b/lib/src/onehz/clinical/illness_cusum.dart @@ -88,6 +88,9 @@ class IllnessDay { /// so this is roughly 4× more specific than the method it used to cite. Say the /// number, not "conservative": red drives a medical-class notification. /// [persistDays] yellow nights required before escalating to red. +/// [maxGapNights] unscored nights (no RHR, or a baseline that abstained) that +/// may sit between two scored nights without breaking the run. A longer gap +/// zeroes the accumulator and every run counter. /// /// Returns a per-night list. Nights before [minBaseline] valid history are /// green with null cusum (we never flag without a baseline — no fabrication). @@ -101,6 +104,7 @@ List illnessCusum( int minBaseline = illnessCusumMinBaseline, double returnZ = 0.5, int recoverDays = 2, + int maxGapNights = 2, }) { final n = rhr.length; final out = []; @@ -174,7 +178,15 @@ List illnessCusum( // guard had just zeroed normalRun so recovery could not clear it either. // A CUSUM is evidence accumulated over consecutive observations; there are // no observations across a gap, so there is no evidence to carry. - if (day[i] - lastScoredDay > 1) { + // + // But a gap of one or two nights must not count as a break. Sleep goes + // undetected and the band comes off precisely while someone is sick, so + // resetting on any single unscored night wiped an episode's evidence at + // the moment it mattered: two +2.5 z nights (cusum 4.0), one night with no + // main sleep, and a third +2.5 night restarted from 0 and read green; on + // alternate-night gaps the alarm could never trip at all. A short gap + // carries the run untouched — it neither advances nor resets it. + if (day[i] - lastScoredDay > 1 + maxGapNights) { yellowRun = 0; normalRun = 0; cusum = 0; diff --git a/lib/src/onehz/clinical/load_trimp.dart b/lib/src/onehz/clinical/load_trimp.dart index e8c7b3d..3c89ca4 100644 --- a/lib/src/onehz/clinical/load_trimp.dart +++ b/lib/src/onehz/clinical/load_trimp.dart @@ -295,8 +295,10 @@ double strainFromNetTrimp(double net) { /// Headline 0–21 strain as a Metric, alongside the raw TRIMP (EST tier). /// -/// LUMP FORM — see [strainScore]. Kept for rescaling a stored TRIMP that has no -/// per-minute series; anything holding the series uses [strainScoreFromSeries]. +/// LUMP FORM — see [strainScore]. The guarded entry point for rescaling a +/// stored TRIMP that has no per-minute series; anything holding the series +/// uses [strainScoreFromSeries]. Calling bare [strainScore] for that rescale +/// skips the non-finite and above-[maxQuietHrr] refusals below. /// /// [trimp] the raw Banister TRIMP for the day/session, [wakeMinutes] the wake /// window it was accumulated over, [quietHrr] this user's quiet-waking level @@ -356,20 +358,39 @@ Metric strainScoreMetric( 'the day\'s own training away', ); } - return Metric( - value: strainScore(trimp, + return _strainEstimate( + strainScore(trimp, wakeMinutes: wakeMinutes, quietHrr: quietHrr, female: female), - confidence: quietSettled ? 0.6 : 0.45, - tier: Tier.estimate, - inputs_used: inputs, - note: quietSettled - ? 'headline 0–21 strain = log map of TRIMP earned above the ' - 'quiet-waking baseline; wrist-HR estimate' - : 'calibrating: quiet-waking level from fewer than ' - '$quietHrrSettledDays days', + inputs: inputs, + quietSettled: quietSettled, + settledNote: 'headline 0–21 strain = log map of TRIMP earned above the ' + 'quiet-waking baseline; wrist-HR estimate', ); } +/// The present-strain envelope shared by [strainScoreMetric] and +/// [strainScoreFromSeries]. One place, because the two are the lump and +/// guarded routes to the same headline: if the calibrating rule changed in one +/// copy only, a backfilled day and a live day with the same quiet level would +/// carry different confidence and the tests of the untouched copy would still +/// pass. +Metric _strainEstimate( + double value, { + required List inputs, + required bool quietSettled, + required String settledNote, +}) => + Metric( + value: value, + confidence: quietSettled ? 0.6 : 0.45, + tier: Tier.estimate, + inputs_used: inputs, + note: quietSettled + ? settledNote + : 'calibrating: quiet-waking level from fewer than ' + '$quietHrrSettledDays days', + ); + /// %HRR at or above which a minute is EXERCISE for the strain floor: ACSM's /// moderate-intensity floor (40 % HRR). Deliberately the same number as /// [maxQuietHrr] and as `Calories.activeHRRFraction`: one boundary between @@ -513,15 +534,12 @@ Metric strainScoreFromSeries( 'HRR in (0, $maxQuietHrr], and at least one worn minute', ); } - return Metric( - value: strainFromNetTrimp(n.net), - confidence: quietSettled ? 0.6 : 0.45, - tier: Tier.estimate, - inputs_used: inputs, - note: quietSettled - ? 'Banister TRIMP above your quiet-waking level, log-mapped to 0–21' - : 'calibrating: quiet-waking level from fewer than ' - '$quietHrrSettledDays days', + return _strainEstimate( + strainFromNetTrimp(n.net), + inputs: inputs, + quietSettled: quietSettled, + settledNote: + 'Banister TRIMP above your quiet-waking level, log-mapped to 0–21', ); } @@ -635,7 +653,7 @@ Metric personalQuietWakingHrr(List priorDailyLevels) { // TRIMP. // ════════════════════════════════════════════════════════════════════════════ -// StrainScorer — Banister TRIMP → 0–100 strain ("Effort"). +// StrainScorer — Banister TRIMP building blocks. // Implementation of published exercise-physiology methods. // // 1. Heart-Rate Reserve (Karvonen): HRR = HRmax − RHR. @@ -644,41 +662,24 @@ Metric personalQuietWakingHrr(List priorDailyLevels) { // from the real timestamps and capped at the stream's median cadence so a // gap is never counted as effort): sample contributes dur × x × y(x), with // y = 0.64·e^(1.92x) (men) / 0.86·e^(1.67x) (women). -// 4. strain = 100 × ln(TRIMP + 1) / ln(D), D = 7201. +// +// There is no strain scale here. The one strain is the 0–21 +// `strainScoreFromSeries` / `strainFromNetTrimp` path above, which weights with +// [StrainScorer.banisterY]. A 0–100 "Effort" map (`trimpToStrain` / +// `StrainScorer.strain` / `trimpStrain`) used to sit on top of these helpers — +// the second scale the file header says is gone. Nothing displayed it, yet its +// existence let a cadence check pass while measuring a number no screen shows. // // References: Karvonen 1957; Banister 1991 (y = 0.64·e^(1.92x) men / // 0.86·e^(1.67x) women); Tanaka 2001 (HRmax = 208 − 0.7·age). // -// NOTE (steps/active-energy floor): strain is PURELY HR-derived -// (Edwards/Banister TRIMP → log map). Steps and active calories are computed as -// SEPARATE, independent daily metrics and are NEVER fused into, nor floor, the -// strain score. +// NOTE (steps/active-energy floor): TRIMP here is PURELY HR-derived. Steps and +// active calories are computed as SEPARATE, independent daily metrics and are +// NEVER fused into, nor floor, the strain score. // ════════════════════════════════════════════════════════════════════════════ -/// StrainScorer constants + the strain pipeline. +/// StrainScorer constants + the Banister TRIMP helpers. class StrainScorer { - /// Minimum HR readings before computing strain on a DENSE stream (≈10 min @1Hz). - static const int minReadings = 600; - - /// Sparse-stream acceptance (#482): a low-cadence strap also qualifies once the - /// HR series SPANS at least [minSpanSeconds] with a small sample floor. - static const int minSparseReadings = 20; - - /// Wall-clock coverage (s) that qualifies a sparse stream (10 min). - static const int minSpanSeconds = 600; - - /// Top of the strain ("Effort") scale (rescaled 21.0 → 100.0). - static const double maxStrain = 100.0; - - /// Logarithmic-map denominator D = 7200 + 1. - /// - /// It was set from the Edwards daily ceiling (top weight 5 sustained 24 h = - /// 7200), and Edwards is gone (MOT-06). Kept unchanged anyway: the Banister - /// ceiling over the same 24 h is 1440 × 1 × y(1) = 6,293, within 15 % of it, - /// and this is a DISPLAY map — rescaling it would rescore every stored day for - /// no gain in truth. - static const double strainDenominator = 7201.0; - /// Fallback per-sample duration (minutes) — 1 s at 1 Hz. static const double fallbackSampleMin = 1.0 / 60.0; @@ -771,13 +772,14 @@ class StrainScorer { /// /// This replaces the old `sampleDurationMinutes`, which read ONE interval /// (the first two timestamps) and applied it to every sample — catastrophic - /// on exactly the sparse/irregular streams [minSparseReadings] admits: 21 - /// samples over 20 min with the first two 1 s apart scored strain 8.08 - /// instead of ~47, and a 1 Hz stream with a 5-min leading gap scored 104. + /// on sparse/irregular streams: 21 samples over 20 min with the first two + /// 1 s apart scored strain 8.08 instead of ~47, and a 1 Hz stream with a + /// 5-min leading gap scored 104. /// EMPTY when the stream's cadence cannot be measured — see /// [sampleCadenceSeconds]. Every duration here is a multiple of that cadence, - /// so an unmeasurable cadence is an unmeasurable effort, and [strain] turns - /// it into an absent metric rather than a small confident one. + /// so an unmeasurable cadence is an unmeasurable effort: a caller must turn + /// the empty list into an absent metric, never hand it to [banisterTRIMP] + /// (which credits [fallbackSampleMin] per sample when given no durations). static List sampleDurationsMinutes(List tsSec) { final n = tsSec.length; if (n == 0) return const []; @@ -806,127 +808,6 @@ class StrainScorer { } return acc; } - - // ── Logarithmic map ───────────────────────────────────────────────────────── - - /// Map accumulated TRIMP onto [0, 100] via 100 × ln(TRIMP+1) / ln(D), 2 dp. - /// TRIMP ≤ 0 → 0; above the D−1 ceiling the score is CLAMPED at [maxStrain] - /// (it used to run off the top of its own documented range: TRIMP 14400 → - /// 107.8, while the sibling [strainScore] clamped correctly). - static double trimpToStrain(double trimp, - {double denominator = strainDenominator}) { - if (trimp <= 0) return 0; - final value = maxStrain * math.log(trimp + 1.0) / math.log(denominator); - final clamped = math.min(maxStrain, math.max(0.0, value)); - return (clamped * 100).roundToDouble() / 100; - } - - // ── TRIMP method ────────────────────────────────────────────────────────────── - - /// Compute strain (0–100) from a time-ordered HR series. APPROXIMATE. - /// - /// [bpm] per-sample HR; [tsSec] their timestamps in SECONDS (same length). - /// Returns null when there isn't enough data to trust the number (fewer than - /// [minReadings] AND less than [minSpanSeconds] coverage), when [maxHR] is - /// missing, or when maxHR ≤ restingHR (invalid HRR). - /// - /// [maxHR] IS REQUIRED-AND-NULLABLE, and null means NO NUMBER (audit MOT-11). - /// It used to fall back to `defaultMaxHR()` = 190, a 220−age ceiling for a 30 - /// y/o applied to whoever's wrist arrived — a fabricated anchor inside the - /// value, caught (if at all) by the caller. Honesty belongs at the source, so - /// that helper is gone too. `AutoWorkoutDetector` keeps its own 190 for the - /// "did you work out?" gate, which is a prompt, not a published number. - static double? strain( - List bpm, - List tsSec, { - required double? maxHR, - double restingHR = defaultRestingHR, - bool female = false, - double denominator = strainDenominator, - }) { - if (maxHR == null) return null; - final effMax = maxHR; - final bool enoughData; - if (bpm.length >= minReadings) { - enoughData = true; - } else if (bpm.length >= minSparseReadings) { - if (tsSec.isEmpty) { - enoughData = false; - } else { - var mn = tsSec[0], mx = tsSec[0]; - for (final t in tsSec) { - if (t < mn) mn = t; - if (t > mx) mx = t; - } - enoughData = (mx - mn) >= minSpanSeconds; - } - } else { - enoughData = false; - } - if (!enoughData || effMax <= restingHR) return null; - - final durations = sampleDurationsMinutes(tsSec); - // No measurable cadence ⇒ no durations ⇒ no effort. Must NOT fall through: - // `banisterTRIMP` credits `fallbackSampleMin` per sample when handed an - // empty list, which is the fabricated 1 s this abstention exists to stop. - if (durations.isEmpty) return null; - final trimp = banisterTRIMP(bpm, restingHR, effMax - restingHR, durations, - female: female); - return trimpToStrain(trimp, denominator: denominator); - } -} - -/// Banister TRIMP strain ("Effort", 0–100) as a Metric, with the honesty -/// envelope: absent when the gates fail (never fabricated). -/// -/// [bpm] per-sample HR (bpm). [tsSec] timestamps (s), same length. [maxHr] / -/// [restingHr] the personal anchors (HRmax resolved by the caller via -/// [StrainScorer.estimateHRmax] / Tanaka). [sex] selects the Banister -/// coefficient. The `method` parameter is gone with the Edwards path (MOT-06). -Metric trimpStrain( - List bpm, - List tsSec, { - double? maxHr, - double? restingHr, - Sex sex = Sex.male, -}) { - const inputs = ['hr_series', 'resting_hr', 'max_hr']; - // used to fall back to StrainScorer.defaultRestingHR/defaultMaxHR when - // these were omitted, which fed a made-up anchor into a confident-looking - // ESTIMATE score with nothing telling anyone it wasn't a real anchor. those - // defaults are for AutoWorkoutDetector's internal gate, not for this public - // honesty-wrapped function - if we don't actually have the person's real - // anchors, we don't have a real number either. - if (maxHr == null || restingHr == null) { - return const Metric.absent( - tier: Tier.estimate, - inputs_used: inputs, - note: 'needs a real maxHr and restingHr, not fabricated defaults', - ); - } - final s = StrainScorer.strain( - bpm, - tsSec, - maxHR: maxHr, - restingHR: restingHr, - female: sex == Sex.female, - ); - if (s == null) { - return const Metric.absent( - tier: Tier.estimate, - inputs_used: inputs, - note: - 'strain needs ≥600 HR samples (or ≥20 spanning ≥600 s) and HRmax>RHR', - ); - } - return Metric( - value: s, - confidence: 0.6, - tier: Tier.estimate, - inputs_used: inputs, - note: 'Banister TRIMP → 0–100 strain ' - '(100·ln(TRIMP+1)/ln(7201)); wrist-HR ESTIMATE, not clinical', - ); } class LoadState { diff --git a/lib/src/onehz/clinical/prsa.dart b/lib/src/onehz/clinical/prsa.dart index e22a3bb..a1fec3a 100644 --- a/lib/src/onehz/clinical/prsa.dart +++ b/lib/src/onehz/clinical/prsa.dart @@ -60,38 +60,62 @@ class PrsaResult { } /// PRSA deceleration capacity. +/// +/// Pass [nnTimesMs] (aligned to [nnMs], e.g. `correctRr(...).nnTimes`) so no +/// anchor window is formed across a hole in the beat record. Metric decelerationCapacity(List nnMs, - {int l = 2, double anchorRatioCap = 0.05}) => - _prsa(nnMs, l: l, deceleration: true, anchorRatioCap: anchorRatioCap); + {int l = 2, double anchorRatioCap = 0.05, List? nnTimesMs}) => + _prsa(nnMs, + l: l, + deceleration: true, + anchorRatioCap: anchorRatioCap, + nnTimesMs: nnTimesMs); -/// PRSA acceleration capacity. +/// PRSA acceleration capacity. See [decelerationCapacity] for [nnTimesMs]. Metric accelerationCapacity(List nnMs, - {int l = 2, double anchorRatioCap = 0.05}) => - _prsa(nnMs, l: l, deceleration: false, anchorRatioCap: anchorRatioCap); + {int l = 2, double anchorRatioCap = 0.05, List? nnTimesMs}) => + _prsa(nnMs, + l: l, + deceleration: false, + anchorRatioCap: anchorRatioCap, + nnTimesMs: nnTimesMs); Metric _prsa( List nnMs, { required int l, required bool deceleration, required double anchorRatioCap, + List? nnTimesMs, }) { final kind = deceleration ? 'DC' : 'AC'; - final inputs = const ['rr_cleaned']; final n = nnMs.length; + final inputs = nnTimesMs == null + ? const ['rr_cleaned'] + : const ['rr_cleaned', 'beat_times']; + // Beat times that do not line up with the intervals cannot place a seam, and + // falling back to the untimed path would anchor across exactly the holes + // the caller passed them to guard against. + if (nnTimesMs != null && nnTimesMs.length != n) { + return Metric.absent( + tier: Tier.relative, + inputs_used: inputs, + note: 'PRSA ($kind): beat times not aligned with intervals', + ); + } if (l < 2) { // Bauer 2006 quantifies PRSA with the Haar contrast at wavelet scale s = 2: // DC = [X(0) + X(1) − X(−1) − X(−2)]/4, which needs TWO profile points on // each side of the anchor. With l < 2 the profile is too short (l = 1 used // to index profile[l−2] = profile[−1] and throw a RangeError). return Metric.absent( - tier: Tier.high, + tier: Tier.relative, inputs_used: inputs, note: 'PRSA ($kind) needs l≥2 for the s=2 Haar contrast', ); } if (n < 2 * l + 4) { return Metric.absent( - tier: Tier.high, + tier: Tier.relative, inputs_used: inputs, note: 'too few beats for PRSA ($kind)', ); @@ -99,8 +123,38 @@ Metric _prsa( // Anchor selection with an artifact-suppression ratio cap: |RR(i)/RR(i-1)-1| // must be ≤ anchorRatioCap so gross jumps (residual artifacts) don't anchor. + // + // Anchors must also sit in a stretch that is contiguous in TIME. correctRr + // drops artifact runs and re-anchors its clock across sensor holes, so two + // beats adjacent in the list can sit either side of a 40-minute gap; one + // such pair used to anchor on its own and publish a DC built from beats an + // hour apart. Same seam test as hrv_time.dart: a pair is contiguous iff the + // elapsed time between the two beat times is the interval itself. + final gapAware = nnTimesMs != null; + // seamBefore[j] = the step from beat j-1 to beat j crosses a hole. A beat + // time that runs backwards is a break in the record too; equal times are + // left contiguous (whole-second source clocks repeat them). + final seamBefore = List.filled(n, false); + if (gapAware) { + for (var j = 1; j < n; j++) { + final dt = nnTimesMs[j] - nnTimesMs[j - 1]; + seamBefore[j] = dt < 0 || dt > nnMs[j] + 0.5; + } + } final anchors = []; for (var i = l; i < n - l; i++) { + if (gapAware) { + // The averaged window is i-l .. i+l-1; every step inside it must be + // contiguous. + var crosses = false; + for (var j = i - l + 1; j <= i + l - 1; j++) { + if (seamBefore[j]) { + crosses = true; + break; + } + } + if (crosses) continue; + } final prev = nnMs[i - 1]; if (prev <= 0) continue; final delta = nnMs[i] - prev; @@ -112,7 +166,7 @@ Metric _prsa( } if (anchors.isEmpty) { return Metric.absent( - tier: Tier.high, + tier: Tier.relative, inputs_used: inputs, note: 'no valid PRSA anchors ($kind)', ); @@ -146,9 +200,13 @@ Metric _prsa( kind: kind, ), confidence: conf, - tier: Tier.high, + tier: Tier.relative, inputs_used: inputs, - note: 'Bauer 2006 PRSA $kind; needs a long overnight RR record for ' - 'mortality-grade stability', + // RELATIVE, not HIGH: the same subject read 7.6-9.9 ms on gen4 and 1.7 ms + // on W5 (header), so the number only trends against the user's own + // baseline on one device family. Callers must not chart it across a + // family change. + note: 'Bauer 2006 PRSA $kind; relative to own baseline on one device ' + 'family only', ); } diff --git a/lib/src/onehz/clinical/readiness_lnrmssd.dart b/lib/src/onehz/clinical/readiness_lnrmssd.dart index 296dcf5..9326651 100644 --- a/lib/src/onehz/clinical/readiness_lnrmssd.dart +++ b/lib/src/onehz/clinical/readiness_lnrmssd.dart @@ -9,18 +9,24 @@ // // HONESTY: needs ≥ minNights of valid history; absent otherwise. PRV. +import 'dart:math' as math; + import '../types.dart'; import '../util.dart'; class ReadinessLnRmssd { final double lnRmssdToday; final double rolling7Mean; + /// The window's MEDIAN — the centre z/SWC/band are taken against. final double rolling7Median; final double cvPct; // CV of lnRMSSD over the window final double? z; // (today - rollingMedian)/robustSD final double? swc; // smallest worthwhile change (0.5×robust SD, Plews-style) - final String band; // 'suppressed' | 'normal' | 'elevated' + /// 'suppressed' | 'normal' | 'elevated': tonight against centre ± spread × + /// [_lnRmssdBandWidening] of the prior-window size (wider than ± 1 spread on + /// short windows, see that function). + final String band; final bool saturationFlag; const ReadinessLnRmssd({ required this.lnRmssdToday, @@ -136,10 +142,30 @@ Metric readinessLnRmssd( // Plews SWC ≈ 0.5 × within-window SD (a small worthwhile change in lnRMSSD). final swc = 0.5 * spread; + // The band is a SINGLE-NIGHT decision, so its threshold has to be a + // single-night one: tonight outside centre ± 1 × spread, Kiviniemi 2007's + // daily reference value (baseline mean ± SD). The SWC above (0.5 × spread) + // is Plews' threshold for a 7-day ROLLING mean, where averaging has already + // cut night-to-night noise by ~1/sqrt(7); applied to one night it sits inside + // ordinary noise and a user whose HRV never changed read 'suppressed' or + // 'elevated' on ~70% of nights. Comparing a rolling mean instead does not + // work against this baseline: the week being tested is six-sevenths the same + // nights the centre comes from, so a real multi-night drop drags the centre + // with it and goes unflagged, and one old bad night keeps an ordinary + // tonight 'suppressed' for the rest of the week. At zero spread any + // departure from the centre is outside the band. + // + // "± 1 SD" only means a ~32% flag rate when the centre and spread are known. + // Here they come from 3-7 prior nights: the median misses the true level and + // the MAD of so few points is biased low and noisy, so a user whose HRV + // never changed still read 'suppressed'/'elevated' on 62% of nights with 3 + // prior nights and 45% with 7. [_lnRmssdBandWidening] restores the intended + // single-night rate for the window size actually used. + final halfWidth = spread * _lnRmssdBandWidening(priorWindow.length); final String band; - if (today < centre - swc) { + if (today < centre - halfWidth) { band = 'suppressed'; - } else if (today > centre + swc) { + } else if (today > centre + halfWidth) { band = 'elevated'; } else { band = 'normal'; @@ -170,3 +196,45 @@ Metric readinessLnRmssd( : 'Plews/Kiviniemi whole-night lnRMSSD readiness; PRV', ); } + +/// Half-width of the single-night band in units of the baseline's robust SD +/// (MAD × 1.4826), by number of prior nights. Each entry is the empirical +/// 68.27th percentile of |tonight − median| / (MAD × 1.4826) when tonight and +/// the prior nights are drawn from one fixed normal distribution (400k trials +/// per n), so a stable user lands outside the band on ~31.7% of nights, +/// Kiviniemi's mean ± SD rate with known parameters. It folds in both the +/// prediction-interval term for a new observation (√(1 + 1/n)) and the small-n +/// bias and noise of the MAD; 3 is above 2 because the MAD of three points is +/// 1.4826 × the smaller non-median deviation and runs especially small. +/// Interpolated log-linearly in n; past the last anchor it decays toward 1 as +/// 1/n, since the MAD then converges to σ. +const List<(int, double)> _lnRmssdBandWideningTable = [ + (2, 2.144), + (3, 2.851), + (4, 1.884), + (5, 1.687), + (6, 1.487), + (7, 1.418), + (8, 1.334), + (10, 1.257), + (14, 1.172), + (20, 1.119), + (28, 1.085), + (45, 1.051), + (90, 1.025), +]; + +double _lnRmssdBandWidening(int n) { + final first = _lnRmssdBandWideningTable.first; + final last = _lnRmssdBandWideningTable.last; + if (n <= first.$1) return first.$2; + if (n >= last.$1) return 1 + (last.$2 - 1) * last.$1 / n; + for (var i = 1; i < _lnRmssdBandWideningTable.length; i++) { + final (n1, f1) = _lnRmssdBandWideningTable[i]; + if (n > n1) continue; + final (n0, f0) = _lnRmssdBandWideningTable[i - 1]; + final t = (math.log(n) - math.log(n0)) / (math.log(n1) - math.log(n0)); + return math.exp(math.log(f0) + t * (math.log(f1) - math.log(f0))); + } + return 1.0; +} diff --git a/lib/src/onehz/clinical/stress_si.dart b/lib/src/onehz/clinical/stress_si.dart index a0fabc1..a9bee0f 100644 --- a/lib/src/onehz/clinical/stress_si.dart +++ b/lib/src/onehz/clinical/stress_si.dart @@ -96,16 +96,22 @@ Metric baevskyStressIndex( } // Window ~5 min of beats (256), 50% overlap; compute SI per window. - const win = 256, step = 128; + // MxDMn grows and AMo shrinks with window length, and the bands below are + // Baevsky's cut-points for ~5-min recordings. A gap-fragmented night yields + // many 30–127 beat segments whose SI runs systematically high (about 2× on + // the same beats), and a median over them reads 'elevated' on a normal + // night. So a window must hold at least [minWin] beats (~4 min at a resting + // HR) to count; shorter fragments are dropped, and if nothing qualifies the + // metric is absent. + const win = 256, step = 128, minWin = 200; final sis = []; final modes = [], amos = [], ranges = []; for (var s = 0; s < segStart.length; s++) { final lo = segStart[s]; final hi = (s + 1 < segStart.length ? segStart[s + 1] : nn.length); final seg = nn.sublist(lo, hi); - for (var start = 0; start + 30 <= seg.length; start += step) { + for (var start = 0; start + minWin <= seg.length; start += step) { final segWin = seg.sublist(start, math.min(start + win, seg.length)); - if (segWin.length < 30) break; final r = _siOfSegment(segWin, minRangeMs); if (r != null) { sis.add(r[0]); @@ -120,7 +126,7 @@ Metric baevskyStressIndex( return const Metric.absent( tier: Tier.estimate, inputs_used: inputs, - note: 'no valid Stress Index window', + note: 'no valid ~5-min Stress Index window', ); } final si = median(sis)!; diff --git a/lib/src/onehz/clinical/vo2max.dart b/lib/src/onehz/clinical/vo2max.dart index ea583df..15401f8 100644 --- a/lib/src/onehz/clinical/vo2max.dart +++ b/lib/src/onehz/clinical/vo2max.dart @@ -75,7 +75,8 @@ Metric vo2maxSubmaxEstimate({ ); } final hrReserve = hrMaxBpm - restingHrBpm; - if (hrReserve < 20) { + // NaN fails every comparison, so `< 20` alone lets a NaN reserve through. + if (!hrReserve.isFinite || hrReserve < 20) { return Metric.absent( tier: Tier.estimate, inputs_used: inputs, @@ -84,11 +85,13 @@ Metric vo2maxSubmaxEstimate({ } final hrr = (avgHrBpm - restingHrBpm) / hrReserve; if (!hrr.isFinite || hrr < vo2maxMinHrr || hrr > vo2maxMaxHrr) { + // round() throws on NaN/Infinity, so a non-finite %HRR must not reach it. + final pct = hrr.isFinite ? '${(hrr * 100).round()}%' : 'non-finite'; return Metric.absent( tier: Tier.estimate, inputs_used: inputs, note: 'bout intensity outside the submaximal %HRR band ' - '(${(hrr * 100).round()}%, need ' + '($pct, need ' '${(vo2maxMinHrr * 100).round()}-${(vo2maxMaxHrr * 100).round()}%)', ); } diff --git a/lib/src/onehz/human/alertness_forecast.dart b/lib/src/onehz/human/alertness_forecast.dart index 7220a78..290373b 100644 --- a/lib/src/onehz/human/alertness_forecast.dart +++ b/lib/src/onehz/human/alertness_forecast.dart @@ -94,10 +94,12 @@ class AlertnessForecast { const double _sUpper = 14.3; // ha — asymptote S recovers toward during sleep const double _sLower = 2.4; // la — asymptote S decays toward while awake const double _tauWake = 1 / 0.0353; // h — 28.33, wake decay (published rate) -const double _tauSleep = 1 / 0.3813; // h — 2.62, sleep recovery (published rate) +const double _tauSleep = + 1 / 0.3813; // h — 2.62, sleep recovery (published rate) const double _cAmplitude = 2.5; // Ca — 24 h circadian amplitude const double _uAmplitude = 0.5; // Ua — 12 h ultradian amplitude const double _uMean = -0.5; // Um — ultradian offset +const double _uTroughLeadHours = 3.0; // U trough lead before the C acrophase const double _wInertia = 5.72; // Wc — sleep-inertia magnitude at wake const double _tauInertia = 0.67; // Wd — h, inertia decay @@ -177,10 +179,11 @@ Metric alertnessForecast({ final s0 = _sUpper - (_sUpper - _sLower) * math.exp(-sleepDurationHours / _tauSleep); - final phase = - (circadianAcrophaseHours != null && circadianAcrophaseHours.isFinite) - ? circadianAcrophaseHours - : _defaultAcrophaseHours; + // One predicate decides both the curve and its disclosure: a non-finite + // phase falls back to the population curve, so the note must say so too. + final fittedPhase = + circadianAcrophaseHours != null && circadianAcrophaseHours.isFinite; + final phase = fittedPhase ? circadianAcrophaseHours : _defaultAcrophaseHours; // THE CURVE ENDS WHEN THE DAY DOES. Until 2026-08-17 it ran a fixed 18 h and // the trough search below could only look at windows that fitted inside that @@ -205,10 +208,16 @@ Metric alertnessForecast({ : _sLower + (s - _sLower) * math.exp(-_stepHours / _tauWake); } final c = _cAmplitude * math.cos(2 * math.pi * (clock - phase) / 24.0); - // Process U — the model's 12 h ultradian harmonic, phase-locked to the same - // acrophase. Missing until 2026-08-17, which flattened the mid-day wobble. - final u = - _uAmplitude * math.cos(2 * math.pi * (clock - phase) / 12.0) + _uMean; + // Process U — the model's 12 h ultradian harmonic. It exists to put the + // post-lunch dip in, so its TROUGH sits 3 h before the circadian acrophase + // (Åkerstedt, Folkard & Portin 2004: U = Um − Ua·cos(2π(t − (p − 3))/12)), + // ~13:48 at the default phase. Until 2026-10-17 it was a peak at the + // acrophase with no offset, which put the dip at 10:48 and a bump at the + // real early-afternoon low. + final u = _uMean - + _uAmplitude * + math.cos( + 2 * math.pi * (clock - (phase - _uTroughLeadHours)) / 12.0); final w = _wInertia * math.exp(-(i * _stepHours) / _tauInertia); raw.add(s + c + u - w); } @@ -269,7 +278,7 @@ Metric alertnessForecast({ tier: Tier.estimate, inputs_used: inputs, note: 'Three-Process Model shape only — NO SCORE. ' - '${circadianAcrophaseHours == null ? _assumedPhaseNote : _fittedPhaseNote}' + '${fittedPhase ? _fittedPhaseNote : _assumedPhaseNote}' 'The low window is where this model bottoms out over your waking day — ' 'usually its last hours. It follows from when you woke and how long you ' 'slept, and is not a dip measured on you. ' diff --git a/lib/src/onehz/human/coaching.dart b/lib/src/onehz/human/coaching.dart index 4db4216..8e00721 100644 --- a/lib/src/onehz/human/coaching.dart +++ b/lib/src/onehz/human/coaching.dart @@ -210,11 +210,15 @@ Metric strainTarget({ targetMin: lo, targetMax: hi, band: band, - rationale: 'Target shaped by recovery and recent load.', + // Without a built-up fitness-fatigue model the window is recovery + // alone; naming load here would claim an input that never entered it. + rationale: hasLoad + ? 'Target shaped by recovery and recent load.' + : 'Target shaped by recovery.', ), confidence: 0.6, tier: Tier.estimate, - inputs_used: const ['recovery', 'load'], + inputs_used: hasLoad ? const ['recovery', 'load'] : const ['recovery'], ); } @@ -619,7 +623,9 @@ class JournalNumericEffect { /// count alone does not explain — currently only /// `need_history:have=N,need=M`, meaning the permutation test's own /// discreteness floor could not clear the multiplicity correction at N days - /// however strong the effect was, and M days would. Null otherwise. + /// however strong the effect was, and M days would; or + /// `unreachable_floor:have=N` when no day count within a season would. + /// Null otherwise. final String? note; const JournalNumericEffect({ @@ -838,12 +844,21 @@ double _binom(int n, int k) { /// field at exactly 8 days CANNOT be published however real the effect. At n = 9 /// it can (1/126 × 4 = 0.032). A test that cannot pass is not a "no", and the /// silence with no reason attached is the part that was wrong. -double _permPFloor(int n, int n1, int b) { +double _permPFloor(int n, int n1, int b) => + math.max(1.0 / (b + 1), _labelPFloor(n, n1)); + +/// The data's own part of [_permPFloor]: the floor with every relabelling +/// enumerated, before any sampling floor is added. +double _labelPFloor(int n, int n1) { final c = _binom(n, n1); final k = (n1 * 2 == n) ? 2.0 : 1.0; - return math.max(1.0 / (b + 1), c <= 0 ? 1.0 : k / c); + return c <= 0 ? 1.0 : k / c; } +/// Upper bound on the family-scaled shuffle count in +/// [journalNumericCorrelations]; enough for a family of 500 tests at 0.05. +const int _maxPermutations = 19999; + /// The smallest number of paired days at which this test could clear [alpha] /// after a family-of-[m] BH correction, holding the observed group balance. /// Null if it cannot get there within a season. @@ -922,8 +937,8 @@ List journalNumericCorrelations({ double? cohensD, int? nWith, int? nWithout, - double? p, - double pFloor, + double Function(int b)? pAt, + double dataFloor, int nSmallerSide, bool floorOk, bool insufficient, @@ -942,8 +957,8 @@ List journalNumericCorrelations({ cohensD: null, nWith: null, nWithout: null, - p: null, - pFloor: 1.0, + pAt: null, + dataFloor: 1.0, nSmallerSide: 0, floorOk: false, insufficient: true, @@ -1027,8 +1042,8 @@ List journalNumericCorrelations({ cohensD: d, nWith: withIt.length, nWithout: without.length, - p: _permTwoSampleP(ys, inGroup, permutations, permutationSeed), - pFloor: _permPFloor(n, withIt.length, permutations), + pAt: (b) => _permTwoSampleP(ys, inGroup, b, permutationSeed), + dataFloor: _labelPFloor(n, withIt.length), nSmallerSide: math.min(withIt.length, without.length), // Both sides exactly constant leaves d undefined; the per-side floor // is already the tag path's zero-spread rule, so a real gap between @@ -1077,8 +1092,8 @@ List journalNumericCorrelations({ // result for CONTINUOUS data applied to tie-heavy ordinal journal fields. // See [_permSpearmanP]. The interval above is unchanged and still what a // screen displays. - final p = - n >= 4 ? _permSpearmanP(xs, ys, permutations, permutationSeed) : null; + final pAt = + n >= 4 ? (int b) => _permSpearmanP(xs, ys, b, permutationSeed) : null; rows.add(( field: field, @@ -1093,11 +1108,11 @@ List journalNumericCorrelations({ cohensD: null, nWith: null, nWithout: null, - p: p, + pAt: pAt, // A rank permutation has n! relabellings — 40,320 already at the n = 8 - // minimum, where 2/n! is far below the sampling floor — so in practice - // only 1/(b+1) binds. Nothing like the binary path's 1/C(n, n₁). - pFloor: math.max(1.0 / (permutations + 1), 2.0 / _factorial(n)), + // minimum — so the data never binds here; only the shuffle count could. + // Nothing like the binary path's 1/C(n, n₁). + dataFloor: 2.0 / _factorial(n), nSmallerSide: n, floorOk: rho.abs() >= minAbsRho, insufficient: false, @@ -1105,25 +1120,39 @@ List journalNumericCorrelations({ } } - final qs = benjaminiHochberg([for (final r in rows) r.p]); // The BH family size, so each row can be asked whether it COULD have passed. - // Computed once from the p's as they stand: dropping the unreachable rows and - // re-running would only shrink m and let more through, which is the wrong - // direction for a gate whose job is to refuse. - final m = rows.where((r) => r.p != null).length; + // Computed once from the rows as they stand: dropping the unreachable rows + // and re-running would only shrink m and let more through, which is the + // wrong direction for a gate whose job is to refuse. + final m = rows.where((r) => r.pAt != null).length; + // The shuffle count scales with the family. 1/(b+1) is a resolution limit of + // the simulation, not a property of the data: with b fixed at 999 it alone + // failed the gate below once m reached 51 (13 fields × 4 outcomes), so one + // more logged field blanked every numeric finding behind a "need more days" + // note that no amount of logging could satisfy, since that term does not + // depend on n. 2m/alpha keeps the sampling floor at half the BH budget, and + // the cap only bounds the cost for a journal with hundreds of fields. + // A non-positive alpha admits nothing, and 2m/0 is not a shuffle count + // (ceil() of infinity throws), so it takes the cap. + final want = fdrAlpha > 0 ? 2 * m / fdrAlpha : double.infinity; + final b = math.max( + permutations, want >= _maxPermutations ? _maxPermutations : want.ceil()); + final ps = [for (final r in rows) r.pAt?.call(b)]; + final qs = benjaminiHochberg(ps); final byField = >{}; for (var i = 0; i < rows.length; i++) { final r = rows[i]; + final p = ps[i]; final q = qs[i]; // STRUCTURALLY UNREACHABLE ⇒ ABSTAIN, and say how many days it would take. // A permutation p cannot go below [_permPFloor]; if even that floor fails // the correction, this test could not have been published however strong // the effect, and reporting it as a computed-but-weak "no" is a verdict we // never actually reached. - final unreachable = r.p != null && r.pFloor * m > fdrAlpha; + final floor = math.max(1.0 / (b + 1), r.dataFloor); + final unreachable = p != null && floor * m > fdrAlpha; final need = unreachable - ? _needForFloor(r.n, r.nSmallerSide, permutations, m, fdrAlpha, - binary: r.binary) + ? _needForFloor(r.n, r.nSmallerSide, b, m, fdrAlpha, binary: r.binary) : null; byField.putIfAbsent(r.field, () => []).add( JournalNumericEffect( @@ -1135,15 +1164,19 @@ List journalNumericCorrelations({ n: r.n, insufficient: r.insufficient || unreachable, meaningful: !unreachable && r.floorOk && q != null && q <= fdrAlpha, + // No N within a season clears the floor (a very lopsided binary + // split): say so rather than promise a day count that would not. note: unreachable - ? 'need_history:have=${r.n},need=${need ?? r.n + 90}' + ? (need != null + ? 'need_history:have=${r.n},need=$need' + : 'unreachable_floor:have=${r.n}') : null, binary: r.binary, delta: r.delta, cohensD: r.cohensD, nWith: r.nWith, nWithout: r.nWithout, - p: r.p, + p: p, q: q, ), ); diff --git a/lib/src/onehz/human/overreaching_conjunction.dart b/lib/src/onehz/human/overreaching_conjunction.dart index 4cde0c3..7c403aa 100644 --- a/lib/src/onehz/human/overreaching_conjunction.dart +++ b/lib/src/onehz/human/overreaching_conjunction.dart @@ -15,6 +15,8 @@ // NO NEW SCORE AND NO NEW NOTIFICATION CLASS. In-app feed only. A push channel // would turn this into exactly the forward alerting MT-10 refuses. +import '../clinical/illness_cusum.dart' + show belowQuantumNote, degenerateBaselineNote; import '../clinical/load_trimp.dart'; import '../foundations/baseline.dart'; import '../types.dart'; @@ -114,13 +116,26 @@ Metric overreachingConjunction({ final rhrBelowQuantum = dispersionBelowQuantum(rhrBaselineWindow, 1.0); final gate = (!rhrBelowQuantum && scale != null && scale > 0) ? 0.5 * scale : null; - if (centre == null || gate == null || !base.sufficient) { + if (centre == null || !base.sufficient) { return Metric.absent( tier: Tier.estimate, inputs_used: inputs, note: needBaselineNote(have: base.nValid, need: minBaseline), ); } + if (gate == null) { + // Enough nights, but no usable dispersion to gate on. This is NOT a cold + // start: more nights of the same whole-bpm values never clear it, so a + // need_baseline note here read as "still collecting" indefinitely. + final sd = stddev(rhrBaselineWindow) ?? 0; + return Metric.absent( + tier: Tier.estimate, + inputs_used: inputs, + note: (rhrBelowQuantum && sd > 0) + ? belowQuantumNote + : degenerateBaselineNote, + ); + } var elevated = 0; for (final v in observed) { diff --git a/lib/src/onehz/human/percentile_of_you.dart b/lib/src/onehz/human/percentile_of_you.dart index 86b54f7..66b9edf 100644 --- a/lib/src/onehz/human/percentile_of_you.dart +++ b/lib/src/onehz/human/percentile_of_you.dart @@ -1,19 +1,18 @@ -// HUMAN LAYER — Percentile-of-you and records. -// Catalog §D "Percentile-of-you + records" [PUB order-statistics]. The streak -// half of that catalog line is deliberately NOT implemented: no streaks in this -// product, so there is nothing here to wire one up from. +// HUMAN LAYER — Percentile-of-you. +// Catalog §D "Percentile-of-you + records" [PUB order-statistics]. Neither the +// streak nor the records half of that catalog line is implemented: there are no +// streaks in this product and no personal-record surface reads one, so an +// unused detector would only be one more thing a baseline change has to keep +// green. // // n-of-1 only: a value's rank is computed against the SAME USER's recent -// history (a robust empirical CDF), never a population leaderboard. Records are -// MDC-gated: a new high/low is only a "record" when it beats the prior extreme -// by more than the metric's minimal detectable change — otherwise it is noise. +// history (a robust empirical CDF), never a population leaderboard. // // HONESTY: within-user percentiles (catalog honesty rule), no validity -// exposure, MDC-gated records, "—" when history is too short. +// exposure, "—" when history is too short. import '../types.dart'; import '../util.dart'; -import '../foundations/baseline.dart'; /// Which direction is "better" for a metric. Required by every function here /// that puts a WORD on a number: a rank is direction-free, a label is not. @@ -99,80 +98,3 @@ String _bandLabel(double pct, Better better) { if (good > 10) return 'worse than usual'; return 'among your worst'; } - -class RecordCheck { - final bool isRecord; // beat the prior extreme by > MDC - final String kind; // 'high' | 'low' | 'none' - final double value; - final double? priorExtreme; - final double? margin; // |value - priorExtreme| - final double? mdc; // gate threshold used - const RecordCheck({ - required this.isRecord, - required this.kind, - required this.value, - required this.priorExtreme, - required this.margin, - required this.mdc, - }); - Map toJson() => { - 'is_record': isRecord, - 'kind': kind, - 'value': round6(value), - if (priorExtreme != null) 'prior_extreme': round6(priorExtreme!), - if (margin != null) 'margin': round6(margin!), - if (mdc != null) 'mdc': round6(mdc!), - }; -} - -/// Is [value] a personal record vs [history], gated by MDC? -/// -/// A record requires the value to beat the prior extreme (max for [Better.higher], -/// min for [Better.lower]) by MORE than the metric's minimal detectable change, -/// so regression-to-mean noise never gets celebrated. When no MDC can be -/// established (degenerate/quantized scale) we refuse to call a record. -Metric personalRecord( - double value, - List history, { - required Better better, - double? typicalError, - int minN = 14, -}) { - const inputs = ['metric_history']; - if (history.length < minN || better == Better.neither) { - return const Metric.absent( - tier: Tier.relative, - inputs_used: inputs, - note: 'records need history and a defined "better" direction', - ); - } - final baseline = robustBaseline(history, minValid: minN); - final gate = mdc(baseline, typicalError: typicalError); - if (gate == null) { - return const Metric.absent( - tier: Tier.relative, - inputs_used: inputs, - note: 'no MDC (degenerate scale) — refusing to claim a record', - ); - } - final prior = better == Better.higher - ? history.reduce((a, b) => a > b ? a : b) - : history.reduce((a, b) => a < b ? a : b); - final beats = better == Better.higher ? value - prior : prior - value; - final isRecord = beats > gate; - final kind = isRecord ? (better == Better.higher ? 'high' : 'low') : 'none'; - return Metric( - value: RecordCheck( - isRecord: isRecord, - kind: kind, - value: value, - priorExtreme: prior, - margin: beats, - mdc: gate, - ), - confidence: (history.length / 60.0).clamp(0.3, 0.9), - tier: Tier.relative, - inputs_used: inputs, - note: 'record only if it beats your prior extreme by > MDC', - ); -} diff --git a/lib/src/onehz/human/sleep_regularity.dart b/lib/src/onehz/human/sleep_regularity.dart index f0f5d21..a790844 100644 --- a/lib/src/onehz/human/sleep_regularity.dart +++ b/lib/src/onehz/human/sleep_regularity.dart @@ -28,6 +28,13 @@ import '../types.dart'; import '../util.dart'; +/// Unconstrained nights needed before [sleepDebt] reports a p75. With one +/// night the "percentile" is that night verbatim, so a single 4.5 h or 11 h +/// Saturday became the user's sleep need — and with it tonight's need, the +/// recommended bedtime and sleep performance. Three is the fewest at which a +/// p75 is not just the maximum or the only value. +const int sleepDebtMinFreeNights = 3; + class SleepDebt { /// p75 of your UNCONSTRAINED nights (h) — "where your longer free nights /// land". NOT Kitamura's optimal sleep duration; see the file header. @@ -55,20 +62,32 @@ class SleepDebt { /// [recentSleepH] recent nightly sleep durations (h), oldest→newest. /// [freeNightSleepH] sleep durations on UNCONSTRAINED nights (no alarm / free /// days). We take their 75th percentile. With no free nights we report habitual -/// sleep but DECLINE to compare (honest: "no free night yet"). +/// sleep but DECLINE to compare (honest: "no free night yet"); the same with +/// fewer than [minFree] free nights, where the percentile is a single night. /// /// NOT Kitamura's OSD — see the file header for what that would require. Metric sleepDebt( List recentSleepH, List freeNightSleepH, { int minRecent = 3, + int minFree = sleepDebtMinFreeNights, }) { const inputs = ['sleep_duration_recent', 'sleep_duration_free']; + // A NaN night is not a night: it must not count toward a minimum or reach + // the median/percentile, where it either sorts to an end or poisons the sum. + recentSleepH = [ + for (final h in recentSleepH) + if (h.isFinite) h + ]; + freeNightSleepH = [ + for (final h in freeNightSleepH) + if (h.isFinite) h + ]; if (recentSleepH.length < minRecent) { - return const Metric.absent( + return Metric.absent( tier: Tier.high, inputs_used: inputs, - note: 'need more recent nights to estimate habitual sleep', + note: 'need_baseline:have=${recentSleepH.length},need=$minRecent', ); } final habitual = median(recentSleepH)!; @@ -82,6 +101,19 @@ Metric sleepDebt( note: 'no free night yet — cannot honestly estimate your sleep need', ); } + // A p75 of one or two nights is just those nights: one short Saturday read + // as "you need 2.5 h less than you sleep", one lie-in as 2.3 h of debt that + // raised tonight's need and bedtime. Below [minFree] we still say a free + // night exists but decline to compare. + if (freeNightSleepH.length < minFree) { + return Metric( + value: SleepDebt(null, habitual, null, true), + confidence: 0.3, + tier: Tier.high, + inputs_used: inputs, + note: 'need_baseline:have=${freeNightSleepH.length},need=$minFree', + ); + } // The upper part of free-night durations: "when unconstrained, this is where // your longer nights land". A percentile of ad-lib nights, NOT a rebound // asymptote — see the file header. diff --git a/lib/src/onehz/motion/steps.dart b/lib/src/onehz/motion/steps.dart index cf8c7a2..ea80abb 100644 --- a/lib/src/onehz/motion/steps.dart +++ b/lib/src/onehz/motion/steps.dart @@ -721,11 +721,18 @@ class DailyMovementEstimate { /// deliberately NO constant fallback: a constant absolute floor is precisely /// the failure this design removes. /// -/// The note counts DAYS, against [personalDynFloorMinDays]. It used to count +/// The note counts DAYS, against [requiredHistoryDays]. It used to count /// them against [personalDynFloorMinMinutes] (2000) while every caller passed a /// day count, so a user three days in was told `have=3,need=2000` — ~1997 more /// of a unit nobody was measuring, for a metric that unlocks in two. /// +/// [requiredHistoryDays] must be the gate the caller actually applies before +/// it stops passing a null floor. It defaults to [enrollmentDaysForFrozenFloor] +/// because a caller that freezes the floor keeps it null until then; reporting +/// against [personalDynFloorMinDays] made days 5-13 read `have=9,need=5`, a +/// met gate on a metric that was still absent. A caller that uses the bare +/// [personalDynFloorFromDailySummaries] median passes that function's minimum. +/// /// Tier is always ESTIMATE. /// /// This function is intentionally independent of the 100 Hz pedometer: it @@ -738,6 +745,7 @@ Metric dailyActiveMinutes( double minSamplesPerMinute = 30, int minBoutMin = 3, int historyDaysAvailable = 0, + int requiredHistoryDays = enrollmentDaysForFrozenFloor, int? expectedMinutes, }) { const inputs = ['dyn_amp_per_min', 'personal_dyn_floor']; @@ -756,10 +764,15 @@ Metric dailyActiveMinutes( return Metric.absent( tier: Tier.estimate, inputs_used: inputs, - note: needBaselineNote( - have: historyDaysAvailable, - need: personalDynFloorMinDays, - ), + // With the history gate already met, a missing or unusable floor is + // not a baseline shortfall; a met `need_baseline` note would read as + // "unlocked" on a metric that is still absent. + note: historyDaysAvailable < requiredHistoryDays + ? needBaselineNote( + have: historyDaysAvailable, + need: requiredHistoryDays, + ) + : 'no usable personal movement floor', ); } // NO CEILING. A `floor x 3` upper bound used to reject "vigorous non-gait" diff --git a/lib/src/onehz/respiration/cvhr_apnea.dart b/lib/src/onehz/respiration/cvhr_apnea.dart index 5eef9b2..fe1ada7 100644 --- a/lib/src/onehz/respiration/cvhr_apnea.dart +++ b/lib/src/onehz/respiration/cvhr_apnea.dart @@ -183,7 +183,15 @@ Metric cvhrApneaScreen( } final perHour = cycleCount / analyzedHours; - final conf = ((1 - artifactFraction) * 0.85).clamp(0.2, 0.85); + // Confidence also scales with how much time was analysed. Any gap-free + // stretch of a minute counts, so a night that dropout left with one clean + // 72 s segment and a single dip read 50 cycles/h at 0.85 and tier HIGH. Full + // weight needs the same observed hours the 30-night screen asks of a night; + // below that the per-hour index is an estimate off a short sample. + final coverage = (analyzedHours / minAnalyzedHoursPerNight).clamp(0.0, 1.0); + final conf = ((1 - artifactFraction) * 0.85 * coverage).clamp(0.01, 0.85); + final tier = + analyzedHours < minAnalyzedHoursPerNight ? Tier.estimate : Tier.high; return Metric( value: CvhrResult( cycleCount: cycleCount, @@ -195,7 +203,7 @@ Metric cvhrApneaScreen( widthQuartilesSec: _quartiles(widths), ), confidence: conf, - tier: Tier.high, + tier: tier, inputs_used: inputs, note: 'CVHR/ACAT (Hayano) apnea SCREEN — NOT a diagnosis, NOT an AHI; ' 'single-night CVHR has substantial night-to-night variability, ' @@ -272,10 +280,15 @@ class CvhrDistribution { /// mean is being read against. final List quartiles; - /// True when the recent nights sit above the retained distribution's own p75. - /// This is the ONLY comparative statement the screen makes, and it is a - /// statement about this user's own history, not about apnea. - final bool aboveOwnUsual; + /// True when the recent nights sit clearly above the OLDER retained nights: + /// above their p75 and by more than their own spread allows for. False when + /// they do not. NULL while there are fewer than + /// [minReferenceNightsForCvhrComparison] older nights: a false there would be + /// read as "stayed inside your range" even for a jump from 5/h to 40/h, and + /// the honest answer is that no comparison exists yet. This is the ONLY + /// comparative statement the screen makes, and it is a statement about this + /// user's own history, not about apnea. + final bool? aboveOwnUsual; const CvhrDistribution({ required this.nightsUsed, @@ -304,14 +317,62 @@ const double minAnalyzedHoursPerNight = 4.0; /// And this many such nights must exist before the screen says anything. const int minNightsForCvhrDistribution = 5; -/// Nights considered, newest-first, once the gates have run. +/// Nights considered, newest-first, once the gates have run. Also the calendar +/// span, in days ending on the newest night, that a night must fall inside. const int cvhrDistributionWindowNights = 30; +/// The comparison needs this many OLDER nights behind the recent ones. It is +/// a separate, higher floor than [minNightsForCvhrDistribution]: at the 5-night +/// distribution floor the reference would be 3 nights, and a stable user then +/// read as "running higher" on about one evaluation in five. +const int minReferenceNightsForCvhrComparison = 6; + +/// Floor on the reference spread, in cycles/h. Identical or quantized nights +/// give MAD = 0, and a zero spread turns any difference at all into a verdict. +/// A 1/h difference in an AHI-scaled index is below any published category +/// boundary (5/15/30). Candidate for personalisation, not a validated constant. +const double minCvhrReferenceSpreadPerHour = 1.0; + +/// One-sided 95% normal quantile for the recent-vs-reference gap. +const double _cvhrComparisonZ = 1.645; + +/// Whether the recent nights sit clearly above the OLDER nights behind them. +/// +/// The recent nights are held out of [reference]: comparing them with a range +/// that contains them makes a single high night its own evidence. Above the +/// reference p75 alone is not enough either — with a handful of nights p75 is +/// close to an order statistic and ordinary night-to-night noise clears it. +/// The recent mean must clear that p75 by z·σ̂·√(1/m + π/(2r)), the standard +/// error of an m-night mean against an r-night median (π/2 is the median's +/// variance inflation for normal data), σ̂ = MAD×1.4826. Measuring the gap +/// from the median instead let a stable user through 5-10% of the time; from +/// p75 it stays at or below what the plain p75 test gave on long histories. +/// +/// Null, not false, when the reference is under the floor: there is nothing +/// to compare against yet, and false would claim the recent nights stayed +/// inside a range nobody measured. +bool? _recentAboveReference({ + required double recentMean, + required int recentNights, + required List reference, +}) { + if (reference.length < minReferenceNightsForCvhrComparison || + recentNights < 1) { + return null; + } + final p75 = percentile(reference, 75)!; + final spread = math.max(mad(reference)!, minCvhrReferenceSpreadPerHour); + final se = + spread * math.sqrt(1 / recentNights + math.pi / (2 * reference.length)); + return recentMean - p75 > _cvhrComparisonZ * se; +} + /// RESP-01 — the 30-night personal CVHR distribution. /// /// [nights] in any order; the newest [cvhrDistributionWindowNights] retained -/// nights are used. Returns ABSENT — never a partial number — until at least -/// [minNightsForCvhrDistribution] nights clear both gates. +/// nights that fall within [cvhrDistributionWindowNights] calendar days of the +/// newest night are used. Returns ABSENT — never a partial number — until at +/// least [minNightsForCvhrDistribution] nights clear both gates. Metric cvhrPersonalDistribution(List nights) { const inputs = ['cvhr_apnea', 'irregular_rhythm_flag']; // Newest first, one row per day (a re-derive can leave two). @@ -322,11 +383,26 @@ Metric cvhrPersonalDistribution(List nights) { final sorted = byDay.values.toList() ..sort((a, b) => b.dayKey.compareTo(a.dayKey)); + // Bounded by calendar days, not just by count: after a wear gap the newest + // 30 rows reach back months, and the "recent" third then mixes the first + // nights back with nights from before the gap. Anchored on the newest night, + // never the wall clock. A label that does not parse is not bounded, which is + // the old count-only behaviour. + final anchorDay = + sorted.isEmpty || DateTime.tryParse(sorted.first.dayKey) == null + ? null + : calendarDays([sorted.first.dayKey])[0]; + var excludedIrregular = 0; var excludedThin = 0; final kept = []; for (final n in sorted) { if (kept.length >= cvhrDistributionWindowNights) break; + if (anchorDay != null && DateTime.tryParse(n.dayKey) != null) { + final day = calendarDays([n.dayKey])[0]; + // Newest-first, so every later night is older still. + if (anchorDay - day >= cvhrDistributionWindowNights) break; + } if (!n.cvhrPerHour.isFinite || !n.analyzedHours.isFinite) continue; // RESP-02 cross-gate FIRST: an AF night is not a thin night, and counting // it as one would misreport why the screen is quiet. @@ -383,7 +459,11 @@ Metric cvhrPersonalDistribution(List nights) { weightedMean: all, recentWeightedMean: rec, quartiles: q, - aboveOwnUsual: rec > q[2], + aboveOwnUsual: _recentAboveReference( + recentMean: rec, + recentNights: recent.length, + reference: [for (final n in kept.skip(recentCount)) n.cvhrPerHour], + ), ), // Relative tier on purpose: this compares the user only with themselves. tier: Tier.relative, @@ -394,7 +474,8 @@ Metric cvhrPersonalDistribution(List nights) { inputs_used: inputs, note: 'CVHR is a CARDIAC SURROGATE, not a breathing measurement: this is ' 'how often the pattern showed up across ${kept.length} of your own ' - 'nights, never an AHI, never a severity, never a per-night value. ' + 'nights in the last $cvhrDistributionWindowNights days, never an AHI, ' + 'never a severity, never a per-night value. ' 'It also fires on atrial fibrillation (those nights are excluded), on ' 'periodic breathing at altitude and on any fragmented night, and it is ' 'blunted by beta-blockers, autonomic neuropathy and diabetes — so a ' diff --git a/lib/src/onehz/respiration/resp_rate.dart b/lib/src/onehz/respiration/resp_rate.dart index 7b5b7dc..6cd3452 100644 --- a/lib/src/onehz/respiration/resp_rate.dart +++ b/lib/src/onehz/respiration/resp_rate.dart @@ -246,6 +246,7 @@ Metric rsaRespRate( final peakHz = []; // the same peaks in Hz, index-aligned final peakPwr = []; // their spectral power, index-aligned var atCeiling = 0; + var atFloor = 0; // The ceilings those sub-windows peaked at, br/min — their own, which // move with heart rate; not the whole input's. var ceilingLo = double.infinity, ceilingHi = 0.0; @@ -304,11 +305,28 @@ Metric rsaRespRate( ceilingHi = math.max(ceilingHi, segHi * 60); continue; } + // The same absence at the BOTTOM edge. A true rate below [rsaLoHz] (slow + // deep-sleep or paced breathing at 7.5–9 br/min) leaves only the skirt of + // its main lobe inside the band, so the in-band maximum sits on or just + // above the floor. Every sub-window then lands near that one value, the + // agreement gate passes, and ~9 br/min went out at confidence 0.8–0.9 for + // a true 8. One resolution width (1/span) was not enough: noise on the + // skirt still put sub-windows at 9.2–10 br/min. Two widths of a full + // [rsaSegmentSec] sub-window (0.4 br/min) withholds the whole skirt; the + // cost is that a true 9.0–9.4 is withheld too, which is an absence rather + // than a wrong number. The guard is capped at that width rather than + // scaling with 1/segSpan: on a short input (a 3-min window gives ~90 s + // sub-windows) 2/segSpan is 1.33 br/min and blanked a true 10 br/min, + // while the floor failure there already sits on the first bin. + if (pk <= rsaLoHz + math.min(2 / segSpan, 2 / rsaSegmentSec)) { + atFloor++; + continue; + } peaks.add(pk * 60.0); peakHz.add(pk); peakPwr.add(_powerAt(ls, pk)); } - final dropped = atCeiling + belowBand + thin + gappy; + final dropped = atCeiling + atFloor + belowBand + thin + gappy; final total = peaks.length + dropped; if (peaks.length < 3) { // The largest counter is the reason; ties go to the earlier entry. @@ -327,6 +345,11 @@ Metric rsaRespRate( '${ceilingHi > ceilingLo ? '–${round6(ceilingHi)}' : ''} ' 'br/min)' ), + ( + atFloor, + '$atFloor of $total sub-windows peaked at the floor of the RSA band ' + '(${round6(rsaLoHz * 60)} br/min) — the true rate may be below it' + ), ( gappy, '$gappy of $total ${round6(segSec)}s sub-windows had less than 80 % ' @@ -358,13 +381,16 @@ Metric rsaRespRate( for (final p in peaks) { if ((p - medBrpm0).abs() <= tolBrpm) within++; } - final consensus = within / peaks.length; + // A floor-pinned sub-window is a vote for "below the band", not an + // abstention: leaving it out let the few skirt sub-windows that cleared the + // floor guard agree among themselves and publish a rate. + final consensus = within / (peaks.length + atFloor); if (consensus < minConsensus) { return Metric.absent( tier: Tier.high, inputs_used: inputs, note: 'HF peak unstable across the window\'s own sub-windows — only ' - '$within of ${peaks.length} ${round6(segSec)}s sub-windows fall ' + '$within of ${peaks.length + atFloor} ${round6(segSec)}s sub-windows fall ' 'within ${round6(tolBrpm)} br/min of the median; withheld', ); } @@ -391,8 +417,7 @@ Metric rsaRespRate( .clamp(0.2, 0.9); return Metric( value: RespEstimate(brpm, bestPeakHz, bestPower, 'rsa', - usableSubwindows: peaks.length, - subwindows: total), + usableSubwindows: peaks.length, subwindows: total), confidence: conf, tier: Tier.high, inputs_used: inputs, diff --git a/lib/src/onehz/types.dart b/lib/src/onehz/types.dart index c0c7f3b..4bee568 100644 --- a/lib/src/onehz/types.dart +++ b/lib/src/onehz/types.dart @@ -44,7 +44,11 @@ class Tier { /// Published but estimate-grade on a wrist 1 Hz signal (e.g. TRIMP). static const String estimate = 'ESTIMATE'; - /// Relative-only — direction/percentile, never an absolute unit. + /// Relative-only — never an absolute unit. Either a direction/percentile, + /// or a number whose scale only holds against the same user's own baseline + /// on one device family (raw-ADC SpO2/temp, PRSA capacity in ms). Such a + /// number may be trended for that user; it must not be compared across + /// people or charted across a device-family change. static const String relative = 'RELATIVE'; static const Set all = {auth, high, estimate, relative}; diff --git a/lib/src/onehz/util.dart b/lib/src/onehz/util.dart index ac69c03..634880d 100644 --- a/lib/src/onehz/util.dart +++ b/lib/src/onehz/util.dart @@ -388,9 +388,6 @@ List freqGrid(double loHz, double hiHz, int n) { return List.generate(n, (i) => loHz + step * i); } -/// Natural log guard: ln of a positive value, else null. -double? safeLn(double x) => x > 0 ? math.log(x) : null; - /// Calendar-day index (days since epoch) for each `yyyy-MM-dd` label. /// /// Trailing windows and "N nights running" counters in this package used to be diff --git a/lib/src/onehz/wellness/anomaly.dart b/lib/src/onehz/wellness/anomaly.dart index b659c53..5a18270 100644 --- a/lib/src/onehz/wellness/anomaly.dart +++ b/lib/src/onehz/wellness/anomaly.dart @@ -48,7 +48,7 @@ const int multivariateAnomalyMinBaseline = 10; class AnomalyDay { final String date; final double? - mahalanobis; // robust Mahalanobis distance (null if no baseline) + mahalanobis; // robust one-sided (illness-direction) Mahalanobis distance; null if no baseline final bool flagged; // distance crossed gate AND persistence satisfied final bool candidate; // distance crossed gate THIS night (pre-persistence) final List drivers; // per-feature signed contribution @@ -206,18 +206,38 @@ List multivariateAnomaly( ]; final cov = _robustCorr(rows, keep, center, scale, ridge); final inv = _invert(cov); - double d2; - if (inv == null) { - // Fall back to identity (uncorrelated) — sum of squared z. - d2 = zc.fold(0.0, (s, v) => s + v * v); - } else { - d2 = 0.0; - for (var a = 0; a < zc.length; a++) { - for (var b = 0; b < zc.length; b++) { - d2 += zc[a] * inv[a][b] * zc[b]; + // ONE-SIDED: d² is a quadratic form, so it is blind to the sign of zc, and + // a run of unusually GOOD nights (RHR well down, HRV well up — recovery + // after illness or a deload, measured against a baseline that still holds + // the bad stretch) cleared the gate exactly like a deterioration and + // surfaced as an illness-family "unusual physiology" flag. Only deviations + // in the oriented illness direction count toward the distance. Drivers + // below keep the signed zc. + // + // Clamping alone does NOT only lower d²: when a baseline pair comes out + // anti-correlated (Pearson on 10-28 rows often does for temp/resp or + // temp/RHR), inv[a][b] > 0, and zeroing the negative member removes a + // negative cross term, so d² RISES — temp +5 with resp -3 crossed a gate + // the two-sided distance stayed under. The [_madInflationTable] was + // measured on the two-sided statistic, so the one-sided distance is taken + // as min(clamped, two-sided): it keeps the one-sided behaviour and stays + // bounded by the statistic the gate is calibrated for. + final zb = [for (final v in zc) v > 0 ? v : 0.0]; + double quad(List z) { + if (inv == null) { + // Fall back to identity (uncorrelated) — sum of squared z. + return z.fold(0.0, (s, v) => s + v * v); + } + var q = 0.0; + for (var a = 0; a < z.length; a++) { + for (var b = 0; b < z.length; b++) { + q += z[a] * inv[a][b] * z[b]; } } + return q; } + + var d2 = math.min(quad(zb), quad(zc)); if (d2 < 0) d2 = 0; // numerical guard final dist = math.sqrt(d2); diff --git a/lib/src/onehz/wellness/changepoint.dart b/lib/src/onehz/wellness/changepoint.dart index eb40b4d..4daec1e 100644 --- a/lib/src/onehz/wellness/changepoint.dart +++ b/lib/src/onehz/wellness/changepoint.dart @@ -59,6 +59,10 @@ class CusumDetection { /// `clinical/illness_cusum.dart` already carries, and the failure mode /// `util.dart`'s [calendarDays] documents. /// +/// A gap resets the accumulator only. A gap longer than [staleAfterDays] also +/// restarts the pre-change regime: the level from before a month off-wrist is +/// not a baseline for the first nights back. +/// /// Omitting [dates] (or passing a mismatched length) keeps the old positional /// behaviour — no caller is made worse off — but it is not the honest reading. /// @@ -86,6 +90,7 @@ List cusumChangePoints( double k = 0.5, double h = 5.0, int minBaseline = 10, + int staleAfterDays = 28, }) { final out = []; if (x.length <= minBaseline) return out; @@ -94,13 +99,21 @@ List cusumChangePoints( var regimeStart = 0; // first index of the current (pre-change) regime var up = 0.0, dn = 0.0; for (var i = 0; i < x.length; i++) { - // A break in the calendar breaks the regime AND the accumulator: the - // pre-change window restarts here and [minBaseline] consecutive days must - // accrue again before anything can be scored. - if (day != null && i > 0 && day[i] - day[i - 1] > 1) { - up = 0; - dn = 0; - regimeStart = i; + // Any break in the calendar breaks the accumulator: there are no + // observations across it, so no evidence to carry. It does NOT break the + // regime unless it is long enough for the baseline to be stale — a wearer + // who charges overnight once a week has a one-night gap every 6-7 rows, and + // restarting the regime on each one meant the baseline never reached + // [minBaseline] and nothing was ever scored. Same split as + // `clinical/illness_cusum.dart`, which resets only its accumulator and keeps + // a 28-day calendar baseline. + if (day != null && i > 0) { + final step = day[i] - day[i - 1]; + if (step > 1) { + up = 0; + dn = 0; + if (step > staleAfterDays) regimeStart = i; + } } final baseline = x.sublist(regimeStart, i); if (baseline.length < minBaseline) continue; diff --git a/lib/src/onehz/wellness/cycle_lengths.dart b/lib/src/onehz/wellness/cycle_lengths.dart index 3bcc2bd..76054d0 100644 --- a/lib/src/onehz/wellness/cycle_lengths.dart +++ b/lib/src/onehz/wellness/cycle_lengths.dart @@ -42,7 +42,14 @@ const int cycleLongIntervalCriterionDays = 60; const int cycleLengthMinCycles = 12; /// Above this, an interval is more likely a missed log than a cycle. -const int cycleLengthImplausibleGapDays = 90; +/// +/// Pinned to the long-interval criterion on purpose: one missed start turns two +/// 28-day cycles into a 56-day one, and two missed starts into 84, so anything +/// past the criterion is exactly where a forgotten log and a real long interval +/// stop being separable. A looser cutoff (this used to be 90) charts 61-90-day +/// holes as cycles that meet the criterion by themselves, and disagrees with +/// the cycle screen, which refuses anything over 60. +const int cycleLengthImplausibleGapDays = cycleLongIntervalCriterionDays; class CycleLengthSeries { /// Days between consecutive logged onsets, oldest first. diff --git a/lib/src/onehz/wellness/temp_circadian.dart b/lib/src/onehz/wellness/temp_circadian.dart index 5b5b91c..34fd260 100644 --- a/lib/src/onehz/wellness/temp_circadian.dart +++ b/lib/src/onehz/wellness/temp_circadian.dart @@ -87,6 +87,7 @@ class _TempCal { /// g of |‖a‖ − 1| above which an epoch is masked. Null = the family sends /// no accelerometer with its temperature (a ring), so nothing is masked. final double? motionGate; + /// One-sided settle band, in [unit]: a sample more than this far BELOW the /// night's median is not settled skin (see [nightlySkinTemp]). Null = this /// family has no measured band, so the settled mean REFUSES for it rather @@ -254,6 +255,14 @@ Metric nightlySkinTemp( /// every import, every raw replay — REFUSES, because a gen4 ADC count and a /// gen5 centi-°C in the same array cannot be told apart afterwards. /// +/// [tzOffsetResolver] maps a unix-second timestamp to the UTC offset (seconds) +/// in effect at that instant. The acrophase is a clock time and the IS profile +/// is an hour-of-day profile, so both must be built on LOCAL time: on raw +/// epoch time the peak comes out as the UTC hour (5.5 h out for a UTC+5:30 +/// user), and a DST change inside the record shifts real behaviour by an hour +/// against the grid. Resolved per timestamp for that reason. When omitted the +/// machine's local offset at each instant is used, as in `segmentSleep`. +/// /// Returns a RELATIVE-tier metric (phase only). Absent if too little data. Metric tempCircadian( List samples, { @@ -261,6 +270,7 @@ Metric tempCircadian( List? accel, double? motionGate, int epochMin = 60, + int Function(int tsSec)? tzOffsetResolver, }) { const inputs = ['skin_temp_adc', 'device_family']; final cal = calibrationFor(_tempCal, deviceFamily); @@ -272,8 +282,14 @@ Metric tempCircadian( ); } final gate = motionGate ?? cal.motionGate; + final int Function(int tsSec) tzAt = tzOffsetResolver ?? + (int tsSec) => DateTime.fromMillisecondsSinceEpoch( + tsSec * 1000, + isUtc: false, + ).timeZoneOffset.inSeconds; // De-mask: drop invalid + high-motion epochs. final ts = []; + final utcTs = []; final adc = []; var deMasked = 0; for (var i = 0; i < samples.length; i++) { @@ -289,16 +305,17 @@ Metric tempCircadian( } } } - ts.add(s.tsMs); + // Local wall-clock ms: everything downstream reads hour-of-day off this. + ts.add(s.tsMs + tzAt((s.tsMs / 1000).floor()) * 1000.0); + utcTs.add(s.tsMs); adc.add(s.adc); } if (adc.length < 4) { return Metric.absent( tier: Tier.relative, - inputs_used: gate == null || accel == null - ? inputs - : [...inputs, 'accel'], + inputs_used: + gate == null || accel == null ? inputs : [...inputs, 'accel'], note: 'too few valid temp epochs for circadian analysis', ); } @@ -308,25 +325,29 @@ Metric tempCircadian( final cos = cosinor(tHours, adc, periodHours: 24); // --- Nonparametric IS/IV/RA/L5/M10 --- - final np = _nonparam(ts, adc, epochMin: epochMin); + final np = _nonparam(utcTs, ts, adc, epochMin: epochMin); + + // This metric is the phase. A fit whose zero-amplitude test did not reject + // has no acrophase, so it is not reported (nor allowed to set confidence); + // the nonparametric half stands on its own. + final fit = cos.value?.acrophaseHours == null ? null : cos.value; - if (cos.value == null && np == null) { + if (fit == null && np == null) { return Metric.absent( tier: Tier.relative, - inputs_used: inputs, + inputs_used: + gate == null || accel == null ? inputs : [...inputs, 'accel'], note: 'circadian fit degenerate', ); } // Confidence: tie to cosinor R² when present, else a modest np-only value. - final conf = cos.value != null ? (cos.value!.r2).clamp(0.1, 0.9) : 0.3; + final conf = fit != null ? fit.r2.clamp(0.1, 0.9) : 0.3; return Metric( - value: TempCircadian(cos.value, np, cal.unit), + value: TempCircadian(fit, np, cal.unit), confidence: conf, tier: Tier.relative, - inputs_used: gate == null || accel == null - ? inputs - : [...inputs, 'accel'], + inputs_used: gate == null || accel == null ? inputs : [...inputs, 'accel'], note: 'RELATIVE skin-temp phase only (no °C/fever/core). Wrist temp is ' 'ANTIPHASE to core; activity-demasked epochs dropped=$deMasked ' '(gate=${gate == null ? 'none, no accel' : '${gate}g'}). Amplitude ' @@ -336,26 +357,36 @@ Metric tempCircadian( } /// Witting/van Someren nonparametric circadian statistics. Bins the series into -/// fixed [epochMin]-minute epochs aligned to wall-clock (hour-of-day from -/// tsMs), averages within epoch, then computes IS/IV/RA/L5/M10. +/// fixed [epochMin]-minute epochs, averages within epoch, then computes IS/IV. +/// +/// The chronology (which epoch follows which, for IV and coverage) runs on +/// [utcMs]; the epoch-of-day each epoch contributes to IS is read off +/// [localMs]. Binning the chronology on local time folds the two occurrences +/// of a repeated fall-back hour into one epoch, dropping a step from IV. CircadianNonparam? _nonparam( - List tsMs, + List utcMs, + List localMs, List adc, { required int epochMin, }) { - if (tsMs.isEmpty) return null; + if (utcMs.isEmpty) return null; final epochsPerDay = (24 * 60) ~/ epochMin; if (epochsPerDay < 5) return null; final epochMs = epochMin * 60 * 1000.0; // Bin into absolute epoch index (continuous timeline), averaging samples. + // The UTC grid is shifted by the first sample's offset so its edges fall on + // local epoch boundaries (offsets change by whole hours). + final refOffMs = localMs[0] - utcMs[0]; final sums = {}; final counts = {}; + final eodOf = {}; var minBin = 1 << 62, maxBin = -(1 << 62); - for (var i = 0; i < tsMs.length; i++) { - final bin = (tsMs[i] / epochMs).floor(); + for (var i = 0; i < utcMs.length; i++) { + final bin = ((utcMs[i] + refOffMs) / epochMs).floor(); sums[bin] = (sums[bin] ?? 0) + adc[i]; counts[bin] = (counts[bin] ?? 0) + 1; + eodOf.putIfAbsent(bin, () => (localMs[i] / epochMs).floor() % epochsPerDay); if (bin < minBin) minBin = bin; if (bin > maxBin) maxBin = bin; } @@ -403,8 +434,7 @@ CircadianNonparam? _nonparam( for (var i = 0; i < series.length; i++) { final v = series[i]; if (v == null) continue; - final eod = (minBin + i) % epochsPerDay; // epoch-of-day - final idx = eod < 0 ? eod + epochsPerDay : eod; + final idx = eodOf[minBin + i]!; // local epoch-of-day phaseSum[idx] += v; phaseN[idx] += 1; } diff --git a/lib/src/onehz/wellness/temp_health.dart b/lib/src/onehz/wellness/temp_health.dart index d3cae0e..61b81c2 100644 --- a/lib/src/onehz/wellness/temp_health.dart +++ b/lib/src/onehz/wellness/temp_health.dart @@ -61,7 +61,8 @@ class TempIllnessDay { /// is the canonical producer. Null = missing night. [luteal] OPTIONAL per-night /// luteal-phase flag (true => suppress). [baselineDays] trailing robust-baseline /// window; [zThresh] flag elevation above this robust z; [persistDays] -/// consecutive elevated nights required. +/// consecutive elevated nights required. [maxGapNights] unscored nights in a +/// row that carry the run instead of breaking it (as in `illnessCusum`). /// /// DO NOT PASS A Z-SCORE. This function standardises internally against its own /// trailing window (`robustBaseline(...).modZ`), so a pre-standardised series @@ -85,6 +86,7 @@ List tempIllnessFlag( double zThresh = 2.0, int persistDays = 2, int minBaseline = tempIllnessMinBaseline, + int maxGapNights = 2, }) { final n = nightlyTemp.length; final out = []; @@ -109,7 +111,8 @@ List tempIllnessFlag( : null; out.add( TempIllnessDay(dates[i], TempFlag.normal, null, lut, 0.0, need: need)); - elevatedRun = 0; + // An unscored night carries the run; the gap check below decides when + // the calendar break is long enough to end it. continue; } final base = robustBaseline(window, minValid: minBaseline); @@ -117,11 +120,16 @@ List tempIllnessFlag( if (zz == null) { // Can't standardize honestly -> stay normal, confidence 0. out.add(TempIllnessDay(dates[i], TempFlag.normal, null, lut, 0.0)); - elevatedRun = 0; continue; } - // "N nights running" means CONSECUTIVE NIGHTS, not consecutive rows. - if (day[i] - lastScoredDay > 1) elevatedRun = 0; + // "N nights running" means CONSECUTIVE NIGHTS, not consecutive rows — but + // a gap of one or two unscored nights is not a break. Skin temp goes + // unscored exactly when someone is ill: the strap gets re-seated (the + // night never settles), sleep goes undetected, or a warm room fails the + // ambient check. Resetting on every such night meant three fever nights in + // four, or any alternate-night pattern, could never reach [persistDays]. + // A short gap neither advances nor resets the run; a long one ends it. + if (day[i] - lastScoredDay > 1 + maxGapNights) elevatedRun = 0; lastScoredDay = day[i]; final elevated = zz >= zThresh; if (elevated) { diff --git a/lib/src/onehz/workout/auto_detect.dart b/lib/src/onehz/workout/auto_detect.dart index 3ad37a8..3c3ccf2 100644 --- a/lib/src/onehz/workout/auto_detect.dart +++ b/lib/src/onehz/workout/auto_detect.dart @@ -89,8 +89,11 @@ class AutoWorkoutDetector { /// you work out?" prompt should ask about. static const double hrrFloorFraction = 0.45; - /// HRmax fallback when the caller passes none (age-predicted is preferred). - static const int defaultMaxHR = 190; + // There is deliberately no HRmax fallback. A population 190 put the %HRR + // gate ~12 bpm too high for a 62-year-old at RHR 62 (120 vs 108 on their + // age line) and too low for anyone whose real ceiling is above 190, and + // nothing in the result said a default was used. No ceiling → no + // suggestions; see [detect] and [autoDetectWorkouts]. /// A suggested WORKOUT is time-bounded. A sustained elevation longer than this /// (observed: ~3 h of steady 131 bpm from an active evening) is lifestyle @@ -189,6 +192,8 @@ class AutoWorkoutDetector { /// /// [hrTs]/[hrBpm] the day's HR samples (parallel, any order; empty → []). /// [restingBpm] nightly resting HR for the day; null → [defaultRestingHR] (60). + /// [maxBpm] the user's HRmax (measured or age-predicted); null → [] (the + /// %HRR gate cannot be placed without it). /// [motion] OPTIONAL continuous motion series for confirmation; null/empty → /// HR-only. [savedSpans] already-saved windows to exclude by overlap. /// [classify] the sport seam — runs on each surviving bout; default @@ -203,7 +208,7 @@ class AutoWorkoutDetector { SportClassifier classify = defaultSportClassifier, }) { final n = hrTs.length; - if (n == 0) return const []; + if (n == 0 || maxBpm == null) return const []; final order = List.generate(n, (i) => i) ..sort((a, b) => hrTs[a].compareTo(hrTs[b])); final ts = [for (final i in order) hrTs[i]]; @@ -211,8 +216,7 @@ class AutoWorkoutDetector { // %HRR floor (primary), never below the RHR+margin safety floor. final rhr = restingBpm ?? defaultRestingHR; - final hrMax = maxBpm ?? defaultMaxHR; - final hrrFloor = (rhr + hrrFloorFraction * (hrMax - rhr)).round(); + final hrrFloor = (rhr + hrrFloorFraction * (maxBpm - rhr)).round(); final floor = math.max(rhr + elevatedMarginBPM, hrrFloor); // --- 1+2+3: grow sustained spans tolerating brief dips --- @@ -303,8 +307,12 @@ class AutoWorkoutDetector { cnt++; } } - meanMotion = cnt == 0 ? 0.0 : sum / cnt; - if (meanMotion < motionConfirmMean) { + // No motion sample inside the window is missing evidence, not a still + // wrist: reading it as 0.0 g failed the gate and dropped the bout, + // and handed the classifier an intensity nobody measured. Judge the + // window on HR alone, as with no motion series at all. + meanMotion = cnt == 0 ? null : sum / cnt; + if (meanMotion != null && meanMotion < motionConfirmMean) { // Low-motion but a genuine exercise-like onset — see [onsetRiseBpm]. final preMean = _meanBpmInRange(ts, bpm, start - onsetLookbackS, start); @@ -371,8 +379,9 @@ class AutoWorkoutDetector { } } -/// Wrap a list of [DetectedWorkout] in the honesty envelope. Always present -/// (an empty list is a valid, honest "no suggested workouts" answer). +/// Wrap a list of [DetectedWorkout] in the honesty envelope. An empty list is a +/// valid, honest "no suggested workouts" answer; a null value means the gate +/// could not be placed because [maxBpm] was absent (`need_input:hr_max`). Metric> autoDetectWorkouts({ required List hrTs, required List hrBpm, @@ -382,6 +391,13 @@ Metric> autoDetectWorkouts({ List savedSpans = const [], SportClassifier classify = defaultSportClassifier, }) { + if (maxBpm == null) { + return const Metric>.absent( + tier: Tier.estimate, + inputs_used: [], + note: 'need_input:hr_max', + ); + } final list = AutoWorkoutDetector.detect( hrTs: hrTs, hrBpm: hrBpm, diff --git a/test/onehz/an_training_human_test.dart b/test/onehz/an_training_human_test.dart index 6a58fa3..2a5279c 100644 --- a/test/onehz/an_training_human_test.dart +++ b/test/onehz/an_training_human_test.dart @@ -33,6 +33,29 @@ void main() { } }); + test( + 'REGRESSION: a non-finite phase draws AND discloses the population curve', + () { + // The curve fell back to the population phase on NaN while the note + // still claimed a phase fitted on the user's own data. + final nan = alertnessForecast( + wakeLocalHour: 7, + sleepDurationHours: 7, + circadianAcrophaseHours: double.nan, + ); + final none = alertnessForecast(wakeLocalHour: 7, sleepDurationHours: 7); + expect(nan.present, isTrue); + expect(nan.note, contains('ASSUMED from the population')); + expect(nan.note, isNot(contains('fitted on your own data'))); + expect(nan.value!.troughLabel, none.value!.troughLabel); + final fitted = alertnessForecast( + wakeLocalHour: 7, + sleepDurationHours: 7, + circadianAcrophaseHours: 16, + ); + expect(fitted.note, contains('fitted on your own data')); + }); + test('the safety refusal is on every output, present or absent', () { final present = alertnessForecast(wakeLocalHour: 7, sleepDurationHours: 7.5); @@ -90,7 +113,23 @@ void main() { if (shape[i] > shape[peak]) peak = i; } final peakClock = 7 + peak * 0.25; - expect(peakClock, closeTo(15.75, 0.5)); + expect(peakClock, closeTo(16.75, 0.5)); + }); + + test('the 12 h ultradian term puts the dip after lunch, not mid-morning', + () { + // Process U's trough sits 3 h before the circadian acrophase (Åkerstedt, + // Folkard & Portin 2004), ~13:48 at the default phase. It used to peak AT + // the acrophase with no offset: a dip at 10:48, a bump at 16:48, so the + // curve climbed steeply through late morning and then sagged in the late + // afternoon while C was still peaking. + final v = + alertnessForecast(wakeLocalHour: 6, sleepDurationHours: 9).value!; + double at(double h) => v.shape[((h - 6) / v.stepHours).round()]; + final lateMorningRise = at(12) - at(10.75); + final lateAfternoonRise = at(16.75) - at(15); + expect(lateAfternoonRise, greaterThan(0)); + expect(lateMorningRise, lessThan(lateAfternoonRise)); }); test('RD-02: the low window moves with sleep length, not with the horizon', @@ -382,6 +421,18 @@ void main() { expect(m.value!.longestIntervalDays, 60); expect(m.value!.maxConsecutiveDifferenceDays, 32); }); + + test('an interval past the long-interval criterion is a hole, not a bar', + () { + // 75 days: under the old 90-day cutoff this was charted and met the + // 60-day criterion on its own, while the app's cycle screen refused the + // same log. One rule, one threshold. + final m = cycleLengthSeries( + onsets([28, 28, 75, 28, 28, 28, 28, 28, 28, 28, 28, 28, 28])); + expect(m.present, isFalse); + expect(m.note, contains('log_has_holes:intervals_over_60d=1')); + expect(cycleLengthImplausibleGapDays, cycleLongIntervalCriterionDays); + }); }); } diff --git a/test/onehz/auto_detect_test.dart b/test/onehz/auto_detect_test.dart new file mode 100644 index 0000000..814c9a5 --- /dev/null +++ b/test/onehz/auto_detect_test.dart @@ -0,0 +1,230 @@ +// AutoWorkoutDetector — synthetic known-answer tests. +// +// RHR 55 / HRmax 190 throughout, so the gate is +// max(55 + 40, round(55 + 0.45·135)) = 116 bpm. +import 'package:test/test.dart'; +import 'package:openstrap_analytics/src/onehz/workout/auto_detect.dart'; +import 'package:openstrap_analytics/src/onehz/workout/sport.dart'; + +const _rhr = 55; +const _max = 190; + +class _Hr { + final ts = []; + final bpm = []; + + /// Constant [v] for seconds [from, to). + void flat(int from, int to, int v) { + for (var t = from; t < to; t++) { + ts.add(t); + bpm.add(v); + } + } + + /// Linear ramp [a] → [b] over seconds [from, to). + void ramp(int from, int to, int a, int b) { + final n = to - from; + for (var t = from; t < to; t++) { + ts.add(t); + bpm.add((a + (b - a) * (t - from) / n).round()); + } + } +} + +List _motion(int from, int to, double v, + {int? holeFrom, int? holeTo}) => + [ + for (var t = from; t < to; t++) + if (holeFrom == null || t < holeFrom || t >= holeTo!) + MotionPoint(t, v), + ]; + +List _detect(_Hr hr, + {List? motion, SportClassifier? classify}) => + AutoWorkoutDetector.detect( + hrTs: hr.ts, + hrBpm: hr.bpm, + restingBpm: _rhr, + maxBpm: _max, + motion: motion, + classify: classify ?? defaultSportClassifier, + ); + +/// 3 h day: rest, a gradual 20-min warm-up 70→150, 30 min at 150, rest. +/// The warm-up is slow enough that the onset bypass cannot fire. +_Hr _gradualBout() => _Hr() + ..flat(0, 3600, 70) + ..ramp(3600, 4800, 70, 150) + ..flat(4800, 6600, 150) + ..flat(6600, 10800, 70); + +/// 3 h day: rest, a sharp step to 150 for 30 min, rest. +_Hr _sharpBout() => _Hr() + ..flat(0, 3600, 70) + ..flat(3600, 5400, 150) + ..flat(5400, 10800, 70); + +void main() { + group('motion series with no sample inside the window', () { + test('is judged on HR alone, like a null motion series', () { + final hr = _gradualBout(); + final hrOnly = _detect(hr); + expect(hrOnly, hasLength(1)); + + // Motion present for the rest of the day, absent across the bout. + final motion = _motion(0, 10800, 0.3, holeFrom: 3600, holeTo: 6600); + final withHole = _detect(hr, motion: motion); + expect(withHole, hasLength(1)); + expect(withHole.single.startSec, hrOnly.single.startSec); + expect(withHole.single.endSec, hrOnly.single.endSec); + }); + + test('hands the classifier no motion features', () { + final seen = []; + final out = _detect( + _sharpBout(), + motion: _motion(0, 10800, 0.3, holeFrom: 3000, holeTo: 6000), + classify: (bout, feats) { + seen.add(feats); + return 'x'; + }, + ); + expect(out, hasLength(1)); + expect(seen, [isNull]); + }); + + test('a measured still wrist still fails the gate without an onset', () { + final out = _detect(_gradualBout(), motion: _motion(0, 10800, 0.01)); + expect(out, isEmpty); + }); + + test('measured motion above the gate confirms and reaches the classifier', + () { + final seen = []; + final out = _detect( + _gradualBout(), + motion: _motion(0, 10800, 0.3), + classify: (bout, feats) { + seen.add(feats); + return 'x'; + }, + ); + expect(out, hasLength(1)); + expect(seen.single!.meanIntensity, closeTo(0.3, 1e-9)); + }); + }); + + group('onset bypass', () { + test('a sharp HR onset survives a low-motion window', () { + final seen = []; + final out = _detect( + _sharpBout(), + motion: _motion(0, 10800, 0.01), + classify: (bout, feats) { + seen.add(feats); + return 'x'; + }, + ); + expect(out, hasLength(1)); + expect(out.single.avgBpm, 150); + expect(seen.single!.meanIntensity, closeTo(0.01, 1e-9)); + }); + + test('abstains when there is no pre-window to measure the onset against', + () { + // Recording starts mid-elevation: nothing before the span. + final hr = _Hr() + ..flat(0, 1800, 150) + ..flat(1800, 5400, 70); + expect(_detect(hr), hasLength(1)); + expect(_detect(hr, motion: _motion(0, 5400, 0.01)), isEmpty); + }); + }); + + group('maxGapS', () { + test('a hole in the record longer than maxGapS splits the span', () { + // Two 15-min efforts with 40 min of no samples between them. + final hr = _Hr() + ..flat(0, 1800, 70) + ..flat(1800, 2700, 150) + ..flat(5100, 6000, 150) + ..flat(6000, 9000, 70); + final out = _detect(hr); + expect(out, hasLength(2)); + expect(out[0].endSec, 2699); + expect(out[1].startSec, 5100); + }); + + test('a hole no longer than maxGapS is bridged', () { + final hr = _Hr() + ..flat(0, 1800, 70) + ..flat(1800, 2700, 150) + ..flat(2760, 3600, 150) + ..flat(3600, 6000, 70); + final out = _detect(hr); + expect(out, hasLength(1)); + expect(out.single.startSec, 1800); + expect(out.single.endSec, 3599); + }); + }); + + group('merge', () { + test('spans less than mergeGapS apart merge into one bout', () { + // 3 min of present, resting HR between two 15-min efforts. + final hr = _Hr() + ..flat(0, 1800, 70) + ..flat(1800, 2700, 150) + ..flat(2700, 2880, 70) + ..flat(2880, 3780, 150) + ..flat(3780, 7000, 70); + final out = _detect(hr); + expect(out, hasLength(1)); + expect(out.single.startSec, 1800); + expect(out.single.endSec, 3779); + }); + + test('spans at least mergeGapS apart stay separate', () { + final hr = _Hr() + ..flat(0, 1800, 70) + ..flat(1800, 2700, 150) + ..flat(2700, 3300, 70) + ..flat(3300, 4200, 150) + ..flat(4200, 7000, 70); + expect(_detect(hr), hasLength(2)); + }); + }); + + group('HRmax', () { + test('absent ceiling abstains instead of assuming 190', () { + final hr = _sharpBout(); + expect(_detect(hr), hasLength(1)); + expect( + AutoWorkoutDetector.detect( + hrTs: hr.ts, hrBpm: hr.bpm, restingBpm: _rhr, maxBpm: null), + isEmpty, + ); + final m = autoDetectWorkouts( + hrTs: hr.ts, hrBpm: hr.bpm, restingBpm: _rhr, maxBpm: null); + expect(m.value, isNull); + expect(m.note, 'need_input:hr_max'); + }); + + test("the gate follows the caller's ceiling", () { + // RHR 62, HRmax 165: gate = max(102, round(62 + 0.45·103)) = 108, so a + // 25-min bout at 115 is a suggestion. A 190 ceiling would put the gate + // at 120 and drop it. + final hr = _Hr() + ..flat(0, 3600, 70) + ..flat(3600, 5100, 115) + ..flat(5100, 9000, 70); + final out = AutoWorkoutDetector.detect( + hrTs: hr.ts, hrBpm: hr.bpm, restingBpm: 62, maxBpm: 165); + expect(out, hasLength(1)); + expect( + AutoWorkoutDetector.detect( + hrTs: hr.ts, hrBpm: hr.bpm, restingBpm: 62, maxBpm: 190), + isEmpty, + ); + }); + }); +} diff --git a/test/onehz/clinical/vo2max_test.dart b/test/onehz/clinical/vo2max_test.dart index 71072d1..caecbb6 100644 --- a/test/onehz/clinical/vo2max_test.dart +++ b/test/onehz/clinical/vo2max_test.dart @@ -94,7 +94,8 @@ void main() { expect(m.note, contains('%HRR band')); }); - test('hr_max too close to resting_hr abstains rather than divide-by-noise', () { + test('hr_max too close to resting_hr abstains rather than divide-by-noise', + () { final m = vo2maxSubmaxEstimate( speedMps: 2.78, avgHrBpm: 90, @@ -116,5 +117,24 @@ void main() { ); expect(m.present, isFalse); }); + test('non-finite heart-rate inputs abstain instead of throwing', () { + for (final (avg, rest, max) in [ + (double.nan, 55.0, 185.0), + (double.infinity, 55.0, 185.0), + (150.0, double.nan, 185.0), + (150.0, 55.0, double.nan), + (150.0, 55.0, double.infinity), + ]) { + final m = vo2maxSubmaxEstimate( + speedMps: 3, + avgHrBpm: avg, + boutDurationSec: 600, + restingHrBpm: rest, + hrMaxBpm: max, + ); + expect(m.present, isFalse, reason: 'avg=$avg rest=$rest max=$max'); + expect(m.value, isNull); + } + }); }); } diff --git a/test/onehz/clinical_test.dart b/test/onehz/clinical_test.dart index 45ad3e7..533ef03 100644 --- a/test/onehz/clinical_test.dart +++ b/test/onehz/clinical_test.dart @@ -1026,6 +1026,96 @@ void main() { // l=2 (the default) still works on the same series. expect(decelerationCapacity(rr, l: 2).present, isTrue); }); + test('REGRESSION: DC is RELATIVE, not HIGH', () { + // The same subject read 7.6-9.9 ms on gen4 and 1.7 ms on W5; publishing + // it HIGH made a strap switch look like an 80% collapse. + final rr = [ + for (var i = 0; i < 400; i++) 1000 + 20 * math.sin(2 * math.pi * i / 20) + ]; + expect(decelerationCapacity(rr).tier, Tier.relative); + expect(accelerationCapacity(rr).tier, Tier.relative); + expect(decelerationCapacity([1000, 1010, 990]).tier, Tier.relative); + }); + test('REGRESSION: no anchor window spans a hole in the beat record', () { + // 600 steady beats, a 40-min hole, 600 steady beats at a slower rate. + // Neither run has a single deceleration; the only one is the pair + // straddling the hole, which used to publish DC = 20 ms on its own. + final rr = [], ts = []; + var t = 1.7e12; + for (var i = 0; i < 600; i++) { + t += 1000; + rr.add(1000); + ts.add((t ~/ 1000) * 1000.0); + } + t += 40 * 60 * 1000; + for (var i = 0; i < 600; i++) { + t += 1040; + rr.add(1040); + ts.add((t ~/ 1000) * 1000.0); + } + final c = correctRr(rr, rrTsMs: ts); + expect(decelerationCapacity(c.nn).present, isTrue, + reason: 'without times the seam still anchors (legacy callers)'); + final dc = decelerationCapacity(c.nn, nnTimesMs: c.nnTimesMs); + expect(dc.present, isFalse); + expect(dc.note, contains('no valid PRSA anchors')); + }); + // Two steady runs whose only deceleration is the pair joining them, so DC + // is present exactly when that pair is allowed to anchor. + (List, List) twoRuns({required double jumpMs}) { + final rr = [], ts = []; + var t = 1.7e12; + for (var i = 0; i < 600; i++) { + t += 1000; + rr.add(1000); + ts.add(t); + } + t += jumpMs; + for (var i = 0; i < 600; i++) { + t += 1040; + rr.add(1040); + ts.add(t); + } + return (rr, ts); + } + + test('REGRESSION: beat times that run backwards are a seam', () { + final (rr, ts) = twoRuns(jumpMs: -60 * 60 * 1000); + expect(decelerationCapacity(rr).present, isTrue); + final dc = decelerationCapacity(rr, nnTimesMs: ts); + expect(dc.present, isFalse); + expect(dc.note, contains('no valid PRSA anchors')); + }); + test('REGRESSION: misaligned beat times abstain instead of going untimed', + () { + final (rr, ts) = twoRuns(jumpMs: 40 * 60 * 1000); + final dc = decelerationCapacity(rr, nnTimesMs: ts.sublist(1)); + expect(dc.present, isFalse); + expect(dc.note, contains('not aligned')); + expect( + accelerationCapacity(rr, nnTimesMs: ts.sublist(1)).present, isFalse); + }); + test('beat times are reported in inputs_used when they are used', () { + final (rr, ts) = twoRuns(jumpMs: 0); + expect(decelerationCapacity(rr).inputs_used, ['rr_cleaned']); + expect(decelerationCapacity(rr, nnTimesMs: ts).inputs_used, + ['rr_cleaned', 'beat_times']); + }); + test('beat times on a contiguous record change nothing', () { + final rr = [ + for (var i = 0; i < 400; i++) 1000 + 20 * math.sin(2 * math.pi * i / 20) + ]; + final ts = []; + var t = 0.0; + for (final r in rr) { + t += r; + ts.add(t); + } + final a = decelerationCapacity(rr); + final b = decelerationCapacity(rr, nnTimesMs: ts); + expect(b.value!.anchors, a.value!.anchors); + expect(b.value!.capacity, a.value!.capacity); + }); }); group('nocturnal RHR + dip', () { @@ -1198,6 +1288,61 @@ void main() { expect(out.last.state, IllnessState.green); expect(out.last.need, belowQuantumNote); }); + group('a short unscored gap carries the run', () { + const base = [ + 58, + 60, + 59, + 61, + 58, + 60, + 59, + 61, + 58, + 60, + 59, + 61, + 58, + 60 + ]; + List datesFor(int n) => [ + for (var i = 0; i < n; i++) + DateTime.utc(2026, 9, 1 + i).toIso8601String().substring(0, 10) + ]; + test( + 'REGRESSION: one night with no RHR does not wipe an accumulating ' + 'episode', () { + final rhr = [...base, 64, 64, null, 64]; + final out = illnessCusum(datesFor(rhr.length), rhr); + final cont = + illnessCusum(datesFor(rhr.length), [...base, 64, 64, 64]); + // The run keeps the evidence of the two nights before the gap. + expect(out[17].cusum!, greaterThan(out[15].cusum!)); + expect(out[17].state, isNot(IllnessState.green)); + expect(cont.last.state, IllnessState.red); + }); + test('REGRESSION: alternate-night gaps can still trip the alarm', () { + final rhr = [...base, 64, null, 64, null, 64, null]; + final out = illnessCusum(datesFor(rhr.length), rhr); + expect(out.map((d) => d.state), contains(IllnessState.red)); + }); + test('a long wear gap still zeroes the accumulator', () { + final rhr = [...base, 64, 64, 64, 64]; + final dates = datesFor(rhr.length); + // 70 nights off-wrist, then one night back at baseline. + final back = DateTime.utc(2026, 9, rhr.length + 70) + .toIso8601String() + .substring(0, 10); + // A baseline window long enough to span the gap, so the night back + // is actually scored rather than abstaining for want of history. + final out = + illnessCusum([...dates, back], [...rhr, 60], baselineDays: 120); + expect(out[rhr.length - 1].state, IllnessState.red); + expect(out.last.cusum, isNotNull); + expect(out.last.cusum!, lessThan(1.0)); + expect(out.last.state, IllnessState.green); + }); + }); }); group('lnRMSSD readiness stack', () { @@ -1268,6 +1413,77 @@ void main() { expect(m.present, isTrue, reason: m.note); expect(m.value!.z, rows.value!.z); }); + test('an ordinary night is not "suppressed" on single-night noise', () { + // tonight sits half a spread under the median: ordinary night-to-night + // noise, nothing physiological. + final hist = [4.05, 4.10, 4.18, 4.20, 4.22, 4.30, 4.35, 4.12]; + final m = readinessLnRmssd(hist); + expect(m.value!.z!, lessThan(-0.5)); + expect(m.value!.band, 'normal'); + }); + test('an old bad night does not label an ordinary tonight', () { + // tonight equals the centre; the 3.568 night six nights back must not + // make it 'suppressed', and the band must agree with z. + final hist = [4.56, 3.568, 4.623, 4.397, 4.662, 4.794, 4.653, 4.62]; + final m = readinessLnRmssd(hist); + expect(m.value!.z!.abs(), lessThan(0.1)); + expect(m.value!.band, 'normal'); + }); + test('a stable user is flagged on ~32% of nights at every window size', + () { + // the band used to be centre ± 1 MAD of the prior 3-7 nights. the MAD of + // so few points runs small, so with no change at all 62% of nights were + // flagged at 3 prior nights and 45% at 7, against the 31.7% a known + // mean ± SD gives. + final rng = math.Random(11); + double gauss() => + math.sqrt(-2 * math.log(1 - rng.nextDouble())) * + math.cos(2 * math.pi * rng.nextDouble()); + for (final prior in [3, 4, 5, 6, 7]) { + var flagged = 0; + const trials = 6000; + for (var t = 0; t < trials; t++) { + final hist = [for (var i = 0; i <= prior; i++) 4.0 + 0.1 * gauss()]; + if (readinessLnRmssd(hist).value!.band != 'normal') flagged++; + } + expect(flagged / trials, inInclusiveRange(0.28, 0.36), + reason: '$prior prior nights'); + } + }); + group('stable user N(4.2, 0.15), 7 prior nights', () { + final rng = math.Random(7); + double gauss() => + math.sqrt(-2 * math.log(1 - rng.nextDouble())) * + math.cos(2 * math.pi * rng.nextDouble()); + double rate(double tonightShiftSd, String band) { + var hit = 0; + const trials = 4000; + for (var t = 0; t < trials; t++) { + final hist = [ + for (var i = 0; i < 7; i++) 4.2 + 0.15 * gauss(), + 4.2 + 0.15 * (tonightShiftSd + gauss()), + ]; + if (readinessLnRmssd(hist).value!.band == band) hit++; + } + return hit / trials; + } + + test('with no change, most nights read "normal"', () { + // against a half-SD band one night at a time, ~70% were flagged. + expect(rate(0, 'normal'), greaterThan(0.5)); + expect(rate(0, 'suppressed'), lessThan(0.27)); + }); + test('a real drop is still caught', () { + // a night 2 SD under the user's own level reads 'suppressed' most of + // the time, and a 3 SD drop nearly always. a band calibrated to a + // ~32% steady-state flag rate on 7 estimated nights catches a 2 SD + // drop ~74% of the time; 84% is only reachable with known parameters. + expect(rate(-2, 'suppressed'), greaterThan(0.7)); + expect(rate(-3, 'suppressed'), greaterThan(0.9)); + // and a drop is never mistaken for the opposite. + expect(rate(-2, 'elevated'), lessThan(0.02)); + }); + }); }); group('cosinor', () { @@ -1292,6 +1508,7 @@ void main() { final y = [for (var i = 0; i < 24; i++) i.isEven ? 1.0 : -1.0]; final m = cosinor(t, y); expect(m.value!.r2, lessThan(0.2)); + expect(m.value!.acrophaseHours, isNull); }); test('REGRESSION: 4 random points never score a confident circadian fit', () { @@ -1309,10 +1526,14 @@ void main() { }); test('REGRESSION: confidence comes from the ADJUSTED R² (3 fitted params)', () { - // 8 points of noise: raw R² is upward-biased, adjusted R² is not. + // A modest real rhythm under noise: raw R² is upward-biased, adjusted + // R² is not, and confidence follows the adjusted value. final rnd = math.Random(4242); - final t = [for (var i = 0; i < 8; i++) i * 3.0]; - final y = [for (var i = 0; i < 8; i++) rnd.nextDouble()]; + final t = [for (var i = 0; i < 24; i++) i * 2.0]; + final y = [ + for (var i = 0; i < 24; i++) + math.cos(2 * math.pi * t[i] / 24) + rnd.nextDouble() + ]; final m = cosinor(t, y); expect(m.present, isTrue); expect(m.value!.r2Adj, lessThan(m.value!.r2)); @@ -1326,6 +1547,94 @@ void main() { expect(good.value!.r2Adj, closeTo(1.0, 1e-6)); expect(good.confidence, 0.95); }); + test('REGRESSION: no acrophase when the zero-amplitude test does not reject', + () { + // A flat or rhythm-free series used to publish an acrophase picked by + // float rounding / noise, which callers then treated as the user's own + // circadian peak. + for (final n in [72, 300]) { + final t = [for (var h = 0; h < n; h++) h.toDouble()]; + final flat = cosinor(t, List.filled(n, 60.0)); + expect(flat.present, isTrue); + expect(flat.value!.acrophaseHours, isNull); + expect(flat.value!.acrophaseRad, isNull); + expect(flat.value!.toJson()['acrophase_hours'], isNull); + expect(flat.value!.amplitudeP, 1.0); + } + final rnd = math.Random(7); + final t = [for (var h = 0; h < 72; h++) h.toDouble()]; + final noise = [ + for (var h = 0; h < 72; h++) 60 + (rnd.nextDouble() * 3 - 1.5) + ]; + final m = cosinor(t, noise); + expect(m.value!.r2Adj, 0); + expect(m.value!.amplitudeP, greaterThan(cosinorAcrophaseAlpha)); + expect(m.value!.acrophaseHours, isNull); + // A weak but real rhythm in noise keeps its acrophase. + final rhythm = [ + for (var h = 0; h < 72; h++) + 60 + + 3 * math.cos(2 * math.pi * (h - 15) / 24) + + (rnd.nextDouble() * 3 - 1.5) + ]; + final r = cosinor(t, rhythm); + expect(r.value!.amplitudeP, lessThan(cosinorAcrophaseAlpha)); + expect(r.value!.acrophaseHours, closeTo(15, 1.0)); + }); + }); + + group('cosinor zero-amplitude gate', () { + // Hourly HR over 3 days (n = 72), built the way the cross-day pipeline + // feeds it. With no 24 h rhythm the fit explains ~0% of the variance, yet + // it used to publish an acrophase at tier HIGH and the alertness forecast + // placed its trough by that noise. + List hours() => [for (var i = 0; i < 72; i++) i.toDouble()]; + + test('REGRESSION: rotating shifts (peaks drift 8 h/day) publish no phase', + () { + final y = [ + for (var i = 0; i < 72; i++) + () { + final day = i ~/ 24; + final peak = [14.0, 22.0, 6.0][day]; + return 70 + 15 * math.cos(2 * math.pi * ((i % 24) - peak) / 24); + }() + ]; + final m = cosinor(hours(), y); + expect(m.value!.acrophaseHours, isNull); + expect(m.value!.acrophaseRad, isNull); + expect(m.value!.amplitudeP, greaterThanOrEqualTo(0.05)); + }); + + test('REGRESSION: irregular noise publishes no phase', () { + var abstained = 0; + for (var seed = 0; seed < 100; seed++) { + final rnd = math.Random(seed); + final y = [ + for (var i = 0; i < 72; i++) 70 + 20 * rnd.nextDouble() + ]; + if (cosinor(hours(), y).value!.acrophaseHours == null) abstained++; + } + // alpha = 0.05: noise passes the gate about 5% of the time, never most. + expect(abstained, greaterThanOrEqualTo(85)); + }); + + test('a constant signal has no rhythm to phase', () { + final m = cosinor(hours(), List.filled(72, 60.0)); + expect(m.value?.acrophaseHours, isNull); + }); + + test('a real but noisy circadian rhythm still fits', () { + final rnd = math.Random(7); + final y = [ + for (var i = 0; i < 72; i++) + 65 + 8 * math.cos(2 * math.pi * ((i % 24) - 16) / 24) + + 6 * (rnd.nextDouble() - 0.5) + ]; + final m = cosinor(hours(), y); + expect(m.present, isTrue); + expect(m.value!.acrophaseHours, closeTo(16, 0.5)); + }); }); group('TRIMP + CTL/ATL/TSB', () { @@ -1520,6 +1829,31 @@ void main() { expect(m.note, contains('PRV not ECG')); }); + test('stage mask is read on the beat clock, not from the first kept beat', + () { + // correctRr keeps the elapsed time of a dropped leading artifact run, so + // the first NN beat can sit minutes past 0 while the caller's mask starts + // at 0. Here the first 10 min were dropped: 600-1200 s is NREM (small + // swing), 1200-2400 s is REM (large swing) and masked out. + final nn = []; + final times = []; + var t = 600000.0; + var i = 0; + while (t < 2400000.0) { + final amp = t < 1200000.0 ? 10.0 : 40.0; + nn.add(1000.0 + amp * math.sin(2 * math.pi * i / 12)); + t += nn.last; + times.add(t); + i++; + } + final mask = List.generate(2500, (s) => s < 1200); + final m = nocturnalRmssd(nn, times, stageMaskPerSec: mask); + expect(m.present, isTrue); + // Only the two NREM windows count (~3.7 ms); reading the mask from the + // first kept beat let two REM windows in and pulled the median to ~9 ms. + expect(m.value!, lessThan(6.0)); + }); + test('absent without enough beats', () { expect(nocturnalRmssd([800, 810], [800, 1610]).present, isFalse); }); @@ -1879,62 +2213,35 @@ void main() { }); }); - group('StrainScorer (Banister TRIMP → 0–100)', () { - test('trimpToStrain pins: 0→0, 7200→~100, monotone, 2dp', () { - expect(StrainScorer.trimpToStrain(0), 0.0); - // ln(7201)/ln(7201)=1 → 100. - expect(StrainScorer.trimpToStrain(7200), closeTo(100.0, 1e-9)); - expect(StrainScorer.trimpToStrain(100) < StrainScorer.trimpToStrain(335), - isTrue); - // Rounded to 2 decimals. - final v = StrainScorer.trimpToStrain(123.456); - expect((v * 100).round() / 100, v); - }); - - test('REGRESSION: trimpToStrain is CLAMPED to maxStrain', () { - // Docstring says "Map accumulated TRIMP onto [0, 100]" but nothing - // clamped: 14400 → 107.8. (The sibling strainScore() always clamped.) - expect(StrainScorer.trimpToStrain(14400), 100.0); - expect(StrainScorer.trimpToStrain(1e9), StrainScorer.maxStrain); - expect(StrainScorer.trimpToStrain(7200), closeTo(100.0, 1e-9)); - // Below the ceiling nothing changed. - expect(StrainScorer.trimpToStrain(335), lessThan(100.0)); - }); - + group('StrainScorer (Banister TRIMP helpers)', () { test( - 'REGRESSION: strain integrates PER-SAMPLE durations, not the first ' + 'REGRESSION: TRIMP integrates PER-SAMPLE durations, not the first ' 'inter-sample gap applied to everything', () { const bpmv = 150.0; - // (a) 21 samples over 20 min whose FIRST two are 1 s apart — exactly the - // sparse stream minSparseReadings admits. sampleDuration was 1 s for all - // 21 samples → strain 8.08 instead of ~47. + double trimp(List ts, int n) => StrainScorer.banisterTRIMP( + List.filled(n, bpmv), + 50, + 140, + StrainScorer.sampleDurationsMinutes(ts)); + // (a) 21 samples over 20 min whose FIRST two are 1 s apart. The old + // single-interval duration credited 1 s to all 21 samples. final tsIrregular = [ 0, 1, for (var i = 1; i < 20; i++) 1 + i * 63.1 ]; - final bpm21 = List.filled(21, bpmv); - final irregular = - StrainScorer.strain(bpm21, tsIrregular, maxHR: 190, restingHR: 50)!; - final uniform = StrainScorer.strain( - bpm21, [for (var i = 0; i < 21; i++) i * 60.0], - maxHR: 190, restingHR: 50)!; - expect(irregular, greaterThan(40.0)); - expect(irregular, closeTo(uniform, 3.0), - reason: 'same HR over the same wall-clock span → similar strain'); + final irregular = trimp(tsIrregular, 21); + final uniform = trimp([for (var i = 0; i < 21; i++) i * 60.0], 21); + expect(irregular, closeTo(uniform, uniform * 0.1), + reason: 'same HR over the same wall-clock span → similar TRIMP'); // (b) The inverse: 700 samples at 1 Hz behind a 300 s leading gap. The - // 5-min first gap became every sample's duration → strain 104.25. + // 5-min first gap used to become every sample's duration. final tsGap = [0, for (var i = 0; i < 699; i++) 300.0 + i]; - final gapped = StrainScorer.strain(List.filled(700, bpmv), tsGap, - maxHR: 190, restingHR: 50)!; - final dense = StrainScorer.strain(List.filled(700, bpmv), - [for (var i = 0; i < 700; i++) i.toDouble()], - maxHR: 190, restingHR: 50)!; - expect(gapped, lessThanOrEqualTo(StrainScorer.maxStrain)); - expect(gapped, closeTo(dense, 1.0), + final gapped = trimp(tsGap, 700); + final dense = trimp([for (var i = 0; i < 700; i++) i.toDouble()], 700); + expect(gapped, closeTo(dense, 1e-9), reason: 'a hole in the stream is not elapsed effort'); - expect(gapped, lessThan(60.0)); // Per-sample durations: gaps capped at the median cadence, tail gets it. final durs = StrainScorer.sampleDurationsMinutes(tsGap); @@ -1943,16 +2250,6 @@ void main() { expect(durs.reduce(math.max), closeTo(1 / 60.0, 1e-12)); }); - test('MOT-11: no HRmax → no number, not a 220−age stand-in', () { - // `strain` used to fall back to defaultMaxHR() = 190 — a 30 y/o's - // ceiling applied to whoever's wrist arrived. The perimeter caught it; - // the source now does. - final bpm = List.filled(700, 150.0); - final ts = [for (var i = 0; i < 700; i++) i.toDouble()]; - expect(StrainScorer.strain(bpm, ts, maxHR: null), isNull); - expect(StrainScorer.strain(bpm, ts, maxHR: 190), isNotNull); - }); - test('Banister monotonic increasing in intensity', () { // Two same-length streams, one strictly higher HR → more Banister TRIMP. final lo = List.filled(30, 100.0); @@ -1982,56 +2279,6 @@ void main() { expect(src2, 'observed'); expect(h2, greaterThan(StrainScorer.tanakaHRmax(30))); }); - - test('gating: too few samples → null; spanning ≥600s with ≥20 → computes', - () { - // 30 samples but spanning only 30s → fails sparse-span gate → null. - final ts30 = [for (var i = 0; i < 30; i++) i.toDouble()]; - expect( - StrainScorer.strain(List.filled(30, 150), ts30, - maxHR: 190, restingHR: 50), - isNull); - // 20 samples spanning 600s → qualifies. - final tsSpan = [for (var i = 0; i < 20; i++) i * 32.0]; // 19*32=608s - final s = StrainScorer.strain(List.filled(20, 150), tsSpan, - maxHR: 190, restingHR: 50); - expect(s, isNotNull); - }); - - test('maxHR ≤ restingHR → null (invalid HRR)', () { - final ts = [for (var i = 0; i < 700; i++) i.toDouble()]; - expect( - StrainScorer.strain(List.filled(700, 100), ts, - maxHR: 50, restingHR: 60), - isNull); - }); - - test( - 'trimpStrain envelope: present/ESTIMATE on enough data, absent otherwise', - () { - final ts = [for (var i = 0; i < 700; i++) i.toDouble()]; - final m = trimpStrain(List.filled(700, 140), ts, - maxHr: 190, restingHr: 50); - expect(m.present, isTrue); - expect(m.tier, 'ESTIMATE'); - expect(m.value!, greaterThan(0)); - final absent = trimpStrain([100, 110], [0, 1], maxHr: 190, restingHr: 50); - expect(absent.present, isFalse); - expect(absent.confidence, 0); - }); - - test( - 'trimpStrain is absent (not fabricated) when maxHr or restingHr is missing', - () { - final ts = [for (var i = 0; i < 700; i++) i.toDouble()]; - final bpm = List.filled(700, 140); - final noMax = trimpStrain(bpm, ts, restingHr: 50); - expect(noMax.present, isFalse); - expect(noMax.confidence, 0); - final noResting = trimpStrain(bpm, ts, maxHr: 190); - expect(noResting.present, isFalse); - expect(noResting.confidence, 0); - }); }); group('baseline-need signals (need_baseline convention)', () { @@ -2163,6 +2410,59 @@ void main() { expect(untimed.present, isTrue); expect(untimed.value!.si, isNot(closeTo(gapAware.value!.si, 1e-6))); }); + + test( + 'REGRESSION: gap fragments far shorter than ~5 min are not scored ' + 'against the 5-min resting bands', () { + // Resting beats around 1000 ms with RSA + slow drift. Over a 5-min + // window MxDMn picks up the drift; over a ~1-min fragment it cannot, + // so a short window reads a narrower range and a taller mode → higher SI. + double rr(int i) => + 1000.0 + + 35.0 * math.sin(2 * math.pi * i / 4.5) + + 50.0 * math.sin(2 * math.pi * i / 300.0) + + 12.0 * math.sin(i * 2.399); + + // A restless night: 60-beat runs separated by 45 s RR dropouts. + List timesWithGaps(List nn, int every, double gapMs) { + var t = 0.0; + return [ + for (var i = 0; i < nn.length; i++) + t += nn[i] + (i > 0 && i % every == 0 ? gapMs : 0.0) + ]; + } + + final frag = [for (var i = 0; i < 3000; i++) rr(i)]; + final fragmented = + baevskyStressIndex(frag, nnTimesMs: timesWithGaps(frag, 60, 45000)); + expect(fragmented.present, isFalse); + expect(fragmented.value, isNull); + + // Long segments mixed with fragments: the fragments must not move the + // median — it equals the SI of the long segments on their own. + final long1 = [for (var i = 0; i < 600; i++) rr(i)]; + final shorts = [for (var i = 600; i < 1200; i++) rr(i)]; + final long2 = [for (var i = 1200; i < 1800; i++) rr(i)]; + final nn = [...long1, ...shorts, ...long2]; + var t = 0.0; + final times = []; + for (var i = 0; i < nn.length; i++) { + final inShorts = i >= 600 && i < 1200; + final gap = (i == 600 || i == 1200 || (inShorts && (i - 600) % 50 == 0)) + ? 45000.0 + : 0.0; + times.add(t += nn[i] + gap); + } + final mixed = baevskyStressIndex(nn, nnTimesMs: times); + final longOnly = baevskyStressIndex([...long1, ...long2], + nnTimesMs: [ + ...[for (var i = 0; i < 600; i++) times[i]], + ...[for (var i = 1200; i < 1800; i++) times[i]], + ]); + expect(mixed.present, isTrue); + expect(longOnly.present, isTrue); + expect(mixed.value!.si, closeTo(longOnly.value!.si, 1e-9)); + }); }); group('cardiac coherence (McCraty & Zayas 2014)', () { diff --git a/test/onehz/coaching_test.dart b/test/onehz/coaching_test.dart index da0ca98..86c6986 100644 --- a/test/onehz/coaching_test.dart +++ b/test/onehz/coaching_test.dart @@ -263,6 +263,25 @@ void main() { expect(m.confidence, closeTo(0.6, 1e-9)); }); + test('REGRESSION: load is named only when the load model contributed', () { + // With no fitness-fatigue history (cold start, or ctlAtlTsb abstained) + // the window comes from recovery alone, so neither the provenance nor + // the rationale may claim recent load shaped it. + final cold = + strainTarget(recovery0to100: 70, ctl: null, atl: null, tsb: null); + expect(cold.inputs_used, ['recovery']); + expect(cold.value!.rationale, isNot(contains('load'))); + + final zero = strainTarget(recovery0to100: 70, ctl: 0, atl: 0, tsb: 0); + expect(zero.inputs_used, ['recovery']); + expect(zero.value!.rationale, isNot(contains('load'))); + + final warm = + strainTarget(recovery0to100: 70, ctl: 300, atl: 300, tsb: 0); + expect(warm.inputs_used, ['recovery', 'load']); + expect(warm.value!.rationale, contains('recent load')); + }); + test('REGRESSION: a recover target is reachable, not below the floor', () { // The bands were sized for a scale the app never produced: "recover 4–8" // sat BELOW what an inactive worn day scored (~13 on the old map), so a diff --git a/test/onehz/human_test.dart b/test/onehz/human_test.dart index 14c648e..38b4041 100644 --- a/test/onehz/human_test.dart +++ b/test/onehz/human_test.dart @@ -98,6 +98,48 @@ void main() { expect(m.value!.debtHours, isNull); expect(m.value!.hasFreeNight, isFalse); }); + test('one or two free nights => no comparison, not a 1-night p75', () { + final recent = [7.5, 7.4, 7.6, 7.5, 7.3, 7.6, 5.0]; + for (final free in [ + [5.0], + [9.8], + [9.8, 5.0], + ]) { + final m = sleepDebt(recent, free); + expect(m.present, isTrue); + expect(m.value!.hasFreeNight, isTrue); + expect(m.value!.osdHours, isNull, reason: '$free'); + expect(m.value!.debtHours, isNull, reason: '$free'); + expect(m.value!.habitualHours, closeTo(7.5, 1e-9)); + } + }); + test('too few free nights => OSD null, need_baseline note', () { + // One 11 h Saturday used to come back verbatim as the p75 and set the + // user's sleep need. + final one = sleepDebt([7.0, 6.8, 7.1, 11.0], const [11.0]); + expect(one.present, isTrue); + expect(one.value!.osdHours, isNull); + expect(one.value!.debtHours, isNull); + expect(one.value!.hasFreeNight, isTrue); + expect(one.value!.habitualHours, closeTo(7.05, 1e-9)); + expect(one.note, 'need_baseline:have=1,need=$sleepDebtMinFreeNights'); + final two = sleepDebt([7.0, 6.8, 7.1, 5.0], const [5.0, 10.5]); + expect(two.value!.osdHours, isNull); + expect(two.note, 'need_baseline:have=2,need=$sleepDebtMinFreeNights'); + }); + test('REGRESSION: a NaN night does not count toward either minimum', () { + final m = sleepDebt([7.0, 6.8, 7.1], const [7.0, 6.5, double.nan]); + expect(m.value!.osdHours, isNull); + expect(m.note, 'need_baseline:have=2,need=$sleepDebtMinFreeNights'); + final r = sleepDebt([7.0, 6.8, double.nan], const [8.0, 8.2, 8.4]); + expect(r.present, isFalse); + expect(r.note, 'need_baseline:have=2,need=3'); + }); + test('too few recent nights => need_baseline note', () { + final m = sleepDebt([7.0, 6.8], const [8.0, 8.2, 8.4]); + expect(m.present, isFalse); + expect(m.note, 'need_baseline:have=2,need=3'); + }); test('free-night rebound => positive debt', () { // Habitual ~6 h, free nights ~8.5 h => debt ≈ 2.5 h. final m = sleepDebt([6.0, 6.1, 5.9, 6.0, 6.2], [8.4, 8.6, 8.5]); @@ -256,19 +298,6 @@ void main() { 'among your highest'); }); - test('record gated by MDC: tiny beat is NOT a record, big beat is', () { - final hist = List.generate(30, (i) => 40.0 + (i % 10)); // spread - final priorMax = hist.reduce((a, b) => a > b ? a : b); - // Barely beats the max => within MDC noise => NOT a record. - final small = - personalRecord(priorMax + 0.01, hist, better: Better.higher); - expect(small.value!.isRecord, isFalse); - // Clearly beats it. - final big = personalRecord(priorMax + 50, hist, better: Better.higher); - expect(big.value!.isRecord, isTrue); - expect(big.value!.kind, 'high'); - }); - test('short history => absent', () { final m = percentileOfYou(10, [1, 2, 3], better: Better.higher); expect(m.present, isFalse); diff --git a/test/onehz/journal_numeric_correlations_test.dart b/test/onehz/journal_numeric_correlations_test.dart index 56388cf..d3fbc9c 100644 --- a/test/onehz/journal_numeric_correlations_test.dart +++ b/test/onehz/journal_numeric_correlations_test.dart @@ -82,6 +82,23 @@ void main() { }); group('journalNumericCorrelations', () { + test('REGRESSION: fdrAlpha 0 admits nothing instead of throwing', () { + final dates = _dates(12); + final caffeine = [1, 2, 3, 4, 5, 6, 1, 2, 3, 4, 5, 6]; + final out = journalNumericCorrelations( + fieldLagDays: const {}, + journal: [ + for (var i = 0; i < dates.length; i++) + JournalNumericDay(dates[i], {'caffeine': caffeine[i]}), + ], + dates: dates, + outcomes: { + 'rmssd': [for (final c in caffeine) 80.0 - 6 * c] + }, + fdrAlpha: 0, + ); + expect(_effect(out, 'caffeine', 'rmssd').meaningful, isFalse); + }); test('finds a strong monotone relationship and signs it correctly', () { final dates = _dates(12); final caffeine = [1, 2, 3, 4, 5, 6, 1, 2, 3, 4, 5, 6]; @@ -636,6 +653,47 @@ void main() { expect(at9.meaningful, isTrue); }); + // The 1/(b+1) sampling floor is a resolution limit of how many shuffles were + // drawn, not a property of the data. With b fixed at 999 it alone failed + // the BH gate once the family reached 51 tests (13 fields × 4 outcomes), so + // one more logged field blanked every numeric insight, with a "need more + // days" note that no amount of logging could ever satisfy. + test('a large family does not make every test structurally unreachable', + () { + List run(int fields) { + const n = 30; + final dates = _dates(n); + final rng = math.Random(7); + final noise = [ + for (var f = 0; f < fields; f++) + [for (var i = 0; i < n; i++) (rng.nextInt(5) + 1).toDouble()], + ]; + return journalNumericCorrelations( + fieldLagDays: const {}, + journal: [ + for (var i = 0; i < n; i++) + JournalNumericDay(dates[i], { + 'f00': (i % 10).toDouble(), + for (var f = 1; f < fields; f++) + 'f${f.toString().padLeft(2, '0')}': noise[f][i], + }), + ], + dates: dates, + outcomes: { + for (final k in ['rmssd', 'rhr', 'readiness', 'efficiency']) + k: [for (var i = 0; i < n; i++) 50.0 + 3.0 * (i % 10) + (i % 2)], + }, + ); + } + + for (final fields in [12, 13, 20]) { + final e = _effect(run(fields), 'f00', 'readiness'); + expect(e.insufficient, isFalse, reason: 'fields=$fields'); + expect(e.note, isNull, reason: 'fields=$fields'); + expect(e.meaningful, isTrue, reason: 'fields=$fields q=${e.q}'); + } + }); + test('empty input is empty output, not a crash', () { expect( journalNumericCorrelations( diff --git a/test/onehz/overreaching_conjunction_test.dart b/test/onehz/overreaching_conjunction_test.dart index b1f457d..2770d5a 100644 --- a/test/onehz/overreaching_conjunction_test.dart +++ b/test/onehz/overreaching_conjunction_test.dart @@ -47,5 +47,32 @@ void main() { expect(m.present, isTrue); expect(m.value!.bothPointSameWay, isTrue); }); + test('a long but coarse baseline is not reported as a cold start', () { + final recent = [65.0, 66.0, 64.0, 67.0, 65.0]; + final flat = overreachingConjunction( + load: load, + rhrRecent: recent, + rhrBaselineWindow: List.filled(28, 58), + ); + expect(flat.present, isFalse); + expect(flat.note, degenerateBaselineNote); + + final alternating = overreachingConjunction( + load: load, + rhrRecent: recent, + rhrBaselineWindow: [for (var i = 0; i < 28; i++) i.isEven ? 58 : 59], + ); + expect(alternating.present, isFalse); + expect(alternating.note, belowQuantumNote); + }); + + test('a short baseline still reports need_baseline', () { + final m = overreachingConjunction( + load: load, + rhrRecent: const [65.0, 65.0, 65.0, 65.0, 65.0], + rhrBaselineWindow: const [54, 56, 58, 60, 55], + ); + expect(m.note, 'need_baseline:have=5,need=14'); + }); }); } diff --git a/test/onehz/readiness_calibration_test.dart b/test/onehz/readiness_calibration_test.dart index 074c988..131a6d0 100644 --- a/test/onehz/readiness_calibration_test.dart +++ b/test/onehz/readiness_calibration_test.dart @@ -174,6 +174,7 @@ void main() { final m = readinessLnRmssd(hist); expect(m.value!.rolling7Mean, closeTo(4.4653, 1e-3)); expect(m.value!.rolling7Median, closeTo(4.623, 1e-9)); + // tonight itself is low: 1.5 robust SDs under the median. expect(m.value!.band, 'suppressed'); // mean/SD said 'normal' expect(m.value!.z!, lessThan(-1)); }); diff --git a/test/onehz/respiration_test.dart b/test/onehz/respiration_test.dart index 3151933..768c602 100644 --- a/test/onehz/respiration_test.dart +++ b/test/onehz/respiration_test.dart @@ -199,7 +199,8 @@ void main() { return (rr: rr, t: t); } - test('RSA: an 8-hour night resolves its time-median rate, not a noise spike', + test( + 'RSA: an 8-hour night resolves its time-median rate, not a noise spike', () { // Measured, same three noise seeds: OLD published 18.48 / 16.94 / 16.06 — // three different answers to the same signal, each passing its own @@ -224,8 +225,7 @@ void main() { final s = rsaNight(11); final corr = correctRr(s.rr, rrTsMs: s.t); final shuffled = [...corr.nn]..shuffle(math.Random(3)); - final m = - rsaRespRate(shuffled, corr.nnTimesMs, artifactFraction: 0); + final m = rsaRespRate(shuffled, corr.nnTimesMs, artifactFraction: 0); expect(m.present, isFalse, reason: 'got ${m.value?.brpm}'); expect(m.confidence, 0); expect(m.note, contains('sub-windows')); @@ -273,8 +273,7 @@ void main() { expect(mg.present, isTrue, reason: 'seed $seed: ${mg.note}'); expect(mg.value!.brpm!, closeTo(mf.value!.brpm!, 1.0)); expect(mg.value!.brpm!, closeTo(15.0, 0.8)); // time-median of 16→14 - expect(mg.value!.usableSubwindows!, - lessThan(mg.value!.subwindows!), + expect(mg.value!.usableSubwindows!, lessThan(mg.value!.subwindows!), reason: 'the holes cost sub-windows, and the value says so'); } }); @@ -321,7 +320,9 @@ void main() { for (var k = 0; k < 4; k++) (k * 0.25 + 0.05, k * 0.25 + 0.25) ]); // The fixture is what it says: four bursts of ~3 min, ~12 min of beats. - final bursts = >[[s.t.first, s.t.first]]; + final bursts = >[ + [s.t.first, s.t.first] + ]; for (final x in s.t.skip(1)) { if (x - bursts.last[1] > 60e3) bursts.add([x, x]); bursts.last[1] = x; @@ -385,10 +386,10 @@ void main() { // 1.4 h at 56 bpm, then 0.6 h at 40 bpm (Nyquist 20 br/min, below the // top of the HF band). The whole input's median clears the band; the // slow stretch's sub-windows do not, and must not be used. - final a = gappyRsaNight(5, hrBpm: 56, hours: 1.4, resp0Brpm: 15, - resp1Brpm: 15); - final b = gappyRsaNight(6, hrBpm: 40, hours: 0.6, resp0Brpm: 15, - resp1Brpm: 15); + final a = + gappyRsaNight(5, hrBpm: 56, hours: 1.4, resp0Brpm: 15, resp1Brpm: 15); + final b = + gappyRsaNight(6, hrBpm: 40, hours: 0.6, resp0Brpm: 15, resp1Brpm: 15); final rr = [...a.rr, ...b.rr]; final t = [...a.t, for (final x in b.t) a.t.last + x]; final c = correctRr(rr, rrTsMs: t); @@ -401,7 +402,8 @@ void main() { reason: 'the 40 bpm stretch holds ~13 sub-windows'); }); - test('RSA: the at-ceiling note names the rejecting sub-windows\' own ' + test( + 'RSA: the at-ceiling note names the rejecting sub-windows\' own ' 'ceilings', () { // A beat-to-beat alternation sits exactly at each sub-window's own // Nyquist, so every sub-window peaks at its ceiling: ~26 br/min in the @@ -433,14 +435,66 @@ void main() { expect(double.parse(nums.group(2) ?? '0'), inInclusiveRange(28.5, 30.0)); }); + test( + 'RSA: slow breathing below the band floor is withheld, not published ' + 'as the floor', () { + // An hour of 7.5–9 br/min RSA (deep sleep / paced breathing) with a + // 0.1 Hz Mayer wave. Only the skirt of the true peak is inside the + // 9 br/min floor, so every sub-window agreed on ~9 br/min and it went + // out at high confidence. + RespEstimate? run(double brpm, + {double spanMs = 3600e3, double artifactFraction = 0, int seed = 1}) { + final rnd = math.Random(seed); + final nn = [], t = []; + var tMs = 0.0; + while (tMs < spanMs) { + final s = tMs / 1000; + final v = 1000 + + 30 * math.sin(2 * math.pi * brpm / 60 * s) + + 10 * math.sin(2 * math.pi * 0.1 * s) + + (rnd.nextDouble() * 20 - 10); + tMs += v; + nn.add(v); + t.add(tMs); + } + final m = rsaRespRate(nn, t, artifactFraction: artifactFraction); + if (m.present) return m.value; + expect(m.note, isNotNull); + return null; + } + + for (final slow in [7.5, 8.0, 8.5, 9.0]) { + final v = run(slow); + expect(v, isNull, reason: 'true $slow br/min reported as ${v?.brpm}'); + } + // In-band rates are unaffected. + expect(run(10.0)?.brpm, closeTo(10.0, 0.3)); + expect(run(12.0)?.brpm, closeTo(12.0, 0.3)); + + // A 3-min window has ~90 s sub-windows. The floor guard must not widen + // with the coarser resolution: a true 10 br/min is still a measurement + // there, while the slow-breathing skirt is still withheld. + for (var seed = 1; seed <= 5; seed++) { + final short10 = + run(10.0, spanMs: 180e3, artifactFraction: 0.15, seed: seed); + expect(short10?.brpm, closeTo(10.0, 0.5), + reason: 'seed $seed: 3-min true 10 br/min was ${short10?.brpm}'); + final short8 = + run(8.0, spanMs: 180e3, artifactFraction: 0.15, seed: seed); + expect(short8, isNull, + reason: 'seed $seed: 3-min true 8 br/min reported as ' + '${short8?.brpm}'); + } + }); + test('RSA: confidence scales with usable sub-windows', () { final s = rsaNight(11); final c = correctRr(s.rr, rrTsMs: s.t); ({List nn, List t}) firstSec(double sec) { // through the first beat PAST [sec], so the span covers all of it - final k = c.nnTimesMs - .indexWhere((t) => t - c.nnTimesMs.first > sec * 1000) + - 1; + final k = + c.nnTimesMs.indexWhere((t) => t - c.nnTimesMs.first > sec * 1000) + + 1; return (nn: c.nn.sublist(0, k), t: c.nnTimesMs.sublist(0, k)); } @@ -603,6 +657,39 @@ void main() { reason: 'pre-fix this fell by a third'); }); + test('CVHR: a few analysed minutes never publish at full confidence', () { + // One 60 s dip in a short clean stretch read ~50 cycles/h at 0.85, tier + // HIGH, because confidence only looked at the artifact fraction. + List times(List nn) { + var t = 0.0; + return [for (final r in nn) t += r]; + } + + final nn = []; + var t = 0.0; + while (t < 72000) { + final s = t / 1000.0; + final hr = s >= 6 && s < 66 + ? 60.0 - 6.0 * math.sin(2 * math.pi * (s - 6) / 60) + : 60.0; + final r = 60000.0 / hr; + nn.add(r); + t += r; + } + final short = cvhrApneaScreen(nn, times(nn), artifactFraction: 0); + expect(short.present, isTrue, reason: short.note); + expect(short.value!.analyzedHours, lessThan(0.05)); + expect(short.confidence, lessThan(0.05)); + expect(short.tier, Tier.estimate); + + // A full night of clean beats keeps the old ceiling. + final night = List.filled(5 * 3600, 1000.0); + final full = cvhrApneaScreen(night, times(night), artifactFraction: 0); + expect(full.present, isTrue, reason: full.note); + expect(full.confidence, closeTo(0.85, 1e-9)); + expect(full.tier, Tier.high); + }); + test('CVHR: no cycles means mean depth/width are ABSENT, not 0', () { // an-motion-7. They were non-nullable and filled with 0, serialising an // absence as a measurement. @@ -730,14 +817,12 @@ void main() { test( 'relativeOdi: a mid-night gap (charging break) does not dilute ' - 'odiPerHour/burdenPct — observed-hours denominator, not raw span', - () { + 'odiPerHour/burdenPct — observed-hours denominator, not raw span', () { // Two clean 300 s segments with 5 dips each, separated by a 2 h gap // (charging break) in ts. Same shape as cvhr_apnea's regression case. const pulsHz = 0.3; - List buildSegment(int startSec, math.Random rnd, - List events, List red, List ir, - List ts) { + List buildSegment(int startSec, math.Random rnd, List events, + List red, List ir, List ts) { for (var s = 0; s < 300; s++) { final t = (startSec + s).toDouble(); final inEvent = events.any((e) => s >= e && s < e + 20); @@ -820,7 +905,8 @@ void main() { 'a desaturation at the segment boundary'); }); - test('relativeOdi: default maxGapSec is a no-op on gap-free input ' + test( + 'relativeOdi: default maxGapSec is a no-op on gap-free input ' '(regression guard)', () { const totalSec = 600; final red = []; @@ -846,8 +932,7 @@ void main() { // span is (totalSec-1) seconds (ts runs 0..totalSec-1 inclusive) — // same first-to-last-span semantics as before this fix, just now // computed as a single observed segment instead of the raw ts delta. - expect(m.value!.analyzedHours, - closeTo((totalSec - 1) / 3600.0, 1e-9)); + expect(m.value!.analyzedHours, closeTo((totalSec - 1) / 3600.0, 1e-9)); }); test('BRV: variable breathing rates -> CV>0 + Theil-Sen slope', () { @@ -1085,6 +1170,99 @@ void _resp01Tests() { expect(rising.toJson().toString(), isNot(contains('ahi'))); }); + test('a wear gap does not pull months-old nights into the window', () { + CvhrNight at(String day, double perHour) => + CvhrNight(dayKey: day, cvhrPerHour: perHour, analyzedHours: 7); + final june = [ + for (var d = 1; d <= 30; d++) + at('2026-06-${d.toString().padLeft(2, '0')}', 20), + ]; + // Four nights back after three months off: too few to say anything, + // rather than 26 June nights standing in for "your usual". + final back = cvhrPersonalDistribution([ + ...june, + for (var d = 1; d <= 4; d++) at('2026-10-0$d', 2), + ]); + expect(back.present, isFalse); + expect(back.note, contains('need_baseline:nights=4/5')); + + // Once there are enough, the screen describes only the new nights. + final settled = cvhrPersonalDistribution([ + ...june, + for (var d = 1; d <= 6; d++) at('2026-10-0$d', 2), + ]); + expect(settled.value!.nightsUsed, 6); + expect(settled.value!.recentWeightedMean, closeTo(2, 1e-9)); + expect(settled.value!.quartiles, everyElement(closeTo(2, 1e-9))); + expect(settled.note, contains('in the last 30 days')); + }); + + test('the window spans 30 calendar days ending on the newest night', () { + final m = cvhrPersonalDistribution([ + CvhrNight(dayKey: '2026-02-01', cvhrPerHour: 99, analyzedHours: 7), + for (var d = 3; d <= 7; d++) + CvhrNight(dayKey: '2026-03-0$d', cvhrPerHour: 10, analyzedHours: 7), + ]); + // 2026-02-01 is 34 days before 2026-03-07. + expect(m.value!.nightsUsed, 5); + expect(m.value!.weightedMean, closeTo(10, 1e-9)); + final edge = cvhrPersonalDistribution([ + CvhrNight(dayKey: '2026-02-06', cvhrPerHour: 40, analyzedHours: 7), + for (var d = 3; d <= 7; d++) + CvhrNight(dayKey: '2026-03-0$d', cvhrPerHour: 10, analyzedHours: 7), + ]); + // 2026-02-06 is 29 days before 2026-03-07: inside. + expect(edge.value!.nightsUsed, 6); + }); + + test('night-to-night noise is not "running higher"', () { + // 0.05/h over the old all-nights p75, recent nights inside the range. + final noise = cvhrPersonalDistribution([ + night(1, 4.8), + night(2, 5.1), + night(3, 4.9), + night(4, 5.3), + night(5, 5.4), + ]); + expect(noise.value!.aboveOwnUsual, isNot(isTrue)); + // A single high newest night is not its own evidence, at the floor or + // a few nights past it. + for (final n in [5, 6, 8, 12]) { + final oneHigh = cvhrPersonalDistribution([ + for (var d = 1; d < n; d++) night(d, 5), + night(n, 9), + ]); + expect(oneHigh.value!.aboveOwnUsual, isNot(isTrue), reason: 'n=$n'); + } + // Recent nights inside an older, noisier spread stay quiet. + final spread = cvhrPersonalDistribution([ + for (var d = 1; d <= 9; d++) night(d, d.isEven ? 4 : 12), + night(10, 13), + night(11, 12), + night(12, 13), + ]); + expect(spread.value!.aboveOwnUsual, isFalse); + }); + + test('too few older nights says nothing comparative', () { + // 8 nights: 3 recent, 5 older — under the reference floor however + // large the jump. + final m = cvhrPersonalDistribution([ + for (var d = 1; d <= 5; d++) night(d, 5), + for (var d = 6; d <= 8; d++) night(d, 40), + ]); + expect(m.present, isTrue); + // null, not false: false would read as "stayed inside your range". + expect(m.value!.aboveOwnUsual, isNull); + expect(m.value!.toJson()['above_own_usual'], isNull); + // one more older night and the same jump reads. + final enough = cvhrPersonalDistribution([ + for (var d = 1; d <= 6; d++) night(d, 5), + for (var d = 7; d <= 9; d++) night(d, 40), + ]); + expect(enough.value!.aboveOwnUsual, isTrue); + }); + test('a re-derived day counts once', () { final m = cvhrPersonalDistribution([ for (var d = 1; d <= 5; d++) night(d, 10), diff --git a/test/onehz/steps_test.dart b/test/onehz/steps_test.dart index d479318..7dabad9 100644 --- a/test/onehz/steps_test.dart +++ b/test/onehz/steps_test.dart @@ -464,8 +464,41 @@ void main() { expect(m.tier, Tier.estimate); // DAYS, not minutes: the only floor builder a storage-bound caller can // use (personalDynFloorFromDailySummaries) is gated on days, and edge - // passes a day count. This used to read have=3,need=2000. - expect(m.note, 'need_baseline:have=3,need=$personalDynFloorMinDays'); + // passes a day count. This used to read have=3,need=2000. The need is + // the FROZEN-floor gate, since that is what keeps the floor null. + expect(m.note, 'need_baseline:have=3,need=$enrollmentDaysForFrozenFloor'); + }); + + test('cold-start note never reads as met while the floor is still null', + () { + // Days 5-13 have enough history for the summary median but not for the + // frozen floor, so the caller still passes null. The note used to say + // have=9,need=5, which reads as an unlocked gate on an absent metric. + for (final have in [5, 9, enrollmentDaysForFrozenFloor - 1]) { + final m = dailyActiveMinutes(day(120, 30), + personalDynFloorG: null, historyDaysAvailable: have); + expect(m.present, isFalse); + expect(m.note, + 'need_baseline:have=$have,need=$enrollmentDaysForFrozenFloor'); + } + // A caller that unlocks on the bare summary median says so. + final early = dailyActiveMinutes(day(120, 30), + personalDynFloorG: null, + historyDaysAvailable: 3, + requiredHistoryDays: personalDynFloorMinDays); + expect(early.note, 'need_baseline:have=3,need=$personalDynFloorMinDays'); + }); + + test( + 'REGRESSION: an unusable floor past the history gate is not a met ' + 'need_baseline note', () { + for (final floor in [double.nan, 0.0, null]) { + final m = dailyActiveMinutes(day(120, 30), + personalDynFloorG: floor, + historyDaysAvailable: enrollmentDaysForFrozenFloor + 6); + expect(m.present, isFalse); + expect(m.note, isNot(startsWith('need_baseline')), reason: '$floor'); + } }); test('an off-wrist HOUR between moving minutes does not make one bout', () { diff --git a/test/onehz/strain_calibration_test.dart b/test/onehz/strain_calibration_test.dart index ed83a3e..5f0e25a 100644 --- a/test/onehz/strain_calibration_test.dart +++ b/test/onehz/strain_calibration_test.dart @@ -248,6 +248,31 @@ void main() { expect(m.tier, Tier.estimate); expect(m.value, greaterThan(14.0)); }); + + test('lump and series routes disclose calibration identically', () { + // The stored-TRIMP rescale and the per-minute scorer report the same + // headline. Same quiet level, same settled state => same confidence and + // the same calibrating note, or a backfilled day reads differently from + // a live one. + const rhr = 50.0, hrMax = 187.0; + final hr = List.filled(960, rhr + 0.3 * (hrMax - rhr)); + final trimp = + banisterTrimp(hr, restingHr: rhr, maxHr: hrMax, sex: Sex.male); + for (final settled in [true, false]) { + final lump = strainScoreMetric(trimp.value, + wakeMinutes: 960, quietHrr: quietWakingHrr, quietSettled: settled); + final series = strainScoreFromSeries(hr, + restingHr: rhr, + maxHr: hrMax, + quietHrr: quietWakingHrr, + sex: Sex.male, + quietSettled: settled); + expect(lump.present && series.present, isTrue); + expect(lump.confidence, series.confidence); + expect(lump.tier, series.tier); + if (!settled) expect(lump.note, series.note); + } + }); }); // MOT-03 shipped the level as an argument; these pin what scoring against it diff --git a/test/onehz/util_test.dart b/test/onehz/util_test.dart index d031213..0b14c64 100644 --- a/test/onehz/util_test.dart +++ b/test/onehz/util_test.dart @@ -265,9 +265,7 @@ void main() { expect(HeartRateZones.timeInZone(stream(600, 20), zoneSet), isNull); }); - test('StrainScorer: 301 s and 600 s produce no durations and no strain', - () { - final bpm = List.filled(40, 150.0); + test('StrainScorer: 301 s and 600 s produce no durations', () { List ts(double cadence) => [for (var i = 0; i < 40; i++) i * cadence]; expect(StrainScorer.medianIntervalSeconds(ts(300)), closeTo(300, 1e-12)); @@ -275,15 +273,7 @@ void main() { expect(StrainScorer.medianIntervalSeconds(ts(600)), isNull); expect(StrainScorer.sampleDurationsMinutes(ts(301)), isEmpty); expect(StrainScorer.sampleDurationsMinutes(ts(600)), isEmpty); - // The abstention has to reach the published number, not stop at the - // helper: `banisterTRIMP` credits `fallbackSampleMin` per sample when - // handed an empty list, which is the fabricated 1 s all over again. - expect( - StrainScorer.strain(bpm, ts(301), maxHR: 190, restingHR: 50), isNull); - expect( - StrainScorer.strain(bpm, ts(600), maxHR: 190, restingHR: 50), isNull); - expect(StrainScorer.strain(bpm, ts(300), maxHR: 190, restingHR: 50), - isNotNull); + expect(StrainScorer.sampleDurationsMinutes(ts(300)), hasLength(40)); }); test('AdvancedSleepStager stages nothing past its own cadence ceiling', () { diff --git a/test/onehz/wellness_test.dart b/test/onehz/wellness_test.dart index d64843c..279ebe1 100644 --- a/test/onehz/wellness_test.dart +++ b/test/onehz/wellness_test.dart @@ -16,6 +16,10 @@ import 'package:test/test.dart'; import 'package:openstrap_analytics/src/onehz/wellness/wellness.dart'; import 'package:openstrap_protocol/openstrap_protocol.dart'; +// The synthetic series below put their phase on the UTC epoch grid; pin the +// local offset to zero so they do not depend on the machine's time zone. +int _utc(int _) => 0; + void main() { // ------------------------------------------------------------------------- // 1. Relative skin-temp circadian — cosinor recovers an injected phase. @@ -33,7 +37,8 @@ void main() { 2000 + 300 * math.cos(2 * math.pi / 24 * (tHours - peakHour)); samples.add(AdcSample(tMs, adc)); } - final m = tempCircadian(samples, deviceFamily: 'gen4', epochMin: 60); + final m = tempCircadian(samples, + deviceFamily: 'gen4', epochMin: 60, tzOffsetResolver: _utc); expect(m.present, isTrue); expect(m.tier, Tier.relative); final fit = m.value!.cosinorFit!; @@ -46,6 +51,93 @@ void main() { expect(np.intradailyVariability, lessThan(0.5)); }); + test('acrophase and IS bins are local clock time, not UTC', () { + // UTC+5:30 user whose wrist temp peaks at 03:00 LOCAL = 21:30 UTC. + const off = 19800; + const peakLocal = 3.0; + final samples = []; + // 7 days @ 1 sample / 5 min, starting at a UTC midnight. + final t0 = DateTime.utc(2026, 3, 1).millisecondsSinceEpoch.toDouble(); + for (var i = 0; i < 7 * 24 * 12; i++) { + final tMs = t0 + i * 5 * 60 * 1000.0; + final localH = ((tMs / 1000 + off) / 3600) % 24; + final adc = + 2000 + 300 * math.cos(2 * math.pi / 24 * (localH - peakLocal)); + samples.add(AdcSample(tMs, adc)); + } + final m = tempCircadian(samples, + deviceFamily: 'gen4', epochMin: 60, tzOffsetResolver: (_) => off); + expect(m.value!.cosinorFit!.acrophaseHours, closeTo(peakLocal, 0.1)); + }); + + test('a DST change inside the record does not blur IS', () { + // Peak at 03:00 local throughout; the offset steps +1 h mid-record. On a + // UTC grid the second half's profile lands an hour late. + final t0 = DateTime.utc(2026, 3, 1).millisecondsSinceEpoch.toDouble(); + final switchMs = t0 + 7 * 24 * 3600 * 1000.0; + int offAt(int s) => s * 1000 < switchMs ? 0 : 3600; + // Square wave on the local clock: hour-aligned plateaus, so an exact + // local grid gives IS == 1 and any 1 h misalignment gives IS < 1. + final samples = []; + for (var i = 0; i < 14 * 24 * 12; i++) { + final tMs = t0 + i * 5 * 60 * 1000.0; + final localH = ((tMs / 1000 + offAt((tMs / 1000).floor())) / 3600) % 24; + samples.add(AdcSample(tMs, localH < 8 ? 2300 : 2000)); + } + final local = tempCircadian(samples, + deviceFamily: 'gen4', epochMin: 60, tzOffsetResolver: offAt); + expect(local.value!.nonparam!.interdailyStability, closeTo(1.0, 1e-9)); + }); + + test('REGRESSION: a repeated fall-back hour stays two epochs', () { + // One sample per UTC hour, offset +1 h then 0 h (fall-back), so local + // 01:00 occurs twice. Each hour carries its own value; IV must be the + // IV of that hourly series, with both 01:00 epochs in it. + final t0 = DateTime.utc(2026, 10, 1).millisecondsSinceEpoch.toDouble(); + final switchMs = t0 + 2 * 24 * 3600 * 1000.0 + 1 * 3600 * 1000.0; + int offAt(int s) => s * 1000 < switchMs ? 3600 : 0; + final vals = []; + final samples = []; + for (var i = 0; i < 4 * 24; i++) { + final tMs = t0 + i * 3600 * 1000.0; + final v = 2000 + 100 * math.sin(i * 1.7) + (i % 5) * 13.0; + vals.add(v); + samples.add(AdcSample(tMs, v)); + } + final m = tempCircadian(samples, + deviceFamily: 'gen4', epochMin: 60, tzOffsetResolver: offAt); + final grand = vals.reduce((a, b) => a + b) / vals.length; + var diffSq = 0.0, varTot = 0.0; + for (var i = 1; i < vals.length; i++) { + diffSq += math.pow(vals[i] - vals[i - 1], 2); + } + for (final v in vals) { + varTot += math.pow(v - grand, 2); + } + final expectedIv = (diffSq / (vals.length - 1)) / (varTot / vals.length); + expect( + m.value!.nonparam!.intradailyVariability, closeTo(expectedIv, 1e-9)); + }); + + test('REGRESSION: a fit with no resolvable phase is not reported', () { + // A pure 12-h rhythm over whole days has no 24-h component, so the + // zero-amplitude test cannot reject and the acrophase is undefined. The + // nonparametric half still has variance to work with. + final samples = []; + for (var i = 0; i < 4 * 24 * 6; i++) { + final tMs = i * 10 * 60 * 1000.0; + final h = tMs / 3.6e6; + samples + .add(AdcSample(tMs, 2000 + 300 * math.cos(2 * math.pi * h / 12))); + } + final m = tempCircadian(samples, + deviceFamily: 'gen4', epochMin: 60, tzOffsetResolver: _utc); + expect(m.present, isTrue); + expect(m.value!.cosinorFit, isNull); + expect(m.value!.nonparam, isNotNull); + expect(m.confidence, 0.3); + }); + test('RD-14: IS is taken about the GRAND mean, not the profile mean', () { // Unbalanced epoch-of-day coverage — the normal case on real data, where // per-day 24-bin coverage runs 14/24 to 24/24 — plus a level drift across @@ -64,7 +156,8 @@ void main() { } } } - final np = tempCircadian(samples, deviceFamily: 'gen4', epochMin: 60) + final np = tempCircadian(samples, + deviceFamily: 'gen4', epochMin: 60, tzOffsetResolver: _utc) .value! .nonparam!; // van Someren's IS puts BOTH variances about the grand mean. Using the @@ -86,13 +179,13 @@ void main() { if (h == 3) continue; for (var k = 0; k < 6; k++) { final tHours = day * 24.0 + h + k / 6.0; - final adc = - 2000 + 150 * math.cos(2 * math.pi / 24 * (tHours - 4)); + final adc = 2000 + 150 * math.cos(2 * math.pi / 24 * (tHours - 4)); samples.add(AdcSample(tHours * 3.6e6, adc)); } } } - final np = tempCircadian(samples, deviceFamily: 'gen4', epochMin: 60) + final np = tempCircadian(samples, + deviceFamily: 'gen4', epochMin: 60, tzOffsetResolver: _utc) .value! .nonparam!; // Correct value with divisor = epochsPerDay; strictly less than the @@ -109,7 +202,8 @@ void main() { samples.add( AdcSample(tMs, 2000 + 300 * math.cos(2 * math.pi / 24 * tHours))); } - final m = tempCircadian(samples, deviceFamily: 'gen4', epochMin: 60); + final m = tempCircadian(samples, + deviceFamily: 'gen4', epochMin: 60, tzOffsetResolver: _utc); final j = m.value!.nonparam!.toJson(); // Not "absent when unobserved" — absent from the shape entirely, because // RA on a median-centred series divides by a quantity that crosses zero. @@ -194,15 +288,16 @@ void main() { // 3 days of an identical ADC reading -> varTot == 0 and diffN's // squared-difference sum == 0. IS/IV are both undefined divisions by // zero here, not "perfectly stable / perfectly regular" measurements. - // The cosinor half can still fit a (degenerate) phase on flat input, - // so the overall Metric stays present — it's specifically `nonparam` - // that must come back null instead of IS=0.0/IV=0.0. + // The cosinor half has no resolvable phase on flat input either, so + // nothing is left to report and the whole metric is absent. final samples = [ for (var i = 0; i < 3 * 24 * 6; i++) AdcSample(i * 10 * 60 * 1000.0, 2000.0), ]; - final m = tempCircadian(samples, deviceFamily: 'gen4', epochMin: 60); - expect(m.value!.nonparam, isNull); + final m = tempCircadian(samples, + deviceFamily: 'gen4', epochMin: 60, tzOffsetResolver: _utc); + expect(m.present, isFalse); + expect(m.value, isNull); }); }); @@ -255,6 +350,33 @@ void main() { expect(out.last.z, isNull); expect(out.last.confidence, 0); }); + + test('a short unscored gap carries the elevated run', () { + // Fever nights with an unscored night between each (strap re-seated, + // ambient-excluded): the run must survive the holes. + final temp = [ + for (var i = 0; i < 14; i++) 2000.0 + (i.isEven ? 3 : -3), + 2060, null, 2065, null, 2070, + ]; + final out = tempIllnessFlag(dates(temp.length), temp, + baselineDays: 21, zThresh: 2.0, persistDays: 2, minBaseline: 7); + expect(out[14].flag, TempFlag.normal); + expect(out[15].z, isNull); + expect(out[16].flag, TempFlag.elevated); + expect(out[18].flag, TempFlag.elevated); + }); + + test('a long unscored gap still ends the elevated run', () { + final temp = [ + for (var i = 0; i < 14; i++) 2000.0 + (i.isEven ? 3 : -3), + 2060, null, null, null, 2065, + ]; + final out = tempIllnessFlag(dates(temp.length), temp, + baselineDays: 21, zThresh: 2.0, persistDays: 2, minBaseline: 7); + // Gap of 4 calendar days > 1 + maxGapNights: a fresh single night. + expect(out.last.z, isNotNull); + expect(out.last.flag, TempFlag.normal); + }); }); // ------------------------------------------------------------------------- @@ -341,6 +463,70 @@ void main() { expect(out[20].candidate, isTrue); expect(out[20].flagged, isFalse, reason: 'one night cannot flag'); }); + + test('a run of unusually GOOD nights is not an illness anomaly', () { + // Recovery after a holiday / deload: RHR well down, HRV well up, resp + // down. Every oriented deviation points the healthy way, so the + // one-sided {RHR↑, HRV↓, temp↑, resp↑} detector must stay quiet. + final r = math.Random(11); + double g() { + final u1 = 1 - r.nextDouble(), u2 = r.nextDouble(); + return math.sqrt(-2 * math.log(u1)) * math.cos(2 * math.pi * u2); + } + + final feats = [ + for (var i = 0; i < 60; i++) + AnomalyFeatures( + rhr: (55 + 1.5 * g()).roundToDouble(), + hrv: 60 + 6 * g(), + resp: 14.5 + 0.4 * g()), + for (var i = 0; i < 3; i++) + const AnomalyFeatures(rhr: 47, hrv: 95, resp: 12.8), + ]; + final out = multivariateAnomaly(dates(feats.length), feats, + baselineDays: 28, minBaseline: 10, persistDays: 2); + for (var i = 60; i < 63; i++) { + expect(out[i].mahalanobis, isNotNull); + expect(out[i].candidate, isFalse, reason: 'night $i improved'); + expect(out[i].flagged, isFalse, reason: 'night $i improved'); + // The drivers still report the signed (healthy) deviation. + expect(out[i].drivers.every((d) => d.contribution < 0), isTrue); + } + }); + + test('clamping a healthy-side z never raises d² above the two-sided one', + () { + // Temp and resp anti-correlated in the baseline, so the inverse + // correlation has a positive off-diagonal. Tonight temp is well up and + // resp well down: zeroing resp's negative z removes a negative cross + // term, which on its own pushes d² past a gate the two-sided distance + // (the one the inflation table was measured on) stays under. + final r = math.Random(7); + double g() { + final u1 = r.nextDouble() + 1e-12, u2 = r.nextDouble(); + return math.sqrt(-2 * math.log(u1)) * math.cos(2 * math.pi * u2); + } + + final feats = []; + for (var k = 0; k < 28; k++) { + final u = g(), w = g(); + feats.add( + AnomalyFeatures(temp: 100 + 2 * u, resp: 15 + (-0.6 * u + 0.8 * w))); + } + for (var k = 28; k < 30; k++) { + final win = feats.sublist(k - 28, k); + final tc = [for (final f in win) f.temp!]; + final rc = [for (final f in win) f.resp!]; + feats.add(AnomalyFeatures( + temp: median(tc)! + 5 * mad(tc)!, + resp: median(rc)! - 3 * mad(rc)!)); + } + final out = multivariateAnomaly(dates(feats.length), feats); + for (var i = 28; i < 30; i++) { + expect(out[i].candidate, isFalse, reason: 'night $i'); + expect(out[i].flagged, isFalse, reason: 'night $i'); + } + }); }); // ------------------------------------------------------------------------- @@ -516,6 +702,29 @@ void main() { expect(cusumChangePoints(x, dates: contiguous, h: 5.0), hasLength(cusumChangePoints(x, h: 5.0).length)); }); + + // A one-night gap every week (the band charging overnight) used to restart + // the regime on every gap, so the baseline never reached minBaseline and a + // sustained shift was never reported with dates passed. + test('a weekly one-night gap keeps the baseline (short gaps)', () { + const noise10 = [0, 1, -1, 2, -2, 1, -1, 0, 2, -2]; + String day(int n) => DateTime.utc(2026, 1, 1) + .add(Duration(days: n)) + .toIso8601String() + .substring(0, 10); + final x = []; + final dates = []; + for (var d = 0; d < 120; d++) { + if (d % 7 == 6) continue; // charging night + x.add((d < 60 ? 55.0 : 63.0) + noise10[d % 10]); + dates.add(day(d)); + } + final dets = cusumChangePoints(x, dates: dates, h: 5.0); + expect(dets, hasLength(1)); + expect(dets.single.direction, 1); + expect( + dates[dets.single.index].compareTo(day(60)), greaterThanOrEqualTo(0)); + }); }); // ------------------------------------------------------------------------- @@ -592,7 +801,20 @@ void main() { // is asserted as a refusal in its own test below. The rescue itself is // unchanged; only a baseline with resolvable dispersion reaches it. final tightBase = [ - 52, 52, 52, 52, 52, 52, 52, 52, 52, 52, 52, 49, 55, 53 + 52, + 52, + 52, + 52, + 52, + 52, + 52, + 52, + 52, + 52, + 52, + 49, + 55, + 53 ]; // minInputs:1 — the subject here is the MAD==0 rescue, not the RD-04 // minimum-inputs gate (which has its own test below). @@ -616,8 +838,9 @@ void main() { // The measured production failure: [58,58,59] bpm + a 52 bpm night => // MAD 0 => mean/SD fallback on an SD of 0.577 bpm => z = -10.97 => // readiness 99.949 at confidence 0.60. Maximal recovery out of 6 bpm. - final m = readinessComposite([rhrInput(52.0, [58, 58, 59])], - minInputs: 1, minWeightSum: 0.0); + final m = readinessComposite([ + rhrInput(52.0, [58, 58, 59]) + ], minInputs: 1, minWeightSum: 0.0); expect(m.present, isFalse, reason: 'was: 99.949 at confidence 0.60'); expect(m.note, contains('need_baseline:have=3,need=14')); }); @@ -630,7 +853,20 @@ void main() { // clamp (a bounded-but-wrong score at confidence 0.60 is harder to catch // than an absurd one). final base = [ - 58, 58, 58, 58, 58, 58, 58, 58, 58, 58, 58, 58, 58, 59 + 58, + 58, + 58, + 58, + 58, + 58, + 58, + 58, + 58, + 58, + 58, + 58, + 58, + 59 ]; final m = readinessComposite([rhrInput(52.0, base)], minInputs: 1, minWeightSum: 0.0); @@ -856,8 +1092,7 @@ void main() { final feats = [ for (var i = 0; i < 12; i++) AnomalyFeatures( - rhr: i == 5 ? 59 : 58, - hrv: 40.0 + (i.isEven ? 1.0 : -1.0)), + rhr: i == 5 ? 59 : 58, hrv: 40.0 + (i.isEven ? 1.0 : -1.0)), const AnomalyFeatures(rhr: 60, hrv: 40.0), // ordinary 2-bpm move ]; final dates = [for (var i = 0; i < feats.length; i++) 'd$i']; @@ -886,7 +1121,20 @@ void main() { // so the baseline still has more than 1 bpm of SD — below that the input // is refused outright rather than scored (A4). final base = [ - 55, 55, 55, 55, 55, 55, 55, 55, 55, 55, 55, 55, 58, 51 + 55, + 55, + 55, + 55, + 55, + 55, + 55, + 55, + 55, + 55, + 55, + 55, + 58, + 51 ]; final m = readinessComposite([rhrInput(60, base)], minInputs: 1, minWeightSum: 0.0); @@ -898,7 +1146,20 @@ void main() { test('a dispersed baseline still discloses robust-z (median+MAD)', () { final base = [ - 50, 52, 54, 56, 58, 60, 62, 51, 53, 55, 57, 59, 61, 63 + 50, + 52, + 54, + 56, + 58, + 60, + 62, + 51, + 53, + 55, + 57, + 59, + 61, + 63 ]; final m = readinessComposite([rhrInput(70, base)], minInputs: 1, minWeightSum: 0.0); @@ -1026,8 +1287,7 @@ void main() { } }); - test('a cold segment 2 °C or more low is trimmed, 25 % of one refuses', - () { + test('a cold segment 2 °C or more low is trimmed, 25 % of one refuses', () { for (final drop in [200.0, 300.0, 500.0]) { final m = nightlySkinTemp( ringNight(coldFrom: 30, coldRecords: 12, coldDrop: drop), @@ -1262,4 +1522,48 @@ void main() { reason: 'a MAD over ten nights is not a known scale'); }); }); + + group('multivariateAnomaly one-sided clamp (regression)', () { + test('an improving feature cannot raise d2 above the two-sided distance', + () { + // Resp negatively coupled to RHR across the baseline, so the robust + // correlation is strongly negative. Tonight RHR is high (illness side) + // and resp is well below its median (the good side). Zeroing the good + // component used to drop the negative cross term and inflate d2 above + // the neutral-resp night and over the gate. + double noise(int i) => + math.sin(i * 1.7) + + 0.6 * math.sin(i * 0.31) + + 0.4 * math.cos(i * 2.9); + final base = [ + for (var i = 0; i < 28; i++) + () { + final rhr = (60 + 2 * noise(i)).roundToDouble(); + return AnomalyFeatures( + rhr: rhr, resp: 15 - 0.25 * (rhr - 60) + 0.2 * noise(i + 50)); + }(), + ]; + final dates = [ + for (var i = 0; i < 29; i++) + DateTime.utc(2026, 9, 1) + .add(Duration(days: i)) + .toIso8601String() + .substring(0, 10) + ]; + double d2(double resp) { + final out = multivariateAnomaly( + dates, [...base, AnomalyFeatures(rhr: 66, resp: resp)]); + final m = out.last.mahalanobis!; + return m * m; + } + + // Neutral = resp exactly at its baseline median (z = 0). + final resps = [for (final f in base) f.resp!]..sort(); + final median = (resps[13] + resps[14]) / 2; + final neutral = d2(median); + final improving = d2(median - 1.5); + expect(improving, lessThan(neutral - 1e-6), + reason: 'slower breathing must not score as more anomalous'); + }); + }); }