PPG-ECG多模态信号预处理与特征提取实战

简介面向脉搏波与心电信号分析需求资源提供可直接运行的预处理与特征提取项目。代码基于nk库完成信号去噪并提取多项关键特征潮波幅值比、重搏波幅值比、收缩面积比、舒张面积比、主波高度、波形周期面积、收缩期面积、舒张期面积以及脉搏传导时间。其中脉搏传导时间由心电信号与脉搏波之间的渡越时间计算可用于评估动脉弹性与血管功能。整套资源共四十三个文件含二十九个文本格式信号数据、八个配置信息、四个备份文件与两个程序脚本压缩包仅四点三一兆下载后即可复现。目前已有七百四十四人浏览学习适合生物医学工程、信号处理方向的学生与科研人员快速入门。数据与代码配套提供便于逐项核对特征计算过程也可作为二次开发与算法优化的起点。1. 先对齐、再清洗、后提特征PPG-ECG预处理到底在解决什么问题同时拿到PPG和ECG两条信号序列最常碰到的不是“特征不够”而是“信号根本没法用”。PPG是光学传感器测的容积脉搏波ECG是电极测的心电活动两者在硬件上就是两套采样系统采样率不一样、时间戳有各自时钟、噪声来源完全不同。你做血压估算、疲劳监测或者心率变异性分析起点都是同一个——把这两路不同形态的信号拉到同一时间基准上洗干净噪声再把能代表心血管状态的指标抽出来。这个过程就是标题里的“预处理及特征提取”。它解决的问题不是让波形好看而是让后续建模用的每一个特征都来自正确的时刻、正确的幅度、稳定的形态。这篇文章面向正在用可穿戴设备数据做算法的工程师和学生我会把时间对齐、滤波、特征提取三条线拆开每个环节给出可直接复现的代码和参数再把你大概率要踩的坑提前指出来。2. 信号对齐多模态特征提取与时序对齐为什么不能直接拼接2.1 两路信号天然不同步采样率、时钟、延迟三处不一致常见可穿戴设备里PPG采样率往往是25Hz、50Hz或125Hz而ECG通道常见是125Hz、250Hz甚至更高。更麻烦的是两路信号来自两个独立ADC片上时钟各自走开机时间和中断延迟不一样导致时间戳不仅存在固定偏移还可能有缓慢漂移。如果把两路原始数据直接逐点拼在一起后续算出来的脉搏波传导时间就是错的。做过血压估算的同学都知道PWTT同一心跳的ECG R波到PPG波峰的时间差通常只有几十到一百多毫秒一次20毫秒的时间戳错位就足以让模型完全失去意义。所以预处理的第一步不是滤波而是把两路信号放到同一时间轴上。这个动作做得好不好直接决定后面特征提取是“在真实时间上计算”还是“在幻觉时间上计算”。2.2 统一采样率用resample_poly做抗混叠重采样我推荐直接用scipy.signal.resample_poly而不是用numpy.interp线性插值。原因很简单线性插值会把PPG的重搏波切迹磨平而resample_poly内部走多相滤波器带抗混叠设计波形形态保留得更好。假设PPG原始采样率50HzECG原始采样率250Hz打算统一到125Hzimport numpy as np from scipy.signal import resample_poly # 原始信号一维数组时间连续单位是样本点 # PPG 从 50Hz 重采样到 125Hzup5, down2 ppg_rs resample_poly(ppg_raw, up5, down2, axis0) # ECG 从 250Hz 重采样到 125Hzup1, down2 ecg_rs resample_poly(ecg_raw, up1, down2, axis0)这里up和down是两个整数重采样比是up/down需要先约分。比如50Hz到125Hz的比例是125/50 5/2所以up5down2。要是原始采样率比目标采样率高比如ECG从250Hz降到125Hz就是up1down2。注意resample_poly默认会做抗混叠滤波因此重采样后的信号带宽受限于新的奈奎斯特频率这是好事但代价是重采样后的PPG如果目标采样率低于80Hz切迹信息会丢失这一点后面避坑章还会说。2.3 互相关校准固定延迟时间戳对齐的最后一步统一采样率之后两路信号还可能有固定的物理延迟比如PPG传感器靠光反射存在几十毫秒的光路和器件响应延迟这个延迟不校掉PWTT会被系统性抬高。常见做法是用互相关找延迟量from scipy.signal import correlate fs 125.0 # 统一后的采样率 seg_len int(10 * fs) # 取10秒干净段做校准 a ppg_rs[:seg_len] b ecg_rs[:seg_len] # 全长互相关找到使两者最对齐的滞后样本数 corr correlate(b, a, modefull) lags np.arange(-len(a) 1, len(b)) delay_samples lags[np.argmax(np.abs(corr))] # 延迟为正表示PPG落后ECG把PPG往前挪 if delay_samples 0: ppg_aligned np.concatenate([ppg_rs[delay_samples:], np.zeros(delay_samples)]) else: ppg_aligned np.concatenate([np.zeros(-delay_samples), ppg_rs[:delay_samples]])这里的逻辑是把ECG当参考时钟计算PPG相对它的滞后样本数再在时间上平移补齐。参数上要注意做互相关一定要选一段波形干净、R波和PPG波峰都明显的片段否则互相关峰值会被运动伪影带偏。另一个细节在校准完固定延迟后如果信号很长还应该分段二次校验因为部分设备的时钟漂移不是常数这时要分段对齐。多模态特征提取与时序对齐的关键就在这里——对不齐时间轴后面所有联合特征都是空中楼阁。对齐做得足够好ECG的每个R波才能准确对应到PPG的同一个心跳周期上时域和频域特征才谈得上可信。3. 去噪与滤波把数据预处理里的数据清洗落在信号层面3.1 噪声拆解运动伪影、基线漂移、工频干扰的三张频谱PPG和ECG的噪声来源差异很大不能共用一套滤波器参数。先看一张我常用的噪声频率对照表噪声类型典型频率范围主要受害信号危害表现工频干扰50Hz/60Hz附近ECG频域特征出现伪峰QRS检测抖动基线漂移0.5Hz以下PPG、ECG幅值特征失真波谷定位漂移运动伪影0.5Hz–10Hz宽带主要是PPG主波和切迹被叠加波淹没肌电干扰20Hz–500HzECGQRS波起点定位不稳定工频干扰和基线漂移都有明确频率区间用带通就能处理。运动伪影最麻烦它和PPG基波所在的0.5–6Hz频段高度重叠单纯滤波无法根除后文会讲怎么用信号质量评估兜底。3.2 一条够用的滤波链中值去基线 带通 工频陷波我一般对不同任务给不同频段但基础滤波链是固定的先去基线再过带通最后按需陷波。用Python写就是import numpy as np from scipy.signal import butter, filtfilt, iirnotch def clean_ppg(x, fs, lo0.5, hi6.0): # 步骤1用移动平均当高通去掉0.5Hz以下基线漂移 win_len int(fs) # 1秒窗口 kernel np.ones(win_len) / win_len baseline np.convolve(x, kernel, modesame) x_nobase x - baseline # 步骤2四阶巴特沃斯带通零相位滤波 b, a butter(4, [lo / (fs / 2), hi / (fs / 2)], btypebandpass) x_clean filtfilt(b, a, x_nobase) return x_clean def clean_ecg(x, fs, lo0.5, hi45.0, notch_freq50.0): # ECG主能量在0.5~45Hz保留QRS形态 b, a butter(2, [lo / (fs / 2), hi / (fs / 2)], btypebandpass) x_band filtfilt(b, a, x) # 50Hz工频陷波Q值取30太大会振铃 if notch_freq: bn, an iirnotch(notch_freq, Q30, fsfs) x_clean filtfilt(bn, an, x_band) return x_clean这里有两个参数非常关键。第一PPG带通上限取6Hz而不是10Hz或更高因为超过6Hz的主要是高频噪声硬要保留只会让切迹识别不稳定。如果你要保留更完整的PPG形态上限可以放宽到8Hz但要接受更多高频毛刺。第二ECG如果后面只做R波检测带通上限取45Hz足矣如果要分析ST段至少保留到100Hz但那需要更高采样率。工频陷波的Q值取30是一个折中——Q太小陷波带宽太宽会吃掉了ECG的频域能量Q太大则时域上产生振铃在QRS波附近制造假峰。使用filtfilt而不是lfilter的原因很简单filtfilt是零相位滤波不产生时间偏移这对后面算PWTT至关重要否则滤波本身就会引入几十毫秒的群延迟让特征集体失真。数据预处理里的数据清洗落到信号上不是把噪声“清零”而是保证噪声不干扰特征提取的决策。所以滤波参数要跟着特征走做HRV的ECG基线要滤干净做AI值的PPG切迹频段要留住。3.3 运动伪影滤波解决不了的部分用信号质量兜底运动伪影是PPG预处理里最玄学的部分。跑过可穿戴数据的人都有这种体验一段信号里滤完波波形看上去平了但峰值检测还是跳出一堆假峰。因为运动伪影频谱和PPG基波重合度高再加一级带通也没用。常见做法有两种一是结合三轴加速度计做自适应滤波比如LMS或NLMS把加速度信号作为参考噪声源二是做好信号质量评估SQI把污染重的片段直接丢弃或降权。我在实际项目里偏向后一种因为它不依赖额外传感器数据逻辑简单可解释。一个最基础但有效的SQI是基于峰值规律性的代码很短from scipy.signal import find_peaks def simple_sqi(signal, fs, min_dist_s0.3): # 先粗峰检测采样率125Hz时最小距离取0.3秒对应200bpm上限 peaks, _ find_peaks(signal, distanceint(min_dist_s * fs)) if len(peaks) 3: return 0.0 # 相邻峰间隔的变异系数越规律质量越高 intervals np.diff(peaks) / fs cv np.std(intervals) / (np.mean(intervals) 1e-9) sqi max(0.0, 1.0 - cv) return sqi这个SQI返回0到1之间的值越接近1表示波形越规律。实战中HRV分析要求SQI不低于0.7低于这个值的窗直接丢弃比硬滤更安全。注意这个指标在心律不齐患者身上会误伤因为RR间期本身就不规律那时候要换成基于模板匹配的SQI。这一步是很多人忽略的预处理不只是滤波还要决定哪些数据“不配进特征提取”。4. 特征提取从时域、频域到形态学特征哪些值得落地4.1 时域特征R波检测、RR间期与心率变异性基础指标特征提取前先把关键事件检测做对。ECG的核心事件是R波PPG的核心事件是主波峰。R波检测我用的是Pan-Tompkins的简化版本差分、平方、滑动平均增强再用find_peaks定位from scipy.signal import find_peaks def detect_r_peaks(ecg, fs, thr_scale0.6): # 一阶差分突出QRS斜率 diff np.diff(ecg) square diff ** 2 # 120ms滑窗平滑叠加相邻心跳的QRS能量 win int(0.12 * fs) ma np.convolve(square, np.ones(win) / win, modesame) # 最小间隔0.25秒对应最高240bpm peaks, _ find_peaks(ma, distanceint(0.25 * fs), heightthr_scale * np.mean(ma)) # 回溯到原始ECG上找精确R波位置窗口取±30ms refined [] for p in peaks: lo max(0, p - int(0.03 * fs)) hi min(len(ecg), p int(0.03 * fs)) refined.append(lo int(np.argmax(np.abs(ecg[lo:hi])))) return np.array(refined)这段代码里滑动平均窗口取120ms是依据QRS波持续时间一般小于120ms最小峰间隔取0.25秒超过240bpm的心率在静息场景下不现实。height阈值用均值乘系数好处是不用人工调绝对门限坏处是噪声大时均值会被抬高所以要配合前面SQI把坏窗先剔掉。R波位置拿到后RR间期、心率、SDNN、RMSSD都是直接计算rr np.diff(r_peaks) / fs hr 60.0 / rr rmssd np.sqrt(np.mean(np.diff(rr) ** 2)) sdnn np.std(rr)RMSSD反映短期迷走神经张力对呼吸和压力敏感SDNN反映整体变异性需要较长窗口才稳定。只给十几秒数据时RMSSD还能用SDNN就会很飘这是应用边界不是代码问题。4.2 频域特征RR间期插值重采样后再做频谱估计频域特征的最大坑是RR间期不是等间隔采样不能直接做FFT。必须先把RR序列插值到均匀时间轴上。我常用的做法是三次样条插值到4Hz再用Welch法估计功率谱提取LF和HF频段的能量比from scipy.interpolate import interp1d from scipy.signal import welch def lf_hf_ratio(rr_times, rr_vals, fs_out4.0): # rr_times: 每个RR间期对应的时间点单位秒 # rr_vals: RR间期长度单位秒 interp interp1d(rr_times, rr_vals, kindcubic, bounds_errorFalse, fill_valueextrapolate) t_uniform np.arange(rr_times[0], rr_times[-1], 1.0 / fs_out) rr_uniform interp(t_uniform) # nperseg取60秒窗口要求数据长度至少两三分钟 freqs, psd welch(rr_uniform, fsfs_out, npersegint(60 * fs_out)) lf psd[(freqs 0.04) (freqs 0.15)].sum() hf psd[(freqs 0.15) (freqs 0.40)].sum() return lf / (hf 1e-9)这里采样率fs_out取4Hz是因为HRV频域能量集中在0.4Hz以下4Hz足够。nperseg取60秒窗口意味着少于2分钟的信号段算LF/HF没什么意义这是频域特征固有的数据量要求。如果你拿30秒窗口硬算频域分辨率不够LF和HF边界会糊在一起。4.3 形态学特征PPG的AI值与ECG-PPG联合的PWTT频域特征适合做状态评估形态学特征则更直接地反映血管和血流动力学。PPG最有价值的形态学特征有两个一个是反射波增强指数AI一个是脉搏波传导时间PWTT。AI的计算依赖于重搏波切迹定位。切迹是主波之后的一个凹陷点位于两个主波峰之间def calc_ai(ppg, fs): # 主波峰距离0.4秒最小突出度取标准差10% peaks, _ find_peaks(ppg, distanceint(0.4 * fs), prominence0.1 * np.std(ppg)) ai_list [] for i in range(len(peaks) - 1): seg ppg[peaks[i]:peaks[i 1]] # 找切迹对段内波形取负找局部极小值 dicrotic, _ find_peaks(-seg, distanceint(0.1 * fs)) if len(dicrotic): # AI 切迹幅值 / 主波幅值 ai seg[dicrotic[-1]] / seg[0] ai_list.append(ai) return np.mean(ai_list)AI和动脉僵硬度相关性较好但前提是切迹清晰。年轻人血管弹性好切迹明显老年人或动脉硬化人群切迹变浅就容易漏检。这种情况下要降低prominence阈值或者改用二阶导数法找拐点但后者对噪声更敏感。另外要强调采样率低于100Hz时切迹经常被重采样磨掉所以第2章统一采样率时不要贪低。PWTT则是PPG和ECG联合特征里最有代表性的一种同一个心拍从ECG R波峰到PPG主波峰的时间间隔。计算它必须依赖第2章做过时间对齐的信号否则测出来的是“设备延迟差”而不是生理量。实战中我用的是最近邻匹配每个ECG R波向后找下一个PPG主波峰时间差就是PWTT。PWTT与血压负相关是不少无袖带血压估计算法的核心特征。特征提取方法不是越复杂越好。时域特征适合短窗和实时场景频域特征适合长窗状态分析形态学特征适合血管功能评估。用的时候先问自己要什么解释而不是把所有特征提一遍扔给模型。5. 避坑指南PPG-ECG预处理与特征提取的五个踩坑点5.1 滤波后PWTT整体偏大几十毫秒现象带通滤波前后PWTT从80毫秒变成了110多毫秒而且每个心拍都稳定偏移。原因使用lfilter或sosfilt做了非零相位滤波ECG和PPG经过滤波器后都产生了群延迟。两者的信号频率成分不同延迟量也不同导致相对时间差被系统性污染。解决统一改用filtfilt零相位滤波不引入时间偏移。如果因为实时性要求必须用因果滤波那就要记录该滤波器的群延迟在特征提取时对每个事件时间戳做补偿。这一条是形态学和联合特征计算里最容易翻车的也是最隐蔽的。5.2 重采样到低采样率后重搏波切迹消失现象把50Hz的PPG直接线性插值到125Hz后AI值算出来普遍偏高切迹检测成功率从80%掉到30%。原因线性插值只能恢复两点之间的直线切迹这种宽度只有十几到几十毫秒的细节在50Hz采样率下本身就没采到几个点插值后也无法创造信息。另一种情况是重采样时目标采样率太低抗混叠滤波器直接把切迹能量滤掉了。解决PPG原始采样率低于100Hz时形态学特征先别急着用。重采样优选resample_poly而非线性插值目标采样率尽量不低于125Hz。如果原始数据只有25HzAI和切迹相关特征就不要做了改用PPG峰峰间期类特征更可靠。5.3 训练验证集划分在预处理之后造成特征泄漏现象在公开数据集上做压力分类验证集AUC达到0.95上线后掉到0.6。重新查看代码发现是先对整段信号做了归一化和滤波再切训练集和测试集。原因全局归一化让训练集和测试集共享了统计信息滤波后的信号边界处还带着前后窗口的信息等价于把测试集的信息提前混进了训练过程。这在时间序列里特别常见叫特征泄漏。解决先按窗口切分原始信号再在每个窗口内独立做滤波、归一化和特征提取。跨窗口的滑动特征也要注意只能用历史窗口的数据来计算当前窗口特征不能用未来窗口。这一点在我做多模态特征提取时反复踩过是数据预处理里最该养成的习惯。5.4 ECG和PPG时间戳漂移导致联合特征不收敛现象HRV特征看起来与临床结果一致但任何用到PPG和ECG联合特征的血压模型训练损失都降不下去。原因HRV只在ECG单模态上计算时间戳错位影响不大而PWTT这类特征对两路信号的相对时间极其敏感。两路信号使用独立时钟长时间采集中时钟漂移会让固定延迟校准失效。解决在校准固定延迟之后按数据段长度分段重复互相关校验。一个粗略经验每5分钟信号做一次延迟估计如果延迟相对于段首变化超过一个样本就分段对齐。不要假设设备出厂时的时间校准一直有效实际数据里什么都会飘。5.5 运动伪影污染被误当成真实生理信号现象采集走路状态时PPG的频域特征里LF频段能量暴涨HRV分析得出“压力极高”的错误结论但ECG频域却完全正常。原因运动伪影与PPG基波频段重叠带通滤波器根本分不开。PPG频域特征在这种情况下反映的不是自主神经活动而是步频及其谐波。解决运动场景下先把PPG的SQI算出来SQI低于阈值的窗口直接丢弃或打上低质量标签或者引入加速度计做自适应滤波。不要把质量差的信号硬送进特征提取器那是拿人工特征去拟合噪声。6. 验证与进阶让特征提取从黑匣子变成可回退的流水线6.1 给特征提取加一道信号质量门控特征提取算法再先进也顶不住输入信号本身就是噪声。我现在的做法是给每段数据先算SQI再做门控决策SQI区间处理策略适用场景0.8–1.0正常提取全部特征静息态HRV、连续血压监测0.5–0.8只保留时域特征禁用频域和切迹特征日常活动、轻度运动低于0.5整窗丢弃或标记缺失剧烈运动、传感器松脱这样做的收益是特征表里不会混入伪影制造的异常值后续模型训练和统计分析的稳定性会明显提升。6.2 用一致性验证代替肉眼检查特征提完了要回答一个问题提得对不对。常见做法是拿你算出的心率序列和图谱设备输出的参考心率做Bland-Altman一致性分析而不是只看一张波形图“感觉差不多”。做HRV的话RR间期序列可以导入标准评价工具或与公开数据集里的标注做误差比对。PWTT特征同样可以拿压力传感器或袖带血压计做趋势一致性检验。特征提取不是画完图就结束一致性验证才证明你提取的指标真正可信。6.3 最小流水线骨架与回退能力最后给你一个最小可运行链路把今天讲的内容串起来def ppg_ecg_pipeline(ppg_raw, ecg_raw, fs_ppg, fs_ecg): # 统一到125Hz up_ppg, down_ppg 125 // fs_ppg, 1 # 实际需化简这里示意 ppg_rs resample_poly(ppg_raw, up_ppg, down_ppg, axis0) ecg_rs resample_poly(ecg_raw, up1, downfs_ecg // 125, axis0) # 互相关固定延迟校准 delay estimate_delay(ppg_rs, ecg_rs, fs125) ppg_rs shift_signal(ppg_rs, delay) # 滤波与SQI门控 ppg_clean clean_ppg(ppg_rs, fs125) ecg_clean clean_ecg(ecg_rs, fs125) if simple_sqi(ppg_clean, 125) 0.7: return None # 不做特征提取返回缺失 # 特征提取 rpeaks detect_r_peaks(ecg_clean, fs125) ppg_peaks detect_ppg_peaks(ppg_clean, fs125) features compute_hrv_features(rpeaks, fs125) features[pwtt] compute_pwtt(rpeaks, ppg_peaks) return features这个骨架的每一段都保留原始输入的副本和日志调试时能回退到任意一步。这是我做信号预处理最看重的习惯特征提取一旦出问题你能用一条命令看出是重采样、滤波、还是对齐出了问题而不是对着特征表瞎猜。这几年来我最大的教训是百分之八十的“模型不准”不是模型的问题而是时间轴没对齐、滤波参数选错、噪声窗口没排除。特征提取与其追求复杂算法不如先保证输入是干净且对齐的真实信号。希望帮到你。本文还有配套的精品资源点击获取