用于睡眠呼吸暂停检测的心电特征稳定性分析*

董孝彤，曲新亮，魏守水△

(1. 山东大学控制科学与工程学院，济南 250061；2.山东第一医科大学第一附属医院，济南 250014)

1 引 言

睡眠呼吸暂停会导致患者白天嗜睡及生活萎靡不振，是高血压，心血管疾病的高风险因素[1-2]。目前，对睡眠呼吸障碍的检测主要是依据多导睡眠仪(PSG)，但其价格昂贵且需医院专业技师操作，无法推广，因此，寻求价格低廉且易于推广的方法一直是医学工程界的共识[3]。睡眠呼吸暂停是指口鼻气流量低于正常的90%，且持续时间超过10 s，低通气是指气流量低于正常的30%且血氧饱和度下降4%[4]。研究表明，睡眠呼吸暂停能引起心率的周期性变化，发生呼吸暂停时心动过缓，呼吸暂停结束时，心动过速且血压升高。呼吸暂停发作时，由于缺氧及高碳酸血症，交感神经活动逐渐增加，在呼吸暂停结束时，交感神经活动最为明显，随着呼吸暂停和呼吸恢复，交感神经活动强度迅速减弱。与健康人相比，即使在白天，睡眠呼吸暂停患者的神经活动水平也很高，而且昼夜心率变异性(HRV)明显减弱[5]。

因此，通过心电信号分析，实现睡眠呼吸暂停事件的检测一直是热门的研究课题。Penzel等[6]在PhysioNet发起的基于心电图的呼吸暂停检测和量化的挑战赛总结中指出，表现最好的算法均使用了HRV的频域特征或EDR(Electrocardiogram-derived Respiration)信号提取的特征。Correa等[7]对3种方法得到的EDR信号进行频域分析，使用阈值法进行分类，实现Se=88.28%，Sp=88.66%。Bsoul等[8]对RR间期及EDR信号进行时、频域分析，使用径向基核支持向量机(RBF-SVM)建立分类模型并根据ROC和F1的表现进行特征选择，实现Acc=86.24%。Manish等[9]使用双正交反对称小波滤波器组将心电信号片段分解为6个小波子带，每个子带信号计算模糊熵和对数能量特征，使用Mann -Whitney U检验进行特征选择，使用35折交叉验证及最小二乘支持向量机(LS-SVM)实现Acc=90.11%，Se=90.87%，Sp=88.88%。梁九兴等[10]使用HRV时频特征及概率神经网络(PNN)构建分类模型，实现Acc=75.97%，Se=82.51%，Sp=76.22%。Sadr等[11]对HRV与EDR信号进行时频分析并计算EDR信号的心肺耦合频谱特征，实现Acc=86.5%。然而，很多研究致力于提高检测的准确率却忽视了特征的稳定性研究，本研究针对睡眠障碍检测中的一些常用心电特征进行稳定性分析，并结合SVM建立睡眠呼吸暂停检测模型。

2 数据与方法

2.1 数据来源

本研究使用的数据库为Physionet的Apnea-ECG数据库。它包含70例睡眠记录，其中，心电信号为V2导联，采样频率为100 Hz，时长约8 h。信号的每分钟均由专家标记是否包含睡眠呼吸暂停和低通气事件，若存在，标记为“A”，不存在，标记为“N”。由于仅训练集的标注公开，本研究使用35例训练集记录，其中30例用作训练，5例用作测试。基于1 min时长实现睡眠呼吸暂停事件的检测。

2.2 数据预处理

基于心电信号提取RR间期信号及EDR信号。使用HRV分析工具箱[12]检测QRS波的位置，将相邻QRS波位置差分，得到RR间期序列，并对RR间期进行校正[13]。为了保留足够的信息量，在分析中排除了校正后累计RR间期时长小于30 s的信号片段。

呼吸过程中，胸腔内的电阻抗会发生变化，当肺部的空气量增加或减少时，电极相对于心脏的位置也会发生变化，因此呼吸运动的信息会反应在ECG中。本研究使用QRS的面积估计EDR信号[7]，图1为真实呼吸信号与EDR信号。

图1 真实呼吸信号与EDR信号

2.3 特征描述

对RR间期序列及EDR信号进行时域、频域、非线性及相关性分析，共提取了45个特征，见表1。频域分析中，使用三次样条插值将RR间期和EDR信号重采样到4 Hz，对重采样后的RR间期与EDR信号进行9层DB4小波变换，计算第4-9层近似系数(频域范围约为0.007～0.5 Hz)的方差(var_RR4- var_RR9、var_EDR4- var_EDR9)[8]。使用快速傅里叶变换(FFT)对重采样后的RR间期与EDR信号估计功率谱密度(PSD)，并分别计算极低频功率VLF(0～0.04 Hz)、低频功率LF(0.04～0.15 Hz)、高频功率HF(0.15～0.4 Hz)及LF/HF。除上述特征外，EDR信号还计算了中心频率(Fc)及主频率(Fm)[7]两个特征。计算公式如下，其中，fi表示频率，N表示PSD的长度。

(1)

(2)

非线性分析中，RR间期的庞加莱图可以反映RR间期的分布情况。其中，SD1为y=-x+2μ(μ为均值)方向上的方差，与高频功率有关。SD2为y=x方向上的方差，与低频及高频功率均相关。两者的比值SD1/SD2反映了长、短期HRV之间的平衡[14]。近似熵(sample entropy，ApEn)[15]是对不稳定时间序列复杂度的一种度量，其思想是检测一个时间序列中产生新子序列的概率。模糊测度熵(fuzzy measure entropy，FMEn)[16]使用模糊函数的隶属度代替ApEn中使用的Heaviside函数作为向量相似判断标准。同时，FMEn利用模糊局部测度熵和模糊全局测度熵反映生理信号中隐含的整体复杂性，弥补了模糊熵只关注局部复杂性的弱点。

相关性分析中，Bianchi等[17]在研究中发现RR间期序列及EDR信号的频谱相关系数在睡眠呼吸暂停事件的识别中表现良好，因此，本研究计算了RR间期序列及EDR信号在极低频段、低频段及高频段的相关系数，分别用REv、REl及REf表示。除此之外，计算了RR间期序列及其延迟(1-5个样本)序列的相关系数(RRcorr1- RRcorr5)[8]。

表1 基于RR间期及EDR信号的特征

2.4 mRMR-RRA稳定特征选择

特征选择可以加快数据挖掘算法的速度，提高学习精度，增强模型的可理解性，但特征选择结果往往因样本差异或噪声等因素而出现不一致。为了获取稳定的特征排序，本研究使用集成特征选择方法进行特征选择[18]。将训练集随机分成5组，分别使用最小冗余最大相关(mRMR)算法[19]进行排序，使用稳健排序聚合(RRA)[20]方法得到最终排序结果，我们将此方法称为mRMR-RRA，整体结构见图2。其中，RRA通过计算每个特征在各次排名间的差异显著性对特征排序，其复杂度与输入维度大小成线性关系，对异常值、噪声和误差具有较强的鲁棒性。mRMR是一种基于相关关系的过滤式特征选择方法，它不仅考虑了特征与标签之间的相关性，而且考虑了特征之间的信息冗余。其中，相关性用互信息(MI)度量，两个离散随机变量X和Y的MI可以定义为:

图2 mRMR-RRA方法的整体结构

Fig.2The overall structure of the mRMR-RRA method

(3)

其中p(x,y)是x和y的联合概率密度函数，p(x)和p(y)分别是x和y的边缘概率密度函数。MI越大，相关性越高。设S是所有特征的集合，Sm-1是已排序特征集合，根据最小冗余最大相关准则，添加到Sm-1中的第m个特征满足以下条件:

(4)

其中，x表示特征，c表示样本标签。

2.5 分类模型及评价指标

使用斯皮尔曼排序相关系数(Ds)[21]对排序结果的稳定性进行评价。设r和r′分别为同一方法在不同数据上得到的排序向量，r(i)和r′(i)是特征i在排序向量r和r′中的位置，c为特征总数，计算方法如下：

(5)

Ds的值在[-1,1]之间，Ds越接近1，排序结果越稳定。对于n次排序结果，计算Ds的平均值：

(6)

使用RBF-SVM建立分类模型，通过10折交叉验证选择最佳特征数量，平台为matlab2017b。为了加快训练速度，对特征进行归一化处理，为了提高分类性能，对分类结果进行中值滤波后处理[22]。使用准确率(Acc)、灵敏性(Se)、特异性(Sp)作为模型的评价指标。

3 结果与讨论

3.1 特征稳定性分析

为了对比mRMR-RRA方法与mRMR方法的稳定性，在训练集中随机抽取数据，每次抽取数据量为总数的80%，共抽取5次，计算每两次间的Ds值和平均Ds值，见表2。从表中可以看出，使用mRMR-RRA方法，特征排序稳定性明显提高。

表2 mRMR- RRA方法和mRMR方法的排序稳定性

图3列出了各个特征在呼吸暂停信号与正常信号间的分布差异及最终排序结果。每个子图中左侧为正常睡眠信号，右侧为呼吸暂停信号。‘R’表示特征来自RR间期，‘E’表示特征来自EDR信号。根据mRMR-RRA排序结果，特征从左到右，从上到下依次排列。排名前10的特征包括：RR间期的频域分析(RRd6，RRd7，LF/HF)，EDR信号的时域分析(var)、频域分析(Fm)及非线性分析(FMEn)，相关性分析中的RRcorr2，RRcorr3，RRcorr4，REh。可以看出，排序靠前的特征在两类间的差异明显大于排序靠后的特征，说明差异性大的特征表现较稳定。虽然有些排序靠后的特征在两类间的差异也很显著，但由于其与排名靠前的特征相关性较高，因此排名靠后，如RRcorr2、RRcorr1等。

图3正常睡眠信号与呼吸暂停信号间特征的分布差异(归一化结果)

Fig.3Differences in the distribution of characteristics between normal sleep signals and apnea signals (Normalized results)

RR间期的频域分析显示，呼吸暂停信号的VLF，LF，LF/HF与正常睡眠信号相比偏高，这与睡眠呼吸暂停发生时，交感神经活性增强的生理学结果吻合。观察RRd9-RRd4的变化可以发现，睡眠呼吸暂停与正常睡眠信号间的差异在低频时更明显。EDR的频域分析与RR间期的相应特征分布相似，说明EDR信号与RR间期信号有较好的相关性。睡眠呼吸暂停信号的Fc和Fm相对正常呼吸信号降低。非线性分析中，RR间期与EDR信号的熵值均显示与ApEn相比，FMEn在两类间的差异更显著，并且由于交感神经兴奋会增加信号的确定度[23]，因此呼吸暂停信号的熵值显著降低。相关性分析中，RR间期及其延迟序列相关性在呼吸暂停信号中显著大于正常睡眠信号，RR间期与EDR在高频段的相关性更能反映两类间的差异。

3.2 特征验证

图4为两种方法10折交叉验证准确率随特征数量的变化情况。可以看出，mRMR-RRA方法整体优于mRMR方法，说明特征的稳定性增加的同时，分类效果也得到了提升。当特征数量为14时，mRMR-RRA方法达到局部最优，此时Acc=83.22%。随着特征数量继续增加，Acc呈小幅度波浪式变化。虽然特征数量大于39时，准确率提高到了86%，但模型的复杂度过大。因此，在最终分类模型的训练中，使用了mRMR-RRA方法得到的前14个特征。

图4 10折交叉验证准确率随特征数量的变化

表3 分类结果

表3为分类模型在测试集上的表现。第一行显示Se=71.30%，远远低于Sp=96.95%，说明假阴性(FN)较高，这在临床中存在非常大的隐患。因此，在训练模型时重新设置代价矩阵，将FN的误分类代价(Cost)设置为假阳性(FP)的两倍，此时得到Se=86.71%，相比之前提升了15.41%，Acc=90.03%，相比之前提升了1.77%。

4 结论

本研究针对心电信号检测睡眠呼吸暂停事件中心电特征的稳定性进行研究。实验结果表明，基于集成稳定特征选择策略，将mRMR特征选择方法与RRA方法结合进行特征稳定性分析，与单独使用mRMR方法相比，特征稳定性及分类准确率均有所提高。ECG指标与EDR指标有较好的相关性，结合ECG信号重构呼吸信号值得进一步研究。本研究给出的特征排序方法与结果，可以为心电信号识别睡眠障碍等相关研究提供一定的参考意义。