Linux 上 wfdb.tar.gz 心电信号分析实战：从解压到 R 波检测

简介这份资源是面向Linux平台心电信号处理与分析的WFDB软件包适合生物医学工程、心电算法研究及信号处理方向的学习者与开发者。WFDB由MIT-BIH心律失常实验室开发支持C、C、Java、Python等多种语言可用于读取、显示、分析心电记录并结合MIT-BIH基准数据库开展算法验证。压缩包共545个文件约1.91MB以C源码、头文件、makefile、hea头文件、dat数据文件及readme说明为主另含sh脚本、png图示与tex文档覆盖编译安装、命令行工具与示例数据。已有213人学习下载。借助其中的rdann、rdsamp、pg等工具与库函数读者可完成滤波去噪、基线漂移校正、PQRST波群检测、心率与RR间期计算及心律失常异常检测快速搭建Linux下的心电分析环境为心脏病诊断与预防研究提供可复用的开发基础。1. 拿到 wfdb.tar.gz 之后Linux 上跑通心电信号分析的第一道坎很多人第一次接触心电信号分析是从 PhysioNet 上拖下来一个wfdb.tar.gz压缩包开始的。解压一看里面是.hea、.dat、.atr这类陌生后缀既不是 CSV 也不是图片用 Excel 打不开用文本编辑器看是乱码。这个包本质上就是 WFDBWaveForm DataBase格式的物理信号文件集合配套的是wfdb这个 Python 工具链。它解决的核心问题很具体把 MIT-BIH、PTB 这类标准心电数据库的原始采样点读进内存转成 numpy 数组再叠加注释文件做 R 波定位、心率变异性分析。适合谁做可穿戴心电算法验证的嵌入式工程师、跑信号处理课程设计的学生、以及需要在 Linux 服务器上批量处理心电记录的运维兼算法角色。如果你手上正好有这个压缩包下面这套流程能让你在两小时内从解压走到画出第一张 RR 间期图。2. 环境搭建与 wfdb 工具链从 tar.gz 到可导入的 Python 模块2.1 为什么优先在 Linux 上处理 WFDB 数据WFDB 的官方参考实现是 C 语言写的Python 包wfdb底层通过 ctypes 调用这些库在 Linux 上的兼容性最稳。Windows 上经常遇到libwfdb.so找不到的问题而 Linux 发行版自带的包管理器能直接补齐依赖。另一个现实原因是心电数据库动辄几个 GBLinux 的文件系统和命令行工具在批量解压、校验、切分上效率明显更高。我一般会在 Ubuntu 22.04 或 Rocky Linux 9 上做这件事两者都验证过。需要提前确认的是 Python 版本wfdb4.x 要求 Python 3.8 以上3.11 也跑得通。2.2 解压与目录结构确认拿到wfdb.tar.gz后不要急着tar -xzf一把梭先看一眼压缩包里的顶层目录名避免解压出一堆散文件污染当前目录。# 先列出压缩包内容确认顶层目录结构 tar -tzf wfdb.tar.gz | head -20 # 确认无误后解压到指定目录保持目录层级 mkdir -p ~/ecg_data tar -xzf wfdb.tar.gz -C ~/ecg_data # 查看解压后的文件类型分布 find ~/ecg_data -type f | sed s/.*\.// | sort | uniq -c第一行tar -tzf只列出不解压是防止压缩包内文件散落的关键习惯。-C参数指定解压目标目录避免在当前工作目录制造混乱。最后那条find配合sed和uniq -c能快速告诉你这个包里到底有几种文件后缀正常应该看到hea、dat、atr三类如果出现xws或edf说明数据来源不止一种后续读取方式要区分。2.3 安装 wfdb 与依赖不要用pip install wfdb直接装完就完事心电分析还需要 numpy、scipy、matplotlib 这三个基础库版本不匹配会在滤波环节报奇怪的错误。# 创建独立虚拟环境避免污染系统 Python python3 -m venv ~/ecg_env source ~/ecg_env/bin/activate # 安装核心依赖指定版本范围避免 API 变动 pip install wfdb4.1,5.0 numpy1.24 scipy1.10 matplotlib3.7 # 验证 wfdb 能否正常加载底层库 python -c import wfdb; print(wfdb.__version__); print(wfdb.rdrecord.__doc__[:80])虚拟环境这一步在 Linux 上尤其重要因为系统自带的 Python 往往被其他工具依赖直接pip install可能触发权限问题或版本冲突。wfdb的版本锁定在 4.x 是因为 3.x 和 4.x 的rdrecord返回对象结构有差异网上很多老教程用的是 3.x 写法照抄会报AttributeError。最后那条验证命令如果打印出版本号和函数文档片段说明底层libwfdb已经就位。2.4 读取一条记录并确认采样参数环境就绪后先用一条记录验证整条链路。假设解压后的目录里有一个名为100的记录MIT-BIH 格式的常见命名。import wfdb import numpy as np # 读取记录path 指向解压目录sampto 限制读取长度便于快速验证 record wfdb.rdrecord( 100, pn_dirNone, sampto5000, physicalTrue ) # 打印关键参数确认采样率和导联数 print(f采样率: {record.fs} Hz) print(f导联名: {record.sig_name}) print(f信号形状: {record.p_signal.shape}) print(f单位: {record.units}) # 取第一导联前 10 个采样点看看数值范围 print(record.p_signal[:10, 0])physicalTrue表示返回的是经过物理单位换算的毫伏值如果设为False则返回原始 ADC 整数。sampto5000只读前 5000 个采样点在调试阶段能省掉大量等待时间。p_signal的 shape 是(采样点数, 导联数)这一点和很多人的直觉相反不是行优先。打印出的数值范围正常应该在 -5 到 5 毫伏之间如果看到几千的整数说明物理换算没生效。3. 信号预处理与 R 波检测把原始采样点变成可分析的心搏序列3.1 基线漂移与工频干扰的处理顺序原始心电信号里混着两类主要噪声基线漂移0.5 Hz 以下的低频漂移和工频干扰50 Hz 或 60 Hz。处理顺序有讲究先做带通滤波把两者一起压掉再做 Notch 滤波专门针对工频。如果反过来Notch 滤波器的瞬态响应会被基线漂移放大反而引入新的伪影。from scipy.signal import butter, filtfilt, iirnotch def bandpass_filter(signal, fs, low0.5, high40.0, order3): 巴特沃斯带通滤波保留 0.5-40 Hz 心电主要能量 nyq fs / 2.0 b, a butter(order, [low / nyq, high / nyq], btypeband) # filtfilt 零相位滤波避免 R 波位置偏移 return filtfilt(b, a, signal) def notch_filter(signal, fs, freq50.0, quality30.0): IIR Notch 滤波压制工频干扰 b, a iirnotch(freq, quality, fs) return filtfilt(b, a, signal) # 对第一导联依次处理 sig record.p_signal[:, 0] sig_bp bandpass_filter(sig, record.fs) sig_clean notch_filter(sig_bp, record.fs, freq50.0)filtfilt而不是lfilter是关键前者做正向和反向两次滤波相位偏移为零R 波位置不会移动。order3是经验值阶数太高会在 QRS 复合波边缘产生振铃。Notch 的quality30决定了陷波带宽值越大带宽越窄对 50 Hz 的压制越精准但对频率漂移越敏感。国内工频是 50 Hz如果数据来自北美设备则要改成 60 Hz这个参数搞错等于没滤。3.2 用 wfdb 自带注释做 R 波定位WFDB 格式的一大优势是注释文件里已经标好了专家标注的 R 波位置不需要自己写检测算法就能拿到基准。但要注意注释的采样点编号和信号采样点是一一对应的。# 读取注释文件 annotation wfdb.rdann(100, atr, sampto5000) # 筛选出正常心搏的注释类型 # N 表示正常搏动L 和 R 表示左右束支阻滞 normal_beats [ (pos, sym) for pos, sym in zip(annotation.sample, annotation.symbol) if sym in (N, L, R) ] # 提取 R 波位置数组 r_peaks np.array([pos for pos, _ in normal_beats]) print(f检测到 {len(r_peaks)} 个正常心搏) # 计算 RR 间期单位毫秒 rr_intervals np.diff(r_peaks) / record.fs * 1000 print(f平均心率: {60000 / np.mean(rr_intervals):.1f} bpm)rdann的第二个参数atr是注释文件的后缀不同数据库可能用qrs或ecg。annotation.symbol是一个字符列表每个字符对应一个心搏类型N是最常见的正常搏动。RR 间期计算用np.diff得到相邻 R 波采样点差值除以采样率再乘 1000 转成毫秒。平均心率用 60000 除以平均 RR 间期这是标准算法但要注意如果 RR 间期方差很大平均值会掩盖心律失常这时候要看 RR 间期的标准差。3.3 没有注释文件时自己写 Pan-Tompkins 检测有些数据库不带注释或者你想验证自己的检测算法就需要实现一个 QRS 检测器。Pan-Tompkins 是经典方案核心是微分、平方、移动窗积分三步。def pan_tompkins(signal, fs): 简化版 Pan-Tompkins QRS 检测 # 1. 带通滤波已在外部完成这里直接微分 diff_sig np.ediff1d(signal, to_begin0) # 2. 平方放大高频成分 squared diff_sig ** 2 # 3. 移动窗积分窗口宽度约 150ms window_size int(0.15 * fs) integrated np.convolve(squared, np.ones(window_size) / window_size, modesame) # 4. 自适应阈值取积分信号的最大值的 0.3 倍作为初始阈值 threshold 0.3 * np.max(integrated) peaks [] for i in range(1, len(integrated) - 1): if integrated[i] threshold and integrated[i] integrated[i-1] and integrated[i] integrated[i1]: # refractory period 200ms 内不重复检测 if not peaks or (i - peaks[-1]) int(0.2 * fs): peaks.append(i) return np.array(peaks)微分用np.ediff1d而不是np.diff前者会在开头补零保持长度一致。移动窗积分用np.convolve配合全一数组除以窗口长度等价于滑动平均。阈值取最大值的 0.3 倍是 Pan-Tompkins 论文里的经验值但实际数据中如果噪声大这个阈值会漏检常见做法是用0.3 * np.percentile(integrated, 99)替代np.max来抗离群点。 refractory period 设为 200ms 是因为生理上两次心搏不可能间隔这么短这个约束能过滤掉 T 波误检。4. 避坑与排查心电信号处理里那些让人怀疑数据的瞬间4.1 读取时报FileNotFoundError但文件明明存在现象wfdb.rdrecord(100)报找不到文件但ls能看到100.hea和100.dat。原因rdrecord的第一个参数是记录名不是文件名。如果当前工作目录不在数据目录下或者记录名带了路径分隔符就会解析失败。解决用os.chdir切到数据目录或者给rdrecord传完整路径但去掉后缀。例如wfdb.rdrecord(/home/user/ecg_data/100)注意不要写成100.hea。4.2 画出来的波形是一条直线现象plt.plot(record.p_signal[:, 0])显示一条平线数值全为零或接近零。原因physicalTrue时如果.hea文件里的增益和基线字段格式不标准物理换算会得到全零。或者读取时channels参数选错了导联索引。解决先用physicalFalse读原始 ADC 值确认信号存在再检查record.adc_gain和record.baseline是否合理。如果原始值正常但物理值为零手动做换算(raw - baseline) / gain。4.3 R 波检测数量远多于实际心搏现象自己写的检测算法报出几百个 R 波但 10 秒数据里正常只有十几次心搏。原因阈值设得太低把 T 波和噪声都算进去了。或者没有做 refractory period 约束。解决把阈值从0.3 * max提高到0.5 * max并强制 200ms 的不应期。更稳妥的做法是先对积分信号做一次 5 点中值滤波把毛刺压掉再检测。4.4 滤波后 R 波幅度变小甚至消失现象带通滤波后波形变平滑但 R 波峰值从 1.5mV 降到 0.3mV。原因带通滤波的高截止频率设得太低比如设成 15 HzQRS 复合波的主要能量在 10-40 Hz 之间被削掉了。解决高截止频率至少设到 40 Hz诊断级心电分析常用 0.05-100 Hz但那样会保留更多噪声。常规分析用 0.5-40 Hz 是平衡点。4.5 注释文件里的采样点和信号对不上现象annotation.sample的最大值超过了record.p_signal的长度。原因读取信号时用了sampto限制长度但读取注释时没加同样的限制导致注释索引越界。解决rdann也要传sampto参数保持和rdrecord一致。或者先读完整信号再截取不要在两处分别限制。5. 批量处理与结果验证把单条记录的经验固化成脚本单条记录跑通只是起点真正干活时面对的是几十上百条记录。我一般会写一个批处理脚本把读取、滤波、R 波检测、RR 间期计算、结果保存串起来并且加一个验证环节——用注释文件的 R 波位置和自己检测的结果做比对算敏感度和阳性预测值。import os import glob import numpy as np import wfdb from scipy.signal import butter, filtfilt def process_record(record_path, fs_targetNone): 处理单条记录返回 RR 间期和检测指标 record wfdb.rdrecord(record_path) annotation wfdb.rdann(record_path, atr) sig record.p_signal[:, 0] fs record.fs # 带通滤波 nyq fs / 2.0 b, a butter(3, [0.5 / nyq, 40.0 / nyq], btypeband) sig_clean filtfilt(b, a, sig) # 用注释作为基准计算 RR 间期 r_peaks annotation.sample rr np.diff(r_peaks) / fs * 1000 # 简单的异常值过滤RR 间期在 300-2000ms 之外视为伪影 valid_rr rr[(rr 300) (rr 2000)] return { record: os.path.basename(record_path), fs: fs, n_beats: len(r_peaks), mean_rr: np.mean(valid_rr) if len(valid_rr) 0 else 0, std_rr: np.std(valid_rr) if len(valid_rr) 0 else 0, hr: 60000 / np.mean(valid_rr) if len(valid_rr) 0 else 0 } # 批量处理目录下所有 .hea 文件 data_dir os.path.expanduser(~/ecg_data) hea_files glob.glob(os.path.join(data_dir, *.hea)) results [] for hea in hea_files: record_name hea[:-4] # 去掉 .hea 后缀 try: res process_record(record_name) results.append(res) print(f{res[record]}: HR{res[hr]:.1f} bpm, beats{res[n_beats]}) except Exception as e: print(f{record_name} 处理失败: {e}) # 汇总统计 if results: hrs [r[hr] for r in results if r[hr] 0] print(f\n共处理 {len(results)} 条记录) print(f平均心率范围: {min(hrs):.1f} - {max(hrs):.1f} bpm)这个脚本里glob.glob匹配所有.hea文件然后用切片去掉后缀得到记录名。try/except包住每条记录的处理避免一条坏数据中断整个批次。RR 间期的 300-2000ms 过滤范围对应心率 30-200 bpm超出这个范围的通常是漏检或误检。最后汇总时只统计心率大于零的记录防止除零错误。验证环节我习惯抽三条记录把注释的 R 波位置和滤波后的信号画在一起肉眼确认对齐。这一步不能省因为有些数据库的注释文件本身就有标注错误盲目相信注释会导致后续分析全盘偏差。画图代码很短import matplotlib.pyplot as plt # 抽一条记录画图验证 record wfdb.rdrecord(os.path.join(data_dir, 100), sampto3000) ann wfdb.rdann(os.path.join(data_dir, 100), atr, sampto3000) fig, ax plt.subplots(figsize(12, 4)) ax.plot(record.p_signal[:, 0], linewidth0.8, labelECG) ax.scatter(ann.sample, record.p_signal[ann.sample, 0], colorred, s20, zorder5, labelAnnotation) ax.set_xlabel(Sample) ax.set_ylabel(mV) ax.legend() plt.tight_layout() plt.savefig(verify_100.png, dpi150)zorder5保证红点画在波形上方s20控制点的大小。保存成 PNG 而不是直接show是因为在服务器上跑没有图形界面保存文件再下载查看更实际。从那以后我每次拿到新的 WFDB 数据包都强制先跑一遍这个验证流程读一条、画一张、对一次注释。这个习惯帮我拦下过好几次采样率标错和导联顺序颠倒的问题。希望帮到你。本文还有配套的精品资源点击获取