首页 / 资讯中心 / 文章详情

Unet+小波变换实现ECG心电信号P波QRS波T波自动分段

Unet+小波变换实现ECG心电信号P波QRS波T波自动分段 ★ FEATURED ARTICLE
简介针对心电图P、QRS、T分量自动识别任务这份Unet-ECG-Segmentation-Wavelet压缩包提供了一套基于Unet架构的完整分割实现。资源采用小波变换将QTDB数据集的心电信号转换到小波域配合WFDB工具对P、QRS、T波段进行标注并借助PyTorchWavelets实现模型训练与测试适合医学信号处理、深度学习分割方向的研究者或开发者参考。包内共27个文件以14个Python脚本包含数据加载、Unet网络定义、训练与测试等、8个编译后的pyc文件及3个Markdown说明文档为主另有Jupyter Notebook演示和License整体仅197KB结构紧凑。目前已有783人学习下载。通过完整代码与预处理流程读者可快速复现实验理解小波域ECG分割思路并在此基础上调整数据集与标注策略用于后续研究或工程实践。1. 一段心电信号里P波、QRS波和T波到底是怎么被切出来的用Unet对ECG做分段识别P波、QRS波和T波这件事听起来像是给医生省眼睛实际做起来比目标检测更磨人。心电信号不像自然图像那样有清晰的边缘P波和T波常常低平、倒置QRS波形态又因人因导联千变万化传统阈值法和模板匹配到了噪声稍大的记录上就失灵。Unet-ECG-Segmentation-Wavelet这种组合方案核心思路是先用小波变换做预处理再让Unet以逐点的方式给每个采样点打上类别标签从而端到端地输出P波、QRS波和T波的分段掩码。这个方向适合正在做心电自动分析、可穿戴设备心律失常筛查、以及心电信号预标注工具链的工程师和研究人员。它能解决的不只是“波群检测”而是把整段信号变成可测量的结构为PR间期、QT间期这类临床参数提供基础。2. 小波去噪与P-QRS-T标签生成训练数据这关怎么过2.1 为什么用db4小波而不是带通滤波器ECG预处理最常被问的问题就是都有带通滤波器了为什么还要用小波带通滤波器的问题是它会同时抹掉有用信号。ECG的有效频率范围大约在0.5Hz到45Hz之间但QRS波群的陡峭上升沿和T波的低频缓变形态都在这个范围内。固定截止频率的带通滤波器往往会把T波削平或者让P波的起始点变得模糊。小波变换的天然优势在于它同时保留时域和频域信息可以在去噪的同时尽量保持波形形态尤其是边界的尖锐程度。db4小波是心电处理里最常见的选择。原因有两个第一db4的支撑长度比较短时域局部性好对QRS波这种短时突变信号的还原能力强第二db4的波形形状和心电信号的尖峰特征有一定相似性分解后能量更集中去噪时不容易把信号分量误判成噪声。实际使用中很多人会直接拿来用但分解层数的选择是需要根据采样率来决定的。常见的做法是用离散小波变换做三层到五层分解对噪声占主导的高频细节分量做软阈值处理再重构信号。这里有一个关键点不要把所有高频层都置零。我曾经见过有人把第一层细节系数全部丢掉结果QRS波的尖端被磨平了后续用R峰检测算法一测检测率直接掉了几个点。正确的做法是只对噪声明显的那几层做阈值收缩保留部分高频信息用于维持波形锐度。2.2 pywt去噪代码与参数调整下面这段代码是我在本地环境里常用的预处理流程。输入是一个numpy数组形式的单导联ECG信号采样率默认设为360Hz这个频率和很多公开心电数据集的采样率一致。import numpy as np import pywt def wavelet_denoise_ecg(signal, sample_rate360, waveletdb4, levelNone): # 根据信号长度自动选择分解层数 if level is None: level min(5, pywt.dwt_max_level(len(signal), wavelet)) # 离散小波分解 coeffs pywt.wavedec(signal, wavelet, levellevel) # 对每一层细节系数做软阈值处理 sigma np.median(np.abs(coeffs[-1])) / 0.6745 threshold sigma * np.sqrt(2 * np.log(len(signal))) coeffs_th [coeffs[0]] for i in range(1, len(coeffs)): coeffs_th.append(pywt.threshold(coeffs[i], threshold, modesoft)) # 重构去噪后的信号 clean_signal pywt.waverec(coeffs_th, wavelet) return clean_signal[:len(signal)]这段代码的逻辑分几步走。先用dwt_max_level算出当前信号长度支持的最大分解层数避免层数过多导致边界失真接着用wavedec做分解其中coeffs[0]是近似系数后面的列表项是各层细节系数。噪声水平估计用经典的中值绝对偏差法计算出统一的threshold后对每一层细节系数做软阈值处理。最后用waverec重构信号返回时裁剪到原始长度。参数调整有两个重点。第一wavelet不要轻易换成db8或sym8这类更长的小波支撑长度越长计算越慢边界处的波形失真越明显。第二modesoft通常比hard效果好因为硬阈值会在信号里留下不连续的毛刺后面对Unet训练不利。如果信号里基线漂移特别严重可以在分解后直接把最后一层近似系数做一次减去均值或高通处理再重构这样比单纯提高阈值更有效。2.3 标注策略基于R峰先验的弱标签生成训练Unet做心电分段最大的困难不是网络结构而是标签。公开数据集里带有P波、QRS波、T波手动标注的数据非常少而且不同医生对P波起点和T波终点的判断并不完全一致。指望拿到大批量精标数据不太现实所以大多数落地项目都采用“自动标注人工抽检”的方式来生产标签。常见做法是先做R峰检测再以每个R峰位置为锚点按照时间窗口向前向后推算P波、QRS波和T波的近似范围。以360Hz采样率为例QRS波可以取R峰前后各0.04秒也就是约15个采样点P波一般取R峰前0.2到0.3秒区域T波从R峰后0.2秒开始持续到R峰后0.4到0.5秒。这样生成的标签是弱标签边界位置不一定精确但作为Unet的训练信号已经足够。下面是生成标签的一个示例逻辑from scipy.signal import find_peaks def generate_weak_labels(ecg_signal, r_peaks, sample_rate360): labels np.zeros(len(ecg_signal), dtypenp.int64) qrs_half_width int(0.04 * sample_rate) p_start_offset int(0.30 * sample_rate) p_end_offset int(0.12 * sample_rate) t_start_offset int(0.15 * sample_rate) t_end_offset int(0.45 * sample_rate) for r_idx in r_peaks: # QRS 标注 labels[max(0, r_idx - qrs_half_width):min(len(ecg_signal), r_idx qrs_half_width)] 2 # P 波标注 p_start max(0, r_idx - p_start_offset) p_end max(0, r_idx - p_end_offset) labels[p_start:p_end] 1 # T 波标注 t_start min(len(ecg_signal), r_idx t_start_offset) t_end min(len(ecg_signal), r_idx t_end_offset) labels[t_start:t_end] 3 return labels需要注意的是这种窗口式标注在心率过快或过慢时会出错。心率80bpm时RR间期约0.75秒P波窗口和T波窗口不会重叠但心率到120bpm以上时RR间期只有0.5秒P波和T波区域会出现重叠后写的标注会覆盖先写的。实际处理中我一般会根据每个心拍的RR间期做动态缩放而不是用固定的时间窗口。2.4 划分数据按患者在划分而不是按心拍ECG分段任务中一个隐蔽但严重的错误就是随机划分心拍。同一个人的不同心拍之间形态高度相似如果同一个人的心拍同时出现在训练集和验证集里模型会“记住”个体特征而不是学习到通用的波形结构导致验证集指标虚高。正确做法是按患者划分数据确保同一个患者的所有数据只出现在一个集合里。操作上先把所有记录按患者ID分组然后按患者ID做随机划分。训练集、验证集、测试集的比例我一般用7:1.5:1.5。如果数据量不足可以采用跨患者交叉验证每次留出一个患者的所有记录做验证。这个细节直接影响模型上线后的真实表现值得在标注环节就规划好。3. 把二维Unet改成一维信号分割网络结构与最小复现代码3.1 从图像分割到一维信号分割要改的三件事Unet最初是为二维图像分割设计的编码器逐层下采样解码器逐层上采样中间用跳跃连接拼接特征。改成ECG分段后不需要大改骨架但有三处必须调整。第一是卷积核和池化全部从二维改成朝向时间轴的一维。图像里用3x3卷积一维信号里用长度为7或9的卷积核。第二是下采样方式。图像分割常使用stride2的卷积或2x2最大池化一维信号用stride2的一维卷积或长度2的最大池化。第三是输入通道数。ECG单导联输入时通道数是1如果使用多导联输入比如12导联则可以把通道数改成12但要注意不同导联间需要做时间对齐。还有一个容易被忽视的点Unet的池化层会让特征图的长度逐层减半如果输入长度不是2的幂次解码器拼接时会出现尺寸对不上的问题。所以输入信号长度一般选256、512、1024、2048或4096这样的值。3.2 1D Unet代码编码器、跳跃连接和输出头下面是一份基于PyTorch的最小实现没有复杂的注意力模块只保留Unet最核心的编码器-解码器-跳跃连接结构。这份代码可以直接跑通适合作为基线。import torch import torch.nn as nn class ConvBlock(nn.Module): def __init__(self, in_ch, out_ch): super().__init__() self.conv nn.Sequential( nn.Conv1d(in_ch, out_ch, kernel_size7, padding3), nn.BatchNorm1d(out_ch), nn.ReLU(inplaceTrue), nn.Conv1d(out_ch, out_ch, kernel_size7, padding3), nn.BatchNorm1d(out_ch), nn.ReLU(inplaceTrue) ) def forward(self, x): return self.conv(x) class EncoderBlock(nn.Module): def __init__(self, in_ch, out_ch): super().__init__() self.conv ConvBlock(in_ch, out_ch) self.pool nn.MaxPool1d(kernel_size2, stride2) def forward(self, x): x self.conv(x) return x, self.pool(x) class DecoderBlock(nn.Module): def __init__(self, in_ch, out_ch): super().__init__() self.up nn.ConvTranspose1d(in_ch, out_ch, kernel_size2, stride2) self.conv ConvBlock(out_ch * 2, out_ch) def forward(self, x, skip): x self.up(x) x torch.cat([x, skip], dim1) return self.conv(x) class UNet1D(nn.Module): def __init__(self, in_channels1, num_classes4): super().__init__() self.enc1 EncoderBlock(in_channels, 32) self.enc2 EncoderBlock(32, 64) self.enc3 EncoderBlock(64, 128) self.enc4 EncoderBlock(128, 256) self.bridge nn.Sequential( nn.Conv1d(256, 512, kernel_size7, padding3), nn.BatchNorm1d(512), nn.ReLU(inplaceTrue) ) self.dec1 DecoderBlock(512, 256) self.dec2 DecoderBlock(256, 128) self.dec3 DecoderBlock(128, 64) self.dec4 DecoderBlock(64, 32) self.head nn.Conv1d(32, num_classes, kernel_size1) def forward(self, x): s1, x self.enc1(x) s2, x self.enc2(x) s3, x self.enc3(x) s4, x self.enc4(x) x self.bridge(x) x self.dec1(x, s4) x self.dec2(x, s3) x self.dec3(x, s2) x self.dec4(x, s1) return self.head(x)这里有两个设计要点。第一卷积核用7而不是3原因是ECG波形在时间尺度上跨度较大较小的卷积核需要堆叠很多层才能获得足够的感受野而7的卷积核在一维信号上计算量并不高。第二每层特征通道数从32开始逐层翻倍到bridge层为512。对大多数ECG数据集来说这个容量足够如果数据量特别大可以把初始通道数调到64。3.3 输入长度选择256还是4096输入长度的选择决定了一个训练样本里包含多少个完整心拍。ECG的RR间期在静息状态下大约为0.6到1.0秒在360Hz采样率下就是216到360个采样点。如果输入长度选256可能只有一个不完整的心拍Unet很难学到P波和T波的完整上下文。我一般选1024或2048这样每个样本能包含3到8个心拍模型既能学到局部波形细节又能感知节奏变化。输入长度还会影响显存占用和推理速度。2048长度配合batch size 32在单张24GB显存的显卡上很轻松。如果部署在边缘设备上可以把长度压到512但要重新验证分段边界质量。长度变化后不需要改网络结构只要把最后一个池化步骤产生的特征图尺寸和对应跳跃连接的尺寸保持一致即可。4. 损失函数与评价指标训练一个不“自欺欺人”的分段模型4.1 损失函数Dice Loss与Focal Loss怎么组合ECG分段本质上是一个逐像素分类问题类别包括背景、P波、QRS波和T波四类。但直接使用交叉熵损失的效果并不好原因在于类别不平衡。一段10秒的心电信号里背景类占比可能超过70%QRS波只占不到5%P波和T波各占10%左右。模型如果全部预测为背景交叉熵损失也不会太高但实际没有任何可用性。通用做法是将Dice Loss和Focal Loss按权重相加。Dice Loss对类别不平衡不敏感它直接优化预测区域和真实区域的重叠度而Focal Loss让模型把注意力集中在难分类的样本上。import torch.nn.functional as F def dice_loss(pred, target, smooth1.0): # pred: [batch, num_classes, length] # target: [batch, length] num_classes pred.shape[1] target_one_hot F.one_hot(target, num_classesnum_classes).permute(0, 2, 1).float() pred_soft F.softmax(pred, dim1) intersection (pred_soft * target_one_hot).sum(dim2) union pred_soft.sum(dim2) target_one_hot.sum(dim2) dice (2.0 * intersection smooth) / (union smooth) return 1.0 - dice.mean()代码里的dice按类别分别计算后取平均这样做可以让QRS波这类小目标不会被背景类淹没。Focal Loss可以直接用PyTorch的实现参数gamma设置为2alpha设置为类别占比的倒数效果比较稳定。组合时我通常把Dice Loss权重设为0.6Focal Loss权重设为0.4。4.2 评价指标像素级F1之外还要看P波起止偏差训练过程中的验证指标如果只看逐像素F1很容易被背景类主导。背景类占比高即使其他三类预测得一团糟F1也可能在0.85以上。所以评估时需要拆开看每一类的F1尤其是P波和T波这两个容易被忽略的类别。更贴近临床使用的评价方式是做“事件级别匹配”。先分别从预测结果和真实标签里找到连续的P波段、QRS波段和T波段再判断预测段与真实段是否有足够的交叠比如IoU大于0.5即认为检测成功。这个方法能直接反映“这一段波有没有被找到”和逐像素打分完全不同。如果Unet把一段T波拆成两段逐像素F1可能不低但事件级别匹配会认为这个T波检测失败。4.3 一组稳定的训练参数基于我跑过的多次实验下面这组参数是稳定且容易复现的起点config { input_length: 2048, batch_size: 32, epochs: 80, lr: 3e-4, optimizer: AdamW, weight_decay: 1e-4, scheduler: cosine_annealing, warmup_epochs: 5, loss_weights: {dice: 0.6, focal: 0.4}, }学习率3e-4配合AdamW在大多数ECG数据集上都不会发散但建议配合warmup使用前5个epoch从1e-5线性升到3e-4。原因是ECG信号的标签存在大量边界噪声网络在初期如果学习率太大很容易被这些噪声标签带偏。cosine_annealing可以让后期的学习率缓慢下降对收敛更加友好。如果显存不够优先减小batch size而不是输入长度因为输入长度影响上下文信息完整性。5. 常见踩坑记录Unet做心电分段的5个典型翻车场景5.1 模型全输出背景类正负样本极端不平衡现象训练了十几个epoch后验证集里P波、QRS波、T波全部预测为背景逐像素准确率却有90%以上。原因背景类占比太高模型找到了一个局部最优解全部预测为背景损失函数值依然很低。这个问题在使用交叉熵损失时尤其常见。解决切换到Dice Loss或Focal Loss并在数据加载时进行更精细的采样。比如每个batch里保证包含至少一个完整心拍而不是随机裁剪任意长度。另外可以给P波、QRS波和T波这几个类别设置更高的权重让损失函数更敏感。5.2 小波去噪后QRS波被“磨平”现象预处理完的信号波形看起来干净了但训练出的模型预测结果里QRS波变宽R峰附近出现了很多细小的碎片预测。原因软阈值去噪时阈值设置过高尤其是对第一层细节系数处理过重导致QRS波的高频成分被削弱。模型学到的QRS波边界因此变得模糊。解决去噪时把第一层细节系数的阈值调低甚至不对第一层做阈值处理只处理第二层及以上的系数。另一个思路是先检测R峰再在R峰附近局部削弱去噪强度以保留波形锐度。5.3 R峰误检导致标签错位现象训练过程中损失值不下降或者说验证集事件匹配率很低但人工检查预测结果时发现预测的P波和T波位置总比真实波形偏了约0.1秒。原因自动标签生成依赖R峰检测但R峰检测在某些噪声上会产生半个心拍左右的偏移。如果某个心拍的R峰被检测到T波位置那么这个心拍的所有标签整体错位训练信号就是错的。加之后续没有做人工修正模型只能学习一个平均偏移。解决在标签生成后加入一致性检查。比如检查RR间期是否在合理范围内P波起点是否在QRS波之前T波终点是否在下一个P波之前。对于有明显矛盾的样本直接丢弃而不是强行保留。5.4 心拍长度剪裁导致P波被截断现象训练数据生成时以R峰为中心裁剪心拍给定长的固定窗口如R峰前后各0.4秒。在高心率样本上P波被截断了一半模型对这类样本的预测结果很差。原因固定窗口宽度在设计时按静息心率设定但没有考虑高心率场景下P波离QRS波更近窗口的前缘不够。解决根据每个样本的瞬时心率动态调整窗口宽度。心率超过100bpm时P波窗口缩短到R峰前0.2秒同时T波窗口短到R峰后0.3秒。更稳定的做法是不裁剪单个心拍而是直接输入固定长度的多心拍片段让模型自己学习心拍之间的上下文。5.5 验证集不按患者划分导致模型“作弊”现象训练时损失很低验证集F1也很高但把模型迁移到另一批新数据上时事件级别匹配率骤降。原因验证集和训练集包含同一个患者的数据模型记住了患者的波形形态和噪声模式。这种“作弊”在数据量小的实验里特别容易发生。解决严格按患者ID划分数据确保一个患者的所有记录只出现在训练集或验证集中的一个。迁移测试时使用一个从未参与过任何训练步骤的外部数据集做最终评估。6. 用R峰锚点修正Unet输出边界让分段结果真正能用于测量6.1 边界为什么模糊以及R峰锚点的作用Unet输出的分割掩码在P波起点和T波终点处往往是“渐变”的因为下采样和上采样过程会给边界带来天然的平滑效应。对于临床测量来说PR间期和QT间期这类参数的精度要求是毫秒级一个采样点都不能差。直接拿Unet的softmax输出去算边界稳定性不够。常见做法是用R峰检测结果作为锚点对Unet输出做几何约束修正。R峰是整个心电信号里最显著的特征绝大多数R峰检测算法在去噪后的信号上都很稳定。有了R峰位置就可以知道每个心拍的周期边界再结合Unet预测出的P波和T波概率分布在真实R峰附近的合理范围内做精确定位。6.2 后处理代码anchored correction下面这段后处理逻辑输入是Unet输出的概率矩阵、R峰位置数组和采样率输出是修正后的分段掩码。import numpy as np def anchor_corrected_mask(prob_map, r_peaks, sample_rate360): prob_map: [num_classes, length], 假设类别顺序为 0背景,1P,2QRS,3T r_peaks: 检测到的R峰索引数组 num_classes, length prob_map.shape mask np.zeros(length, dtypenp.int64) for r in r_peaks: # QRS 范围直接在R峰附近取 qrs_start max(0, int(r - 0.04 * sample_rate)) qrs_end min(length, int(r 0.04 * sample_rate)) # 在窗口内取概率最大的类别 qrs_window np.argmax(prob_map[:, qrs_start:qrs_end], axis0) mask[qrs_start:qrs_end] qrs_window # P 波在前面的搜索窗口内找非背景最可能的点 p_win_start max(0, int(r - 0.36 * sample_rate)) p_win_end max(0, int(r - 0.10 * sample_rate)) p_probs prob_map[1, p_win_start:p_win_end] if p_probs.max() 0.5: p_peak p_win_start np.argmax(p_probs) mask[p_win_start:p_win_end] 1 # T 波在后面的搜索窗口内修正 t_win_start min(length, int(r 0.12 * sample_rate)) t_win_end min(length, int(r 0.48 * sample_rate)) t_probs prob_map[3, t_win_start:t_win_end] if t_probs.max() 0.5: t_peak t_win_start np.argmax(t_probs) mask[t_win_start:t_win_end] 3 return mask这段代码的核心逻辑不是简单取最大概率而是用R峰位置把搜索范围限制在生理合理的窗口内再在窗口内寻找P波和T波的概率峰值。这样一来即使Unet预测在远处出现了一个虚假的高概率段也不会被纳入修正结果因为搜索窗口限制了候选区域。阈值0.5可以按实际数据调整如果模型预测偏保守阈值可以降低到0.4如果误检多则提到0.6。6.3 落地验证的指标与最终建议最后一步验证不要只看平均指标要把每个信号单独跑一遍输出每个P波、QRS波、T波的起点和终点然后和医生标注做时间差分析。我习惯用三个数字判断是否达标P波起点偏差小于20毫秒的比例、T波终点偏差小于30毫秒的比例、QRS波边界重叠率大于80%的比例。只有这三个指标同时达到要求分段结果才敢拿去做QT间期分析。我做这个项目的最大教训是网络结构其实不是瓶颈数据标签质量和后处理逻辑才是决定成败的关键。别急着上注意力机制和Transformer先把手里的预处理流程跑稳再逐段检查Unet的输出你会发现问题大多出在数据侧。希望这篇笔记能帮你在做ECG分段实验时少走几步弯路。本文还有配套的精品资源点击获取
阅读完成 · 觉得有帮助?
咨询建站