【基于生成对抗网络的QRS波检测算法分析7100字（论文）】

基于生成对抗网络的QRS波检测算法分析目录TOC\o"1-3"\h\u7099基于生成对抗网络的QRS波检测算法分析 1258391.1引言 142361.2检测算法理论模型 1326301.2.1数据介绍 118351.2.2算法理论基础 2137271.2.3总体框架 5206631.3实验平台与评价指标 6262441.1.1实验平台 683371.1.2评价指标 6117761.4心电图数据预处理 723311.4.1心电图去噪 728161.4.2数据集0-1变换及划分 9131511.5算法框架 10146021.5.1算法流程图 10207081.5.2生成网络的设计 12243771.5.2判别网络的设计 143191.5.3损失函数的设计 1420631.6实验结果及分析 15289651.6.1网络训练结果 15124811.6.2QRS波检测结果 16263721.6.3心率计算 17296701.6.4结果分析 181.1引言QRS波检测是处理ECG信号的关键步骤，对于后续心率(HR)计算和异常情况分析具有重要意义。当前，静态心电信号的QRS波检测精度越来越高，但在受噪声严重干扰或心电信号形态异常时，QRS波检测和HR计算精度仍难以达到医学检测标准。主要问题在于两个方面：ECG记录来自可穿戴或动态心电信号采集设备时，会出现伪影、噪声等干扰；各种心律不齐(比如LBBB、RBBB、PVC、窦性心动过速、窦性心电过缓等)患者的心电信号形态会出现明显的个体变异性的情况，此时传统QRS波检测算法往往无法获得高精度检测结果。针对动态心电信号，本章将运用生成对抗网络对ICBEB数据库的动态心电信号进行高精度的R波检测，结合高精度R波检测值实现心率的计算，为医学心电信号研究提供参考方案。1.2检测算法理论模型1.2.1数据介绍本章数据来自第八届国际生物医学工程与生物技术会议(ICBEB)所提供的心电数据库。ICBEB数据库中，原始数据是从心血管疾病患者中收集的2000组单导联心电信号记录，如表1.1所示。每组记录持续10s，采样率为500Hz，共计5000个点。考虑到ICBEB数据库官方的要求[54]，根据实际应用对每组原始心电信号数据截去前500ms和后500ms的记录，同时对2000组原始心电信号的R波波峰位置进行标注处理，将处理之后的数据作为本章算法的训练数据集[63]。表1.1ICBEB数据库具体情况Table1.1TheparametersofICBEB心电信号单导联信号时长/s10采样频率/Hz500数据维度(2000,5000)数据总量100000001.2.2算法理论基础Pix2Pix网络于2017年提出[64],主要创新在于有效结合了DCGAN与CGAN的原理，目标优化仍为CGAN的损失函数。CGAN的损失函数[65]如下所示:： (1.1)表示判别网络对真实样本和标签的评价结果，代表结果的概率。表示判别网络对真实数据和生成网络生成标签的评价结果，表示结果反面的概率。在公式(1.1)的基础上，一般采用L1或L2范数对损失函数进行优化，提升模型的适应度。L1范数方式如下所示： (1.2)L2范数方程为： (1.3)经过实验证明，针对大部分情况，以L1优化对抗损失函数时，模型的过拟合现象有所改善，只有某些特殊情况，L2优化效果才强于L1。因此，将式(1.1)和(1.2)结合，可以得到以L1优化的损失函数计算公式： (1.4)上式中，为超参量，当时表明不使用L1损失函数，Pix2pix的损失函数一般都是在CGAN的基础上以L1范数方程进行优化，计算公式如(1.4)所示[66]。在神经网络之中，若只存在线性转换，网络只能学习简单的映射关系，且输出还要进行线性转换的组合。激活函数作用就是进行非线性转换，以确保网络能够进行复杂映射关系的学习。生成对抗网络中包含了生成网络和判别网络，因此也需要激活函数。一般常用的激活函数有如下几种：①Sigmoid激活函数：该函数定义见(1.5)，图1.1为Sigmoid函数的示意图。 (1.5)从图中可以看到当z大于0会逐渐趋近于1，当小于0时会逐渐趋近0，因此Sigmoid函数经常用于二分类。但是它的缺点为：一是该函数的均值为0.5，不等于0。这样的情况会导致网络在训练过程中，神经网络就会出现全正或全负的反馈；二是对Sigmoid函数进行求导，如下所示： (1.6)该函数在z无穷大或非常小的时候，函数的导数会接近0，很容易出现梯度弥散的问题，因此当进行网络训练时，容易出现神经网络不能收敛训练失败的情况。图1.1Sigmoid函数的示意图Fig.1.1ThediagramofSigmoidfunction ②Tanh激活函数(双切正切函)：该函数为双曲函数中的一种，该函数的范围是[-1,1]，如图1.2所示，函数表达式如下所示：图1.2Tanh函数示意图Fig.1.2ThediagramofTanhfunction (1.7)从图中可以看出该函数的均值为零，从而优化了Sigmoid均值为0.5达不到0的情况，这就导致了一般情况下Tanh函数应用比Sigmoid函数更广泛。对Tanh函数进行求导，从这里可以看到当y非常大或非常小的时候，g(y)’会等于0，这样会导致权重更新很缓慢。Tanh函数的特点是能够很好的区分差别很大的特征，并且在过程中不会扩大这个特征效果。③ReLU激活函数：为了缓解神经网络梯度消失情况，Nair等人于2010首次运用了ReLU激活函数[67]，ReLU函数被称为非饱和修正线性函数，如图1.3所示，定义如式(1.8)。 (1.8)当输入k大于0时，输出不为0，但当输入k小于等于0时，输出结果为0，这样可以有效削弱输出与输入之间的依赖性，从而避免了过拟合现象的发生。但ReLU函数在训练过程中某些神经元可能会无法更新，从而丢失了数据的多样化。在具体应用中可以通过设置合理的学习率减少这种情况的发生概率。图1.3ReLU函数Fig.1.3ThediagramofReLUfunction为了克服ReLU激活函数的缺点，出现了LeakyReLU激活函数。在k小于0时，，其中是一个不为零的小量，这样输出不为0，改进了神经元不可逆死亡的问题。具体定义如式(1.9)所示。 (1.9)1.2.3总体框架本文算法流程如图1.4所示，可分为如下三个模块：(1)数据预处理模块。针对ICBEB数据库中ECG中存在的肌电干扰、工频干扰和基线漂移等噪声，采用滤波的方法进行滤除。(2)数据重整与网络分析模块。本文所采用的Pix2Pix网络结构是端到端输入和输出，即输入数据和输出数据维度需保持一致。ICBEB数据库中输入原始心电信号长度为5000，而输出结果为R波波峰位置，数据长度远小于5000，因此需对数据进行重整，进行0-1转换；随后，利用Pix2Pix网络对训练数据进行训练，并对剩余数据进行测试。(3)QRS波检测模块。将经过Pix2Pix网络输出的数据进行阈值设置和最优选择算法，得到R波波峰位置。图1.4QRS波检测算法流程图Fig.1.4FlowchartofQRSwavedetectionalgorithm1.3实验平台与评价指标1.1.1实验平台本文的实验都是在配备Intel(R)Core(TM)i7-8700k、1.7GHzCPU、16GB内存，NVIDIATitanXP显卡的服务器上进行，编程环境使用Python1.6，深度学习框架使用的是Tensorflow1.14。1.1.2评价指标为了验证算法的准确性，引入了FP(假阳性)、FN(假阴性)和TP指标。FP是指将其他干扰误认为是R波波峰位置的数量，即误检；FN是指遗漏的R波波峰位置的数量，即漏检。规定输出的R波波峰位置在标注的R波波峰位置左右75ms之内，即两个输出位置与标注位置坐标之差d≤37.5认为是正确的检测[68]。TP表示正确的检测数量。评分规则如式(1.10)和(1.11)所示： (1.10) (1.11)其中QRSscore表示具体每组心电信号QRS检测分数；QRSacc表示QRS波检测最终分数。每组持续时间10s的心电信号取5.5s-9.5s的心电信号计算心率[54]，具体评分规则如下所示： (1.12) (1.13)上式中HRref表示计算得到的心率，HRtest表示ICBEB数据库给出的参考心率，HRscore表示具体每组心电信号HR检测分数。引入灵敏度(Se)和阳性预测(+P)两个指标用于下文的比较。指标定义如下： (1.14) (1.15)1.4心电图数据预处理1.4.1心电图去噪心电信号是非线性的生理信号，由于心电信号本身的特点，容易受到多种噪声的干扰：环境噪声、采集设备产生的噪声、测试者运动产生的噪声，如测试者随意的呼吸动作或习惯性的皮肤摩擦等肢体活动等。总体而言，心电信号噪声主要包含三类：分别是肌电干扰、工频干扰、基线漂移[69]。(1)肌电干扰：当测试者在肌肉发生兴奋及收缩时会出现微弱电势，这个电势会使ECG波形幅值被提高。肌电干扰通常持续时间较短且频谱特性近似于白噪声。(2)工频干扰：心电采集设备在供电过程中会产生电磁干扰，且难以避免。由于中国通用电网电压频率是50Hz，所以工频干扰常常被认为是由测量设备导线引入的50Hz或50Hz整倍数噪声，且工频干扰的幅值一般是50%左右的心电信号波峰值。(3)基线干扰：基线干扰是一种低频干扰噪声，一般情况下主要来自测试者的呼吸和肢体活动，以及测量设备电极的材质问题。因为基线干扰噪声频率范围很窄，所以它容易混淆在心电信号中，常常导致心电信号失真严重，出现心电波形ST段呈现压低或抬高等情况，误导医生做出错误诊断。针对上述三类噪声，本文采用频域滤波的方法加以抑制。通常心电信号频率位于0.01Hz~100Hz[70]，一般集中于0.25Hz~40Hz。但肌电信号频率分布较宽，通常为40Hz至数万赫兹之间[71]，相对而言，心电信号位于低频段，因此利用低通滤波器滤除肌电干扰信号是一种可行的方案。工频干扰通常由50Hz正弦波及其谐波分量组成[72]，因此，采用窗函数法设计陷波器，即用低通滤波器叠加一个高通滤波器加以滤除。基线漂移噪声频率通常位于0.05Hz~2.0Hz之间，接近心电信号的频率[73]，则采用IIR数字滤波器加以抑制。三种滤波方法如图1.5所示。图1.6为ICBEB数据库ECG去噪对比图，(a)为原始ECG信号，(b)为去噪之后的ECG信号，显然，信号质量得到了一定程度的提升，而信号中的波峰波谷的幅度和位置滤波后并未发生明显改变。图1.5心电图去噪Fig.1.5ThedenoisingofECG(a)(b)图1.6ICBEB数据库去噪对比图。(a)原始心电信号；(b)去噪之后的心电信号Tab.1.6ThecomparisonchartofdenoisingforICBEB.(a)ThesignaloforiginalECG;(b)ThesignalofECGafterdenoising1.4.2数据集0-1变换及划分图1.70-1转换过程Fig.1.7Thetransformationprocessof0-1ICBEB数据库心电信号是连续且每组维度均是(1,5000)，其中1表示组数，5000表示每组的点数。而标注的心电信号中每组数据维度基本不一致，导致网络模型的输入输出不一致，无法进行训练。为了解决该问题，本文将滤波之后的心电信号进行0-1转换，实现维度的统一。0-1变换具体转换过程如图1.7所示，图中*表示ICBEB数据库标注好的R波波峰位置，1×n维(n为可变，且远小于5000)表示一组R波波峰位置的点数有n个。为了映射到1×5000维，将n个R波波峰位置转换成n个1，其余位置用0填充。为了解决维度扩大，但样本数据过少的问题，采用加权的方式，给予标注位置更大的权重，加强模型的训练准确度。具体操作为：将n个1扩大20倍即n个R波波峰位置转换成n个20。然后根据交叉验证的方式对数据集进行划分，在2000组数据中随机选择30%的数据为测试数据集，70%的数据为训练数据集。1.5算法框架1.5.1算法流程图图1.7Pix2Pix网络流程图Fig.1.7ThenetworkflowchartofPix2Pix算法流程图如图1.7所示。(1)将经过预处理后的2000组心电数据作为输入数据x，按照(1,5000,1)的维度输入至生成网络G中，生成一个R波波峰位置的预测数据：G(x)。(2)将输入数据x和G(x)在通道上进行数据融合然后作为判别网络D的输入。同时，基于通道维度大小，将ICBEB库中已标注的R波波峰位置数据(即真实标注数据)Y和x同样进行通道的数据融合，作为判别网络的第二个输入，由此构成判别网络D的总体损失值。判别网络D的目标是区分虚假(为“0”)和真实(为“1”)的输入，生成网络G的目标是使生成的R波波峰标注数据蒙蔽判别网络的鉴别能力。(3)当自定义目标准确率达到98%(阈值)以上，网络训练结束，保存训练好的网络模型。(4)调用训练好的Pix2Pix网络，对R波波峰进行检测，然后再利用最优选择算法对检测的R波波峰位置进行优化，得到最优结果，最后利用评价指标对结果进行分析。1.5.2生成网络的设计Pix2Pix网络中生成网络结构如图1.8(b)所示，以U-Net结构为基础，包含下采样(反卷积层)和上采样(卷积层)。与1.8(a)结构相比，U-Net结构采用了跳跃连接，实现了编码层到解码层数据之间的语义融合，加强了梯度反向传播的能力。本文生成网络的U-Net结构采用6层的卷积和反卷积层。(a)(b)图1.8生成网络结构：(a)Encoder-Decoder结构；(b)U-Net结构Fig.1.8Thestructureofthegeneratenetwork：(a)ThestructureofEncoder-Decoder;(b)ThestructureofU-Net如图1.9为生成网络模型的架构图，包括输入层、卷积层、激活层等模块，并加入了批量标准化处理，从而有效防止梯度消失的问题。具体参数如表1.2所示，表中conv2d表示下采样卷积，deconv2d和UpSampling2D是上采样卷积。(1,5000,1)表示(矩阵宽度，矩阵列数，通道)，第1层中输出(1,2500,64)表示(矩阵宽度，矩阵列数，通道)。padding采用“same”方式对当前输入矩阵进行0填充。第1-6层为下采样，激活层采用LeakyReLU。7-12层是上采样，激活层采用ReLU。13层是输出层，激活层采用Sigmoid。第1-6层卷积核数量分别是64、128、256、512、512、512，第7-11层卷积核数量为512、512、256、128、64。图1.9生成网络模型架构图Fig.1.9Thearchitecturediagramofgeneratenetwork表1.2生成网络模型参数Table1.2Theseparametersofgeneratenetwork层数输入输出1(1,5000,1)+conv2d(1,2500,64)2(1,2500,64)+conv2d(1,2500,128)3(1,2500,128)+conv2d(1,625,256)4(1,625,256)+conv2d+BN(1,125,512)5(1,125,512)+conv2d+BN(1,25,512)6(1,25,512)+conv2d+BN(1,5,512)7(1,5,512)+deconv2d(1,25,1024)8(1,25,1024)+deconv2d(1,125,1024)9(1,125,1024)+deconv2d(1,625,512)10(1,624,512)+deconv2d(1,1250,256)11(1,1250,256)+deconv2d(1,2500,128)12(1,2500,128)+UpSampling2D(1,5000,128)13(1,5000,128)+conv2d(1,5000,1)1.5.2判别网络的设计为了更详细的了解数据样本的局部情况，判别网络模型采用“PatchGAN”结构。判别网络将数据按照一定大小进行划分，进而形成多个区域块，判断每个区域的真假值。如图1.10是判别网络模型架构图，图中输入的数据主要有两个：一是条件信息(1,5000,1)和真实标注数据Y(1,5000,1)的融合数据，二是条件信息和生成网络生成的标注数据G(x)的融合数据，数据融合后的大小均为(1,5000,2)。图1.10判别网络模型架构Fig.1.10Thearchitecturediagramofdiscriminantnetwork表1.3是判别网络模型的具体参数。表中(1，5000，2)表示(矩阵宽度，矩阵列数，通道)。1-4层卷积核数量分别为64、128、256、256，第5层是输出层。表1.3判别网络模型参数Table1.3Theseparametersofdiscriminantnetwork层数输入输出1(1,5000,2)+conv2d(1,2500,64)2(1,2500,64)+conv2d(1,1250,128)3(1,1250,128)+conv2d(1,625,256)4(1,625,256)+conv2d(1,125,256)5(1,125,256)+conv2d(1,125,1)1.5.3损失函数的设计生成对抗网络训练过程容易出现训练不稳定、梯度消失或爆炸的问题，本文主要以Wassertein距离评价网络的具体损失，使训练尽快达到稳定。Wassertein距离又称为EM距离(Earth-Mover)，定义如下所示： (1.16)上式中，inf表示下界的最大值。是kr和kg所有可能组合的联合分布的集合。1.6实验结果及分析1.6.1网络训练结果(b)(a)(b)(a)(c)(d)1.11模型训练结果图Fig.1.11Theresultpictureofmodeltraining图1.11中蓝色的直线表示ICBEB数据库中标注的R波位置，红色的直线表示网络训练结束所生成的R波位置，纵坐标value表示权重，7-1000表示第7个周期里面的第1000次训练的模型检测效果。图中，标注的R波位置是固定的，但是网络训练结束生成的R波位置在标注位置附近有些位置是多余的，因此需要设定阈值和运用最优选择算法进行筛选。筛选规则如下所示： (1.17) (1.18)首先对纵坐标设定阈值，规则如公式(1.17)。当纵坐标大于0.25时，能够最大程度去掉标注R波位置附近多余的位置。然后进行横坐标的筛选，运用最优选择算法将连续且很多条直线取x值的最大值和最小值，计算规则如公式(1.18)，只有一条直线直接保存x值。1.6.2QRS波检测结果表1.4是600组测试集QRS波检测结果，表中P表示标准R波波峰位置个数，总计有8163个R波波峰；TP表示准确的检测R波波峰位置个数，该个数总计有7730个；FP表示将其他干扰误认为是R波波峰位置个数，总计有187个；FN表示标准R波波峰位置检测遗漏的个数，总计有246个遗漏的R波波峰；Accuracy表示正确率，正确率运用TP/P×100%计算得到，R波波峰位置检测准确率达到94.69%。表1.4ICBEB数据库QRS波检测结果Table1.4QRSwavedetectionresultwithICBEBdatabasenumberP/pcsTP/pcsFP/pcsFN/pcsAccuracy/%data_0000115141091.33data_00002161600100data_00003161600100data_00004161600100………………data_01998972077.78data_01999972077.78data_0200012110191.67总计8163773018724694.69图1.12是本文算法检测的R波波峰位置图，检测的R波波峰位置范围是250-4750。从图中可以看到(a)和(c)图完全检测到了R波波峰位置，但是(b)图中显示多检测了R波波峰位置，而(d)中显示有未被检测到的R波波峰位置。结果会出现多检和漏检两种情况，但这些多检和漏检的情况占总体检测数据的比例不大(如表1.4所示)，故本算法具有较强的鲁棒性。(a)(b)(c)(d)图1.12本算法检测R波波峰位置的结果图Fig.1.12TheresultofthisalgorithmdetectingthepositionabouttheRwavepeak1.6.3心率计算心率是正常人在安静状态下每分钟的心跳次数，一般情况下心率计算有三种方法，分别是查表法、公式法和目测法[74]。常以自然数60除以RR间期计算心率，但是RR间期是通过相邻R波波峰位置进行计算得到的，因此，心率的计算需要以R波波峰标注位置为前提。通常选取的R波波峰的范围为5.5-9.5s(即2750-4750),心率的计算公式如(1.19)和(1.20)所示： (1.19) (1.20)上式中，lij表示第i组第j+1个R波波峰位置与第j个R波波峰位置的差距。表示lij的均值，fs=500Hz。针对R波波峰多检的问题，本文设置了相应的评价体系，如公式(1.13)所示。具体解释为：若当前多检的R波波峰位置与相邻R波波峰位置距离相差位置处于[10,95]区间，则表示当前检测的R波波峰位置无效，应不与记录；反之，则表示当前检测的R波波峰位置有效。通过对总体数据进行计算，当相差距离为95时，心率计算得分取得最优结果。表1.5HR计算结果Table1.5TheresultsofHRcalculation编号HRrefHRtestdata_00001104.0104.0data_00002101.0101.0data_00003101.0101.0data_00004104.0104.0………data_0199862.061.0data_0199961.062.0data_0200081.081.0总计6006001.6.4结果分析为定量评价本文算法的有效性，我们将识别结果输入评分系统进行评分。评分标准如式(1.10)和式(1.11)所示。同时利用经典的PanandTompkins(P&T)算法和基于CNN的R波识别算法进行对比，最终得分结果如表1.6所示。表1.6各类算法比较Table1.6ComparisonofAccuracywithdifferentalgorithms算法得到评分本算法0.9913P&T算法[54]0.3345CNN算法[75]0.7958由表1.6所得，在相同的评价体系下，本文算法对于R波波峰位置的检测准确率高达99.13%，在动态心电信号的QRS波检测方面明显优于P&T算法和CNN算法。故本文设计的生成对抗网络算法可以有效检测动态心电信号的R波波峰位置，提高了心率计算的准确率。(a)(b)(c)(d)(e)(f)图1.13本文算法、P&T算法和CNN算法结果对比图。(a)(b)本文算法结果；(c)(d)P&T算法结果；(e)(f)CNN算法结果Fig.1.13Comparisonchartsoftheresultsofthisalgorithm,P&TandBasedonCNNalgorithm.(a)(b)Theresultsofthisalgorithm;(c)(d)TheresultsofP&T;(e)(f)TheresultsofCNN.图1.13给出了本文算法、P&T算法和基于CNN算法对两组心电数据(编号data_00004和编号data_00