基于Harris的纸质心电图波形提取方法

杨佳倩　孙连山　杨晓立

摘 要：ECG是心电图机从体表记录心脏每一周期所产生的电活动变化图形。纸质心电图广泛应用于病案管理、取证等方面，针对临床诊断中纸质心电图不便于保存的问题，文章采用Harris角点检测算法与阈值法相结合，提取得到更加精准的纸质心电图的波形。

关键词：纸质心电图；Harris角点；心电图波形；阈值法

心电图是临床检测患者的心脏问题的重要依据。将过去医院产生的心电图纸转变为便于传输、存储和管理的数字化对象，有利于学者开展基于大数据的医学研究，增加临床诊断的效率和可靠性。

图像处理是心电图纸数字化的主流方法。常见的方法包括预测搜索算法、阈值变换法、Snake模型法、基于K-means算法的波形提取方法等。这些现有算法精确度不高且实验耗时，亟待改进。本文借鉴现有算法的优点，提出阈值法与Harris角点检测算法相结合的方式提取心电波形，提高纸质心电图数字化的精度。

1 心电图预处理

对纸质心电图进行数字化的第一步是将其扫描为心电图图像，并对其进行必要的预处理，作为后续提取心电图波形的基础。心电图像预处理包括缩小图片、倾斜检测及矫正以及导联范围截取等。

1.1 缩小图片

由于测试图片分辨率太高，图片过大，所以采用opencv自带函数resize中的最邻近插值法，按一定比例将图像缩小，以便提高角点检测的准确性。

1.2 心电图倾斜检测及矫正

扫描心电图获得的图像通常是倾斜的，有必要采用Hough变换对其进行倾斜矫正。通过Hough变换，就可以将图像空间中直线的检测问题转化为参数空间中对点的检测问题。

1.3 心电图导联范围自动定位

心电图由于其背景与网格特征差异明显，先用阈值法将网格去除，再进行水平投影。针对不同的导联采用水平投影法自动定位，将每行中像素点最多的前12行的水平像素都提取出来，按照顺序取平均值，即为分割图像的宽。自动定位波形范围的目的是保证在界面输入导联即可直接得到波形结果和数字信号。对于图1所示的格式纸质心电图，由于此类型的心电图包括很多字符，自动定位很难准确分割，因此采用手动截图的方法进行波形范围选取。图1中黑框部分为手动框选的范围，将其保存为待处理图片。

2 基于Harris的纸质心电图波形提取方法

本文提出基于Harris算法和阈值法相结合的纸质心电图波形提取方法。首先，用阈值法将图像的背景网格全部去除。其次，采用Harris算法对原始图像检测这些拐点，再将这些拐点与阈值法提取的断点重新排列，连接起来，拟合心电图波形。

图1 格式心电图手动框选

2.1 阈值法去网格

首先采用阈值法，阈值法是一种基于区域的图像分割技术。它适用于目标和背景占据不同灰度级范围的图像，不仅可以极大地压缩数据量，而且也大大简化了分析和处理步骤，因此在很多情况下，是进行图像分析、特征提取与模式识别之前的必要的图像预处理过程。

图像阈值化的目的是要按照灰度级，对像素集合进行一个划分，得到的每个子集形成一个与现实景物相对应的区域，各个区域内部具有一致的属性，而相邻区域不具有这种一致属性。这样的划分可以通过从灰度级出发选取一个或多个阈值来实现。

为按照图像灰度的特性将图像分成背景和目标两部分，由于背景网格与心电波形灰度相近，通过最佳阈值的选取，对原图进行处理，将背景网格全部去除，该阈值的选取需要通过实验来确定，对于给定的图像通过迭代法确定最佳阈值。

式（1）中：T0代表初始阈值，Zmax和Zmin分别表示最大灰度值和最小灰度值。式（2）中：小于T0的部分为TA，大于T0的部分为TB，将T1作为新的阈值代替T0，重复以上过程如此迭代，直至TK收敛，即TK+1=TK。阈值分割之后的心电波形是不连续的散点，保存这些像素点。

2.2 基于Harris算法提取波形

角点是图像的一个重要的局部特征，它保留了图像中重要的特征信息。在计算机视觉中，对于角点的定义有不同的表述。目前的角点检测可分为两类，基于图像边缘的方法和基于图像灰度的方法，前者往往需要对图像边缘进行编码，这很大程度上依赖于图像边缘的提取，具有相当大的难度和计算量，因此这种方法的使用范围很小，基于图像灰度的方法主要通过曲率及梯度来检测角点，但是可能会出现漏检测角点和伪角点等问题。

本实验主要采用Harris角点检测算子对原始波形进行处理，其具有计算简单、提取点的特征均匀且合理、稳定等特点。波形提取只要可以准确识别波峰波谷等像素点即可，所以这里采用Harris算子，设定角点数量，将心电波形的波峰波谷及拐点检测出来，并将角点的像素保存下来。通过定义角点响应函数公式（3）：

式（3）中的CRF可以通过式（4）和式（5）求出，其中K取0.04，式（4）表示矩陣M的行列式，式（5）表示主对角线数之和，CRF为大数值整数时则为角点，CRF为大数值负数时则为边缘，绝对值CRF为小数时则为平坦区，Harris角点检测步骤如下。

Step1：利用Sobel计算像素点在x，y方向的梯度值Ix和Iy，计算出Ix2，Iy2，IxIy。

Step2：得到局部自相关矩阵

M =

对图像进行高斯滤波，得到新的M，本实验高斯窗口选9×9。

Step3：计算原图像中每个像素的角点响应函数CRF的值，对每个像素点，利用一个以它为中心的9×9窗口进行判断，如果中心像素点的角点响应值是这个窗口的极大值，就提取中心像素为角点。

Step4：设定阈值，CRF值大于阈值的点作为最终选取的角点。

对于原图用Harris角点检测算法将角点全部识别出来，关于阈值的确定，取决于波形提取图片的导联范围，对于非最后一行心电图的某一导联，角点数取30～50，即可将波峰波谷的角点全部识别出来，关于角点数的选择，需要通过欧式距离进行结果评估，确定最佳角点数。对于边缘有时会将背景网格误检出来，需要将边缘的误选角点去除，再将角点的像素保存下来。

最后将阈值法和Harris角点检测算法保存的像素进行从小到大排列，连接这些像素点，进行心电波形复原，并与原心电波形进行对比。基本全部覆盖原始波形。

3 实验及结果评估

为验证本文算法的有效性，对测试图像进行实验。本实验环境为Matlab 2014和Visual Studio 2015，并配置Opencv2.4.6，操作系统为64位Windows10，CPU为2.70 GHz，安装内存为4 GB。

实验数据来自北京大学人民医院心内科，一共包括100张纸质心电图报告单，其中有55张图1格式的纸质心电图，其走纸速度为25 mm/s，定标电压为10 mm/MV，横线每格0.04 s，纵线每格0.1 MV。

本文选取5组测试图片，分别为标准肢体导联两组，加压单极肢体导联和胸导联两组。对于测试图片，角点检测的阈值设定分别为30，35，40，45，50，先通过欧式距离评估结果确定角点最佳阈值，再和其他现有算法提取的结果进行对比。评估时首先将提取到心电波形的像素根据以下方法转化为数字信号。

式（6）中：Pn表示心电图纸每个格子包括像素点的个数，t表示每一小格代表的时间0.04s，v表示心电图机的走纸速度为25 mm/s，DPI表示分辨率，可通过式（7）求得，式（7）中：Ph表示水平像素点的个数，Pv表示垂直像素点的个数，size表示电脑屏幕尺寸，本实验的DPI结果为300，实验结果需要将原始波形通过插值法降采样为实验结果的像素数，通过计算可以得出Pn为11.8。式（8）中：Fu表示每个像素点与电压的比例，一维信号为每个竖直像素点与Fu的乘积，心电图纸上每个小正方形代表的信号电强度为0.1 MV，每个像素点与电压的比例Fu可由公式（8）得出。计算结果为8.47 μV/pixel，将每像素的竖直坐标与Fu相乘即可得到一维数字信号。

使用非线性插值技术将原始心电信号重新采样缩放到波形数字化所需要的采样率和振幅分辨率，将提取的一维信号跟降采样后的一维信号进行对比，将欧式距离和余弦距离相似性度量法相结合，进行评估实验结果。欧氏距离（Euclid Distance，ED）是一个通常采用的距离定义，它是在m维空间中两个点之间的真实距离。在二维和三维空间中的欧式距离的就是两点之间的距离，本实验只对一维信号进行评估，其公式是：

式（9）中x为原始波形的竖直坐标转化后的数字信号，欧氏距离可看作信号的相似程度，距离越近就越相似，准确率就越高。单独使用欧氏距离度量并不足以说明结果准确性，需要引入余弦距离进行再次评估。余弦距离是衡量空间内向量的夹角，在本实验场景中，对于每一个波形即n维空间内的向量，通过余弦距离评估两个波形在形状上的相似度，夹角余弦越大表示两个向量的夹角越小，两条波形在形状上越相似，式（10）表示余弦距离的公式：

（10）

通过欧式距离对波形提取结果进行评估找出最佳角点阈值。对于5组测试图片用本文方法与其他算法进行实验，通过欧氏距离和余弦距离相似性度量（CD）相結合方式进行结果评估，根据结果对比可见，本文方法与其他实验方法相比，其欧式距离最小，余弦相似度最接近1，准确率均高于其他方法，心电波形提取的准确率得到明显提高。

4 结语

纸质心电图由纵横交错的背景网格、心电波形曲线、导联字符及校准电压组成。纸质心电图波形提取是通过图像处理技术将心电图波形提取出来并以数字信号形式保存的过程。但因存在很多技术问题，导致提取结果精度较低，实际应用并不广泛。所以本文提出了结合阈值法及Harris角点检测算法的方法，提取出高精度的波形，使患者信息保存的更加完整，更完全。但提取结果还是与原纸质心电图有细微差别，医学中极小的差别也会为临床诊断带来影响，所以本文采用的算法的精度还有待提高。