陈茜 欧阳绳武 马新宇 谢泉
关键词: 心电信号; 医疗诊断; 数字滤波算法; 整形IIR滤波器; 多阶反馈; 原始信号
中图分类号: TN911.72?34 文献标识码: A 文章编号: 1004?373X(2019)04?0045?04
Digital filtering algorithm optimization for ECG signal acquisition
CHEN Qian, OUYANG Shengwu, MA Xinyu, XIE Quan
(College of Big Data and Information Engineering, Guizhou University, Guiyang 550025, China)
Abstract: There exist various random noises in ECG signals as the actual environmental condition of ECG signal detection is not ideal. These noises make ECG signals fuzzy and seriously affect the reliability of medical diagnosis. The current problem that needs to be optimized is how to better suppress the noises that seriously affect the ECG signal detection waveforms. Though the filtering algorithm in this respect has got a certain achievements, the detail characteristics of ECG signals are not handled carefully enough. Therefore, a design method of a multi?order feedback shaping IIR filter is researched in this paper on the basis of analyzing the IIR digital filtering and IIR direct design method from the digital filtering perspective, so as to improve the feature point definition of ECG signals. The simulation results show that the multi?order feedback shaping IIR filter can retain the characteristics of original signals more obviously and clearly, which is beneficial to the latter diagnosis.
Keywords: ECG signal; medical diagnosis; digital filtering algorithm; shaping IIR filter; multi?order feedback; original signal
当今社会正迅速发展,同时对人们工作要求越来越高,大家的压力也越来越大,导致到了一定年纪后心脏病发作比例越来越高[1]。同时追求更高生活质量的人类越发看重医疗健康,但是医院的繁琐医疗程序与高成本的医疗检查并不能满足大部分患者的需求,所以促使了便携性家庭形式的医疗产品迅速发展,家庭式的心电监护仪器就是其中一种。如何从低成本、高效率、高准确度上考虑医疗仪器是现在需要研究的问题[2?3]。在临床治疗分析病患病因时,心电信号(ECG)的监测与波形分析占据很重要的位置。ECG包含了大量的心脏疾病信息,能帮助了解心脏临床诊断,包括快速定位功能与状况、帮助诊测和判断心血管疾病、提高对治疗方法的有效性判断。
心电滤波常用的是IIR滤波器和FIR滤波器[4]。IIR滤波器的特点是在同阶数的算法中有较好的滤波能力,但对于滤波器的相位特性没有很好的控制效果,在自适应IIR滤波器的研究中发现,在反馈回路时因为极点有可能散落到稳定区域外,造成滤波器的不稳定[5]。FIR滤波器的首要优点是因为没有设计反馈回路,不存在不稳定的问题;而且,FIR滤波器有准确的线性相位特性,同时幅度特性任意选取。由此可见,相比IIR滤波器,FIR滤波在稳定和线性相位特性有着显著优势。其他优点还包含:线性的设计结构、对硬件要求低、滤波器过渡过程具有有限区间。不过FIR滤波器阶次很高,因此造成了延迟过长。IIR与FIR复合自适应滤波器克服了各自的缺点,能够在低配置硬件条件下得到比较真实可视化的心电图形,但是在一些心电信号细节特征上却没能做到很好的优化[6]。
心电信号细节特征主要包括T,P,Q,R,S波峰细节特征。本文从数字滤波角度出发,研究了一种多阶反馈式整形IIR滤波器设计方法,提高心电信号特征点清晰度。
递归是IIR滤波器的特性,其输出结果不仅仅受到输入数据的影响,而且还与之前的输出结果有关联,它以FIR一种经典的线性滤波器为基础[7]。IIR滤波器具有递归结构,在相同的滤波器结构中, IIR滤波器比FIR滤波器具有更好的滤波效果。由于IIR滤波器的设计比FIR滤波器的设计更加复杂,不仅要考虑信号的衰减系数,而且还要考虑IIR滤波器带有反馈环。通过滤波器稳定理论可以发现它是不稳定的。整个滤波器的稳定性至关重要,以此目的来约束极点位置:首先必须是在单位圆的域内,然后由相应的极点和零点的设计相互抵消[8]。虽然软件滤波器在理论上可以达到任意衰减系数,但在工程实现中还必须考虑到计算的复杂度,处理器的字长。在信号处理中采用传递函数描述一个信号的转换过程,其公式为:
[H(z)=b0+b1z-1+…+bNz-Na0+a1z-1+…+aMz-M] (1)
因为上述IIR数字滤波器的计算量相对比较大,故而不能满足对实时输出心电波形高要求的环境。在20世纪末提出一个由约束零极点位置设计一个函数的变量整形滤波器的传递函数系数。由于系数为整数,所以不需要浮点运算,使计算大大降低,这个设计可满足本文的心电图数字滤波器[9]。
整数型数字滤波器的设计方法是,根据需要滤除噪声频率来计算单位圆上每个位置的频率,在这些位置上的布零和重叠,并与频率极点配置的频率,让极点零点互相偏移[10?11]。在该域中,每个信号被压缩到单位圆的边缘,并用边缘上的每个点表示对应的数字角频率。每个点具有零点的距离相乘作为分子,再与多个极点距离的分母相比,就得到频率在传递函数的过渡系数,这也是一个递归滤波器结构[12]。
2.1 IIR整形低通滤波器的设计实现
z域中的数字角频率:
[Ω=1fad·fs·2πn, n=0,1,2,…] (2)
式中:[fad]为采样频率;[fs]为当前频率。需要滤除的频率为[fc=50] Hz,频率采样为[fad=480] Hz时,对应的[Ωc=15πn],其中[n=0,1,2,…]。
为了最大幅度[fad=480 Hz] 消除50 Hz噪声。在z平面的单位圆中放置n个零点(n=0,1,2,…)。但这时又有一个问题,即0 Hz与50 Hz在(1,0)处重合。为了避免频率重合,在(1,0)处放置一个极点与零点相互抵消。同时为了使传递函数稳定且传递函数的分子为整数形式,又必须在原点(0,0)处放置极点。
从图1可以看到小圈是零点、小叉为极点,单位圆上显示的是50 Hz的数字角频率0.2 πn,n=0,1,2,…,9,所以均匀分布了10个零点。由于(1,0)这个点是0 Hz的数字角速度所在位置,而且0 Hz是通带频率。为了不滤除0 Hz,在(1,0)处放置极点与零点相互抵消。最后为了使得整个传递函数稳定,在单位圆内部的原点处放置9个极点。这样传递函数可以写成:
[H(jΩ)=(ej0-ejΩ)(ej15π-ejΩ)(ej25π-ejΩ)…(ej35π-ejΩ)(ej0-ejΩ)ej9Ω =ej10Ω-1ej10Ω-ej9Ω] (3)
令[z=ejΩ],得:
[H(z)=z10-1z10-z9=1-z-101-z-1=Y(N)X(N)] (4)
式(5)为二阶低通数字滤波器,增加滤波器阶数的目的是增加低通滤波器滚降特性,其传递函数为:
[H(z)=(1-z-10)2(1-z-1)2=1-2z-10+z-201-2z-1+z-2=Y(N)X(N)] (5)
得到差分方程为:
[y(n)=2y(n-1)-y(n-2)+x(n)- 2x(n-10)+x(n-20)] (6)
通过差分方程实现低通滤波器Java源代码,如下:
private static int OyLowPassFilter1(int data) {
int y0;
x2[n2 + 21] = data;
x2[n2] = x2[n2 + 21];
y0 = (y12 << 1) ? y22 + x2[n2] ? (x2[n2 + 10] << 1) + x2[n2 + 20];
y22 = y12;
y12 = y0;
y0 /= 100;
if (??n2 < 0)
n2 = 20;
return y0;
}
private static int OyLowPassFilter2(int data) {
int y0;
x1[n1 + 11] = data;
x1[n1] = x1[n1 + 11];
y0 = y11 + x1[n1] ? x1[n1 + 10];
y11 = y0;
y0 = y0 / 40;
if (??n1 < 0)
n1 = 10;
return y0;
}
2.2 IIR整形高通滤波器的设计实现
首先设计一个低通滤波器,因为截止频率为2 Hz,需要滤除的频率为[fc=2] Hz,采样频率为[fad=480] Hz时,对应的[Ωc≈1120πn],其中[n=0,1,2,…] 。
根据整数型滤波器设计原则:
[Hlp(z)=Y(z)X(z)=1-z-2401-z-1] (7)
利用全通滤波器与低通滤波器作差设计高通滤波器可得[10]:
[Hhp(z)=P(z)X(z)=z-120-Hlp(z)240] (8)
而高通濾波器的传递函数化为差分方程:
[y(n)=y(n-120)-y(n-1)+x(n)-y(n-240)240] (9)
实现高通滤波器Java源代码,如下:
private static int OyHighPassFilter(int data) {
int y0, Z;
x[n] = x[n + 241] = data;
y0 = y1 + x[n] ? x[n + 240];
y1 = y0;
if (??n < 0) {
n = 240;
}
Z = (x[n + 70] ? (y0 / 240)) ? 250;
// if(Z<0) Z=0;
return Z;
}
3 多阶反馈式整形IIR滤波器的设计
整形IIR滤波器处理后的ECG信号会产生一些失真,改进整形IIR滤波器设计会得到更完整的心电信号。整形IIR滤波器会滤除很多信号细节特征(如QRS特征点),所以为了保留更多信号细节特征,提出一种多阶反馈式整形IIR滤波器[4]如图2所示。
反馈式整形IIR滤波设计步骤:
1) 心电信号In(n)通过IIR整形滤波器滤波得到信号Out1(n);
2) 将原始心电信号减去滤波后信号得到差值信号Sub(n):
[Sub(n)=In(n)-Out1(n)] (10)
3) 将差值信号Sub(n)继续进行IIR整形滤波得到的信号Out2(n),加到第一次滤波后的信号之上得到滤波信号Out(n):
[Out(n)=Out2(n)+Out1(n)] (11)
4) 将原始心电信号减去Out(n)信号得到差值信号Sub1(n):
[Sub1(n)=In(n)-Out(n)] (12)
5) 将差值信号Sub1(n)再进行IIR整形滤波得到信号Out3(n),加到Out(n)信号之上得到滤波信号Out4(n):
[Out4(n)=Out3(n)+Out(n)] (13)
6) 以此循环。
以图3结构进行心电信号处理,结果如图4所示。
从MIT/BIH标准心律失常数据库中取一组标准的ECG真实测试数据,用来验证所提出的多阶反馈式整形滤波去噪性能。图4中a)为有工频干扰的原始心电信号In(n);图4中b)是经过未改进低通滤波后的心电信号Out1(n);图4中e)是经过改进低通滤波后的心电信号Out(n),以及其他滤波过程。
结合之前的分析,从图4中可以看出,Out2(n)是被整形滤波器滤出过多的部分,导致了Out1(n)某些细节特征点消失,变得平滑,不利于医生做详细观察,因此Out1(n)加上Out2(n)这部分才是还原到较为真实的心电波形。所设计的低通滤波器可以很好地滤除掉工频干扰,并且心电图信号很好地保留了细节特征。多阶反馈式整形IIR滤波器使细节可以达到更好的效果。
本文从数字滤波角度出发,基于分析IIR数字滤波和IIR直接设计方法,研究一种多阶反馈式整形IIR滤波器设计方法,以提高心电信号特征点清晰度。由仿真结果可以看出,多阶反馈式整形IIR滤波器能够更加明显清晰地保留原始信号特征,有利于后期诊断。
参考文献
[1] 蒙俊甫.基于嵌入式的医疗监护PDA无线终端设计[D].成都:成都理工大学,2010.
MENG Junfu. The design of wireless medical PDA based on embedded system [D]. Chengdu: Chengdu University of Technology, 2010.
[2] 谷蓉.心电信号处理算法与应用的研究[D].南京:南京邮电大学,2014.
GU Rong. Research of ECG processing algorithm and its application [D]. Nanjing: Nanjing University of Posts and Telecommunications, 2014.
[3] 石进平.心电图机计量指标的临床意义[J].大众标准化,2009(z2):41?43.
SHI Jinping. Clinical significance of measuring index of electrocardiograph [J]. Popular standardization, 2009(S2): 41?43.
[4] 欧阳绳武.便携式移动医疗心电诊断系统研究与设计[D].贵阳:贵州大学,2016.
OUYANG Shengwu. Research and design of portable mobile medical ECG diagnosis system [D]. Guiyang: Guizhou University, 2016.
[5] 郑蕾.基于小波变换的心音信号分析方法的研究[D].兰州:兰州理工大学,2010.
ZHENG Lei. Research of heart sound signal analysis method based on wavelet transform [D]. Lanzhou: Lanzhou University of Technology, 2010.
[6] 彭云.基于EEMD消除ECG基线漂移的新算法[J].信息系统工程,2014(8):128?129.
PENG Yun. A new algorithm based on EEMD to eliminate ECG baseline drift [J]. China CIO news, 2014(8): 128?129.
[7] 宋春丽.怎样识读MIT?BIH中的心电信号[J].科技资讯,2010(9):27.
SONG Chunli. How to read the ECG signal in MIT?BIH [J]. Science & technology information, 2010(9): 27.
[8] 汪学东,向晋涛,熊军.利用时间RR间期散点图及逆向技术分析和诊断房性并行心律[J].中国心脏起搏与心电生理杂志,2013,27(1):79?82.
WANG Xuedong, XIANG Jintao, XIONG Jun. Analysis and diagnosis of atrial concurrent arrhythmia using time RR interval diagram and reverse technique [J].Chinese journal of cardiac pacing and electrophysiology, 2013, 27(1): 79?82.
[9] KANG W S, YUN S, CHO K. ECG denoise method based on wavelet function learning [C]// Proceedings of IEEE sensors. [S.l.]: IEEE, 2012: 1?4.
[10] 毛玲,张国敏,孙即祥.心电图ST段形态分析方法研究[J].信号处理,2009,25(9):1360?1365.
MAO Ling, ZHANG Guomin, SUN Jixiang. Research on shape analysis of ST segments in ECG signal [J]. Signal processing, 2009, 25(9): 1360?1365.
[11] 宋晋忠,严洪,姚宇华,等.基于心电图ST?T段的心肌缺血检测方法研究进展[J].航天医学与医学工程,2011,24(2):146?150.
SONG Jinzhong, YAN Hong, YAO Yuhua, et al. Research progress in detecting methods for myocardial ischemia based on electrocardiogram ST?T complex [J]. Space medicine & medical engineering, 2011, 24(2): 146?150.
[12] MUKHERJEE A, GHOSH K K. An efficient wavelet analysis for ECG signal processing [C]// Proceedings of International Conference on Informatics, Electronics & Vision. Dhaka: IEEE, 2012: 411?415.