宋盟春 熊念
1 广东省医疗器械质量监督检验所 (广州 510663)
2 暨南大学信息科学技术学院 (广州 510632)
胎儿监护及围产期监护是产科中一项重要的监护措施。常用的胎儿监护方法有:超声多普勒监护、胎心音监护和胎心电监护。其中,胎儿心电监护因其从母体腹部提取胎儿心电信号,对母体和胎儿无创,特别是能同时观察胎儿心电图,其优点明显,是目前产科常用的诊断方法。
虽然腹部心电信号提取方便,但是,母体腹部信号复杂,信噪比低,特别是一般母体心电信号比胎儿心电信号大10~20 倍,且胎儿心电约有10%~30%与母体心电重合,形成复杂的QRS 波群,并且胎儿心电的频谱与母体心电的频谱重叠。这些因素给胎儿心电的检测带来了比较大的困难。解决这个问题常用方法有:匹配滤波法、自适应滤波法等[1~3]。本系统采用匹配滤波法,在MATLAB6.5 平台上实现[4~5]。
心电信号的频带响应一般只有0.05~100Hz,从母体腹腔提取的胎儿心电信号比较弱,而且含有很多噪声。
从母体腹腔体表所采集的综合心电信号由纯净的心电信号和随机的噪声信号两部分组成,可表示为:
其中S(t)表示纯净的心电信号(包括母亲心电信号和胎儿心电信号);N(t)表示随机的噪声信号,该躁声信号可表示如下:
其中N1(t) 为50Hz 工频干扰信号,N2(t)为高频干扰信号(>100Hz),N3(t)为低频信号(<0.05Hz),N4(t)表示基线漂移。所以X(t)可表示为:首先设计一个通带为100HZ 的低通数字滤波器,去除高频干扰信号N2(t),再使用50Hz,0.05HZ和100HZ 的数字陷波器,最后得到是有基线漂移的母亲胎儿心电信号,可表示如下:利用信号叠加平均法得到漂移基线N4(t),与(4)式相减得到只含有母亲心电和胎儿心电的信号S(t),其可以表示为:
S(t) 的 频 率 范 围 是0.05~100Hz,SM(t) 为母体的心电信号,SF(t)为胎儿的心电信号,其中SM(t)的频率范围主要集中在0.05~40Hz,而SF(t)的频率范围主要集中在0.05~70Hz[6]。可以使用匹配滤波的方法从S(t)中提取出模板,将S(t)与提取的模板相减得到胎儿心电信号SF(t)。
提取母体腹腔胎儿心电信号的程序设计主要由9 部分组成,其程序设计流程如图1。
信号预处理主要功能是滤除工频干扰和抑制基线漂移。信号预处理可以分为滤除高频信号,陷波器设计(主要是滤除工频干扰),抑制基线漂移4 部分。
(1)滤除高频信号
心电信号的频带响应只有0.05~100Hz,因此,需要设计一个低通滤波器将大于100Hz 的信号滤除。这里使用巴特沃思数字低通滤波器,其幅频响应图如下图2 所示。滤波后的信号为S1(t)。
图1. 信号处理流程
图2. 低通滤波幅频响应
图3. 50Hz陷波幅频响应
图4. 滤除高频信号后的波形
图5. 高通滤波幅频响应
图6. 母亲胎儿心电信号
(2)消除工频干扰
采集的信号中常混入50Hz 的工频干扰,必须加以抑制。这里使用50Hz 的数字陷波的方法将其滤除。选择使用有限冲击响应数字滤波器(FIR),调用Matlab 中函数fir2 和基于Parks-McClellan 算法的函数remez 来设计滤波器参数。其幅频响应图如图3 所示。滤除频率大于100Hz信号,消除工频干扰后的波形如图3,由于高频分量的消减,波形畸变较小,但有明显的基线漂移。设消除工频干扰的信号为S2(t)。
(3)抑制基线漂移
基线漂移的原因有多种,其中由于运动造成体表电极与皮肤之间细微的滑动和硬件电子元件热效应造成参数变化是主要原因。造成的基线漂移大部分是低频的信号,因此我们可以利用高通滤波器加以滤除。这里我们先进行高通滤波,滤除频率低于0.05HZ 的低频信号,再使用信号叠加平均法,最终抑制基线漂移。使用有限冲击响应数字滤波器(FIR),调用Matlab 中函数fir2 来设计高通滤波器参数,其幅频响应如图5 所示。
再利用公式
得到基线Baseline(n), 用S2(t)信号减去基线Baseline(n)得到只剩下母亲和胎儿心电的信号S3(t),如图6 所示。从图中可以看出基线已经被消除了,波形比较平稳。
(1) R 峰检测
R 峰值是指母亲PQRST 中的最大幅值点,在正常情况下整个周期中该点的幅值最大,因此,我们先设定一个参考值peak=180,R 峰值应该在大于peak 的范围。
将S2(t)信号值逐一跟peak 进行比较,当S2(t)>peak,同时又满足:
则将S2(t)设为R 峰,再将其与S2(t)后面连续15 个值进行比较,若S2(t)值为最大则确认S2(t)为R 峰。因为所处理的临床信号的采样频率为429,为了避免由于出现较大的干扰而影响峰值的判断,在找出一个R 峰后,在接下来的250个点内不再进行峰值检测判断。这样既可以避免峰值的重复检测,同时也可以提高检测效率,更易于进行实时处理。
(2) 模板的提取
在确定了R 峰后,就可以采取动态加权平均的方法来得到每一时刻的母亲心电图的模板,考虑到采样频率为429,由经验又可以判断出比较明显的QRST 波群范围为40 点左右,本系统对信号的处理每个周期使用了100 个点的数据(R峰前40 个数据点,R 峰后35 个数据点),在这个范围内求取母亲心电的模板信号。在找到每一个R 峰值后可以根据公式:
取得。这个模板信号随着检测过程的时间推移而动态改变。该方法提取的母亲心电信号模板如图7 所示。
在S3(t)信号中减去MODEstandard 便得到了胎儿心电信号。波形如图8 所示。
心率的计算公式为:
其中,N 为在某一段时间内得到的R 峰个数,M 是这段时间内的所有数据点的个数。
采用与寻找母亲心电R 峰同样的方法确认胎儿心电的R 峰值。而前面母亲心电的R 峰值已经找出,利用公式(10)可以求得母亲和胎儿心率,并显示出来,如图9 所示。
经处理后的心电信号波形可以动态显示出来,这里采用扫屏和滚屏两种显示方式,同时显示出母亲和胎儿混合信号,母亲心电模板信号,胎儿心电信号以及母亲和胎儿心率。另外还设计了用于信号分析的控件,如放大、缩小、暂停等。完整的操作和显示界面如图10 所示。
图7. 母亲心电信号模板
图8. 胎儿心电信号
图9. 母亲和胎儿心率图
图10. 操作和显示界面
本系统采用MATLAB 处理平台利用匹配滤 波方式对胎儿心电信号进行检测提取,运算效率高,处理速度快,编程简单,尤其是使用系统所带的信号处理工具中的滤波器设计,大大简化了信号调理的处理过程,节省了开发时间和系统搭建成本,具有较好的经济效益。同时,在胎儿心电信号波群的检测上还有待完善。
[1] 尹泽明 丁春利.精通MATLAB 6 清华大学出版社,2002
[2] 周辉 懂正宏。数字信号处理基础及MATLAB 实现。北京希望电子出版社,2006
[3] 王美华. 胎儿心电图的检测与分析. 桂林电子工业学报,1998;18(2)
[4] 周礼 梁 厚 叶大田.关于胎儿心电信号检测和处理的新方法研究.上海生物医学工程
[5] 莫玮等. 微机化胎儿心电图仪的研制. 中国医疗器械杂志,1997;21(1)
[6] Steven L. Horner and William M. Holls III Real-time signal processing of an intrauterine catheter’s fetal electrocardiogram Digital Signal Processing 13 (2003) 368–395