CR: feat/hr-estimator-guard (sleepypod/core#840)

modules/common/biometrics/manifest.v1.json

@@ -1,6 +1,6 @@
 {
   "id": "sleepypod-biometric-bundle",
-  "version": "1.0.0",
+  "version": "1.1.0",
   "schemaVersion": 1,
   "canonicalization": "finite-binary64-json-v1",
   "sampleRateHz": 500,
@@ -30,11 +30,20 @@
             "default": 5,
             "unit": "windows",
             "description": "Number of accepted heart-rate estimates retained by the tracker."
+          },
+          "minimumPeakScore": {
+            "type": "number",
+            "minimum": 0.1,
+            "maximum": 0.8,
+            "default": 0.4,
+            "unit": "score",
+            "description": "Minimum subharmonic-summation autocorrelation score for a heart-rate estimate. Weaker windows report no heart rate instead of a guess; raising it trades coverage for accuracy."
           }
         },
         "required": [
           "trackingToleranceBpm",
-          "trackingHistoryLength"
+          "trackingHistoryLength",
+          "minimumPeakScore"
         ],
         "additionalProperties": false
       },
@@ -375,7 +384,8 @@
   "defaults": {
     "heart_rate": {
       "trackingToleranceBpm": 15,
-      "trackingHistoryLength": 5
+      "trackingHistoryLength": 5,
+      "minimumPeakScore": 0.4
     },
     "breathing_rate": {
       "intervalLowerRatio": 0.5,
@@ -624,7 +634,7 @@
     }
   },
   "sourceHashes": {
-    "modules/common/biometrics/signal.py": "fff73469815673c8daebba9e2c9e7bcf97708939ad7a5639a52f0a53c9d75c18",
+    "modules/common/biometrics/signal.py": "0b1447c1162245232c053aeeb5cb3e75e7da3392c38a43c465f28596f084b32f",
     "modules/piezo-processor/main.py": "75db2ea2f005fd537d994b306332c1abf97c0c28abb496a8d3af03816ccf8550",
     "src/lib/sleep-stages.ts": "b1800e49ce8525d314333da6a29a63131a002d160bb294fa97b24df6e655e228",
     "modules/piezo-processor/runtime_adapter.py": "0a22a2fa05aae6569bfc4b53fba80d75e89f702f021ccb61f82ce64e284c5bb8",

modules/common/biometrics/signal.py

@@ -27,6 +27,7 @@ class ProcessingParameters:
     breathing_interval_upper_ratio: float
     hrv_harmonic_tolerance: float
     hrv_minimum_windows: int
+    heart_rate_minimum_peak_score: float
 
     def __post_init__(self):
         resolve_config([
@@ -36,6 +37,7 @@ def __post_init__(self):
             {"path": "breathing_rate.intervalUpperRatio", "value": self.breathing_interval_upper_ratio},
             {"path": "hrv_index.harmonicTolerance", "value": self.hrv_harmonic_tolerance},
             {"path": "hrv_index.minimumWindows", "value": self.hrv_minimum_windows},
+            {"path": "heart_rate.minimumPeakScore", "value": self.heart_rate_minimum_peak_score},
         ])
 
     @classmethod
@@ -46,7 +48,8 @@ def from_config(cls, config):
                    config["breathing_rate"]["intervalLowerRatio"],
                    config["breathing_rate"]["intervalUpperRatio"],
                    config["hrv_index"]["harmonicTolerance"],
-                   config["hrv_index"]["minimumWindows"])
+                   config["hrv_index"]["minimumWindows"],
+                   config["heart_rate"]["minimumPeakScore"])
 
 
 DEFAULT_PARAMETERS = ProcessingParameters.from_config(manifest()["defaults"])
@@ -111,7 +114,8 @@ def _compute_autocorr(filtered: np.ndarray, fs: float) -> Optional[np.ndarray]:
 
 def subharmonic_summation_hr(samples: np.ndarray, fs: float = SAMPLE_RATE,
                              bpm_range: tuple = (45, 120),
-                             n_harmonics: int = 3):
+                             n_harmonics: int = 3,
+                             min_score: float = DEFAULT_PARAMETERS.heart_rate_minimum_peak_score):
     """Subharmonic summation for robust fundamental HR detection.
 
     Only scores autocorrelation PEAKS (not every lag) to avoid edge effects
@@ -120,7 +124,10 @@ def subharmonic_summation_hr(samples: np.ndarray, fs: float = SAMPLE_RATE,
 
     Hermes 1988 (JASA), adapted for BCG by Bruser et al. 2011.
 
-    Returns (heart_rate_bpm, score) or (None, 0.0) on failure.
+    Returns (heart_rate_bpm, score) or (None, 0.0) on failure. A best score
+    below min_score also returns (None, score): on a weak heartbeat the tallest
+    autocorrelation peak is usually not the heart (the true rate is often not a
+    candidate at all), so reporting nothing beats reporting a confident guess.
     """
     filtered = _bandpass(samples, HR_BAND[0], HR_BAND[1], fs)
     acr = _compute_autocorr(filtered, fs)
@@ -161,7 +168,7 @@ def subharmonic_summation_hr(samples: np.ndarray, fs: float = SAMPLE_RATE,
             best_score = score
             best_lag = lag
 
-    if best_lag is None or best_score < 0.1:
+    if best_lag is None or best_score < min_score:
         return None, best_score
 
     hr = 60.0 * fs / best_lag
@@ -483,7 +490,7 @@ def estimate(self, timestamp_s: float, hr_samples: np.ndarray,
         # Match live warmup exactly; never pad unavailable history with zeroes.
         if len(hr_samples) < int(10 * fs):
             return WindowEstimate(timestamp_s, None, None, None, None, 0.0, ("hr_warmup",))
-        hr_raw, hr_score = subharmonic_summation_hr(hr_samples, fs)
+        hr_raw, hr_score = subharmonic_summation_hr(hr_samples, fs, min_score=self.parameters.heart_rate_minimum_peak_score)
         hr = self.hr_tracker.update(hr_raw, hr_score)
         br = compute_breathing_rate(breathing_samples, fs, self.parameters)
         hrv = (compute_hrv(hrv_samples, fs, self.parameters)

modules/piezo-processor/test_shared_signal.py

@@ -77,3 +77,26 @@ def test_explicit_event_time_and_warmup_never_consult_wall_clock():
         estimator.estimate(20, data, data, data, fs=250)
     with pytest.raises(ValueError, match="finite"):
         estimator.estimate(float("nan"), data, data, data)
+
+
+def test_minimum_peak_score_abstains_on_weak_windows_only():
+    # A heartbeat buried in noise still yields a tallest peak, just a weak one.
+    weak = samples()[-15000:] + np.random.default_rng(7).normal(0, 2.5, 15000)
+    hr, score = core.subharmonic_summation_hr(weak)
+    assert hr is not None and 0.1 < score < 0.8
+    assert core.subharmonic_summation_hr(weak, min_score=score - 1e-6)[0] == pytest.approx(hr)
+    assert core.subharmonic_summation_hr(weak, min_score=score + 1e-6) == (None, pytest.approx(score))
+    strict = core.WindowEstimator(parameters("heart_rate.minimumPeakScore", score + 1e-6))
+    estimate = strict.estimate(1700000000, weak, samples()[-30000:], samples())
+    assert estimate.heart_rate is None and "no_hr" in estimate.rejection_reasons
+    # 0.1 is the previous fixed floor; the shipped default is now 0.4.
+    legacy = core.WindowEstimator(parameters("heart_rate.minimumPeakScore", 0.1))
+    assert legacy.estimate(1700000000, weak, samples()[-30000:], samples()).heart_rate == pytest.approx(hr)
+    assert core.DEFAULT_PARAMETERS.heart_rate_minimum_peak_score == 0.4
+    assert core.subharmonic_summation_hr(weak)[0] == (hr if score >= 0.4 else None)
+
+
+@pytest.mark.parametrize("value", [0.05, 0.85])
+def test_minimum_peak_score_range_is_validated(value):
+    with pytest.raises(ValueError):
+        parameters("heart_rate.minimumPeakScore", value)

scripts/biometrics-local-fit/README.md

@@ -68,8 +68,14 @@ Outputs:
 - `out/<candidate>/<split>.json`: per-night HR, breathing and HRV error, plus
   the full evaluator report. A metric without paired reference data reads
   `no_paired_reference`, never zero error.
-- `out/ranking.json` and `out/ranking.md`: ranked by the tuning split's
-  equal-night mean HR MAE difference from defaults, with validation alongside.
+- `out/ranking.json` and `out/ranking.md`: ranked by tuning HR error on each
+  candidate's own readings, among candidates that keep at least `--min-coverage`
+  (default 75%) of the defaults' Watch coverage; validation and holdout are shown
+  alongside. Each split also shows the common-support difference against
+  defaults. That difference reads 0 for a change that only withholds readings
+  (such as `heart_rate.minimumPeakScore`), which is why it is not the ranking key.
+- `--save-readings` also writes every replay's derived readings to
+  `out/readings/<night>-<config>.json` for error analysis.
 - `out/run.json`: timing, worker count and the Python build identity.
 
 ## 4. Check parity, then confirm on the pod

scripts/biometrics-local-fit/cli.ts

@@ -8,7 +8,7 @@ import { compareWithPod, parseLocalReport, rankCandidates, rankingMarkdown } fro
 
 const usage = `Usage:
   verify BUNDLE
-  fit BUNDLE CANDIDATES.json OUTPUT_DIR [--jobs N] [--splits tuning,validation] [--include-holdout] [--protocol FILE] [--python PATH]
+  fit BUNDLE CANDIDATES.json OUTPUT_DIR [--jobs N] [--splits tuning,validation] [--include-holdout] [--save-readings] [--min-coverage 0.75] [--protocol FILE] [--python PATH]
   sweep BUNDLE GRID.json|default OUTPUT_DIR [same options as fit]
   parity OUTPUT_DIR POD_RESULT.json [POD_RESULT.json ...]
 BUNDLE is the tar from GET /api/biometric-tuning/datasets/{id}/export, or a directory extracted from it.`
@@ -17,7 +17,7 @@ function options(args: string[]) {
   const values: Record<string, string> = {}
   const positional: string[] = []
   for (let index = 0; index < args.length; index++) {
-    if (args[index] === '--include-holdout') values['include-holdout'] = 'true'
+    if (args[index] === '--include-holdout' || args[index] === '--save-readings') values[args[index].slice(2)] = 'true'
     else if (args[index].startsWith('--')) values[args[index].slice(2)] = args[++index] ?? ''
     else positional.push(args[index])
   }
@@ -43,6 +43,7 @@ async function fit(bundlePath: string, candidates: Candidate[], outputDirectory:
   const splits = (flags.splits?.split(',') ?? [...LOCAL_SPLITS, ...(includeHoldout ? ['holdout'] : [])]) as LocalSplit[]
   const started = performance.now()
   const result = await fitCandidates({ bundle, candidates, outputDirectory, splits, includeHoldout,
+    readingsDirectory: flags['save-readings'] === 'true' ? join(outputDirectory, 'readings') : undefined,
     jobs: flags.jobs ? Number(flags.jobs) : undefined,
     protocol: flags.protocol ? await readJson(flags.protocol) as never : undefined,
     python: flags.python ? defaultPythonLocation(flags.python) : undefined,
@@ -51,7 +52,7 @@ async function fit(bundlePath: string, candidates: Candidate[], outputDirectory:
     await mkdir(join(outputDirectory, report.candidate.id), { recursive: true })
     await writeFile(join(outputDirectory, report.candidate.id, `${report.split}.json`), `${JSON.stringify(report)}\n`)
   }
-  const ranking = rankCandidates(result.reports)
+  const ranking = rankCandidates(result.reports, flags['min-coverage'] ? Number(flags['min-coverage']) : undefined)
   await writeFile(join(outputDirectory, 'ranking.json'), `${JSON.stringify(ranking, null, 2)}\n`)
   await writeFile(join(outputDirectory, 'ranking.md'), rankingMarkdown(ranking))
   const totalWallMs = performance.now() - started

scripts/biometrics-local-fit/fit.ts

@@ -96,6 +96,8 @@ export interface FitOptions {
   splits?: readonly LocalSplit[]
   /** Score holdout nights too, once, for a candidate that is already chosen. */
   includeHoldout?: boolean
+  /** Also write each replay's derived readings here, one JSON file per night and configuration. */
+  readingsDirectory?: string
   jobs?: number
   protocol?: z.input<typeof experimentProtocolSchema>
   python?: PythonLocation
@@ -145,6 +147,11 @@ export async function fitCandidates(options: FitOptions) {
       { version: 'retained-waveform-work-v1', input, configs: { [snapshot.side]: resolved, [otherSideSide]: defaults.resolved } },
       { timeoutMs: 3600000, expectedBuildHash: build.hash })
     outputs.set(`${item.cohortId}:${item.hash}`, { cohortId: item.cohortId, configHash: item.hash, ...output })
+    if (options.readingsDirectory) {
+      await mkdir(options.readingsDirectory, { recursive: true })
+      await writeFile(join(options.readingsDirectory, `${item.cohortId}-${item.hash.slice(0, 12)}.json`),
+        `${JSON.stringify({ cohortId: item.cohortId, configHash: item.hash, side: snapshot.side, readings: output.readings })}\n`)
+    }
     options.onProgress?.(++done, items.length, { cohortId: item.cohortId,
       candidateIds: candidates.filter(row => row.resolvedHash === item.hash).map(row => row.id), wallMs: performance.now() - itemStarted })
   })

scripts/biometrics-local-fit/report.test.ts

@@ -24,17 +24,18 @@ afterAll(() => {
 })
 
 /** Score a candidate with fixed synthetic readings (references are 62 bpm). */
-function localReport(id: string, overrides: Array<{ path: string, value: unknown }>, split: 'tuning' | 'validation' | 'holdout', heartRates: { default: number, candidate: number }) {
+function localReport(id: string, overrides: Array<{ path: string, value: unknown }>, split: 'tuning' | 'validation' | 'holdout', heartRates: { default: number | Array<number | null>, candidate: number | Array<number | null> }) {
   const defaults = { manifestHash, resolvedHash: configHash(defaultConfig), resolved: defaultConfig }
   const candidate = resolveConfig(overrides as never)
   const configs = [{ label: 'default' as const, revisionId: null, ...defaults },
     { label: 'candidate' as const, revisionId: null, manifestHash, resolvedHash: candidate.hash, resolved: candidate.resolved }]
   const night = exportNights.find(row => row.split === split) ?? exportNights[0]
   const artifact = fixture.snapshot.version === 'waveform-dataset-v1' ? fixture.snapshot.waveformInputs.find(row => row.cohortId === night.cohortId)?.input.artifacts[0].artifact.sha256 ?? '' : ''
   const outputs = [...new Map(configs.map(config => [config.resolvedHash, config])).values()].map((config) => {
-    const heartRate = config.label === 'default' ? heartRates.default : heartRates.candidate
-    const readings = Array.from({ length: 10 }, (_, index) => ({ side: 'left' as const, timestampMs: night.startMs + index * 60000, heartRate, hrvIndex: 30,
-      hrvMetric: 'window_ibi_variability_v1' as const, breathingRate: 14, qualityScore: 1, flags: null, hrRaw: heartRate, sourceSide: 'left' as const,
+    const chosen = config.label === 'default' ? heartRates.default : heartRates.candidate
+    const rate = (index: number) => Array.isArray(chosen) ? chosen[index] : chosen
+    const readings = Array.from({ length: 10 }, (_, index) => ({ side: 'left' as const, timestampMs: night.startMs + index * 60000, heartRate: rate(index), hrvIndex: 30,
+      hrvMetric: 'window_ibi_variability_v1' as const, breathingRate: 14, qualityScore: 1, flags: null, hrRaw: rate(index), sourceSide: 'left' as const,
       configHash: config.resolvedHash, limitations: [] }))
     return { cohortId: night.cohortId, configHash: config.resolvedHash, readings, run: { version: 'retained-waveform-run-v1' as const, build_hash: pythonBuildHash,
       decoded_records: 10, raw_records: 10, counts: { records: 10, accepted_frames: 10, empty_channel_frames: 0, pump_gated_frames: 0, preroll_outputs: 0, scored_outputs: readings.length },
@@ -66,31 +67,54 @@ describe('local fit reports', () => {
     expect(() => parseLocalReport({ ...report, label: 'pod' })).toThrow()
   })
 
-  it('ranks by tuning HR error difference and keeps validation alongside', () => {
+  it('ranks by tuning HR error on own readings and keeps validation alongside', () => {
     const better = [{ path: 'heart_rate.trackingToleranceBpm', value: 10 }], worse = [{ path: 'heart_rate.trackingToleranceBpm', value: 25 }]
     const reports = [
       localReport('worse', worse, 'tuning', { default: 60, candidate: 58 }), localReport('worse', worse, 'validation', { default: 60, candidate: 59 }),
       localReport('better', better, 'tuning', { default: 60, candidate: 61.5 }), localReport('better', better, 'validation', { default: 60, candidate: 60 }),
       localReport('defaults', [], 'tuning', { default: 60, candidate: 60 }), localReport('defaults', [], 'validation', { default: 60, candidate: 60 }),
     ]
     const ranking = rankCandidates(reports)
-    expect(ranking).toMatchObject({ label: 'exploratory-local', evidence: 'not_activation_evidence' })
+    expect(ranking).toMatchObject({ label: 'exploratory-local', evidence: 'not_activation_evidence', minimumCoverage: 0.75 })
     expect(ranking.rows.map(row => [row.rank, row.id])).toEqual([[1, 'better'], [2, 'defaults'], [3, 'worse']])
-    expect(ranking.rows[0].heartRate).toEqual({ tuning: { maeDifference: -1.5, mae: 0.5, scoredNights: 1, coverageReducedNights: 0 },
-      validation: { maeDifference: 0, mae: 2, scoredNights: 1, coverageReducedNights: 0 } })
+    expect(ranking.rows[0].heartRate).toEqual({
+      tuning: { maeDifference: -1.5, mae: 0.5, ownMae: 0.5, coverage: 1, scoredNights: 1, coverageReducedNights: 0 },
+      validation: { maeDifference: 0, mae: 2, ownMae: 2, coverage: 1, scoredNights: 1, coverageReducedNights: 0 } })
+    expect(rankingMarkdown(ranking)).toContain('| 1 | better | 0.500 | -1.500 | 100% | 2.000 | 0.000 | 100% | n/a | n/a |')
     const tied = rankCandidates([localReport('b-two', [...better, { path: 'heart_rate.trackingHistoryLength', value: 8 }], 'tuning', { default: 60, candidate: 61 }),
       localReport('a-one', [...better, { path: 'hrv_index.minimumWindows', value: 10 }], 'tuning', { default: 60, candidate: 61 })])
     expect(tied.rows.map(row => [row.id, row.changedDefaults])).toEqual([['a-one', 1], ['b-two', 2]])
-    expect(rankingMarkdown(ranking)).toContain('| 1 | better | 0.500 | -1.500 | 2.000 | 0.000 | 0/0 | n/a | n/a |')
+  })
+
+  it('credits abstaining on bad windows but not below the coverage floor', () => {
+    const gate = (value: number) => [{ path: 'heart_rate.minimumPeakScore', value }]
+    // References sit on readings 0, 2, 4, 6 and 8 at 62 bpm. Defaults are right
+    // except one wild reading; each gate only removes readings, never changes them.
+    const defaults = [62, 60, 62, 60, 62, 60, 62, 60, 90, 60]
+    const dropsWild = defaults.map((value, index) => index === 8 ? null : value)
+    const dropsMost = defaults.map((value, index) => index >= 2 && index % 2 === 0 ? null : value)
+    const ranking = rankCandidates([
+      localReport('defaults', [], 'tuning', { default: defaults, candidate: defaults }),
+      localReport('gate-0.45', gate(0.45), 'tuning', { default: defaults, candidate: dropsWild }),
+      localReport('gate-0.6', gate(0.6), 'tuning', { default: defaults, candidate: dropsMost }),
+    ])
+    const rows = Object.fromEntries(ranking.rows.map(row => [row.id, row]))
+    expect(rows.defaults.heartRate.tuning).toMatchObject({ ownMae: 5.6, coverage: 1 })
+    expect(rows['gate-0.45'].heartRate.tuning).toMatchObject({ ownMae: 0, coverage: 0.8, maeDifference: 0 })
+    expect(rows['gate-0.6'].heartRate.tuning).toMatchObject({ ownMae: 0, coverage: 0.2, maeDifference: 0 })
+    // Common support scores both gates as no change; own error with the coverage floor does not.
+    expect(ranking.rows.map(row => row.id)).toEqual(['gate-0.45', 'defaults', 'gate-0.6'])
+    expect(rankCandidates([localReport('gate-0.6', gate(0.6), 'tuning', { default: defaults, candidate: dropsMost }),
+      localReport('defaults', [], 'tuning', { default: defaults, candidate: defaults })], 0.1).rows[0].id).toBe('gate-0.6')
   })
 
   it('shows holdout error alongside only when holdout nights were scored', () => {
     const overrides = [{ path: 'heart_rate.trackingToleranceBpm', value: 10 }]
     const without = rankingMarkdown(rankCandidates([localReport('tol10', overrides, 'tuning', { default: 60, candidate: 61 })]))
     expect(without).not.toContain('holdout')
     const ranking = rankCandidates([localReport('tol10', overrides, 'tuning', { default: 60, candidate: 61 }), localReport('tol10', overrides, 'holdout', { default: 60, candidate: 62 })])
-    expect(ranking.rows[0].heartRate.holdout).toEqual({ maeDifference: -2, mae: 0, scoredNights: 1, coverageReducedNights: 0 })
-    expect(rankingMarkdown(ranking)).toContain('| 1 | tol10 | 1.000 | -1.000 | n/a | n/a | 0.000 | -2.000 |')
+    expect(ranking.rows[0].heartRate.holdout).toEqual({ maeDifference: -2, mae: 0, ownMae: 0, coverage: 1, scoredNights: 1, coverageReducedNights: 0 })
+    expect(rankingMarkdown(ranking)).toContain('| 1 | tol10 | 1.000 | -1.000 | 100% | n/a | n/a | n/a | 0.000 | -2.000 | 100% |')
   })
 
   it('compares local default HR error with pod results within 0.1 bpm', () => {