1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132 133 134 135 136 137 138 139 140 141 142 143 144 145 146 147 148 149 150 151 152 153 154 155 156 157 158 159 160 161 162 163 164 165 166 167 168 169 170 171 172 173 174 175 176 177 178 179 180 181 182 183 184 185 186 187 188 189 190 191 192 193 194 195 196 197 198 199 200 201 202 203 204 205 206 207 208 209
| """ EDA + HRV 驾驶员压力监测算法
核心方法: 1. EDA 检测急性压力事件 2. HRV 检测慢性压力水平 3. 多模态融合判定压力等级 """
import numpy as np from typing import Tuple, Dict
def extract_eda_features(eda_signal: np.ndarray, fps: int = 4) -> Dict[str, float]: """ 提取 EDA 特征 Args: eda_signal: EDA 信号序列, shape=(N,), 单位:μS(微西门子) fps: 采样率(Hz) Returns: features: { 'scl': 皮肤电导水平(基线), 'scr_count': SCR 事件次数, 'scr_amplitude': SCR 幅度, 'scr_frequency': SCR 频率(次/分钟) } 驾驶员压力应用: 急性压力事件: - SCR 频率 >10 次/分钟 → 急性压力高 - SCR 幅度 >0.5μS → 强烈压力反应 """ from scipy.signal import butter, filtfilt nyquist = fps / 2 low_cutoff = 0.05 / nyquist b, a = butter(2, low_cutoff, btype='low') scl = filtfilt(b, a, eda_signal) high_cutoff = 0.05 / nyquist b, a = butter(2, high_cutoff, btype='high') scr = filtfilt(b, a, eda_signal) threshold = 0.05 scr_peaks = detect_scr_peaks(scr, threshold) scr_count = len(scr_peaks) scr_amplitude = np.max(scr) if scr_count > 0 else 0 duration_minutes = len(eda_signal) / fps / 60 scr_frequency = scr_count / duration_minutes if duration_minutes > 0 else 0 return { 'scl': np.mean(scl), 'scr_count': scr_count, 'scr_amplitude': scr_amplitude, 'scr_frequency': scr_frequency }
def detect_scr_peaks(scr_signal: np.ndarray, threshold: float) -> list: """ 检测 SCR 峰值 Args: scr_signal: SCR 信号 threshold: 峰值阈值 Returns: peaks: 峰值索引列表 """ peaks = [] for i in range(1, len(scr_signal) - 1): if scr_signal[i] > scr_signal[i-1] and scr_signal[i] > scr_signal[i+1]: if scr_signal[i] > threshold: peaks.append(i) return peaks
def calculate_hrv_stress(rri_sequence: np.ndarray) -> float: """ 计算 HRV 压力指数 Args: rri_sequence: R-R 间期序列, 单位:ms Returns: stress_index: 压力指数(0-100) 驾驶员压力应用: HRV 降低 → 压力增加 RMSSD <20ms → 高压力 """ diff_rri = np.diff(rri_sequence) rmssd = np.sqrt(np.mean(diff_rri ** 2)) stress_index = max(0, 100 - rmssd * 2) return stress_index
def fusion_stress_assessment( eda_features: Dict[str, float], hrv_stress: float, weights: Dict[str, float] = None ) -> Tuple[float, str]: """ EDA + HRV 融合压力评估 Args: eda_features: EDA 特征 hrv_stress: HRV 压力指数 weights: 权重配置 Returns: stress_score: 综合压力分数(0-100) stress_level: 压力等级(low/medium/high/extreme) Garmin CIRQA 方法: 权重配置: - EDA SCR 频率:40%(急性压力) - HRV 压力指数:40%(慢性压力) - EDA SCL:20%(基线压力) """ if weights is None: weights = { 'scr_frequency': 0.4, 'hrv': 0.4, 'scl': 0.2 } scr_frequency_normalized = min(eda_features['scr_frequency'] / 20.0 * 100, 100) scl_normalized = min((eda_features['scl'] - 5) / 15.0 * 100, 100) stress_score = ( scr_frequency_normalized * weights['scr_frequency'] + hrv_stress * weights['hrv'] + scl_normalized * weights['scl'] ) if stress_score >= 80: stress_level = "extreme" elif stress_score >= 60: stress_level = "high" elif stress_score >= 40: stress_level = "medium" else: stress_level = "low" return stress_score, stress_level
if __name__ == "__main__": fps_eda = 4 duration = 60 n_samples_eda = fps_eda * duration t_eda = np.linspace(0, duration, n_samples_eda) scl_baseline = 10.0 scr_events = np.zeros(n_samples_eda) for i in range(15): event_time = i * 4 event_idx = int(event_time * fps_eda) if event_idx < n_samples_eda: scr_events[event_idx:event_idx+20] = np.sin(np.linspace(0, np.pi, 20)) * 0.5 eda_signal = scl_baseline + scr_events + np.random.randn(n_samples_eda) * 0.02 eda_features = extract_eda_features(eda_signal, fps_eda) rri_sequence = np.random.normal(700, 10, 100) hrv_stress = calculate_hrv_stress(rri_sequence) stress_score, stress_level = fusion_stress_assessment(eda_features, hrv_stress) print("="*60) print("驾驶员压力监测测试(EDA + HRV 融合)") print("="*60) print(f"EDA SCL(基线): {eda_features['scl']:.2f} μS") print(f"EDA SCR 频率: {eda_features['scr_frequency']:.1f} 次/分钟") print(f"HRV 压力指数: {hrv_stress:.1f}") print(f"\n综合压力分数: {stress_score:.1f}") print(f"压力等级: {stress_level}")
|