李志强 林永武 李晓东
(1.佛山市质量计量监督检测中心 广东省佛山市 528200 2.佛山科学技术学院 广东省佛山市 528225)
数字心电图机是医学上常用来记录和显示心脏活动时产生的生理电信号的重要设备。随着人们生活水平的提高,心血管方面的疾病逐渐增多,心电图作为一种心血管疾病诊断的最有效手段在临床医疗中日益广泛应用。数字心电图机利用ECG 电极从人体表记录心电活动变化形成心电图形,医生通过分析心电波形的变化作为心血管疾病的辅助诊疗手段。为了保证心电图机的准确可靠工作,国家规定需要定期对心电图机进行检测。根据国家检定规程(JJG 1041-2008)[1]的要求,检测时通过对心电图机各电极馈入标准心电信号,并将心电图机输出的各导联的信号打印在心电图记录纸上,检测人员对打印的心电图进行测量和分析,然后给出检测结果。
但是,人工对打印的心电图进行检测是一项非常繁重和复杂的工作。需要测量数据涉及I,II,III,aVR,aVF,aVL,V1 ~V6 等12个导联[1],需要测量的数据项包括数百项之多。例如,其中一个导联的幅度-时间参数测量结果就包括幅度测量误差A1 ~A10,时间间隔测量误差T1 ~T11 等二十多项数据。由于测量和计算量都很大,给检测工作人员带来了繁重的工作压力,工作效率很低。虽然在实际中,人们经常只测量部分导联的数据来减少工作量,但工作量仍然巨大,并且人为因素也可能带来测量结果的不可靠。
随着计算机图像处理和图像测量技术的发展和应用,利用计算机技术对打印的心电图图像进行自动测量分析是可以实现的。本文以下将讨论基于图像的ECG 波形检测分析系统开发的关键技术问题。
基于图像的ECG 波形检测分析系统的核心工作流程如图1 所示。利用计算机图像处理技术实现对打印的心电图像进行分析测量,首先需要将打印的心电图扫描转换为数字ECG 图像,然后通过图像处理和分析技术对心电图像进行处理、分析及测量[10-12],再根据检定规程要求测量出相关数据项并给出检测结果和检测报告。
通过工作流程分析,基于图像的ECG 波形检测分析系统的总体模块结构如图2 所示。其中系统核心功能模块包括:图像处理模块,图像标定模块,ECG特征波定位和检测模块,数据处理分析模块,报告生成等核心模块。所涉及的关键技术包括:ECG 图像处理和标定技术,ECG 特征波形提取、定位和测量技术等。
图1:ECG 波形检测分析系统的主要工作流程
图2:ECG 波形检测分析系统的主要模块
图3:ECG 波形分割结果
在基于图像的ECG 波形检测分析系统,其核心任务就是实现从ECG 图像中提取ECG 波形,并对波形进行分析测量。图像处理模块实现ECG 波形提取,图像参数标定模块解决图像基础参数测量问题,特征波定位模块实现各类特征波形的快速准确定位,数据测量模块主要解决波形特征点的幅度、时间等数据的测量。数据处理和分析模块主要根据检定规程完成相关数据的计算,报告生成模块主要实现检测报告的生成和打印功能。其他功能主要包括:检测数据查询功能,系统数据配置功能,用户管理功能和报告审核功能等。因此,在本系统中,需要解决的最核心技术就是ECG 波形提取和特征波形的准确定位等技术。
图4:ECG 波形细化后的结果
图5:提取的网格图像及投影直方图
图6:ECG 波形定位及检测结果
在ECG 图像中ECG 波形是以图像的形式存在,ECG 波形的提取就是要将ECG 波形图像从图像中定位并分割出来。为准确地定位分割ECG 图形,首先需要对图像进行预处理,通过图像预处理实现对ECG 图像进行无失真的去噪增强处理(一般的滤波可能会损失图像的部分边沿信息),然后实现ECG 图像分割提取。考虑到坐标纸的底色及打印波形的颜色,为了能准确可靠地从CEG图形从图像中分割出来,需要将图像从RGB 模型空间转换到HSV空间,在HSV 空间进行ECG 波形分割。RGB 颜色空间到HSV 颜色空间的转换方法如下:
由于HSV 模型能很好地反映彩色图像的色彩、饱和度及亮度关系,对彩色图像的分割效果较好。此外,为提高图像边缘的二值化精度,必要时需要根据相关模型确定图像的亚像素边缘[7,8]。本系统为了简化,根据原图像的分辨率大小,将原图像放大2-5 倍,增加的像素采用双线性插值方式计算,通过分析实际效果可满足测量需要。为避免亚像素边缘不平滑,避免波形出现断裂,在此过程中增加了形态学闭运算。ECG 分割效果如图3 所示。
在ECG 波形图像分割提取后,实际上不是单像素宽度。因此为定位的准确,需要对ECG 波形图像进行中心骨架提取即图像细化,为了确保细化时不影响波形的定位精度采用对称的细化模板。从而保证特征波形位置不改变(限于篇幅细化算法不在此讨论)。细化后结果如图4 所示。
对ECG 图像的参数标定通常需要标定板来实现,但本系统的ECG 图像是打印在坐标纸上,坐标纸上的坐标线成为可用的标定依据。为了准确标定图像参数,我们需要提取ECG 图像中的坐标网格并进行细化处理。在实际中,ECG 图像通过扫描得到,图像一般不必进行几何矫正(在必要时进行倾斜矫正)。图像中坐标线构成的每一个小格为1mm,通过连续测量N 个小格的图像宽度,就可以计算每mm的像素个数,这将成为后续ECG波形测量的基础。图5 中(a)为提取的ECG 网格图像;(b)为细化后的网格图像;(c)为其网格投影直方图。设连续记录的直方图峰值个数为m,其对应的像素个数为k。则图像中像素密度为:
在ECG 波形中,ECG 特征波包括QRS 波,T 波,P 波等。尽管实际ECG 波形比较复杂,但各特征波形比较特征明显,特别是幅度最大的R 波。此外,ECG 中各种特征波形在时间上存在固定的顺延关系,因此波形定位时,只需先定位某些特征波(如R 波或QRS 波群),就可以快速定位其他特征波形。
目前,关于QRS 波和R 波的定位有很多比较成熟的算法[4-6],但这些算法主要是针对心电图机输出的ECG 波形信号而言的,不能直接运用于ECG 波形图像。尽管可以将提取和细化后的ECG 波形图像转化为一维的波形,但转化时波形的幅度值会出现非单值的情形,需要做特殊的处理。由于本系统是在心电图检测环境下,实际波形受环境干扰比较小,且R 波具有幅度最大的特点,因此,以此可快速定位R 波或QRS 波群。在定位R 波或QRS 波群后,就可以进行分段快速模板匹配算法,定位其他特征波形。
模板匹配定位算法[9]的最大特点在于匹配定位准确,适合于ECG 特征波形的准确定位。分段模板匹配算法是将心电图不同特征波图像作为模板进行分段匹配搜索。由于模板比较大,搜索过程中计算量很大,匹配速度较慢,因此,预先快速定位R 波峰值位置后,可以限定匹配的范围,从而大大节省模板的匹配计算量。
表1:测试结果与人工检测结果对比
在图像模板匹配算法中,设选定的特征图模板可表示为一个二维矩阵:匹配时以此模版为基础,通过对模板进行伸缩变换,得到一系列模版图像然后用这些模版分别于被检测图像区域进行匹配计算。匹配时,将模板在ECG 图像上平移滑动,分别计算其误差值或相关度。最后找到误差值最小的位置即为匹配的特征波形位置。误差的计算采用下列公式:
考虑到ECG 图形的特点,在计算模板与子图对象的匹配度时,实际只需要匹配心电图曲线上的点,因此,二维图像匹配问题可转化为一维图形的匹配度的计算问题。设一维模板为:
此时,伸缩模板可表示为:
其中,a 为幅度伸缩因子,b 为尺度伸缩因子。此时匹配误差的计算如下:
当E(i)取得最小值时,i 的值即为特征波形匹配的位置。
采用本文技术和算法我们开发了基于图像的ECG 波形检测分析系统,图6 为幅度-时间参数测量模块的检测结果。实验表明本系统能成功地检测出ECG 各导联的相关参数。为了验证算法结果的准确性和可行性,本系统检测的数据委托ECG 专业分析机构进行了检测数据和结果的对比分析(如表1 所示)。比对结果表明:本系统定位测量的数据结果与人工测量数据结果一致。
从检测结果来看,采用本文所讨论的技术或算法得到的检测数据与专业人员测试的结果一致。可见,采用以上研究技术开发的基于图像的ECG 波形检测分析系统可以实现了对ECG 各波形的准确定位,测量精度也达到了国家检定规程的要求。因此采用本系统可以在很大程度上代替相关检测人员完成相关的检测工作,使相关人员从繁重复杂的工作中解脱出来,并避免了人为因素导致的检测不准确的风险。