### 使用 Pan-Tompkins 算法实现 ECG R 波峰定位 为了在 PyCharm 中使用 Python 实现 ECG 的 R 波峰定位功能，可以采用经典的 Pan-Tompkins 算法。该算法通过一系列滤波器来增强 QRS 复合波并抑制其他成分，从而有效地检测到 R 峰。 #### 导入必要的库 首先需要安装 `wfdb` 和 `scipy` 库用于读取数据和处理信号： ```bash pip install wfdb scipy numpy matplotlib ``` 接着导入所需的模块： ```python import numpy as np import matplotlib.pyplot as plt import wfdb from scipy.signal import find_peaks, butter, lfilter ``` #### 加载心电图数据 加载 MIT-BIH 心律失常数据库中的记录文件，并提取所需通道的数据： ```python record = wfdb.rdrecord('mit-bih-arrhythmia-database-1.0.0/106', sampfrom=0, sampto=1000, physical=True, channels=[0]) ECG = record.p_signal[:, 0] fs = record.fs time = np.arange(len(ECG)) / fs plt.plot(time, ECG) plt.title("原始心电信号") plt.xlabel("时间 (秒)") plt.ylabel("振幅") plt.show() ``` #### 设计带通滤波器 设计一个带通滤波器去除低频基线漂移以及高频噪声干扰： ```python def bandpass_filter(signal, lowcut, highcut, fs, order=5): nyquist = 0.5 * fs low = lowcut / nyquist high = highcut / nyquist b, a = butter(order, [low, high], btype='band') y = lfilter(b, a, signal) return y filtered_ecg = bandpass_filter(ECG, 5, 15, fs) plt.figure(figsize=(8,4)) plt.plot(time, filtered_ecg) plt.title("经过带通滤波后的信号") plt.xlabel("时间 (秒)") plt.ylabel("振幅") plt.show() ``` #### 计算导数 计算一阶差分近似求解微分方程，突出QRS波群特征： ```python diff_ecg = np.diff(filtered_ecg) plt.figure(figsize=(8,4)) plt.plot(np.linspace(0,time[-2],len(diff_ecg)), diff_ecg) plt.title("一阶差分结果") plt.xlabel("时间 (秒)") plt.ylabel("变化率") plt.show() ``` #### 平方运算和平滑化处理 对上述得到的结果做平方操作再进行移动平均平滑化处理： ```python square_ecg = diff_ecg ** 2 window_size = int(round(fs*0.12)) # 设置窗口大小约为采样频率的十分之一周期长度 smoothed_square_ecg = np.convolve(square_ecg, np.ones(window_size)/float(window_size), mode="same") plt.figure(figsize=(8,4)) plt.plot(time[:-1], smoothed_square_ecg[:len(time)-1]) plt.title("平方后并经卷积平滑化的信号") plt.xlabel("时间 (秒)") plt.ylabel("强度") plt.show() ``` #### 阈值判定与峰值查找 设定合适的阈值条件筛选出可能存在的R波位置；利用SciPy工具包内的find_peaks函数自动寻找局部极大值点作为候选R波的位置： ```python threshold = max(smoothed_square_ecg)*0.75 peaks,_ = find_peaks(smoothed_square_ecg, height=threshold,distance=int(round(fs*0.3))) r_peaks_time = peaks/fs print(f"找到 {len(peaks)} 个疑似R波:") for i,t in enumerate(r_peaks_time): print(f"{i}: 时间={t:.3f} 秒 ") fig, ax = plt.subplots(figsize=(10,6)) ax.plot(time, ECG,label="原始信号") ax.scatter(r_peaks_time, ECG[peaks], color='red', label="预测R波", zorder=5,s=100,alpha=.9) ax.legend(loc='upper right') plt.title("最终识别出来的R波标记") plt.xlabel("时间 (秒)") plt.ylabel("电压(mV)") plt.grid(True) plt.tight_layout() plt.show() ``` 此方法能够较为准确地标记出大部分正常情况下的R波位置[^1]。