在全球的死亡案例中，心血管疾病（CVD）是主要的致死原因之一。心音分类识别在心血管疾病的早期发现中起着关键作用。正常心音和异常心音之间的区别并不明显，本文为提升心音分类模型的准确度，提出一种基于双谱分析的心音特征提取方法，并将其与卷积神经网络（CNN）结合，对心音进行分类。该算法能够有效地利用双谱分析来抑制高斯噪声，而且不需要准确分割心音信号就能提取其特征，同时结合了卷积神经网络的强大分类性能，从而实现对心音的准确分类。根据实验结果显示，在相同的数据和实验条件下，本文提出的算法在准确率、灵敏度和特异性方面分别达到了0.910、0.884和0.940。与其他心音分类算法相比，本文算法提升明显，并具有较强的鲁棒性和泛化能力，因此有望应用于先心病的辅助检测。
引用本文： 彭利勇, 全海燕. 基于双谱特征提取和卷积神经网络的心音分类算法. 生物医学工程学杂志, 2024, 41(5): 977-985, 994. doi: 10.7507/1001-5515.202310016 复制
0 引言
根据世界卫生组织（World Health Organization，WHO）的数据显示，心血管疾病（cardiovascular disease，CVD）是全球死亡率上升的主要原因之一，它严重影响人们的预期寿命[1]。心音信号包含了大量心血管系统的信息，心音信号的分析对于诊断和预防心血管疾病有重要意义[2]。自动心音分析检测系统，是预防和诊断心血管疾病的重要手段。在早期，心音分析完全依赖于经验丰富的医生以听诊的方式进行，然而这种方法存在主观性较大的问题，不够客观、准确。如今，随着计算机辅助心音分析技术的快速发展，越来越多的自动化心音分析方法受到关注。这些技术的出现为心音分析带来了新的可能性，有望提高诊断准确性并降低主观误差。
目前，心音分类识别是自动心音分析技术的一大挑战。在心音分类任务中的心音信号特征提取和分类识别算法是国内外许多研究者长期关注的领域。Bozkurt等[3]对心音信号的特征进行了分类，分为时域、频域以及时域和频域相结合的特征；时域特征包括第一心音、第二心音以及收缩期和舒张期，而频域特征包括梅尔频率倒谱系数（Mel-frequency cepstrum coefficient，MFCC）滤波器的特定系数能量和梅尔谱图特征，其方法达到了81.5%的分类精度。Chen等[4]对心音信号进行短时傅里叶变换（short-time Fourier transform，STFT），将得到的特征谱图作为卷积神经网络（convolutional neural networks，CNN）的输入，然后对心音进行分类识别，该方法最终在测试集上得到了95.49%的准确率（accuracy，Acc）。Asmare等[5]从心音样本中提取了对数梅尔频率谱系数（log Melfrequency spectral coefficients，MFSC）分量作为心音特征图，结合CNN进行分类，采用170个心音样本，得到了96.10%的准确率。文献[4]和文献[5]都取得了较高的分类准确率，可是它们所使用的心音数据样本量较少，其准确率的可信度较低，很难在实际应用中推广。此外，文献[5]使用MFSC方法可能会引入冗余特征，从而降低分类准确率。Li等[6]使用MFSC作为特征与深度CNN（deep CNN，DCNN）相结合的分类模型，实现了89.6%的无权平均召回，其灵敏度（sensitivity，Se）和特异性（specific，Sp）分别达到了89.5%和89.7%。Ryu等[7]使用汉明窗滤波器去除原始心音信号中的不相关噪声，然后将固定长度的心音贴片送入四层一维CNN进行分类。Rubin等[8]将原心音信号转换为MFCC的二维时频特征，然后通过与CNN相结合的方法实现心音的分类，最终得到了84.80%的准确率。Nilanon等[9]使用心音信号的功率谱密度（power spectral density，PSD）特征谱，然后结合CNN实现心音分类，最终得到88.2%的准确率。文献[8]通过隐藏半马尔可夫模型采用第一心音开始的3 s心音片段作为输入，而文献[9]采用固定长度的滑动片段作为输入。
上述方法在心音的特征提取和分类识别方面都取得了很好的效果，但还存在许多的问题，如模型复杂、使用数据样本较少、准确度低和提取后的特征包含信息较少等。如何设计一种更轻量型、分类效果好且能包含更多心音信息的特征提取方式，对心音分类任务来说具有十分重要的意义[10]。当前，大多数心音信号分析方法主要采用时域和频域分析，通常结合了时频域特征。这些方法主要依赖于二阶统计量作为数学工具，但忽略了有用的高阶信息，导致其结果难以令人满意。为了解决这些问题，需要采用更高阶的方法来更准确地表征信号特征。高阶谱分析（higher order spectral analysis，HOSA），是信号处理学科的前沿性研究方向，是处理非线性、非高斯信号的有力工具。高阶谱分析通过考虑更高阶概率结构来描述随机信号，弥补了仅依赖于二阶统计量的方法的不足。这种分析方法不仅能够反映信号的能量分布，还能保留信号的相位信息。此外，从理论上来看，高阶谱具备完全抑制高斯噪声的能力，因此在降噪方面具有卓越表现。本文采用高阶谱中阶数最低的三阶谱来提取心音信号的特征。一般称三阶谱为双谱（bispectrum），双谱相对于更高阶的谱分析来说计算量最小，但它仍然具备高阶谱的特性，能在保留心音高阶特征信息的同时也具有很好的降噪性能。
综上所述，为简化分类过程，在尽可能保留心音信号更多信息的同时提升分类的准确率，本文提出一种基于双谱分析特征提取与CNN相结合的心音信号分类识别方法，通过对心音信号进行双谱分析得到心音双谱特征图，然后利用CNN对双谱特征图进行分类，并在公开心音数据集上进行测试，以期取得较好的分类效果。本文所提算法提供了一种将传统心音信号处理算法与深度学习算法有机结合的新思路，这将有助于推动心脏信号处理算法和心脏疾病辅助检测技术的进一步探索和发展。
1 方法分析
1.1 算法框架
通常心音分类算法大致分为心音预处理（降噪、分割等）、心音特征提取和心音分类识别等步骤。本文算法首先对心音进行分帧加窗处理，将每一帧作为一个样本心音。然后，对这些样本心音进行双谱分析，生成心音双谱特征图，从而实现将一维信息转换成二维谱图形式[11]。接下来，对双谱特征图进行非线性归一化处理，以增强心音双谱特征的细节，并避免部分特征图损失关键信息。最终，利用CNN实现对心音信号的分类识别。本文算法框图如图1所示。
图1
心音分类模型框图
Figure1.
Block diagram of the heart sound classification model
1.2 数据来源
本文所用数据集为公开数据集，来自于正常/异常心音记录分类：2016年物理网/心脏病学计算挑战赛（classification of normal/abnormal heart sound recordings: the physionet/computing in cardiology challenge 2016）[12]。该数据集由6个子数据集（trainning_a～ trainning_f）组成，数据集中包含来自健康受试者和各种心脏病（主要为：心脏瓣膜缺陷和冠状动脉疾病）患者的3 240份原始心音信号。这些信号是由各种电子听诊器采集，信号长度从5～120 s不等，如表1所示，每个心音样本均被标记为异常或正常，录音都被重新取样到2 000 Hz，并以.wav格式保存。这个心音数据库采集了来自人体4个不同部位的心音信号，包括主动脉区、肺心病区、三尖瓣区和二尖瓣区。正常心音源自于健康人群，而异常心音则来自于确诊患心脏疾病的患者。需要指出的是，异常心音记录没有进一步的具体分类。这一数据集存在类别不平衡问题，正常心音有2 575条，而异常心音仅有665条[13]。
表1
实验数据集
Table1.
Experimental data set
1.3 数据预处理
对心音信号分类的主要目的是检测信号中是否存在异常而不是进一步去识别它，所以无需对心音进行精确分割[14]。连续第一心音（或第二心音）之间的时间间隔构成一个心动周期。正常心音收缩的持续时间范围为300～400 ms，而舒张期的持续时间范围在500～600 ms，典型的心动周期通常持续时间范围在800～1 000 ms[15]。心音信号的周期性特点对于心音分类识别研究具有十分重要的意义。
本文首先去除每个心音样本最开始1 s和最后1 s，消除仪器接入和撤离的干扰，然后用五阶巴特沃斯滤波器对每个心音样本进行滤波，滤波器的通频带为20～800 Hz，消除低频伪影和高频噪声的影响，结果如图2所示。
图2
降噪前后对比
Figure2.
Before and after noise reduction
由于CNN对输入的要求，尽量保证每个输入样本的形状相同。本文对每个数据分帧加窗处理，帧长为4 s，帧移为2 s，对每一帧的样本加汉明窗并保存，这样得到处理后的心音样本帧序列，正常样本19 964个，异常样本6 109个。
1.4 预特征提取
对于音频信号分类来说，特征提取是分类任务中的关键。由于心音数据为时域的高维信息，直接用来识别非常困难，特征提取通过变换的方法将高维信息进行维数压缩。一般经典谱估计法在特征提取中最为常用。而本文所使用的双谱分析考虑了音频信号的非线性相互关系，这意味着它可以捕获到音频信号中更复杂的频域结构和相互作用。而心音信号包含许多不同频率成分，使用双谱分析可以更好地捕获心音信号中的频率交互和非线性特征，包括心音的共振结构和特定频率组合，这有助于提取更丰富的心音特征，以区分不同心脏病症或不同的心音状态。此外，作为高阶谱的一种，双谱完整地保留了信号的幅度、频率和相位等信息，具有时移不变性、尺度变化性和相位保持性等特性，并能够抑制高斯有色噪声，所以被广泛应用于信号分析和处理中[16]。
1.4.1 双谱定义
信号的双谱（即三阶谱）是一种阶数最低的高阶谱，相对于其他高阶谱应用较为广泛，理论和方法比较成熟[17]。一般满足短时平稳的信号，其三阶累积量的二阶傅里叶变换被称做双谱[18]，其定义式如式(1)所示：
 |
其中，ω1、ω2分别代表两个不同的频率分量，且|ω1| ≤ π，|ω2| ≤ π，|ω1+ω2| ≤ π；τ1、τ2是任意时间延迟；C3x(τ1, τ2)为三阶累积量且满足绝对可和；Bx(∙)是双谱特征矩阵；j为虚数单位。
双谱具有对称性，只要知道双谱的主域信息，如图3阴影部分所示，就可以根据对称性得出双谱的全部信息[19]。如图3所示，本文根据后续处理分类需求，取ω1 ≥ 0、ω2 ≥ 0的部分区域来进行研究。
图3
双谱对称区域图
Figure3.
Bispectral symmetry map
1.4.2 双谱计算方法
双谱的计算方法通常有参数法和非参数模型法两种，本文采用非参数模型法计算心音信号的双谱。具体步骤如下：
步骤1：将信号{x(0)，x(1)，…，x(N − 1)}数据分为K段，每段包含M个样本，允许两段相邻数据之间有重叠，本文取50%的重复率，对每段数据都减去该段数据的均值，第k段记为xk(n)，其中，k = 1，2，···，K；n = 0，1，···，M − 1。
步骤2：进一步计算每段的傅里叶系数，如式(2)所示：
 |
其中，Xk(λ)为第k段的傅里叶系数；λ = 0，1，2，···，M/2。
步骤3：通常，设M = M1 × N0，此处M1为正整数，且与任意变量L1的关系为：M1 = 2L1 + 1；N0为总频率样本数。计算傅里叶系数的三重相关，如式(3)所示：
 |
其中，
(·)表示第k段信号的傅里叶系数的三重相关， i1与 i2表示两不同的任意频移。Δ0 = fs/N0表示在双谱区域沿水平和垂直方向的两个频率采样点之间的间隔，fs为信号采样率；λ1和λ2表示任意两不同的频率分量，且 0 ≤ λ2 ≤ λ1，λ1 + λ2 ≤ fs/2。
步骤4：由K段双谱估计的均值可以求出给定样本信号的双谱，如式(4)所示：
 |
其中， 
(·)表示信号的双谱特征矩阵，ω1 = ((2πfs)/N0)λ1，ω2 = ((2πfs)/N0)λ2。
1.4.3 双谱的选择
由1.4.2小节里双谱的计算可以看出，传统非参数直接计算法计算双谱是将每个输入看作多个段，最后计算其均值，以满足信号做双谱计算必须为短时平稳的要求。而对于心音信号而言，其病理性特点并非均匀分布于每一个心音周期内，且部分心音的重要信息分布具有不规律性，这种直接计算均值的方式会使得最后提取的特征信息量有所减少。根据异常心音的这一特点，本文需对1.4.2小节中步骤4进行实验改进，以改进均值计算的缺点，从而确定最优的双谱特征计算结果。
1.4.4 获取双谱特征图
本文将预处理后的每4 s的样本帧进行1 024个点的离散傅里叶变换，可视化结果如图4所示。
图4
心音样本离散傅里叶变换结果图
Figure4.
DFT results of heart sound samples
图4中，样本心音信号的有效信息主要频带在0～200 Hz之间。如式(2)～式(4)所示，对样本做双谱分析后得到一个1 024 × 1 024的二维复数矩阵。
取复数矩阵的实部构成特征矩阵，其数值代表了信号中两个频率分量之间的相关性大小。数值越大，表示两个频率分量之间的相关性越强。结合图4心音的频谱特征和1.4.1小节中双谱的对称性特点截取特征图，最终得到大小为102×102的心音双谱特征图，如图5所示，为正常心音（左）和异常心音（右）样本经特征提取后得到的心音双谱特征图。
图5
正常（左）和异常（右）心音双谱特征图
Figure5.
Bispectral characteristics of normal (left) and abnormal (right) heart sounds
1.4.5 双谱特征图增强处理
图像增强处理涉及将原始图像进行转换，以获得一个通用标准图像的过程。对特征图像增强的好处在于它能够提升后续算法中特征图的适用性[20]。本文中，由于得到的双谱特征图中不同的灰度都可能包含心音的有用信息，因此为避免使用CNN对后续图像特征再提取和分类计算时，由于灰度分布不均导致忽略部分低灰度区域有用信息，本文在分类前，将得到的二维双谱图进行特征图像非线性归一化处理。本文采用对数变换法对双谱特征图进行处理，对数变换法公式如式(5)所示：
 |
其中，r为输入，S为输出，v为底数，c为倍数因子，本文中取：c = 255/log(1 + max(r))。根据输入图像中的最大像素值，本文使用式（5）中c的取值来动态调整，以便将输出限定在标准的0～255范围内，以利于保持图像的视觉质量和对比度。
由上述对数变换法变换处理后，二维和三维特征图像变换前后的对比如图6所示。
图6
对数变换前后对比
Figure6.
Before and after logarithmic transformation
当引入了对数变换法增强特征图后，部分低灰度区域的特征被增强，使得双谱特征图的特征更加丰富，避免了在后续CNN中进行卷积运算的时候将有用信息忽略，从而使得心音的特征更明显，以便于提升分类的准确度，增强分类模型的性能。
1.5 分类模型
近年来，深度学习发展迅速，许多研究者用深度学习构建分类模型取得很不错的效果。CNN是深度学习的主要研究方向之一，其在语音识别、自然语言处理、图像识别等方面的应用较为成熟[21]。CNN独特的结构和图像处理方式使得它非常适合于图像分类问题[6]。本文将使用CNN对心音双谱特征进行特征提取和分类。
如表2所示，本文所搭建的CNN结构，由三个卷积层、三个池化层、两个全连接层和一个输出层组成，激活函数选用线性整流函数（rectified linear unit，ReLU），池化选择最大池化方式。首先，展开特征图为一个向量，通过一个含128个神经元的全连接层传输到输出层，使用输出层抑制过拟合，其值设置为0.5。输出层只有一个神经元，而本文研究的是二分类任务，所以只需要一个输出来表示样本属于正类别的概率。输出层，采用S型生长曲线（sigmoid）作为激活函数，计算公式如式（6）所示：
表2
CNN网络模型架构
Table2.
CNN network model architecture
 |
其中，σ(∙)是sigmoid激活函数，x为输入的特征图的值，函数将输入x映射到一个0～1之间的概率值。通常，将0.5作为阈值，当sigmoid函数的输出大于0.5时，模型将样本分类为正类别；当输出小于0.5时，模型将样本分类为负类别。
在训练过程中，使用二元交叉熵损失函数来衡量模型的预测与真实标签之间的差异。由于本文所有数据集存在数据类不平衡的问题，因此引入类别权重。通过对损失函数（以符号Loss表示）中的类赋予不同的权重，以使模型更加关注数据量较少的类，从而达到平衡数据的作用，如式(7)所示：
 |
其中，y是二元标签，取0或者1；p(yi)是输出属y标签的概率，通常在0～1之间；N为输出个数。本文设正负数据类权重分别为α和β，可根据数据类不平衡程度对其值进行相应调整。本文正样本19 964个，负样本6 109个，于是将α的值设置为1/4，β设为3/4。
2 实验分析
2.1 实验环境与数据说明
本研究软件环境：算法主体编程软件为Pycharm2023专业版（JetBrains Inc.，捷克），编程语言为Python3.9（Python Software Foundation Inc.，美国），所使用CNN分类网络模型在TensorFlow2.0（Google Inc.，美国）框架中实现。
本研究硬件环境：中央处理器i7-13700H CPU @2.4G Hz（Intel Inc.，美国），独立显卡（NVIDIA Inc.，美国）。
本研究所有实验均使用同一数据集且在同一台设备进行。
2.2 评估指标
为了与目前最先进的基于CNN的心音分类算法模型进行比较，本文引入下述评价指标对模型性能进行评价：准确率、灵敏度、特异性。准确率，是指模型在整个测试数据集上正确分类的样本数量与总样本数量之比；它衡量了模型正确分类样本的能力，但在不平衡数据集中可能不够敏感。灵敏度，衡量了模型对于正例样本的识别能力；它表示在所有实际正例中，模型正确识别为正例的比例。特异性，衡量了模型对于负例样本的识别能力；它表示在所有实际负例中，模型正确识别为负例的比例。高特异性，意味着模型能够有效地排除负例，降低了误报的可能性[22]。
2.3 实验结果分析
为证明本文所提分类模型的有效性，分别进行了4组实验。实验一：对前文1.4.2小节心音双谱特征计算中步骤4改进，对比三种计算方法对分类指标的影响，得出更适合的双谱特征图。实验二：对比式(5)中不同对数底数值v对实验结果的影响，从而得到较为理想的v的取值。实验三：通过实验对比不同特征提取方法下心音分类模型的评价指标，证明本文所提基于高阶统计量的心音信号双谱特征提取方法的有效性。实验四：将本文所提心音分类模型与当前已有的心音分类模型做对比，通过各种模型评价指标证实本文所提基于双谱特征提取与CNN相结合的心音分类算法的有效性和优越性。
实验一：如式(4)所示，对心音样本进行双谱分析特征提取计算，将每一个小片段的三重相关结果求和，然后求取平均值。为解决1.4.3小节中所提问题，本文采取对比实验计算三种双谱特征结果，确定更适合的心音双谱特征，分别为均值（mean，M）、方差（variance，Var）和标准差（standard deviation，SD）。该实验所使用的分类网络模型结构如表2所示，实验所使用数据集均为表1所示数据，均使用了底数为10的对数变换法对特征图像增强，实验结果如表3所示。
表3
不同双谱特征结果对比
Table3.
Comparison of results for different bispectral features
表3中，加粗的数据代表同一评价指标下实验结果最好的数据。由实验结果可以看出，基于标准差计算的心音双谱特征更适合作为心音分类模型的特征，这可能是因为病理性的特征分布不均导致的。许多算法在心音预处理过程中仅进行了简单的心音分割，导致某些片段所含的病理性特征较少，而噪声较多。如果对这些异常特征较少的样本片段直接采用传统的双谱计算方法进行均值处理，所得的双谱特征图与正常心音的双谱特征相差无几，从而增加了分类的难度，导致异常样本的识别率较低。这在实验结果中表现为特异性较低，且与实验事实相符。
实验二：不同的特征图增强效果直接影响分类的结果，本文采用非线性图像增强方法中对数变换法对双谱特征图进行增强，其中对数的底数变化对于特征图的增强效果有差异，如图7所示。
图7
不同底数对数变换函数
Figure7.
Logarithmic transformation function with different bases
对数变换法对特征图低灰度值有拉伸效果，避免了后续计算中部分特征的丢失。不同的底数v的值对特征图的增强效果不同，其中v = (r + 1)c/r表示没有进行对数变换增强的特征图，以对比当没有进行图像增强时的模型分类性能。本文实验了5种不同底数，分别是：(r + 1)c/r、8、10、15、20，对心音分类模型的分类准确率的影响。实验所使用的分类网络模型结构如表2所示，5个实验所使用数据均为表1所示数据集，且均使用了基于标准差计算的双谱特征图作为心音特征谱图。实验结果如表4所示。
表4
不同底数对分类结果的影响
Table4.
Effect of different base numbers on classification results
从实验结果来看，没有进行图像增强的心音特征图分类效果最差；而当v值约为10时，模型分类效果最好。如图7所示，当v = 10时，其拉伸曲线与不进行对数变换的虚线相交于点（255，255），这表示v = 10时的对数变换计算对像素值为255附近的像素拉伸或压缩不明显。对特征图进行非线性变换的目的在于增强灰度较低的部分，使得其不至于被后续识别网络所忽略，从而保证特征图的信息量不减少；与此同时，也应尽量不对高灰度区域的值进行过度的改变，从而保证特征图的高对比度；这也是选择非线性方法的一个重要原因。
实验三：本研究对比了另外5种特征提取方式与CNN结合的心音分类实验，所有分类网络模型均采用表1中数据，且均使用了表2所示的CNN结构，其中基于双谱特征的分类实验所使用的双谱特征是基于标准差的分类实验结果，且对双谱特征进行了对数变换增强，其底数为10。所有分类模型评价指标的对比实验结果，如表5所示。其中，一维+CNN，表示直接使用一维心音序列结合CNN对心音信号进行分类识别。
表5
不同特征提取方式分类性能对比
Table5.
Comparison of classification performance of different feature extraction methods
从表5的实验结果来看，基于高阶统计量的心音双谱特征相比于其他基于一二阶心音特征的分类模型来说，分类的准确率更高，这与双谱特征包含有更为丰富的心音特征有关。心音信号通常包含多个频率成分，每个频率成分都对应于心脏的某种生理过程。这些频率成分之间可能存在复杂的相互作用和耦合关系，例如：心室颤动、心室扑动等心律失常疾病。而双谱可以用来分析这些不同频率成分之间的相互作用，以识别各种疾病更具体的特征。
实验四：为证明本文所提基于双谱特征提取和CNN的心音分类算法的有效性，将其与其他4种当前较为先进的心音分类算法做性能比较[8-9, 23-24]，评价指标除前文2.2小节中所述准确率、灵敏度、特异性以外，还增加了模型参量、是否精确分割等两个指标，对各种心音分类模型进行综合评估。5组实验均使用表1所示数据集，且在同一设备进行。实验结果如表6所示，其中“—”代表无数据可查。
表6
不同心音分类模型分类性能对比
Table6.
Comparison of classification performance of different heart sound classification models
从表6的实验结果可以看出，本文所提出的基于双谱特征和CNN相结合的心音分类算法模型性能达到了较为先进的水平，文献[8]虽然也具有很好的准确率和特异性，但是其模型参量较大，而本文的算法模型参量仅为0.04 × 106个，且不需要对心音数据进行精确分割。
综上所述，本文所提心音分类算法模型综合性能优于其他算法，这说明基于高阶统计量的心音特征提取方法与CNN结合的心音分类算法模型在提高心音分类模型的综合性能方面更具优势。
3 结论
在心音分类任务中，心音的病理及生理特征不仅仅包含在幅度和频率中。为更有效地提取心音特征，本文针对低阶统计量（功率谱等）提取心音特征不足的局限，提出采用基于高阶统计量中的双谱对心音信号进行特征提取。高阶谱具有很强的消噪能力，理论上能完全抑制高斯噪声。就双谱的物理意义而言，体现了信号在频谱中的斜度，而以高斯噪声为例，其在频谱中斜度为零，因此双谱可直接对比心音信号中非线性、非对称性的特征成分。本文还采用对数变换法对心音特征图像进行增强，将特征图标准化的同时也能突出更多的信号特征，这有利于CNN对特征图像的分类处理，从而增强模型分类效果。本文实验结果表明，基于双谱特征提取和CNN相结合的方式对心音进行分类，效果显著提升，准确率达到0.91，且模型参量低，无需对原始信号进行精确分割，综合性能优于其他同类算法，有望应用于心脏疾病的辅助检查等方面。但，本文研究也存在一定局限，比如本文只实现了正常心音和异常心音的分类，对于代表具体心脏疾病的异常心音信号并没有做详细的分类，这将成为本课题组今后的研究重点。
重要声明
利益冲突声明：本文全体作者均声明不存在利益冲突。
作者贡献声明：彭利勇负责实验构思与设计、数据收集、数据分析和论文写作；全海燕负责实验指导、论文写作指导、数据审查和论文整体审查。
0 引言
根据世界卫生组织（World Health Organization，WHO）的数据显示，心血管疾病（cardiovascular disease，CVD）是全球死亡率上升的主要原因之一，它严重影响人们的预期寿命[1]。心音信号包含了大量心血管系统的信息，心音信号的分析对于诊断和预防心血管疾病有重要意义[2]。自动心音分析检测系统，是预防和诊断心血管疾病的重要手段。在早期，心音分析完全依赖于经验丰富的医生以听诊的方式进行，然而这种方法存在主观性较大的问题，不够客观、准确。如今，随着计算机辅助心音分析技术的快速发展，越来越多的自动化心音分析方法受到关注。这些技术的出现为心音分析带来了新的可能性，有望提高诊断准确性并降低主观误差。
目前，心音分类识别是自动心音分析技术的一大挑战。在心音分类任务中的心音信号特征提取和分类识别算法是国内外许多研究者长期关注的领域。Bozkurt等[3]对心音信号的特征进行了分类，分为时域、频域以及时域和频域相结合的特征；时域特征包括第一心音、第二心音以及收缩期和舒张期，而频域特征包括梅尔频率倒谱系数（Mel-frequency cepstrum coefficient，MFCC）滤波器的特定系数能量和梅尔谱图特征，其方法达到了81.5%的分类精度。Chen等[4]对心音信号进行短时傅里叶变换（short-time Fourier transform，STFT），将得到的特征谱图作为卷积神经网络（convolutional neural networks，CNN）的输入，然后对心音进行分类识别，该方法最终在测试集上得到了95.49%的准确率（accuracy，Acc）。Asmare等[5]从心音样本中提取了对数梅尔频率谱系数（log Melfrequency spectral coefficients，MFSC）分量作为心音特征图，结合CNN进行分类，采用170个心音样本，得到了96.10%的准确率。文献[4]和文献[5]都取得了较高的分类准确率，可是它们所使用的心音数据样本量较少，其准确率的可信度较低，很难在实际应用中推广。此外，文献[5]使用MFSC方法可能会引入冗余特征，从而降低分类准确率。Li等[6]使用MFSC作为特征与深度CNN（deep CNN，DCNN）相结合的分类模型，实现了89.6%的无权平均召回，其灵敏度（sensitivity，Se）和特异性（specific，Sp）分别达到了89.5%和89.7%。Ryu等[7]使用汉明窗滤波器去除原始心音信号中的不相关噪声，然后将固定长度的心音贴片送入四层一维CNN进行分类。Rubin等[8]将原心音信号转换为MFCC的二维时频特征，然后通过与CNN相结合的方法实现心音的分类，最终得到了84.80%的准确率。Nilanon等[9]使用心音信号的功率谱密度（power spectral density，PSD）特征谱，然后结合CNN实现心音分类，最终得到88.2%的准确率。文献[8]通过隐藏半马尔可夫模型采用第一心音开始的3 s心音片段作为输入，而文献[9]采用固定长度的滑动片段作为输入。
上述方法在心音的特征提取和分类识别方面都取得了很好的效果，但还存在许多的问题，如模型复杂、使用数据样本较少、准确度低和提取后的特征包含信息较少等。如何设计一种更轻量型、分类效果好且能包含更多心音信息的特征提取方式，对心音分类任务来说具有十分重要的意义[10]。当前，大多数心音信号分析方法主要采用时域和频域分析，通常结合了时频域特征。这些方法主要依赖于二阶统计量作为数学工具，但忽略了有用的高阶信息，导致其结果难以令人满意。为了解决这些问题，需要采用更高阶的方法来更准确地表征信号特征。高阶谱分析（higher order spectral analysis，HOSA），是信号处理学科的前沿性研究方向，是处理非线性、非高斯信号的有力工具。高阶谱分析通过考虑更高阶概率结构来描述随机信号，弥补了仅依赖于二阶统计量的方法的不足。这种分析方法不仅能够反映信号的能量分布，还能保留信号的相位信息。此外，从理论上来看，高阶谱具备完全抑制高斯噪声的能力，因此在降噪方面具有卓越表现。本文采用高阶谱中阶数最低的三阶谱来提取心音信号的特征。一般称三阶谱为双谱（bispectrum），双谱相对于更高阶的谱分析来说计算量最小，但它仍然具备高阶谱的特性，能在保留心音高阶特征信息的同时也具有很好的降噪性能。
综上所述，为简化分类过程，在尽可能保留心音信号更多信息的同时提升分类的准确率，本文提出一种基于双谱分析特征提取与CNN相结合的心音信号分类识别方法，通过对心音信号进行双谱分析得到心音双谱特征图，然后利用CNN对双谱特征图进行分类，并在公开心音数据集上进行测试，以期取得较好的分类效果。本文所提算法提供了一种将传统心音信号处理算法与深度学习算法有机结合的新思路，这将有助于推动心脏信号处理算法和心脏疾病辅助检测技术的进一步探索和发展。
1 方法分析
1.1 算法框架
通常心音分类算法大致分为心音预处理（降噪、分割等）、心音特征提取和心音分类识别等步骤。本文算法首先对心音进行分帧加窗处理，将每一帧作为一个样本心音。然后，对这些样本心音进行双谱分析，生成心音双谱特征图，从而实现将一维信息转换成二维谱图形式[11]。接下来，对双谱特征图进行非线性归一化处理，以增强心音双谱特征的细节，并避免部分特征图损失关键信息。最终，利用CNN实现对心音信号的分类识别。本文算法框图如图1所示。
图1
心音分类模型框图
Figure1.
Block diagram of the heart sound classification model
1.2 数据来源
本文所用数据集为公开数据集，来自于正常/异常心音记录分类：2016年物理网/心脏病学计算挑战赛（classification of normal/abnormal heart sound recordings: the physionet/computing in cardiology challenge 2016）[12]。该数据集由6个子数据集（trainning_a～ trainning_f）组成，数据集中包含来自健康受试者和各种心脏病（主要为：心脏瓣膜缺陷和冠状动脉疾病）患者的3 240份原始心音信号。这些信号是由各种电子听诊器采集，信号长度从5～120 s不等，如表1所示，每个心音样本均被标记为异常或正常，录音都被重新取样到2 000 Hz，并以.wav格式保存。这个心音数据库采集了来自人体4个不同部位的心音信号，包括主动脉区、肺心病区、三尖瓣区和二尖瓣区。正常心音源自于健康人群，而异常心音则来自于确诊患心脏疾病的患者。需要指出的是，异常心音记录没有进一步的具体分类。这一数据集存在类别不平衡问题，正常心音有2 575条，而异常心音仅有665条[13]。
表1
实验数据集
Table1.
Experimental data set
1.3 数据预处理
对心音信号分类的主要目的是检测信号中是否存在异常而不是进一步去识别它，所以无需对心音进行精确分割[14]。连续第一心音（或第二心音）之间的时间间隔构成一个心动周期。正常心音收缩的持续时间范围为300～400 ms，而舒张期的持续时间范围在500～600 ms，典型的心动周期通常持续时间范围在800～1 000 ms[15]。心音信号的周期性特点对于心音分类识别研究具有十分重要的意义。
本文首先去除每个心音样本最开始1 s和最后1 s，消除仪器接入和撤离的干扰，然后用五阶巴特沃斯滤波器对每个心音样本进行滤波，滤波器的通频带为20～800 Hz，消除低频伪影和高频噪声的影响，结果如图2所示。
图2
降噪前后对比
Figure2.
Before and after noise reduction
由于CNN对输入的要求，尽量保证每个输入样本的形状相同。本文对每个数据分帧加窗处理，帧长为4 s，帧移为2 s，对每一帧的样本加汉明窗并保存，这样得到处理后的心音样本帧序列，正常样本19 964个，异常样本6 109个。
1.4 预特征提取
对于音频信号分类来说，特征提取是分类任务中的关键。由于心音数据为时域的高维信息，直接用来识别非常困难，特征提取通过变换的方法将高维信息进行维数压缩。一般经典谱估计法在特征提取中最为常用。而本文所使用的双谱分析考虑了音频信号的非线性相互关系，这意味着它可以捕获到音频信号中更复杂的频域结构和相互作用。而心音信号包含许多不同频率成分，使用双谱分析可以更好地捕获心音信号中的频率交互和非线性特征，包括心音的共振结构和特定频率组合，这有助于提取更丰富的心音特征，以区分不同心脏病症或不同的心音状态。此外，作为高阶谱的一种，双谱完整地保留了信号的幅度、频率和相位等信息，具有时移不变性、尺度变化性和相位保持性等特性，并能够抑制高斯有色噪声，所以被广泛应用于信号分析和处理中[16]。
1.4.1 双谱定义
信号的双谱（即三阶谱）是一种阶数最低的高阶谱，相对于其他高阶谱应用较为广泛，理论和方法比较成熟[17]。一般满足短时平稳的信号，其三阶累积量的二阶傅里叶变换被称做双谱[18]，其定义式如式(1)所示：
|
其中，ω1、ω2分别代表两个不同的频率分量，且|ω1| ≤ π，|ω2| ≤ π，|ω1+ω2| ≤ π；τ1、τ2是任意时间延迟；C3x(τ1, τ2)为三阶累积量且满足绝对可和；Bx(∙)是双谱特征矩阵；j为虚数单位。
双谱具有对称性，只要知道双谱的主域信息，如图3阴影部分所示，就可以根据对称性得出双谱的全部信息[19]。如图3所示，本文根据后续处理分类需求，取ω1 ≥ 0、ω2 ≥ 0的部分区域来进行研究。
图3
双谱对称区域图
Figure3.
Bispectral symmetry map
1.4.2 双谱计算方法
双谱的计算方法通常有参数法和非参数模型法两种，本文采用非参数模型法计算心音信号的双谱。具体步骤如下：
步骤1：将信号{x(0)，x(1)，…，x(N − 1)}数据分为K段，每段包含M个样本，允许两段相邻数据之间有重叠，本文取50%的重复率，对每段数据都减去该段数据的均值，第k段记为xk(n)，其中，k = 1，2，···，K；n = 0，1，···，M − 1。
步骤2：进一步计算每段的傅里叶系数，如式(2)所示：
|
其中，Xk(λ)为第k段的傅里叶系数；λ = 0，1，2，···，M/2。
步骤3：通常，设M = M1 × N0，此处M1为正整数，且与任意变量L1的关系为：M1 = 2L1 + 1；N0为总频率样本数。计算傅里叶系数的三重相关，如式(3)所示：
|
其中，(·)表示第k段信号的傅里叶系数的三重相关， i1与 i2表示两不同的任意频移。Δ0 = fs/N0表示在双谱区域沿水平和垂直方向的两个频率采样点之间的间隔，fs为信号采样率；λ1和λ2表示任意两不同的频率分量，且 0 ≤ λ2 ≤ λ1，λ1 + λ2 ≤ fs/2。
步骤4：由K段双谱估计的均值可以求出给定样本信号的双谱，如式(4)所示：
|
其中， (·)表示信号的双谱特征矩阵，ω1 = ((2πfs)/N0)λ1，ω2 = ((2πfs)/N0)λ2。
1.4.3 双谱的选择
由1.4.2小节里双谱的计算可以看出，传统非参数直接计算法计算双谱是将每个输入看作多个段，最后计算其均值，以满足信号做双谱计算必须为短时平稳的要求。而对于心音信号而言，其病理性特点并非均匀分布于每一个心音周期内，且部分心音的重要信息分布具有不规律性，这种直接计算均值的方式会使得最后提取的特征信息量有所减少。根据异常心音的这一特点，本文需对1.4.2小节中步骤4进行实验改进，以改进均值计算的缺点，从而确定最优的双谱特征计算结果。
1.4.4 获取双谱特征图
本文将预处理后的每4 s的样本帧进行1 024个点的离散傅里叶变换，可视化结果如图4所示。
图4
心音样本离散傅里叶变换结果图
Figure4.
DFT results of heart sound samples
图4中，样本心音信号的有效信息主要频带在0～200 Hz之间。如式(2)～式(4)所示，对样本做双谱分析后得到一个1 024 × 1 024的二维复数矩阵。
取复数矩阵的实部构成特征矩阵，其数值代表了信号中两个频率分量之间的相关性大小。数值越大，表示两个频率分量之间的相关性越强。结合图4心音的频谱特征和1.4.1小节中双谱的对称性特点截取特征图，最终得到大小为102×102的心音双谱特征图，如图5所示，为正常心音（左）和异常心音（右）样本经特征提取后得到的心音双谱特征图。
图5
正常（左）和异常（右）心音双谱特征图
Figure5.
Bispectral characteristics of normal (left) and abnormal (right) heart sounds
1.4.5 双谱特征图增强处理
图像增强处理涉及将原始图像进行转换，以获得一个通用标准图像的过程。对特征图像增强的好处在于它能够提升后续算法中特征图的适用性[20]。本文中，由于得到的双谱特征图中不同的灰度都可能包含心音的有用信息，因此为避免使用CNN对后续图像特征再提取和分类计算时，由于灰度分布不均导致忽略部分低灰度区域有用信息，本文在分类前，将得到的二维双谱图进行特征图像非线性归一化处理。本文采用对数变换法对双谱特征图进行处理，对数变换法公式如式(5)所示：
|
其中，r为输入，S为输出，v为底数，c为倍数因子，本文中取：c = 255/log(1 + max(r))。根据输入图像中的最大像素值，本文使用式（5）中c的取值来动态调整，以便将输出限定在标准的0～255范围内，以利于保持图像的视觉质量和对比度。
由上述对数变换法变换处理后，二维和三维特征图像变换前后的对比如图6所示。
图6
对数变换前后对比
Figure6.
Before and after logarithmic transformation
当引入了对数变换法增强特征图后，部分低灰度区域的特征被增强，使得双谱特征图的特征更加丰富，避免了在后续CNN中进行卷积运算的时候将有用信息忽略，从而使得心音的特征更明显，以便于提升分类的准确度，增强分类模型的性能。
1.5 分类模型
近年来，深度学习发展迅速，许多研究者用深度学习构建分类模型取得很不错的效果。CNN是深度学习的主要研究方向之一，其在语音识别、自然语言处理、图像识别等方面的应用较为成熟[21]。CNN独特的结构和图像处理方式使得它非常适合于图像分类问题[6]。本文将使用CNN对心音双谱特征进行特征提取和分类。
如表2所示，本文所搭建的CNN结构，由三个卷积层、三个池化层、两个全连接层和一个输出层组成，激活函数选用线性整流函数（rectified linear unit，ReLU），池化选择最大池化方式。首先，展开特征图为一个向量，通过一个含128个神经元的全连接层传输到输出层，使用输出层抑制过拟合，其值设置为0.5。输出层只有一个神经元，而本文研究的是二分类任务，所以只需要一个输出来表示样本属于正类别的概率。输出层，采用S型生长曲线（sigmoid）作为激活函数，计算公式如式（6）所示：
表2
CNN网络模型架构
Table2.
CNN network model architecture
|
其中，σ(∙)是sigmoid激活函数，x为输入的特征图的值，函数将输入x映射到一个0～1之间的概率值。通常，将0.5作为阈值，当sigmoid函数的输出大于0.5时，模型将样本分类为正类别；当输出小于0.5时，模型将样本分类为负类别。
在训练过程中，使用二元交叉熵损失函数来衡量模型的预测与真实标签之间的差异。由于本文所有数据集存在数据类不平衡的问题，因此引入类别权重。通过对损失函数（以符号Loss表示）中的类赋予不同的权重，以使模型更加关注数据量较少的类，从而达到平衡数据的作用，如式(7)所示：
|
其中，y是二元标签，取0或者1；p(yi)是输出属y标签的概率，通常在0～1之间；N为输出个数。本文设正负数据类权重分别为α和β，可根据数据类不平衡程度对其值进行相应调整。本文正样本19 964个，负样本6 109个，于是将α的值设置为1/4，β设为3/4。
2 实验分析
2.1 实验环境与数据说明
本研究软件环境：算法主体编程软件为Pycharm2023专业版（JetBrains Inc.，捷克），编程语言为Python3.9（Python Software Foundation Inc.，美国），所使用CNN分类网络模型在TensorFlow2.0（Google Inc.，美国）框架中实现。
本研究硬件环境：中央处理器i7-13700H CPU @2.4G Hz（Intel Inc.，美国），独立显卡（NVIDIA Inc.，美国）。
本研究所有实验均使用同一数据集且在同一台设备进行。
2.2 评估指标
为了与目前最先进的基于CNN的心音分类算法模型进行比较，本文引入下述评价指标对模型性能进行评价：准确率、灵敏度、特异性。准确率，是指模型在整个测试数据集上正确分类的样本数量与总样本数量之比；它衡量了模型正确分类样本的能力，但在不平衡数据集中可能不够敏感。灵敏度，衡量了模型对于正例样本的识别能力；它表示在所有实际正例中，模型正确识别为正例的比例。特异性，衡量了模型对于负例样本的识别能力；它表示在所有实际负例中，模型正确识别为负例的比例。高特异性，意味着模型能够有效地排除负例，降低了误报的可能性[22]。
2.3 实验结果分析
为证明本文所提分类模型的有效性，分别进行了4组实验。实验一：对前文1.4.2小节心音双谱特征计算中步骤4改进，对比三种计算方法对分类指标的影响，得出更适合的双谱特征图。实验二：对比式(5)中不同对数底数值v对实验结果的影响，从而得到较为理想的v的取值。实验三：通过实验对比不同特征提取方法下心音分类模型的评价指标，证明本文所提基于高阶统计量的心音信号双谱特征提取方法的有效性。实验四：将本文所提心音分类模型与当前已有的心音分类模型做对比，通过各种模型评价指标证实本文所提基于双谱特征提取与CNN相结合的心音分类算法的有效性和优越性。
实验一：如式(4)所示，对心音样本进行双谱分析特征提取计算，将每一个小片段的三重相关结果求和，然后求取平均值。为解决1.4.3小节中所提问题，本文采取对比实验计算三种双谱特征结果，确定更适合的心音双谱特征，分别为均值（mean，M）、方差（variance，Var）和标准差（standard deviation，SD）。该实验所使用的分类网络模型结构如表2所示，实验所使用数据集均为表1所示数据，均使用了底数为10的对数变换法对特征图像增强，实验结果如表3所示。
表3
不同双谱特征结果对比
Table3.
Comparison of results for different bispectral features
表3中，加粗的数据代表同一评价指标下实验结果最好的数据。由实验结果可以看出，基于标准差计算的心音双谱特征更适合作为心音分类模型的特征，这可能是因为病理性的特征分布不均导致的。许多算法在心音预处理过程中仅进行了简单的心音分割，导致某些片段所含的病理性特征较少，而噪声较多。如果对这些异常特征较少的样本片段直接采用传统的双谱计算方法进行均值处理，所得的双谱特征图与正常心音的双谱特征相差无几，从而增加了分类的难度，导致异常样本的识别率较低。这在实验结果中表现为特异性较低，且与实验事实相符。
实验二：不同的特征图增强效果直接影响分类的结果，本文采用非线性图像增强方法中对数变换法对双谱特征图进行增强，其中对数的底数变化对于特征图的增强效果有差异，如图7所示。
图7
不同底数对数变换函数
Figure7.
Logarithmic transformation function with different bases
对数变换法对特征图低灰度值有拉伸效果，避免了后续计算中部分特征的丢失。不同的底数v的值对特征图的增强效果不同，其中v = (r + 1)c/r表示没有进行对数变换增强的特征图，以对比当没有进行图像增强时的模型分类性能。本文实验了5种不同底数，分别是：(r + 1)c/r、8、10、15、20，对心音分类模型的分类准确率的影响。实验所使用的分类网络模型结构如表2所示，5个实验所使用数据均为表1所示数据集，且均使用了基于标准差计算的双谱特征图作为心音特征谱图。实验结果如表4所示。
表4
不同底数对分类结果的影响
Table4.
Effect of different base numbers on classification results
从实验结果来看，没有进行图像增强的心音特征图分类效果最差；而当v值约为10时，模型分类效果最好。如图7所示，当v = 10时，其拉伸曲线与不进行对数变换的虚线相交于点（255，255），这表示v = 10时的对数变换计算对像素值为255附近的像素拉伸或压缩不明显。对特征图进行非线性变换的目的在于增强灰度较低的部分，使得其不至于被后续识别网络所忽略，从而保证特征图的信息量不减少；与此同时，也应尽量不对高灰度区域的值进行过度的改变，从而保证特征图的高对比度；这也是选择非线性方法的一个重要原因。
实验三：本研究对比了另外5种特征提取方式与CNN结合的心音分类实验，所有分类网络模型均采用表1中数据，且均使用了表2所示的CNN结构，其中基于双谱特征的分类实验所使用的双谱特征是基于标准差的分类实验结果，且对双谱特征进行了对数变换增强，其底数为10。所有分类模型评价指标的对比实验结果，如表5所示。其中，一维+CNN，表示直接使用一维心音序列结合CNN对心音信号进行分类识别。
表5
不同特征提取方式分类性能对比
Table5.
Comparison of classification performance of different feature extraction methods
从表5的实验结果来看，基于高阶统计量的心音双谱特征相比于其他基于一二阶心音特征的分类模型来说，分类的准确率更高，这与双谱特征包含有更为丰富的心音特征有关。心音信号通常包含多个频率成分，每个频率成分都对应于心脏的某种生理过程。这些频率成分之间可能存在复杂的相互作用和耦合关系，例如：心室颤动、心室扑动等心律失常疾病。而双谱可以用来分析这些不同频率成分之间的相互作用，以识别各种疾病更具体的特征。
实验四：为证明本文所提基于双谱特征提取和CNN的心音分类算法的有效性，将其与其他4种当前较为先进的心音分类算法做性能比较[8-9, 23-24]，评价指标除前文2.2小节中所述准确率、灵敏度、特异性以外，还增加了模型参量、是否精确分割等两个指标，对各种心音分类模型进行综合评估。5组实验均使用表1所示数据集，且在同一设备进行。实验结果如表6所示，其中“—”代表无数据可查。
表6
不同心音分类模型分类性能对比
Table6.
Comparison of classification performance of different heart sound classification models
从表6的实验结果可以看出，本文所提出的基于双谱特征和CNN相结合的心音分类算法模型性能达到了较为先进的水平，文献[8]虽然也具有很好的准确率和特异性，但是其模型参量较大，而本文的算法模型参量仅为0.04 × 106个，且不需要对心音数据进行精确分割。
综上所述，本文所提心音分类算法模型综合性能优于其他算法，这说明基于高阶统计量的心音特征提取方法与CNN结合的心音分类算法模型在提高心音分类模型的综合性能方面更具优势。
3 结论
在心音分类任务中，心音的病理及生理特征不仅仅包含在幅度和频率中。为更有效地提取心音特征，本文针对低阶统计量（功率谱等）提取心音特征不足的局限，提出采用基于高阶统计量中的双谱对心音信号进行特征提取。高阶谱具有很强的消噪能力，理论上能完全抑制高斯噪声。就双谱的物理意义而言，体现了信号在频谱中的斜度，而以高斯噪声为例，其在频谱中斜度为零，因此双谱可直接对比心音信号中非线性、非对称性的特征成分。本文还采用对数变换法对心音特征图像进行增强，将特征图标准化的同时也能突出更多的信号特征，这有利于CNN对特征图像的分类处理，从而增强模型分类效果。本文实验结果表明，基于双谱特征提取和CNN相结合的方式对心音进行分类，效果显著提升，准确率达到0.91，且模型参量低，无需对原始信号进行精确分割，综合性能优于其他同类算法，有望应用于心脏疾病的辅助检查等方面。但，本文研究也存在一定局限，比如本文只实现了正常心音和异常心音的分类，对于代表具体心脏疾病的异常心音信号并没有做详细的分类，这将成为本课题组今后的研究重点。
重要声明
利益冲突声明：本文全体作者均声明不存在利益冲突。
作者贡献声明：彭利勇负责实验构思与设计、数据收集、数据分析和论文写作；全海燕负责实验指导、论文写作指导、数据审查和论文整体审查。
