CN103954443B - 自适应总体平均经验模式分解eemd协助噪声大小确定方法 - Google Patents

Abstract

本发明涉及一种自适应总体平均经验模式分解EEMD协助噪声大小确定方法。本发明包括如下步骤：(1)初始化噪声加入最大值e_max及网格数N，每个网格i(i＝1,2…N)对应噪声大小为e_i＝i/N*e_max；(2)对信号x(t)依次加入幅值大小为e_i的高斯白噪声n_i(t)；(3)求取信号x_i(t)的极大值点所在位置Max(j)、极小值点所在位置Min(k)，得到新序列a_i(j)和b_i(j)；计算新序列a_i、b_i的标准差Sa_i、Sb_i，得到噪声分布特性Q_i与e_i的对应曲线关系；(4)Q_i随e_i迅速下降的点对应的e_r，即为自适应总体平均经验模式分解EEMD协助噪声大小的最优值。本发明减少了人为选择参数的盲目性。

Description

自适应总体平均经验模式分解EEMD协助噪声大小确定方法

技术领域

本发明涉及一种自适应总体平均经验模式分解EEMD协助噪声大小确定方法，属于机械设备故障诊断领域，可以自适应地选择EEMD方法中协助噪声幅值大小最优值，实现机械设备早期故障的有效诊断。

背景技术

随着现代工业和科学技术的飞速发展，工业已成为衡量一个国家科学技术发展的标准。机械设备是工业发展的载体，为工业发展提供关键技术，在国民经济中发挥着越来越重要的作用。同时，机电设备也越来越朝着大型化、复杂化、精密化发展，设备的功能越来越多，性能指标越来越高，其组成与结构越来越复杂，这样势必会使得故障出现的概率大大增加。

由于机电设备工况的复杂多样化，机械故障也越来越复杂，设备的故障特征往往是非平稳、非线性的，传统的机械设备故障诊断方法往往只能针对平稳信号才能得到较好的分解结果，对非平稳、非线性信号束手无策。经验模式分解EMD是针对非线性、非平稳信号而提出的一种基于信号局部极值点的处理方法。但是当信号极值点分布不均时，EMD分解结果往往出现模式混淆现象，在其基础上提出了总体平均经验模式分解EEMD，如公开日为2011年11月23日，公开号为CN102254103A的中国专利中，公开了一种自适应总体平均经验模式分解EEMD筛选次数确定方法，以及公开日为2011年11月23日，公开号为CN102254103A的中国专利中，公开了一种自适应总体平均经验模式分解EEMD筛选次数确定方法。通过对信号加入高斯白噪声，改善信号极值点分布，减少模式混淆不足，已被成功应用于转子、齿轮、轴承等旋转机械故障诊断中。

然而，EEMD分解结果较大程度的依赖于高斯白噪声幅值大小，目前，噪声幅值大小较大程度的依赖于人为选择。为了减少参数人为选择导致的分解结果盲目性，从极值点角度出发，考虑加入噪声后信号的极值点分布情况，在不改变原始信号极值点分布的情况下，最大程度的改善信号极值点，使其分布最为均匀。考虑到当加入的高斯白噪声幅值较大时，噪声振动剧烈，破坏了信号高频成分的极值点分布，EEMD分解结果会产生多余成分，出现模式混淆；当加入的高斯白噪声幅值较小时，噪声振动微弱，无法改善信号低频成分的极值点分布，EEMD分解结果与EMD类似，无法克服模式混淆不足。因此，加入的噪声大小应在未破坏原始信号高频成分极值点分布特性基础上尽可能改善低频成分极值点分布。

发明内容

本发明的目的在于克服现有技术中存在的上述不足，而提供一种减少了参数人为选择带来盲目性的自适应总体平均经验模式分解EEMD协助噪声大小确定方法。

本发明解决上述问题所采用的技术方案是：该自适应总体平均经验模式分解EEMD协助噪声大小确定方法，其特征在于：包括如下步骤：

(1)初始化噪声加入最大值e_max及网格数N，每个网格i(i＝1,2…N)对应噪声大小为e_i＝i/N*e_max；

(2)对信号x(t)依次加入幅值大小为e_i的高斯白噪声n_i(t)，得到x_i(t)＝x(t)+n_i(t)；

其中n_i(t)＝e_i*m_i(t)，m_i(t)是均值为0，且标准差为1的高斯白噪声序列；

(3)求取信号x_i(t)的极大值点所在位置Max(j)、极小值点所在位置Min(k)，得到新序列a_i(j)＝Max(j+1)-Max(j)，b_i(j)＝Min(k+1)-Min(k)；计算新序列a_i、b_i的标准差Sa_i、Sb_i，得到构造信号x_i(t)的噪声分布特性得到Q_i与e_i的对应曲线关系；

(4)根据Q_i与e_i的对应曲线关系，Q_i随e_i迅速下降的点对应的e_r，即为自适应总体平均经验模式分解EEMD协助噪声大小的最优值。

作为优选，本发明所述步骤(4)中：

①对Q_i与e_i用最小二乘法进行拟合得到P_i，P_i的最小值为P_min，提取所有P_y≥1.2*P_min的y值构成集合Y；

②e_r根据所有满足P_i≥P_i+L(i∈Y,i+L∈Y)所确定的最小i值i_min对应的噪声大小e_r＝i_min/N*e_max。

本发明与现有技术相比，具有以下优点和效果：核心是实现了EEMD方法协助噪声大小的自适应选取，减少了参数人为选择带来的盲目性。以极值点分布特性为评价参数，既保证了噪声幅值不会过小导致的无法改善信号极值点分布，又避免了噪声幅值不会过大导致的破坏信号高频成分的极值点分布。

附图说明

图1为自适应总体平均经验模式分解EEMD协助噪声大小确定方法的流程图。

图2(a)为仿真信号及其各个组成部分，(b)为仿真信号极值点与噪声大小的对应关系，(c)为仿真信号极值点与噪声大小的对应关系的最小二乘拟合曲线，(d)为自适应EEMD分解结果，(e)为EMD分解结果。

图3(a)某转子实际振动时域信号，(b)是该转子实际振动信号的频谱，(c)是该实际信号极值点与噪声大小的对应关系，(d)是该实际信号极值点与噪声大小的对应关系的最小二乘拟合曲线，(e)为实际信号的自适应EEMD分解结果，(f)为实际信号的EMD分解结果。

具体实施方式

下面结合附图并通过实施例对本发明作进一步的详细说明，以下实施例是对本发明的解释而本发明并不局限于以下实施例。

实施例。

参见图1至图3，本实施例中的自适应总体平均经验模式分解EEMD协助噪声大小确定方法包括如下步骤：

(1)初始化噪声加入最大值e_max及网格数N，每个网格i(i＝1,2…N)对应噪声大小为e_i＝i/N*e_max。

(2)对信号x(t)依次加入幅值大小为e_i的高斯白噪声n_i(t)，得到x_i(t)＝x(t)+n_i(t)；

其中n_i(t)＝e_i*m_i(t)，m_i(t)为均值为0，标准差为1的高斯白噪声序列。

(3)求取信号x_i(t)的极大值点所在位置Max(j)、极小值点所在位置Min(k)，得到新序列a_i(j)＝Max(j+1)-Max(j)，b_i(j)＝Min(k+1)-Min(k)；计算新序列a_i、b_i的标准差Sa_i、Sb_i，得到构造信号x_i(t)的噪声分布特性得到Q_i与e_i的对应曲线关系。

(4)根据Q_i与e_i的对应曲线关系，对Q_i与e_i用最小二乘法进行拟合得到P_i，P_i的最小值为P_min，提取所有P_y≥1.2*P_min的y值构成集合Y。

(5)e_r根据所有满足P_i≥P_i+L(i∈Y,i+L∈Y)所确定的最小i值i_min对应的噪声大小e_r＝i_min/N*e_max。

根据上述发明内容和图1的自适应选择总体平均经验模式分解(EEMD)噪声大小技术的流程图，首先初始化噪声最大值及网格数，得到各网格对应的噪声大小，然后依次将各网格对应幅值大小的高斯白噪声加入信号中，评估噪声协助信号的极值点分布特性，得到极值点分布与噪声大小的对应关系，最后，对Q_i与e_i用最小二乘法进行拟合得到P_i，P_i的最小值为P_min，提取所有P_y≥1.2*P_min的y值构成集合Y；e_r根据所有满足P_i≥P_i+L(i∈Y,i+L∈Y)所确定的最小i值i_min对应的噪声大小e_r＝i_min/N*e_max。

按照上面选择的最优噪声大小对信号进行EEMD分解，得到的分解结果就是自适应求得噪声大小的总体平均经验模式分解(EEMD)结果。

为了证明上述发明内容的有效性，先仿真一组信号，仿真信号s如图2中(a)所示，由冲击信号c₁、调制信号c₂、谐波c₃组成。对合成信号s利用本发明提出的自适应EEMD方法分解。初始化噪声最大值e＝0.5，网格数1000，得到信号的极值点分布与噪声大小的曲线如图2(b)所示，最小二乘拟合得到光滑曲线如图2(c)。从图2(c)中可以看出，当加入的噪声较小时(e_i<0.008)，噪声振动微弱，极值点分布的离散性随噪声的增加而增大，这是由于加入的噪声振动微弱，仅仅改善了信号中某些低频成分的极值点分布，导致极值点分布呈现较大的离散型；当噪声大小进一步增大，达到e_i＝0.008时，极值点分布特性F_i达到最大值，信号的极值点离散型最大。随后，当加入的噪声较大时(e_i>0.008)，极值点分布的不均匀性随噪声的增加而减小，最后趋于稳定，主要是由于随着噪声较大，振动剧烈，剧烈振动的噪声破坏了高频成分的振动特性，噪声协助信号的极值点全部是由加入噪声引起。因此，根据图2(c)可以确定，最优噪声大小为0.02，将确定大小的高斯白噪声加入仿真信号s进行EEMD分解，分解得到的IMF如图2(d)所示。可以看出，分解得到的前三个IMF分别是冲击、调制、谐波分量，依次对应原始信号中的各组成成分，每个IMF的幅值与构成原始信号的各成分几乎相等，各IMF之间不存在模式混淆现象。对仿真信号进行EMD分解，分解结果如图2(e)所示。分解得到的第一个IMF中，既包括了冲击成分，也含有部分调制成分，很明显分解结果发生了模式混淆现象。此外，第三个IMF相比原始信号的谐波成分，存在严重失真。

同时将上述方法应用在实际数据的分析中。图3(a)为用位移传感器测得的某炼油厂重催机组振动信号，图3(b)为其对应的频谱，频谱中的三个主要频率为f1＝25.35Hz，f2＝97.65Hz与f3＝193.36Hz，分别对应齿轮箱低速轴的转频、机组转频和二次谐波。因此从信号的时域波形及频谱中均无法判别振动过大的原因。利用自适应EEMD对信号进行分析，初始化噪声最大值e＝0.5，网格数1000，得到极值点分布特性与噪声大小的对应关系如图3(c)，最小二乘拟合得到光滑曲线如图3(d)，从图3(d)中可以看出，随着加入噪声的不断增加，极值点的分布不均匀性逐渐减少。在噪声大小为0.02时，信号极值点分布不均匀特性趋于平稳，之后，随着噪声的不断增大，信号极值点分布不均匀性急剧下降，最后趋于平稳，这是由于信号的极值点特性已被加入的噪声破坏。因此，最佳噪声大小为0.02。将确定大小的高斯白噪声加入实际信号进行EEMD分解，分解得到的前三个IMF及趋势项如图3(e)所示。从图中可以看出，IMF1、IMF2、IMF3分别为分解得到的前三个分量，R为趋势项。IMF2对应机组的工频振动，频率为f2＝97.65Hz，IMF3代表齿轮箱低速轴的转频，频率为f1＝25.35Hz。仔细观察IMF1，其中含有周期为0.031s的微弱冲击成分，接近机组旋转周期的3倍(1/97.65Hz＝0.01024s)，表明冲击信号的频率是1/3倍的机组工频，表明轴与轴瓦之间产生了早期碰摩故障。为了进行对比，利用EMD方法对信号进行分解，分解得到的前三个IMF及趋势项如图3(f)所示，可以看出，IMF1中周期性冲击不明显，且IMF2与IMF3已出现了模式混淆，无法判别故障原因。

通过以上仿真信号和实际信号的分解结果对比，可以得到本发明的自适应选取协助噪声大小的EEMD方法能够在一定程度上减少模式混叠，能够分解出物理意义明确的本征模式函数，分解结果比较好、操作方便简单，说明本发明的自适应总体平均经验模式分解EEMD协助噪声大小确定方法可以更好的实现机械设备的故障诊断。

虽然本发明已以实施例公开如上，但其并非用以限定本发明的保护范围，任何熟悉该项技术的技术人员，在不脱离本发明的构思和范围内所作的更动与润饰，均应属于本发明的保护范围。

Claims (1)

1.一种自适应总体平均经验模式分解EEMD协助噪声大小确定方法，其特征在于：包括如下步骤：

(1)初始化噪声加入最大值e_max及网格数N，每个网格i，i＝1,2…N，对应噪声大小为e_i＝i/N*e_max；

(2)对信号x(t)依次加入幅值大小为e_i的高斯白噪声n_i(t)，得到x_i(t)＝x(t)+n_i(t)；

其中n_i(t)＝e_i*m_i(t)，m_i(t)是均值为0，且标准差为1的高斯白噪声序列；

(3)求取信号x_i(t)的极大值点所在位置Max(j)、极小值点所在位置Min(k)，得到新序列a_i(j)＝Max(j+1)-Max(j)，b_i(j)＝Min(k+1)-Min(k)；计算新序列a_i、b_i的标准差Sa_i、Sb_i，得到构造信号x_i(t)的噪声分布特性得到Q_i与e_i的对应曲线关系；

(4)根据Q_i与e_i的对应曲线关系，Q_i随e_i迅速下降的点对应的e_r，即为自适应总体平均经验模式分解EEMD协助噪声大小的最优值；

所述步骤(4)中：

①对Q_i与e_i用最小二乘法进行拟合得到P_i，P_i的最小值为P_min，提取所有P_y≥1.2*P_min的y值构成集合Y；

②e_r根据所有满足P_i≥P_i+L，i∈Y,i+L∈Y，所确定的最小i值i_min对应的噪声大小e_r＝i_min/N*e_max。