汪 梅，王 将，李远成，董立红，马 天，李铭禹

(1.西安科技大学 计算机科学与技术学院，陕西 西安 710054；2.西安科技大学 电气与控制工程学院，陕西 西安 710054)

0 引 言

矿工的不安全情绪可能会导致矿难的发生，识别矿工情绪能够判断其当前的情绪状态，为煤矿的安全生产减少人为因素的安全威胁[1]。脑电信号以其客观性和不易隐藏性的特点，被广泛应用于情绪识别领域[2-3]。但脑电信号由于受到在采集过程中外界环境干扰和采集设备等因素的影响，而被引入噪声导致信号质量下降[4]。这些噪声有可能影响信号的主要特征，对于后续的特征分析和情绪识别有一定影响[5-7]。因此，在利用脑电信号识别矿工情绪的过程中，对于噪声的滤除就显得尤为必要。

经典脑电滤波方法主要包括傅里叶分解或小波分解及重构。傅里叶分解由于其构造函数为周期性的正、余弦波，导致其对非周期性或局部特征较明显的信号处理效果较差[8]。小波分解具有良好的时频分析能力，可以很好地分辨信号的突变部分[9-11]。近年来，经验模态分解及其改进方法在脑电滤波算法研究中日渐增多[12-14]。相比于小波分解，经验模态分解不需要预先设定母小波和分解层次就可以自适应地对非线性、非平稳信号进行分解处理，但此方法缺乏良好的数学理论，并且对采样数据和噪声都很敏感[15-18]。为此，变分模态分解(variation mode decomposition，VMD)以其完备的数学理论支持和较好的噪声鲁棒性，已广泛应用在生物电信号滤波领域。

KAUR等将VMD算法分别与离散小波变换和小波包变换结合起来对脑电信号进行滤波，发现基于小波包变换的VMD法性能更优[19]。DORA等利用VMD提取脑电信号中的眼电分量，并结合回归的方法获得干净的脑电信号[20]。XIAO等通过结合VMD和小波阈值的方法对肌电信号进行滤波，实验滤波效果要优于单一的小波阈值法和经验模态分解法[21]。卢莉蓉等则利用VMD和小波软阈值方法来去除心电信号中肌电的干扰[22]。但是VMD的模态数和带宽分别取决于分解个数K和惩罚因子α的预设值，而且K和α值的大小都会影响滤波效果。

针对上述VMD在对脑电信号滤波时，VMD的分解效果会受到参数K和α选取的影响，导致信号滤波效果差的问题，提出一种优化变分参数与改进小波软阈值重构滤波算法。该算法结合了VMD能提供有效频率划分和小波阈值时频分析能力强的优点，改进了由K和α选取不当造成的影响和传统小波阈值存在逼近程度较差或平滑性不足问题。文中创新点包括以下4个方面：第1，给出乌燕鸥参数优化的VMD算法；第2，扩展相关系数差值比的方法来区分有效分量和含噪分量；第3，利用改进的小波软阈值处理含噪分量得到去噪分量；第4，提出优化变分参数与改进小波软阈值重构的滤波算法。最后，对模拟信号和情绪脑电信号进行滤波处理，检验所提方法的滤波性能。

1 乌燕鸥参数优化变分模态分解

1.1 变分模态分解算法

VMD是一种自适应的信号处理方法，与经验模态分解法相比，VMD提供了更有效的频率划分，可以更好地避免模态混叠等问题[23-24]。非平稳多分量的信号f(t)通过VMD可以得到K个具有特定中心频率ωk的模态分量uk。为了确定uk和ωk，需要构造一个约束变分问题，见式(1)。

(1)

式中 {uk}，{ωk}分别为模态分量集合和中心频率集合。

引入拉格朗日乘数λ求解上式，可得到

(2)

(3)

(4)

(5)

式中n为迭代次数；α为惩罚因子；τ为保真系数；ε为收敛精度判据。

1.2 乌燕鸥优化算法

1.2.1 迁移行为

迁移过程中乌燕鸥应该满足以下3个条件。

条件1：避免碰撞

通过引入参数SA来计算乌燕鸥个体的新位置C(z)，以避免相邻个体之间的碰撞

C(z)=SA×P(z)

(6)

式中C(z)为不与其他乌燕鸥发生碰撞的新位置；P(z)为乌燕鸥个体当前位置；z为当前迭代次数。SA为避免碰撞的变量参数，计算公式如下

SA=Cf-(z×(Cf/Maxiterations)

(7)

z=0,1,2,…,Maxiterations

(8)

式中Cf为控制变量，使得SA的值从Cf线性递减到0，设置Cf值为2，Maxiterations为最大迭代次数。

条件2：确定相对间距

在避免碰撞之后，当前个体会向最优个体靠拢，以获得更好的位置

M(z)=CB×(P*(z)-P(z))

(9)

式中M(z)为当前个体位置P(z)与最优个体位置P*(z)的相对间距。CB是一个随机变量，能够使个体更好地进行全局搜索。CB的计算公式为

CB=0.5×Rand

(10)

式中Rand是介于[0，1]的随机数。

条件3：向最优个体靠拢

确定无碰撞位置C(z)和相对间距M(z)后，下面开始向最优个体位置靠拢，以到达新位置D(z)。

D(z)=C(z)+M(z)

(11)

式中D(z)为乌燕鸥迁移行为后的位置。

1.2.2 攻击行为

在乌燕鸥定位猎物后，它们会改变自身的飞行速度和角度对猎物进行攻击，从而在空中产生螺旋状运动轨迹。其在x′，y′和z′平面的运动行为的数学模型为

x′=Radius×sin(i)

(12)

y′=Radius×cos(i)

(13)

z′=Radius×i

(14)

Radius=u×eiv

(15)

式中Radius为螺旋运动的半径；i为在[0，2π]范围内的螺旋角；e为自然对数的底数；u和v为常量，设置u和v的值均为1。

乌燕鸥攻击行为后位置的计算公式为

P(z+1)=(D(z)×(x′+y′+z′))×P*(z)

(16)

式中P(z+1)为乌燕鸥攻击行为后的位置。

1.3 乌燕鸥参数优化变分模态分解

在VMD中，参数K和α的选取将会影响最终的分解效果。若K值过大，会造成分解过度，导致相邻模态中心频率的间距较近；而如果K值太小，会造成分解不足，导致分解失去了实际意义[25]。类似的，如果α选取不恰当，也会出现模态混叠的问题。因此，合理选择K和α的值对于信号分解结果的准确性尤为重要。

由于乌燕鸥优化算法(sooty tern optimization algorithm，STOA)具有寻优能力强、精度高等特点[26]。利用STOA优化VMD以确定最佳参数组合[K，α]，并将文献[27]所提的包络熵作为适应度函数，将VMD参数的优化过程转化为利用STOA寻求最小包络熵值的过程。通过VMD将原始信号分解为K个模态分量，如果模态分量中包含的噪声分量越多，其与原始信号相关的特征信息就越不明显，则包络熵越大；反之，包络熵越小。图1为乌燕鸥参数优化变分模态分解(STOA-VMD)算法流程。

图1 STOA-VMD算法流程Fig.1 Flow of STOA-VMD algorithm

STOA-VMD算法的具体步骤如下。

步骤1：初始化STOA参数，并随机形成一个数量为N的乌燕鸥种群，设置迭代次数为z，寻优维数为2，则乌燕鸥个体位置可表示为：P(z)={P1(z)，P2(z)}，其中，P1(z)=K，P2(z)=α。

步骤2：以当前乌燕鸥个体位置对应的K和α作为VMD的输入参数对信号进行分解。

步骤3：计算不同位置的乌燕鸥个体对应的包络熵值，更新得到当前最小包络熵值。

步骤4：判断是否达到迭代终止条件，如果达到，寻优停止，输出最优的K和α，否则，令z=z+1，并更新种群位置，返回步骤2继续进行迭代优化。

2 模态判断和改进小波软阈值算法

2.1 相关系数差值比的模态判断

通过STOA-VMD将原始信号分解为K个模态分量，需要从中判断并重构合适的分量来对信号进行滤波。相关系数可以度量原始信号与模态分量之间的相似性，相关系数值越大则相似性越强，反之则越弱。文中分别将与原始信号相关系数较大和较小的模态分量称为有效分量和含噪分量。文献[28]利用固定的相关系数值来区分这2类分量，难以对不同信号作出灵活的调整，具有一定的局限性。

文中提出一种相关系数差值比的判断条件，根据相邻模态与原始信号相关系数的差值比来确定有效分量和含噪分量的临界点，判断条件如下

(17)

(18)

式中Rm为原始信号f(t)与第m个分量um(t)之间的相关系数，E[·]和D[·]分别为数学期望和方差。

根据临界点判断条件找出有效分量和含噪分量之间的临界点um(t)。由于模态分量按照频率从低到高排列，而有效分量的频带大多分布在低频部分。因此，对临界点um(t)之前的有效分量保留并利用所提改进小波阈值处理其余含噪分量。

2.2 改进小波软阈值的分解重构算法

根据小波阈值滤波的基本原理可知[29]，小波基、分解尺度的确定和阈值函数的构造都会影响其滤波效果。软阈值和硬阈值是2种常用的阈值函数，但软阈值法会造成重构信号与原信号之间逼近程度较差的缺点；硬阈值法的不连续则会引起重构信号平滑性不足的问题。

(19)

为了直观地反映形状调节因子β的作用效果，选取T=1,β=0.5,2,5,10,20,30，比较改进阈值函数与软、硬阈值函数的特性，如图2所示。通过调整β的值，改进阈值函数可以在软、硬阈值之间变动，这也使得它在实际应用上更为灵活。随着|wz|的增加，改进阈值函数越来越接近于硬阈值函数，这就弥补了软阈值法的缺点。同时，改进阈值函数具有高阶可导性，可以克服振荡，提高信号的平滑度。

图2 阈值函数特性比较Fig.2 Comparison of threshold function characteristics

3 优化变分与改进软阈值重构算法

在上述理论基础上，提出优化变分与改进软阈值重构算法(STOA-VMD-IWS)，如图3所示。

图3 STOA-VMD-IWS算法流程Fig.3 Flow of STOA-VMD-IWS algorithm

首先，通过STOA-VMD算法找到最佳参数组合[K，α]，利用优化的VMD分解原始信号得到一组模态分量{u1,u2,…,uK}。然后，利用相关系数差值比的判断条件来区分有效分量和含噪分量，并用IWS算法对含噪分量进行滤波。最后，重构有效分量和去噪分量，实现信号滤波。

4 实验与结果分析

4.1 参数优化变分模态实验与结果分析

实验选取公开脑电数据集DEAP中一段预处理后的脑电信号(electroencephalogram，EEG)作为研究对象，采样率为128 Hz，采样个数为600。在其中分别加入强度为-10～10 dB且间隔为5 dB的高斯白噪声，得到5种含噪EEG信号。

下面以含有10 dB高斯白噪声的EEG信号为例，利用STOA-VMD算法寻找最佳参数组合[K，α]。设置乌燕鸥种群数为30，最大迭代次数Maxiterations为10，K的范围为[2，10]，且K为整数，α的范围为[500，4 000]。图4为STOA-VMD算法的寻优过程。

图4 STOA-VMD算法的寻优过程Fig.4 Optimization process of STOA-VMD algorithm

从图4可以看出，当迭代数为5时，最小包络熵值为0.969 2，此时[K，α]值为[10，3 348]。因此，以K=10，α=3 348作为VMD的输入参数对EEG信号进行分解，结果如图5所示。

图5 STOA-VMD算法的分解结果 Fig.5 Decomposition results of STOA-VMD algorithm

可以看出，u1～u10的中心频率之间相互独立，且并未发生模态混叠的现象，证明了STOA-VMD算法的分解结果较好。为了验证STOA参数寻优的准确性，将原始EEG信号频谱中对应的主频率值与u1～u10的中心频率值进行比较，结果见表1。

表1 EEG的主频率和u1～u10的中心频率比较Table 1 Comparison of the main frequency of EEG and the center frequency of u1～u10

从表1可知，u1～u10的中心频率值基本对应于EEG信号的主频率值，能很好地反映EEG信号的频率特性。此外，如图5所示，u1～u10之间没有混叠，证明STOA-VMD算法参数优化的有效性。

4.2 相关系数差值比实验与结果分析

在得到模态分量u1～u10后，计算了u1～u10和原始EEG信号之间的相关系数，并得到相邻模态分量相关系数差的绝对值见表2。

表2 相邻模态分量相关系数差的绝对值Table 2 Absolute value of the correlation coefficient difference between adjacent modal components

利用公式(17)的模态判断条件计算相关系数差值比，可以得到：|R6-R5|/|R5-R4|=0.076 5，|R5-R4|/|R4-R3|=9.136 7，则u5为有效EEG分量和含噪EEG分量的临界点。

4.3 改进小波阈值滤波实验与结果分析

根据相关系数差值比的模态判断条件，得到有效EEG分量和含噪EEG分量的临界点为u5，将临界点前的模态分量作为有效EEG分量保留，对其余含噪EEG分量利用IWS算法处理后与有效EEG分量重构得到滤波后EEG信号，如图6所示。

图6 信号滤波效果Fig.6 Filtering effects of the signal

经所提STOA-VMD-IWS算法滤波后噪声得到了明显地抑制。为了验证所提算法的优势，文中将其与传统VMD算法、STOA-VMD和小波硬阈值滤波算法(STOA-VMD-HWTF)、STOA-VMD和小波软阈值滤波算法(STOA-VMD-SWTF)进行了对比实验，并以信噪比(signal-to-noise ratio，SNR)和均方根误差(root mean square error，RMSE)作为评价指标，对比结果见表3。

由表3可知，传统VMD算法的SNR较低，RMSE较高，其滤波效果较差。在传统VMD算法的基础上，尽管STOA-VMD-HWTF算法和STOA-VMD-SWTF算法的滤波效果得到了一定的改善，但所提STOA-VMD-IWS算法在SNR和RMSE方面效果更好，脑电滤波效果更佳。

表3 EEG信号滤波评价指标对比Table 3 Comparison of EEG signals filtering evaluation indexes

4.4 EEG信号STOA-VMD-IWS滤波实验分析

为了进一步验证所提滤波算法对实测EEG信号的有效性，根据情绪二维模型理论，实验选取煤矿井下工作环境相关的3类情绪图片，即积极、消极和中性，在E-Prime 2.0系统上对被试进行情绪诱发，并利用Neuroscan公司的NuAmps设备完成对被试EEG信号的采集。

下面以在消极情绪图片诱发下得到的一段EEG信号为例，图7为STOA-VMD-IWS算法对该EEG信号的滤波效果。

图7 4种EEG信号的滤波效果Fig.7 Filtering effects of 4 EEG signals

其中，图7(a)的EEG信号分别来自FP1,C3,O1和T4电极，可以看出EEG信号具有大量的毛刺，说明其含有大量的高频噪声。图7(b)是利用所提STOA-VMD-IWS算法对这4个电极处EEG信号进行滤波处理后的信号。经所提STOA-VMD-IWS算法处理后，高频噪声得到了有效抑制，也保留了EEG信号的大部分特征，可以得到其幅值变化及波动情况。

由于实验无法预先得到“干净”的EEG信号，导致无法利用SNR和RMSE来评价所提算法的滤波效果，故文中利用噪声抑制比(noise suppression ratio，NSR)来评价所提算法和传统VMD算法对4种EEG信号的滤波效果，计算结果见表4。NSR越大，说明滤波效果越好，其计算公式如下

表4 4种EEG信号的NSR计算结果Table 4 NSR calculation results of 4 EEG signals

(20)

从表4可知，2种滤波算法对C3和O1的EEG信号滤波效果要优于其他2种EEG信号，且所提算法对4种EEG信号的滤波效果均优于传统VMD算法，证明所提算法在脑电滤波方面的有效性。

5 结 论

1)给出变分模态分解参数的乌燕鸥优化方法，即利用乌燕鸥优化算法对变分模态分量个数和惩罚因子进行优化，得到参数优化的VMD算法，解决了VMD有效分解问题。

2)扩展相关系数差值比判别方法，用于判别VMD分解后的有效分量和含噪分量，以便对含噪分量做进一步处理，解决了模态分量选取问题。

3)对小波软阈值函数进行改进，用于对含噪分量的小波分解与重构，实现对含噪分量的去噪，解决了硬阈值不连续引起重构信号平滑性不足和软阈值重构信号的逼近度较差问题。

4)提出优化变分参数与改进小波软阈值重构的滤波算法，对模拟信号和情绪脑电信号进行滤波处理，信噪比最大提高了4.568 3 dB，平均提高了3.284 7 dB；均方根误差最大降低了0.169 1，平均降低了0.069 5，解决了传统VMD算法对脑电滤波效果较差问题，有助于提高脑电的信噪比和矿工情绪识别的准确率，对减少因矿工不安全情绪导致的矿井安全事故的发生具有重要意义。