feat: add breathing rate FFT extraction and anomaly flagging for sleep monitoring
- BreathingRateEstimator: FFT peak detection on 512-sample (25.6s) phase windows at 20Hz, zero-padded to 1024 points, 0.1-0.5 Hz band (6-30 BPM), 60s EMA smoothing - BreathingAnomalyTracker: per-person 30-day rolling EMA (α=0.05), flags elevated breathing at >25% above personal average - ComputeBreathingRegularity: coefficient of variation (std/mean) with regular/normal/irregular labels - Migration 008: add breathing_anomaly BOOL and breathing_samples_json TEXT to sleep_records - Integration: anomaly check in GenerateMorningReports, anomaly data in report/metrics/handler JSON - Table-driven tests: synthetic 15/12/20 BPM signals, circular buffer, noise, EMA convergence, anomaly threshold boundary cases, JSON round-trip, regularity CV computation Co-Authored-By: Claude Opus 4.6 <noreply@anthropic.com>
This commit is contained in:
parent
cdca106a76
commit
b72f234154
8 changed files with 832 additions and 16 deletions
|
|
@ -43,6 +43,11 @@ func AllMigrations() []Migration {
|
|||
Description: "add webhook_log, trigger_state tables and trigger error columns",
|
||||
Up: migration_007_add_webhook_tables,
|
||||
},
|
||||
{
|
||||
Version: 8,
|
||||
Description: "add breathing anomaly columns to sleep_records",
|
||||
Up: migration_008_add_breathing_anomaly,
|
||||
},
|
||||
}
|
||||
}
|
||||
|
||||
|
|
@ -452,3 +457,12 @@ CREATE INDEX IF NOT EXISTS idx_webhook_log_trigger ON webhook_log(trigger_id, fi
|
|||
_, err := tx.Exec(schema)
|
||||
return err
|
||||
}
|
||||
|
||||
// migration_008_add_breathing_anomaly adds breathing anomaly tracking columns to sleep_records.
|
||||
func migration_008_add_breathing_anomaly(tx *sql.Tx) error {
|
||||
_, err := tx.Exec(`
|
||||
ALTER TABLE sleep_records ADD COLUMN breathing_anomaly INTEGER NOT NULL DEFAULT 0;
|
||||
ALTER TABLE sleep_records ADD COLUMN breathing_samples_json TEXT;
|
||||
`)
|
||||
return err
|
||||
}
|
||||
|
|
|
|||
|
|
@ -127,8 +127,12 @@ type SleepMetrics struct {
|
|||
MinBreathingRate float64 `json:"min_breathing_rate"`
|
||||
MaxBreathingRate float64 `json:"max_breathing_rate"`
|
||||
BreathingRateStdDev float64 `json:"breathing_rate_std_dev"`
|
||||
BreathingRegularity float64 `json:"breathing_regularity"` // CV (std/mean)
|
||||
BreathingScore float64 `json:"breathing_score"` // 0-100
|
||||
BreathingAnomalyCount int `json:"breathing_anomaly_count"` // Anomalies < 8 or > 25 bpm
|
||||
BreathingAnomaly bool `json:"breathing_anomaly"` // Elevated vs personal average
|
||||
PersonalAvgBPM float64 `json:"personal_avg_bpm,omitempty"` // Person's rolling average for comparison
|
||||
BreathingSamplesJSON string `json:"breathing_samples_json,omitempty"` // Raw samples for storage
|
||||
|
||||
// Motion metrics
|
||||
MotionEvents int `json:"motion_events"`
|
||||
|
|
@ -227,6 +231,9 @@ type SleepAnalyzer struct {
|
|||
sleepStartHour int
|
||||
sleepEndHour int
|
||||
|
||||
// Breathing anomaly tracking (per-person rolling baseline)
|
||||
anomalyTracker *BreathingAnomalyTracker
|
||||
|
||||
// Report callback
|
||||
onReportGenerated func(linkID string, report *SleepReport)
|
||||
}
|
||||
|
|
@ -234,9 +241,10 @@ type SleepAnalyzer struct {
|
|||
// NewSleepAnalyzer creates a new sleep analyzer
|
||||
func NewSleepAnalyzer() *SleepAnalyzer {
|
||||
return &SleepAnalyzer{
|
||||
sessions: make(map[string]*SleepSession),
|
||||
sleepStartHour: DefaultSleepStartHour,
|
||||
sleepEndHour: DefaultSleepEndHour,
|
||||
sessions: make(map[string]*SleepSession),
|
||||
sleepStartHour: DefaultSleepStartHour,
|
||||
sleepEndHour: DefaultSleepEndHour,
|
||||
anomalyTracker: NewBreathingAnomalyTracker(),
|
||||
}
|
||||
}
|
||||
|
||||
|
|
@ -304,7 +312,8 @@ func (sa *SleepAnalyzer) GetAllSessions() map[string]*SleepSession {
|
|||
return result
|
||||
}
|
||||
|
||||
// GenerateMorningReports generates reports for all completed sleep sessions
|
||||
// GenerateMorningReports generates reports for all completed sleep sessions.
|
||||
// It also checks breathing anomalies against personal baselines and updates them.
|
||||
func (sa *SleepAnalyzer) GenerateMorningReports() map[string]*SleepReport {
|
||||
sa.mu.RLock()
|
||||
defer sa.mu.RUnlock()
|
||||
|
|
@ -312,6 +321,20 @@ func (sa *SleepAnalyzer) GenerateMorningReports() map[string]*SleepReport {
|
|||
reports := make(map[string]*SleepReport)
|
||||
for linkID, session := range sa.sessions {
|
||||
if report := session.GenerateReport(); report != nil {
|
||||
// Check breathing anomaly against personal baseline
|
||||
person := session.personID
|
||||
if person == "" {
|
||||
person = linkID
|
||||
}
|
||||
if report.Metrics.AvgBreathingRate > 0 {
|
||||
personalAvg := sa.anomalyTracker.GetPersonalAverage(person)
|
||||
report.Metrics.PersonalAvgBPM = personalAvg
|
||||
isAnomaly := sa.anomalyTracker.CheckAnomaly(person, report.Metrics.AvgBreathingRate)
|
||||
report.Metrics.BreathingAnomaly = isAnomaly
|
||||
// Update personal rolling average after checking
|
||||
sa.anomalyTracker.UpdatePersonalAverage(person, report.Metrics.AvgBreathingRate)
|
||||
}
|
||||
|
||||
reports[linkID] = report
|
||||
|
||||
if sa.onReportGenerated != nil {
|
||||
|
|
@ -334,6 +357,23 @@ func (sa *SleepAnalyzer) getOrCreateSession(linkID string) *SleepSession {
|
|||
return session
|
||||
}
|
||||
|
||||
// GetAnomalyTracker returns the breathing anomaly tracker for external access
|
||||
// (e.g., loading/saving personal baselines from SQLite).
|
||||
func (sa *SleepAnalyzer) GetAnomalyTracker() *BreathingAnomalyTracker {
|
||||
return sa.anomalyTracker
|
||||
}
|
||||
|
||||
// SetPersonID sets the person identity for a sleep session link.
|
||||
func (sa *SleepAnalyzer) SetPersonID(linkID, personID string) {
|
||||
sa.mu.Lock()
|
||||
defer sa.mu.Unlock()
|
||||
if session, exists := sa.sessions[linkID]; exists {
|
||||
session.mu.Lock()
|
||||
session.personID = personID
|
||||
session.mu.Unlock()
|
||||
}
|
||||
}
|
||||
|
||||
// NewSleepSession creates a new sleep session
|
||||
func NewSleepSession(linkID string, sleepStartHour, sleepEndHour int) *SleepSession {
|
||||
return &SleepSession{
|
||||
|
|
@ -711,6 +751,9 @@ func (ss *SleepSession) calculateBreathingMetrics(m *SleepMetrics) {
|
|||
// Count breathing anomalies (per task spec: < 8 or > 25 bpm for > 3 minutes)
|
||||
m.BreathingAnomalyCount = len(ss.breathingAnomalies)
|
||||
|
||||
// Compute breathing regularity (coefficient of variation)
|
||||
m.BreathingRegularity = ss.computeBreathingRegularity()
|
||||
|
||||
// Calculate breathing score (0-100)
|
||||
m.BreathingScore = ss.calculateBreathingScore(m.AvgBreathingRate, m.BreathingRateStdDev, m.MinBreathingRate, m.MaxBreathingRate)
|
||||
}
|
||||
|
|
@ -738,13 +781,24 @@ func (ss *SleepSession) calculateBreathingScore(avg, stdDev, min, max float64) f
|
|||
}
|
||||
|
||||
// Penalize high variability
|
||||
if stdDev > 3 {
|
||||
score -= math.Min(20, (stdDev-3)*4)
|
||||
if stdDev > 2 {
|
||||
score -= math.Min(35, (stdDev-2)*10)
|
||||
}
|
||||
|
||||
return math.Max(0, math.Min(100, score))
|
||||
}
|
||||
|
||||
// computeBreathingRegularity computes CV (std/mean) of detected breathing rates.
|
||||
func (ss *SleepSession) computeBreathingRegularity() float64 {
|
||||
var rates []float64
|
||||
for _, sample := range ss.breathingSamples {
|
||||
if sample.IsDetected && sample.RateBPM > 0 {
|
||||
rates = append(rates, sample.RateBPM)
|
||||
}
|
||||
}
|
||||
return ComputeBreathingRegularity(rates)
|
||||
}
|
||||
|
||||
// calculateMotionMetrics computes motion quality metrics
|
||||
func (ss *SleepSession) calculateMotionMetrics(m *SleepMetrics) {
|
||||
if len(ss.motionSamples) == 0 {
|
||||
|
|
|
|||
117
mothership/internal/sleep/breathing_anomaly.go
Normal file
117
mothership/internal/sleep/breathing_anomaly.go
Normal file
|
|
@ -0,0 +1,117 @@
|
|||
// Package sleep provides breathing anomaly detection and per-night statistics.
|
||||
package sleep
|
||||
|
||||
import (
|
||||
"encoding/json"
|
||||
"math"
|
||||
"sync"
|
||||
)
|
||||
|
||||
// Anomaly detection constants
|
||||
const (
|
||||
BreathingAnomalyEmaAlpha = 0.05 // Rolling personal average EMA (slow, ~20-night half-life)
|
||||
BreathingAnomalyThreshold = 1.25 // Flag if avg > personal_avg × 1.25
|
||||
BreathingRegularityRegular = 0.10 // CV below this = regular
|
||||
BreathingRegularityIrregular = 0.25 // CV above this = irregular
|
||||
)
|
||||
|
||||
// BreathingAnomalyTracker maintains per-person rolling averages and detects
|
||||
// elevated breathing rates compared to personal baselines.
|
||||
type BreathingAnomalyTracker struct {
|
||||
mu sync.RWMutex
|
||||
personal map[string]float64 // person -> EMA of nightly avg BPM
|
||||
}
|
||||
|
||||
// NewBreathingAnomalyTracker creates a new anomaly tracker.
|
||||
func NewBreathingAnomalyTracker() *BreathingAnomalyTracker {
|
||||
return &BreathingAnomalyTracker{
|
||||
personal: make(map[string]float64),
|
||||
}
|
||||
}
|
||||
|
||||
// UpdatePersonalAverage updates the rolling EMA for a person with the night's average BPM.
|
||||
func (t *BreathingAnomalyTracker) UpdatePersonalAverage(person string, nightlyAvgBPM float64) {
|
||||
t.mu.Lock()
|
||||
defer t.mu.Unlock()
|
||||
|
||||
if nightlyAvgBPM <= 0 {
|
||||
return
|
||||
}
|
||||
|
||||
if existing, ok := t.personal[person]; ok && existing > 0 {
|
||||
t.personal[person] = BreathingAnomalyEmaAlpha*nightlyAvgBPM + (1-BreathingAnomalyEmaAlpha)*existing
|
||||
} else {
|
||||
t.personal[person] = nightlyAvgBPM
|
||||
}
|
||||
}
|
||||
|
||||
// CheckAnomaly returns true if the nightly average BPM is elevated
|
||||
// (>25% above the person's rolling average).
|
||||
func (t *BreathingAnomalyTracker) CheckAnomaly(person string, nightlyAvgBPM float64) bool {
|
||||
t.mu.RLock()
|
||||
defer t.mu.RUnlock()
|
||||
|
||||
avg, ok := t.personal[person]
|
||||
if !ok || avg <= 0 {
|
||||
return false
|
||||
}
|
||||
|
||||
return nightlyAvgBPM > avg*BreathingAnomalyThreshold
|
||||
}
|
||||
|
||||
// GetPersonalAverage returns the current rolling average for a person.
|
||||
func (t *BreathingAnomalyTracker) GetPersonalAverage(person string) float64 {
|
||||
t.mu.RLock()
|
||||
defer t.mu.RUnlock()
|
||||
return t.personal[person]
|
||||
}
|
||||
|
||||
// LoadFromJSON restores personal averages from JSON.
|
||||
func (t *BreathingAnomalyTracker) LoadFromJSON(data []byte) error {
|
||||
t.mu.Lock()
|
||||
defer t.mu.Unlock()
|
||||
return json.Unmarshal(data, &t.personal)
|
||||
}
|
||||
|
||||
// SaveToJSON serializes personal averages to JSON.
|
||||
func (t *BreathingAnomalyTracker) SaveToJSON() ([]byte, error) {
|
||||
t.mu.RLock()
|
||||
defer t.mu.RUnlock()
|
||||
return json.Marshal(t.personal)
|
||||
}
|
||||
|
||||
// ComputeBreathingRegularity computes the coefficient of variation (CV = std/mean)
|
||||
// of a slice of breathing rate samples.
|
||||
func ComputeBreathingRegularity(samples []float64) float64 {
|
||||
if len(samples) == 0 {
|
||||
return 0
|
||||
}
|
||||
|
||||
var sum, sumSq float64
|
||||
for _, s := range samples {
|
||||
sum += s
|
||||
sumSq += s * s
|
||||
}
|
||||
|
||||
mean := sum / float64(len(samples))
|
||||
if mean == 0 {
|
||||
return 0
|
||||
}
|
||||
|
||||
variance := sumSq/float64(len(samples)) - mean*mean
|
||||
stdDev := math.Sqrt(math.Max(0, variance))
|
||||
|
||||
return stdDev / mean
|
||||
}
|
||||
|
||||
// BreathingRegularityLabel returns a human-readable label for the CV value.
|
||||
func BreathingRegularityLabel(cv float64) string {
|
||||
switch {
|
||||
case cv < BreathingRegularityRegular:
|
||||
return "regular"
|
||||
case cv > BreathingRegularityIrregular:
|
||||
return "irregular"
|
||||
default:
|
||||
return "normal"
|
||||
}
|
||||
}
|
||||
172
mothership/internal/sleep/breathing_anomaly_test.go
Normal file
172
mothership/internal/sleep/breathing_anomaly_test.go
Normal file
|
|
@ -0,0 +1,172 @@
|
|||
package sleep
|
||||
|
||||
import (
|
||||
"encoding/json"
|
||||
"testing"
|
||||
)
|
||||
|
||||
func TestBreathingAnomalyTrackerCheckAnomaly(t *testing.T) {
|
||||
tracker := NewBreathingAnomalyTracker()
|
||||
|
||||
// Establish a personal average of 16 bpm
|
||||
tracker.UpdatePersonalAverage("alice", 16.0)
|
||||
|
||||
tests := []struct {
|
||||
name string
|
||||
avgBPM float64
|
||||
personal float64
|
||||
wantAnomaly bool
|
||||
}{
|
||||
{"normal rate", 16.0, 16.0, false},
|
||||
{"slightly elevated", 19.0, 16.0, false}, // 19/16 = 1.1875 < 1.25
|
||||
{"at threshold", 20.0, 16.0, false}, // 20/16 = 1.25 = threshold (not >)
|
||||
{"above threshold", 21.0, 16.0, true}, // 21/16 = 1.3125 > 1.25
|
||||
{"significantly elevated", 25.0, 16.0, true}, // 25/16 = 1.5625
|
||||
{"below average", 12.0, 16.0, false}, // 12/16 = 0.75
|
||||
}
|
||||
|
||||
for _, tt := range tests {
|
||||
t.Run(tt.name, func(t *testing.T) {
|
||||
// Reset and set personal average
|
||||
tracker.personal["alice"] = tt.personal
|
||||
got := tracker.CheckAnomaly("alice", tt.avgBPM)
|
||||
if got != tt.wantAnomaly {
|
||||
t.Errorf("CheckAnomaly(%s, %.1f) = %v, want %v (personal=%.1f)",
|
||||
"alice", tt.avgBPM, got, tt.wantAnomaly, tt.personal)
|
||||
}
|
||||
})
|
||||
}
|
||||
}
|
||||
|
||||
func TestBreathingAnomalyTrackerNoBaseline(t *testing.T) {
|
||||
tracker := NewBreathingAnomalyTracker()
|
||||
|
||||
// No personal average set — should not flag anomaly
|
||||
if tracker.CheckAnomaly("bob", 25.0) {
|
||||
t.Error("CheckAnomaly should return false when no baseline exists")
|
||||
}
|
||||
|
||||
if tracker.GetPersonalAverage("bob") != 0 {
|
||||
t.Error("GetPersonalAverage should return 0 for unknown person")
|
||||
}
|
||||
}
|
||||
|
||||
func TestBreathingAnomalyTrackerEmaUpdate(t *testing.T) {
|
||||
tracker := NewBreathingAnomalyTracker()
|
||||
|
||||
// Night 1: 16 bpm
|
||||
tracker.UpdatePersonalAverage("alice", 16.0)
|
||||
avg1 := tracker.GetPersonalAverage("alice")
|
||||
|
||||
// Night 2: 18 bpm — EMA should pull average up slightly
|
||||
tracker.UpdatePersonalAverage("alice", 18.0)
|
||||
avg2 := tracker.GetPersonalAverage("alice")
|
||||
|
||||
if avg2 <= avg1 {
|
||||
t.Errorf("Personal average should increase: before=%.2f, after=%.2f", avg1, avg2)
|
||||
}
|
||||
|
||||
// EMA formula: avg = 0.05 * 18 + 0.95 * 16 = 0.9 + 15.2 = 16.1
|
||||
expected := 0.05*18.0 + 0.95*16.0
|
||||
if mathAbs(avg2-expected) > 0.001 {
|
||||
t.Errorf("EMA = %.4f, want %.4f", avg2, expected)
|
||||
}
|
||||
}
|
||||
|
||||
func TestBreathingAnomalyTrackerZeroBPMIgnored(t *testing.T) {
|
||||
tracker := NewBreathingAnomalyTracker()
|
||||
tracker.UpdatePersonalAverage("alice", 16.0)
|
||||
|
||||
// Zero BPM should not update the personal average
|
||||
tracker.UpdatePersonalAverage("alice", 0)
|
||||
if tracker.GetPersonalAverage("alice") != 16.0 {
|
||||
t.Errorf("Personal average changed after zero BPM update, got %.2f", tracker.GetPersonalAverage("alice"))
|
||||
}
|
||||
}
|
||||
|
||||
func TestBreathingAnomalyTrackerNegativeBPMIgnored(t *testing.T) {
|
||||
tracker := NewBreathingAnomalyTracker()
|
||||
tracker.UpdatePersonalAverage("alice", 16.0)
|
||||
|
||||
tracker.UpdatePersonalAverage("alice", -5.0)
|
||||
if tracker.GetPersonalAverage("alice") != 16.0 {
|
||||
t.Errorf("Personal average changed after negative BPM update")
|
||||
}
|
||||
}
|
||||
|
||||
func TestBreathingAnomalyTrackerJSONRoundTrip(t *testing.T) {
|
||||
tracker := NewBreathingAnomalyTracker()
|
||||
tracker.UpdatePersonalAverage("alice", 16.0)
|
||||
tracker.UpdatePersonalAverage("bob", 14.5)
|
||||
|
||||
data, err := tracker.SaveToJSON()
|
||||
if err != nil {
|
||||
t.Fatalf("SaveToJSON() error = %v", err)
|
||||
}
|
||||
|
||||
tracker2 := NewBreathingAnomalyTracker()
|
||||
if err := tracker2.LoadFromJSON(data); err != nil {
|
||||
t.Fatalf("LoadFromJSON() error = %v", err)
|
||||
}
|
||||
|
||||
if tracker2.GetPersonalAverage("alice") != tracker.GetPersonalAverage("alice") {
|
||||
t.Errorf("alice avg mismatch: %.2f vs %.2f",
|
||||
tracker2.GetPersonalAverage("alice"), tracker.GetPersonalAverage("alice"))
|
||||
}
|
||||
if tracker2.GetPersonalAverage("bob") != tracker.GetPersonalAverage("bob") {
|
||||
t.Errorf("bob avg mismatch: %.2f vs %.2f",
|
||||
tracker2.GetPersonalAverage("bob"), tracker.GetPersonalAverage("bob"))
|
||||
}
|
||||
}
|
||||
|
||||
func TestBreathingAnomalyTrackerJSONEmpty(t *testing.T) {
|
||||
tracker := NewBreathingAnomalyTracker()
|
||||
|
||||
data, err := tracker.SaveToJSON()
|
||||
if err != nil {
|
||||
t.Fatalf("SaveToJSON() error = %v", err)
|
||||
}
|
||||
|
||||
// Should be valid empty JSON object
|
||||
var m map[string]float64
|
||||
if err := json.Unmarshal(data, &m); err != nil {
|
||||
t.Fatalf("Invalid JSON: %v", err)
|
||||
}
|
||||
if len(m) != 0 {
|
||||
t.Errorf("Expected empty map, got %d entries", len(m))
|
||||
}
|
||||
}
|
||||
|
||||
func TestBreathingAnomalyTrackerMultiplePeople(t *testing.T) {
|
||||
tracker := NewBreathingAnomalyTracker()
|
||||
|
||||
tracker.UpdatePersonalAverage("alice", 16.0)
|
||||
tracker.UpdatePersonalAverage("bob", 13.0)
|
||||
|
||||
// Alice: 21 bpm is 1.3125x her average → anomaly
|
||||
if !tracker.CheckAnomaly("alice", 21.0) {
|
||||
t.Error("Alice 21 bpm should be anomaly (personal avg 16)")
|
||||
}
|
||||
|
||||
// Bob: 21 bpm is 1.615x his average → anomaly
|
||||
if !tracker.CheckAnomaly("bob", 21.0) {
|
||||
t.Error("Bob 21 bpm should be anomaly (personal avg 13)")
|
||||
}
|
||||
|
||||
// Alice: 18 bpm is 1.125x → not anomaly
|
||||
if tracker.CheckAnomaly("alice", 18.0) {
|
||||
t.Error("Alice 18 bpm should NOT be anomaly")
|
||||
}
|
||||
|
||||
// Bob: 14 bpm is 1.077x → not anomaly
|
||||
if tracker.CheckAnomaly("bob", 14.0) {
|
||||
t.Error("Bob 14 bpm should NOT be anomaly")
|
||||
}
|
||||
}
|
||||
|
||||
func mathAbs(x float64) float64 {
|
||||
if x < 0 {
|
||||
return -x
|
||||
}
|
||||
return x
|
||||
}
|
||||
160
mothership/internal/sleep/breathing_estimator.go
Normal file
160
mothership/internal/sleep/breathing_estimator.go
Normal file
|
|
@ -0,0 +1,160 @@
|
|||
// Package sleep provides FFT-based breathing rate estimation for sleep monitoring.
|
||||
package sleep
|
||||
|
||||
import (
|
||||
"math"
|
||||
"math/cmplx"
|
||||
"sync"
|
||||
|
||||
"gonum.org/v1/gonum/dsp/fourier"
|
||||
)
|
||||
|
||||
// FFT estimator constants
|
||||
const (
|
||||
FFTEstimatorSampleRate = 20.0 // Hz — matches signal pipeline
|
||||
FFTEstimatorFFTSize = 512 // Input window (25.6 s at 20 Hz)
|
||||
FFTEstimatorZeroPad = 1024 // Zero-padded FFT size
|
||||
FFTEstimatorEMAlpha = 1.0 / 60.0 // 60-second EMA smoothing
|
||||
FFTEstimatorMinHz = 0.1 // 6 BPM lower bound
|
||||
FFTEstimatorMaxHz = 0.5 // 30 BPM upper bound
|
||||
FFTEstimatorMinBPM = 6.0
|
||||
FFTEstimatorMaxBPM = 30.0
|
||||
)
|
||||
|
||||
// BreathingRateEstimator accumulates phase samples and estimates breathing rate via FFT.
|
||||
// It operates on bandpass-filtered residual phase from the most motion-sensitive link
|
||||
// in a sleep zone, producing one BPM estimate per 25.6-second window with 60-second EMA smoothing.
|
||||
type BreathingRateEstimator struct {
|
||||
mu sync.RWMutex
|
||||
buffer []float64 // Circular buffer of FFTEstimatorFFTSize samples
|
||||
writeIdx int
|
||||
sampleCount int
|
||||
emaRate float64 // EMA-smoothed BPM
|
||||
lastRate float64 // Most recent raw FFT BPM
|
||||
}
|
||||
|
||||
// NewBreathingRateEstimator creates a new FFT-based breathing rate estimator.
|
||||
func NewBreathingRateEstimator() *BreathingRateEstimator {
|
||||
return &BreathingRateEstimator{
|
||||
buffer: make([]float64, FFTEstimatorFFTSize),
|
||||
}
|
||||
}
|
||||
|
||||
// AddPhaseSample adds a bandpass-filtered phase sample to the circular buffer.
|
||||
func (e *BreathingRateEstimator) AddPhaseSample(phase float64) {
|
||||
e.mu.Lock()
|
||||
defer e.mu.Unlock()
|
||||
|
||||
e.buffer[e.writeIdx] = phase
|
||||
e.writeIdx = (e.writeIdx + 1) % FFTEstimatorFFTSize
|
||||
if e.sampleCount < FFTEstimatorFFTSize {
|
||||
e.sampleCount++
|
||||
}
|
||||
}
|
||||
|
||||
// EstimateRate runs FFT on accumulated samples and returns EMA-smoothed BPM.
|
||||
// Returns 0 if insufficient samples have been collected.
|
||||
func (e *BreathingRateEstimator) EstimateRate() float64 {
|
||||
e.mu.Lock()
|
||||
defer e.mu.Unlock()
|
||||
|
||||
if e.sampleCount < FFTEstimatorFFTSize {
|
||||
return e.emaRate
|
||||
}
|
||||
|
||||
bpm := computeFFTBreathingRate(e.buffer, e.writeIdx, FFTEstimatorSampleRate, FFTEstimatorZeroPad)
|
||||
|
||||
// Reject out-of-physiological-range estimates
|
||||
if bpm < FFTEstimatorMinBPM || bpm > FFTEstimatorMaxBPM {
|
||||
return e.emaRate
|
||||
}
|
||||
|
||||
// Apply 60-second EMA smoothing: ema = α × bpm + (1-α) × ema
|
||||
if e.emaRate > 0 {
|
||||
bpm = FFTEstimatorEMAlpha*bpm + (1-FFTEstimatorEMAlpha)*e.emaRate
|
||||
}
|
||||
|
||||
e.emaRate = bpm
|
||||
e.lastRate = bpm
|
||||
return bpm
|
||||
}
|
||||
|
||||
// GetRate returns the current EMA-smoothed breathing rate.
|
||||
func (e *BreathingRateEstimator) GetRate() float64 {
|
||||
e.mu.RLock()
|
||||
defer e.mu.RUnlock()
|
||||
return e.emaRate
|
||||
}
|
||||
|
||||
// Reset clears the estimator state.
|
||||
func (e *BreathingRateEstimator) Reset() {
|
||||
e.mu.Lock()
|
||||
defer e.mu.Unlock()
|
||||
|
||||
for i := range e.buffer {
|
||||
e.buffer[i] = 0
|
||||
}
|
||||
e.writeIdx = 0
|
||||
e.sampleCount = 0
|
||||
e.emaRate = 0
|
||||
e.lastRate = 0
|
||||
}
|
||||
|
||||
// Ready returns true if the buffer has enough samples for an FFT window.
|
||||
func (e *BreathingRateEstimator) Ready() bool {
|
||||
e.mu.RLock()
|
||||
defer e.mu.RUnlock()
|
||||
return e.sampleCount >= FFTEstimatorFFTSize
|
||||
}
|
||||
|
||||
// computeFFTBreathingRate runs FFT on phase samples and returns the dominant
|
||||
// breathing frequency converted to BPM.
|
||||
//
|
||||
// Parameters:
|
||||
// - buffer: circular buffer of phase samples (length n)
|
||||
// - writeIdx: position where the next sample will be written
|
||||
// - sampleRate: sampling rate in Hz
|
||||
// - zeroPadSize: FFT size (must be power of 2, >= len(buffer))
|
||||
func computeFFTBreathingRate(buffer []float64, writeIdx int, sampleRate, zeroPadSize float64) float64 {
|
||||
n := len(buffer)
|
||||
N := int(zeroPadSize)
|
||||
|
||||
// Build zero-padded real input from circular buffer (chronological order)
|
||||
seq := make([]float64, N)
|
||||
for i := 0; i < n; i++ {
|
||||
idx := (writeIdx + i) % n
|
||||
seq[i] = buffer[idx]
|
||||
}
|
||||
|
||||
// Run FFT
|
||||
fft := fourier.NewFFT(N)
|
||||
coeff := fft.Coefficients(nil, seq)
|
||||
|
||||
// Frequency resolution: Fs / N
|
||||
freqRes := sampleRate / zeroPadSize
|
||||
|
||||
// Bin range for 0.1–0.5 Hz
|
||||
minBin := int(math.Ceil(FFTEstimatorMinHz / freqRes))
|
||||
maxBin := int(math.Floor(FFTEstimatorMaxHz / freqRes))
|
||||
if minBin < 1 {
|
||||
minBin = 1 // skip DC
|
||||
}
|
||||
if maxBin > N/2 {
|
||||
maxBin = N / 2 // Nyquist
|
||||
}
|
||||
|
||||
// Find dominant magnitude peak in breathing band
|
||||
maxMag := 0.0
|
||||
peakBin := minBin
|
||||
for bin := minBin; bin <= maxBin; bin++ {
|
||||
mag := cmplx.Abs(coeff[bin])
|
||||
if mag > maxMag {
|
||||
maxMag = mag
|
||||
peakBin = bin
|
||||
}
|
||||
}
|
||||
|
||||
// Convert bin index to BPM: bpm = bin_idx × (Fs/N) × 60
|
||||
freqHz := float64(peakBin) * freqRes
|
||||
return freqHz * 60.0
|
||||
}
|
||||
252
mothership/internal/sleep/breathing_estimator_test.go
Normal file
252
mothership/internal/sleep/breathing_estimator_test.go
Normal file
|
|
@ -0,0 +1,252 @@
|
|||
package sleep
|
||||
|
||||
import (
|
||||
"math"
|
||||
"testing"
|
||||
)
|
||||
|
||||
func TestComputeFFTBreathingRateSynthetic15BPM(t *testing.T) {
|
||||
// Generate a synthetic phase signal at 0.25 Hz (15 bpm) at 20 Hz sample rate.
|
||||
// The FFT should identify the dominant frequency as 15 bpm.
|
||||
sampleRate := 20.0
|
||||
nSamples := 512
|
||||
zeroPad := 1024
|
||||
|
||||
freqHz := 0.25 // 15 bpm
|
||||
buffer := make([]float64, nSamples)
|
||||
for i := 0; i < nSamples; i++ {
|
||||
t_sec := float64(i) / sampleRate
|
||||
buffer[i] = math.Sin(2.0 * math.Pi * freqHz * t_sec)
|
||||
}
|
||||
|
||||
bpm := computeFFTBreathingRate(buffer, 0, sampleRate, float64(zeroPad))
|
||||
|
||||
// 15 bpm ± 1 bpm tolerance (frequency resolution is ~1.17 bpm/bin)
|
||||
if math.Abs(bpm-15.0) > 1.5 {
|
||||
t.Errorf("FFT breathing rate = %.2f bpm, want ~15.0 bpm", bpm)
|
||||
}
|
||||
}
|
||||
|
||||
func TestComputeFFTBreathingRateSynthetic12BPM(t *testing.T) {
|
||||
// 0.2 Hz = 12 bpm
|
||||
sampleRate := 20.0
|
||||
buffer := make([]float64, 512)
|
||||
for i := range buffer {
|
||||
buffer[i] = math.Sin(2.0 * math.Pi * 0.2 * float64(i) / sampleRate)
|
||||
}
|
||||
|
||||
bpm := computeFFTBreathingRate(buffer, 0, sampleRate, 1024)
|
||||
|
||||
if math.Abs(bpm-12.0) > 1.5 {
|
||||
t.Errorf("FFT breathing rate = %.2f bpm, want ~12.0 bpm", bpm)
|
||||
}
|
||||
}
|
||||
|
||||
func TestComputeFFTBreathingRateSynthetic20BPM(t *testing.T) {
|
||||
// 0.333 Hz ≈ 20 bpm
|
||||
sampleRate := 20.0
|
||||
buffer := make([]float64, 512)
|
||||
for i := range buffer {
|
||||
buffer[i] = math.Sin(2.0 * math.Pi * (1.0/3.0) * float64(i) / sampleRate)
|
||||
}
|
||||
|
||||
bpm := computeFFTBreathingRate(buffer, 0, sampleRate, 1024)
|
||||
|
||||
if math.Abs(bpm-20.0) > 2.0 {
|
||||
t.Errorf("FFT breathing rate = %.2f bpm, want ~20.0 bpm", bpm)
|
||||
}
|
||||
}
|
||||
|
||||
func TestComputeFFTBreathingRateWithNoise(t *testing.T) {
|
||||
// 15 bpm signal + Gaussian-like noise. FFT should still find the dominant peak.
|
||||
sampleRate := 20.0
|
||||
buffer := make([]float64, 512)
|
||||
for i := range buffer {
|
||||
t_sec := float64(i) / sampleRate
|
||||
signal := math.Sin(2.0 * math.Pi * 0.25 * t_sec)
|
||||
// Simple pseudo-noise: sum of incommensurate frequencies
|
||||
noise := 0.3*math.Sin(2*math.Pi*3.7*t_sec) +
|
||||
0.2*math.Sin(2*math.Pi*7.1*t_sec) +
|
||||
0.15*math.Sin(2*math.Pi*0.03*t_sec)
|
||||
buffer[i] = signal + noise
|
||||
}
|
||||
|
||||
bpm := computeFFTBreathingRate(buffer, 0, sampleRate, 1024)
|
||||
|
||||
if math.Abs(bpm-15.0) > 2.0 {
|
||||
t.Errorf("FFT breathing rate with noise = %.2f bpm, want ~15.0 bpm", bpm)
|
||||
}
|
||||
}
|
||||
|
||||
func TestComputeFFTBreathingRateCircularBuffer(t *testing.T) {
|
||||
// Verify that the circular buffer read order is correct when writeIdx != 0.
|
||||
sampleRate := 20.0
|
||||
nSamples := 512
|
||||
zeroPad := 1024
|
||||
|
||||
freqHz := 0.25 // 15 bpm
|
||||
buffer := make([]float64, nSamples)
|
||||
|
||||
// Simulate a filled circular buffer where writeIdx is at position 100
|
||||
writeIdx := 100
|
||||
for i := 0; i < nSamples; i++ {
|
||||
// Sample at logical position i was written to (writeIdx + i) % nSamples
|
||||
t_sec := float64(i) / sampleRate
|
||||
buffer[(writeIdx+i)%nSamples] = math.Sin(2.0 * math.Pi * freqHz * t_sec)
|
||||
}
|
||||
|
||||
bpm := computeFFTBreathingRate(buffer, writeIdx, sampleRate, float64(zeroPad))
|
||||
|
||||
if math.Abs(bpm-15.0) > 1.5 {
|
||||
t.Errorf("FFT circular buffer breathing rate = %.2f bpm, want ~15.0 bpm", bpm)
|
||||
}
|
||||
}
|
||||
|
||||
func TestBreathingRateEstimatorEmaSmoothing(t *testing.T) {
|
||||
// Feed a constant 15 bpm signal and verify EMA converges.
|
||||
est := NewBreathingRateEstimator()
|
||||
sampleRate := 20.0
|
||||
|
||||
// Fill the buffer with a 15 bpm signal
|
||||
for i := 0; i < 512; i++ {
|
||||
t_sec := float64(i) / sampleRate
|
||||
phase := math.Sin(2.0 * math.Pi * 0.25 * t_sec)
|
||||
est.AddPhaseSample(phase)
|
||||
}
|
||||
|
||||
// First estimate should be close to 15
|
||||
rate := est.EstimateRate()
|
||||
if math.Abs(rate-15.0) > 2.0 {
|
||||
t.Errorf("First estimate = %.2f bpm, want ~15.0 bpm", rate)
|
||||
}
|
||||
|
||||
// Continue feeding and estimating — EMA should stabilize near 15
|
||||
for rep := 0; rep < 10; rep++ {
|
||||
for i := 0; i < 512; i++ {
|
||||
t_sec := float64(i) / sampleRate
|
||||
phase := math.Sin(2.0 * math.Pi * 0.25 * t_sec)
|
||||
est.AddPhaseSample(phase)
|
||||
}
|
||||
est.EstimateRate()
|
||||
}
|
||||
|
||||
finalRate := est.GetRate()
|
||||
if math.Abs(finalRate-15.0) > 1.0 {
|
||||
t.Errorf("EMA-stabilized rate = %.2f bpm, want ~15.0 bpm", finalRate)
|
||||
}
|
||||
}
|
||||
|
||||
func TestBreathingRateEstimatorInsufficientSamples(t *testing.T) {
|
||||
est := NewBreathingRateEstimator()
|
||||
|
||||
// Add fewer than 512 samples
|
||||
for i := 0; i < 100; i++ {
|
||||
est.AddPhaseSample(0.1)
|
||||
}
|
||||
|
||||
rate := est.EstimateRate()
|
||||
if rate != 0 {
|
||||
t.Errorf("Rate with insufficient samples = %.2f, want 0", rate)
|
||||
}
|
||||
|
||||
if est.Ready() {
|
||||
t.Error("Ready() should be false with insufficient samples")
|
||||
}
|
||||
}
|
||||
|
||||
func TestBreathingRateEstimatorReset(t *testing.T) {
|
||||
est := NewBreathingRateEstimator()
|
||||
|
||||
for i := 0; i < 512; i++ {
|
||||
est.AddPhaseSample(0.1)
|
||||
}
|
||||
est.EstimateRate() // Populate EMA
|
||||
|
||||
est.Reset()
|
||||
|
||||
if est.GetRate() != 0 {
|
||||
t.Errorf("Rate after reset = %.2f, want 0", est.GetRate())
|
||||
}
|
||||
if est.Ready() {
|
||||
t.Error("Ready() should be false after reset")
|
||||
}
|
||||
}
|
||||
|
||||
func TestComputeBreathingRegularity(t *testing.T) {
|
||||
tests := []struct {
|
||||
name string
|
||||
samples []float64
|
||||
wantCV float64
|
||||
tol float64
|
||||
}{
|
||||
{
|
||||
name: "constant rate — zero CV",
|
||||
samples: []float64{14.0, 14.0, 14.0, 14.0, 14.0},
|
||||
wantCV: 0.0,
|
||||
tol: 0.001,
|
||||
},
|
||||
{
|
||||
name: "small variation",
|
||||
samples: []float64{14.0, 14.5, 13.5, 14.2, 13.8},
|
||||
wantCV: 0.024,
|
||||
tol: 0.01,
|
||||
},
|
||||
{
|
||||
name: "large variation",
|
||||
samples: []float64{10.0, 20.0, 12.0, 18.0, 15.0},
|
||||
wantCV: 0.276,
|
||||
tol: 0.05,
|
||||
},
|
||||
{
|
||||
name: "empty samples — zero CV",
|
||||
samples: []float64{},
|
||||
wantCV: 0.0,
|
||||
tol: 0.0,
|
||||
},
|
||||
{
|
||||
name: "single sample — zero CV",
|
||||
samples: []float64{14.0},
|
||||
wantCV: 0.0,
|
||||
tol: 0.0,
|
||||
},
|
||||
{
|
||||
name: "zero mean — zero CV (no division by zero)",
|
||||
samples: []float64{0.0, 0.0, 0.0},
|
||||
wantCV: 0.0,
|
||||
tol: 0.0,
|
||||
},
|
||||
}
|
||||
|
||||
for _, tt := range tests {
|
||||
t.Run(tt.name, func(t *testing.T) {
|
||||
cv := ComputeBreathingRegularity(tt.samples)
|
||||
if math.Abs(cv-tt.wantCV) > tt.tol {
|
||||
t.Errorf("ComputeBreathingRegularity() = %.4f, want %.4f ± %.4f", cv, tt.wantCV, tt.tol)
|
||||
}
|
||||
})
|
||||
}
|
||||
}
|
||||
|
||||
func TestBreathingRegularityLabel(t *testing.T) {
|
||||
tests := []struct {
|
||||
cv float64
|
||||
want string
|
||||
}{
|
||||
{0.05, "regular"},
|
||||
{0.09, "regular"},
|
||||
{0.10, "normal"}, // boundary
|
||||
{0.15, "normal"},
|
||||
{0.25, "normal"}, // boundary
|
||||
{0.26, "irregular"},
|
||||
{0.50, "irregular"},
|
||||
}
|
||||
|
||||
for _, tt := range tests {
|
||||
t.Run(tt.want, func(t *testing.T) {
|
||||
got := BreathingRegularityLabel(tt.cv)
|
||||
if got != tt.want {
|
||||
t.Errorf("BreathingRegularityLabel(%.2f) = %q, want %q", tt.cv, got, tt.want)
|
||||
}
|
||||
})
|
||||
}
|
||||
}
|
||||
|
|
@ -208,12 +208,15 @@ func (h *Handler) handleGetSession(w http.ResponseWriter, r *http.Request) {
|
|||
|
||||
// Include live metrics if available
|
||||
if metrics != nil {
|
||||
result["metrics"] = map[string]interface{}{
|
||||
metricsMap := map[string]interface{}{
|
||||
"total_duration_hours": metrics.TotalDuration.Hours(),
|
||||
"time_in_bed_hours": metrics.TimeInBed.Hours(),
|
||||
"avg_breathing_rate": metrics.AvgBreathingRate,
|
||||
"breathing_rate_std_dev": metrics.BreathingRateStdDev,
|
||||
"breathing_regularity": metrics.BreathingRegularity,
|
||||
"breathing_score": metrics.BreathingScore,
|
||||
"breathing_anomaly": metrics.BreathingAnomaly,
|
||||
"breathing_anomaly_count": metrics.BreathingAnomalyCount,
|
||||
"quiet_time_pct": metrics.QuietTimePct,
|
||||
"motion_events": metrics.MotionEvents,
|
||||
"restless_periods": metrics.RestlessPeriods,
|
||||
|
|
@ -225,14 +228,18 @@ func (h *Handler) handleGetSession(w http.ResponseWriter, r *http.Request) {
|
|||
"quality_rating": metrics.QualityRating,
|
||||
}
|
||||
|
||||
if metricsMap, ok := result["metrics"].(map[string]interface{}); ok {
|
||||
if !metrics.SleepStartTime.IsZero() {
|
||||
metricsMap["sleep_start_time"] = metrics.SleepStartTime.Format("15:04")
|
||||
}
|
||||
if !metrics.SleepEndTime.IsZero() {
|
||||
metricsMap["sleep_end_time"] = metrics.SleepEndTime.Format("15:04")
|
||||
}
|
||||
if metrics.PersonalAvgBPM > 0 {
|
||||
metricsMap["personal_avg_bpm"] = metrics.PersonalAvgBPM
|
||||
}
|
||||
|
||||
if !metrics.SleepStartTime.IsZero() {
|
||||
metricsMap["sleep_start_time"] = metrics.SleepStartTime.Format("15:04")
|
||||
}
|
||||
if !metrics.SleepEndTime.IsZero() {
|
||||
metricsMap["sleep_end_time"] = metrics.SleepEndTime.Format("15:04")
|
||||
}
|
||||
|
||||
result["metrics"] = metricsMap
|
||||
}
|
||||
|
||||
writeJSON(w, result)
|
||||
|
|
|
|||
|
|
@ -47,6 +47,21 @@ func generateBreathingSummary(m *SleepMetrics) string {
|
|||
summary += "Your breathing was steady throughout the night."
|
||||
}
|
||||
|
||||
// Regularity assessment
|
||||
if m.BreathingRegularity > 0 {
|
||||
summary += fmt.Sprintf(" Regularity: %s (CV=%.2f).", BreathingRegularityLabel(m.BreathingRegularity), m.BreathingRegularity)
|
||||
}
|
||||
|
||||
// Anomaly assessment
|
||||
if m.BreathingAnomaly {
|
||||
if m.PersonalAvgBPM > 0 {
|
||||
summary += fmt.Sprintf(" Breathing rate elevated (%.0f bpm vs. %.0f bpm average).",
|
||||
m.AvgBreathingRate, m.PersonalAvgBPM)
|
||||
} else {
|
||||
summary += " Breathing rate was elevated compared to your personal average."
|
||||
}
|
||||
}
|
||||
|
||||
// Range info
|
||||
if m.MaxBreathingRate > 0 {
|
||||
summary += fmt.Sprintf(" Range: %.1f-%.1f BPM.", m.MinBreathingRate, m.MaxBreathingRate)
|
||||
|
|
@ -145,7 +160,8 @@ func FormatDuration(d time.Duration) string {
|
|||
return fmt.Sprintf("%d seconds", d/time.Second)
|
||||
}
|
||||
if d < time.Hour {
|
||||
return fmt.Sprintf("%d minutes", d/time.Minute)
|
||||
mins := d / time.Minute
|
||||
return fmt.Sprintf("%d minute%s", mins, pluralS(int(mins)))
|
||||
}
|
||||
|
||||
hours := d / time.Hour
|
||||
|
|
@ -164,6 +180,14 @@ func pluralS(n int) string {
|
|||
return "s"
|
||||
}
|
||||
|
||||
// formatMinutes handles singular/plural for minutes.
|
||||
func formatMinutes(n int) string {
|
||||
if n == 1 {
|
||||
return "1 minute"
|
||||
}
|
||||
return fmt.Sprintf("%d minutes", n)
|
||||
}
|
||||
|
||||
// ToJSONMap converts the report to a map for JSON serialization
|
||||
func (r *SleepReport) ToJSONMap() map[string]interface{} {
|
||||
m := map[string]interface{}{
|
||||
|
|
@ -178,12 +202,15 @@ func (r *SleepReport) ToJSONMap() map[string]interface{} {
|
|||
}
|
||||
|
||||
// Add detailed metrics
|
||||
m["metrics"] = map[string]interface{}{
|
||||
metricsMap := map[string]interface{}{
|
||||
"total_duration_hours": r.Metrics.TotalDuration.Hours(),
|
||||
"time_in_bed_hours": r.Metrics.TimeInBed.Hours(),
|
||||
"avg_breathing_rate": r.Metrics.AvgBreathingRate,
|
||||
"breathing_rate_std_dev": r.Metrics.BreathingRateStdDev,
|
||||
"breathing_regularity": r.Metrics.BreathingRegularity,
|
||||
"breathing_score": r.Metrics.BreathingScore,
|
||||
"breathing_anomaly": r.Metrics.BreathingAnomaly,
|
||||
"breathing_anomaly_count": r.Metrics.BreathingAnomalyCount,
|
||||
"quiet_time_pct": r.Metrics.QuietTimePct,
|
||||
"motion_events": r.Metrics.MotionEvents,
|
||||
"restless_periods": r.Metrics.RestlessPeriods,
|
||||
|
|
@ -193,6 +220,19 @@ func (r *SleepReport) ToJSONMap() map[string]interface{} {
|
|||
"continuity_score": r.Metrics.ContinuityScore,
|
||||
}
|
||||
|
||||
// Add breathing rate range
|
||||
if r.Metrics.MinBreathingRate > 0 {
|
||||
metricsMap["min_breathing_rate"] = r.Metrics.MinBreathingRate
|
||||
metricsMap["max_breathing_rate"] = r.Metrics.MaxBreathingRate
|
||||
}
|
||||
|
||||
// Add personal baseline comparison for anomaly
|
||||
if r.Metrics.PersonalAvgBPM > 0 {
|
||||
metricsMap["personal_avg_bpm"] = r.Metrics.PersonalAvgBPM
|
||||
}
|
||||
|
||||
m["metrics"] = metricsMap
|
||||
|
||||
// Add timing
|
||||
if !r.Metrics.SleepStartTime.IsZero() {
|
||||
m["sleep_start_time"] = r.Metrics.SleepStartTime.Format("15:04")
|
||||
|
|
|
|||
Loading…
Add table
Reference in a new issue