快速自适应经验模态分解方法的基本原理及其性能评估

周　义， 李鸿光

(上海交通大学 机械系统与振动国家重点实验室，上海　200240)



快速自适应经验模态分解方法的基本原理及其性能评估

周义， 李鸿光

(上海交通大学 机械系统与振动国家重点实验室，上海200240)

摘要：经验模态分解是一种有效的信号分解方法，尤其是针对非平稳非线性信号。然而，随着研究的深入，学者们发现该方法中存在着诸多弊端。根据Bhuiyan的研究，提出了一种针对一维信号的快速自适应经验模态分解方法。通过大量的数值仿真，证明这种方法不但能克服传统方法的弊端、得到高质量的分解结果，还能大幅度地提高计算效率。

关键词：经验模态分解；快速自适应经验模态分解；数值仿真

传统的信号处理方法大致可以分为两类：时域统计分析和傅里叶变换。然而，这些方法都基于一个共同的假设，处理对象是平稳线性信号。然而，工程实际中遇到的都是非平稳、非线性的瞬态信号，因此如果仍使用传统方法，往往会得到错误的结果。为此，学者们提出了许多针对非平稳非线性信号的时频分析方法，如小波变换。但是小波变换需要人为地选择基函数，只有在基函数的形状与信号性质时，才能产生高幅值小波系数，这种非自适应性极大地限制了这种方法的发展。

Huang[1]提出了一种自适应信号处理方法，经验模态分解(Empirical Mode Decomposition，EMD)。该方法基于信号本身的局部特征时间尺度，能够将一个非平稳非线性信号分解成为一系列完备的、近正交的分量，称为“固有模态函数”(Intrinsic Mode Function，IMF)。IMF作为EMD方法的基函数，是嵌入在原信号中的固有振荡模态。它不需要提前人为决定，而是由信号本身在分解过程中决定。因此，经验模态分解是从根本上脱离傅里叶变换，而从信号本身进行分析处理，具有完全的自适应性、无监督性。

经过十几年的发展，EMD的理论研究逐渐深入，在一维信号处理方面的应用更为广泛。Huang等主要建立了EMD的基本框架，分析了基本原理，比较了EMD算法和小波以及其他信号分析方法的区别。其他学者在理论框架的基础上，成功地将EMD应用在地震物理学[2]、生物医学[3]、设备故障诊断[4-5]等领域。

然而，随着研究的深入，学者发现EMD存在诸多弊端：①终止准则：Rilling等[6]将终止准则中的SD，变为双阈值准则。②端点效应：郭明威等[7]使用基于时间尺度的LS-SVM端点延拓方法，蔡艳平等[8]使用Lyapunov指数预测模型进行边界延拓，时培明等[4]提出一种波形特征匹配延拓与余弦窗函数相结合的方法。③模态混叠：Wu等[9]提出集合经验模态分解，消除模态混叠现象。④运算效率：Rilling等[10]又把单变量EMD扩展到双变量EMD，减少计算量；Damerval等[11]提出基于Delauney三角和分段三次多项式插值方法提高EMD的分解速度。

尽管如此，还没有一种算法，可以综合考虑上述四个问题，只是针对某一个问题展开讨论。因此，本文根据学者Bhuiyan的研究[12]，提出了一种旨在提高效率的新算法——快速自适应经验模态分解(Fast and Adaptive Empirical Mode Decomposition，FAEMD)，同时指出此算法能在不同程度上优于目前现有的改进方法，能较好的同时解决上述几个问题。

1传统经验模态分解

经验模态分解，作为时频分析方法的基础，可以将一个非平稳非线性信号分解成多个被称为“固有模态函数”的振荡信号和趋势项，所有IMF按频率从高到底排列，最后一项为趋势项。

固有模态函数是满足单一振荡模式以及单分量信号物理解释的一类信号。单分量信号是相对于多分量信号而言的，只有单分量信号的瞬时频率才有物理意义。Huang定义的固有模态函数必须满足以下两个条件[1]：①在整个数据范围内，极值点和过零点的数目必须相等或最多相差一个；②在任何时刻，由极大值点定义的上包络线和极小值定义的下包络线的平均值为零。

Huang认为，任何信号，无论是平稳还是非平稳信号，都能通过不断迭代分解出满足条件的固有模态函数，这个迭代过程即为“筛分(Sifting)”。筛分具有两个目的：一是消除模态波形的叠加；二是使波形轮廓更加对称。每一轮筛分过程都将分离出1个固有模态函数，然后将余量继续执行筛分操作，当全部固有模态函数提取出来后，就剩下一个余量或趋势项，原始信号x(t)可表示为：

(1)

式中：ci(t)表示第i个IMF分量，r(t)表示筛选到最后剩下的趋势项。

2快速自适应经验模态分解

2.1基本原理

Bhuiyan[12]提出了快速自适应二维经验模态分解，完全突破了传统经验模态分解的实现模式。基于此，我们做了大量改进工作，提出了针对一维信号的快速自适应经验模态分解方法。通过使用顺序统计滤器替代样条插值函数获取上下包络线，并在保证分解质量的前提下，将控制迭代次数，不仅降低了程序的计算时间，结果在许多指标上也比传统EMD出色。

顺序统计滤波器(Order Statistics Filter, OSF)是一种非线性滤波器，最早用于图像处理，增强图像的平滑度。其思路是：以极值点间距的最小值作为滤波器宽度对信号进行顺序统计滤波，求上下包络；再用均值滤波器做平滑处理。

假设原信号含有n个极大值点，OSF统计每个极大值点与相邻极大值点的最小距离，并存入数组dadj-max中，显然该数组有n个数据；同理，对于极小值点，也能得到一个数组dadj-min。可以通过如下四种准则确定OSF窗口宽度w：

w=d1=min{min{dadj-max},min{dadj-min}}

(2)

w=d2=max{min{dadj-max},min{dadj-min}}

(3)

w=d3=min{max{dadj-max},max{dadj-min}}

(4)

w=d4=max{max{dadj-max},max{dadj-min}}

(5)

式中：max代表取数组中所有元素的最大值；min代表取数组中所有元素的最小值。分别定义式中四种准则为1型、2型、3型和4型，四个所得值之间的关系为：d1

确定OSF窗口宽度后，用宽度为w、平行于时间坐标轴、相切于极值点的线段，初步绘制包络线，形成上下包络线：

(6)

(7)

式中：i,j代表提取第i阶IMF中的第j次迭代；uij(t)为上包络，lij(t)为下包络；hmaxij(t)代表极大值谱，hminij(t)代表极小值谱；zt是窗口宽度为w的线段。

在源信号上(下)方绘制了若干个线段，但这些线段之间有的没有接壤，有的接壤却不是平滑过渡，因此还需要进行均值平滑滤波。这种平滑滤波在方程所得结果上直接进行，推荐使用与次序统计滤波滤波器同样的窗口宽度，具体操作如下所示：

(8)

(9)

式中：w′是均值平滑滤波器的窗口宽度，令w′=w。这一步平滑滤波操作是一种算术平均滤波，使数据内的局部变量更加平滑。

此外，还要设置具体的单次筛分迭代次数，一般1~3次就可以保证精度。其余步骤与传统EMD一致。

以仿真信号为例，图1为仿真信号通过FAEMD在获取包络均值的过程。显然通过顺序统计滤波器得到的包络均值存在许多一阶不可导点。所以我们在获得包络均值后，对包络均值用算数平均的方式做平滑处理，获得图2所示光滑曲线。

图1　仿真信号在获取包络均值的过程Fig.1 The process obtaining mean envelop for simulated signal

图2　均值平滑滤波处理过程Fig.2 The process of average smoothing filtering

2.2窗口宽度设定

根据大量的仿真计算可知，式(2)~(5)中：1型窗口尺寸最小，4型窗口尺寸最大；使用3型和4型提取第i+1阶IMF，窗口尺寸会大于第i阶IMF；使用1型和2型提取第i+1阶IMF，窗口尺寸有时可能并不会大于第i阶IMF，这种现象有悖于经验模态分解的基本原则，即后阶IMF一定是含有频率较低的成分。因此，有必要额外添加一个操作，令后面的窗口尺寸大于前面的，例如可以令当前窗口尺寸为前一个窗口尺寸的1.5倍。尽管这个操作可能没有必要，但它可以充分保证经验模态分解的基本原则。

一般情况下，建议使用统一的OSF窗口宽度绘制上下包络。如果分别确定上下OSF窗口宽度，可能会遇到如下问题：①源信号中极小值点足够多而极大值点只有一个，会导致数组dadj-max为空，无法绘制上包络面，却没有满足终止准则，不能终止筛分迭代的过程。②对于同一阶IMF的提取，使用相同宽度的OSF，有助于从源信号提取尺度相近的成分；而如果宽度不一致，得到的IMF分量中尺度信息可能相差很远，这就是 “模式混叠”现象。

2.3小结

综上所述，FAEMD与原始EMD的最大不同之处在于：

(1)传统EMD方法均使用插值技术绘制包络曲线，插值法本身就是一个复杂的迭代过程，需要大量的计算资源；而FAEMD突破传统的插值法，使用估计法绘制包络曲线，这种方法尤其适用于含噪信号，可以快速绘制包络线。

(2)EMD方法采用筛分的方法提取IMF分量，每一次筛分都需要几十次甚至上百次的循环迭代，极度耗时；而FAEMD在保证分解质量的前提下，把单次筛分的迭代次数降低到5次以内，甚至是1次，大大提高了计算速度。

3方法的性能对比实验

虽然EMD发展了10余年，但仍有许多问题需要解决：①终止准则;②端点效应;③模式混叠;④运算效率。

3.1终止准则

过多的筛分会导致IMF分量变成纯粹的频率调制信号，而幅度趋于恒定。因而存在筛分过程的终止准则问题。为了保证IMF在幅值和频率上都有物理定义，必须定义一个筛分过程的停止准则。这个准则可以通过限制标准差的大小来实现。标准差SD通过两个连续的处理结果来计算得出：

(10)

式中：T代表信号长度，hij(t)代表提取第i阶IMF过程中的第j次迭代所得信号，hi(j-1)(t)代表提取第i阶IMF过程中的第(j-1)次迭代所得信号。SD称为“筛分阀值”，一般取0.2。如果计算所得SD若小于此值，筛分过程停止，认为为第一阶IMF；否则，还要继续循环迭代，直到满足该条件为止。

①课程内容按照学科知识点的逻辑顺序组织。为了促进学生主动思维能力的发展，教师设定每个知识单元学习目标和完成标准，指引学生跟随任务完成学习。如在教学三相异步电动机的正反转控制时，教师事先提出预习思考题：异步电动机的转动原理是什么？如何实现异步电动机的正转和反转？使学生带着学习目标学习并思考。

终止准则会对分解质量和分解速度产生极大的影响：如果SD值过大，虽然分解速度很快，但分解精度变粗，每个IMF会含有过宽的频率成分，而不再是局部窄带的单分量信号；反之，如果SD值过小，虽然能够保证分解精度，但运算时间过长。对于不同工况下的工程信号，找到其合适的SD值是比较困难的，往往需要多次使用不同的SD值，人工地寻找合理的分解结果，但这种“试探”的过程，与小波分解中基函数的寻找相类似，都破坏了方法的自适应性。

而本文提出的FAEMD方法，使用顺序统计滤波器绘制包络线，窗口宽度w的产生是基于信号极值点。经过大量的实验可知，在w尺度下，无须迭代或较少次迭代即可获得一个IMF，直接避免了过度筛分导致的幅值同化。

3.2端点效应

在使用EMD处理信号时，端点效应是常遇到的一个问题，特别是在处理长度较短的信号时，端点效应会成为影响分析结果精度的原因之一。

EMD的一个关键步骤是根据极值点作三次样条插值获得包络线。除非信号两端均为极值点，否则不能简单地认为端点就是极值点，这将导致插值产生的包络线存在拟合误差。由于EMD的筛分过程需要反复求取信号的包络线，而端点处极值的不确定性导致每次求包络时都存在拟合误差。随着误差不断累积，第一阶IMF的端点处就会存在较大误差。而第二阶IMF的分解是建立在原信号减去第一阶IMF的残余项的基础上进行的。所以第一阶IMF中的误差会传递给第二阶IMF，导致第二阶IMF有更大的误差。以此类推，随着分解的进行，误差会不断的积累并从端点向内部传递。严重时会使分解后的IMF失去意义。

构造一个简单的信号函数如下：

(11)

使用传统的EMD方法处理仿真信号，分解结果如图4所示。可以看到，分解结果，不能如实反映原信号的两个单分量成分，而是出现许多不能解其物理含义的分量。这就是因为，边界效应的影响不但在边界处，还会随着分解的进行，深入到信号内部，最终导致分解失败。

图3　仿真信号及其真实分量Fig.3 The simulated signal and its real components

图4　仿真信号的传统EMD分解结果Fig.4 The decomposition results by applying traditional EMD to the simulated signal in Fig. 3

因此，使用本文提出的FAEMD方法，处理图3所示仿真信号，所得结果如图5所示。从结果可以看出，FAEMD很好地从原信号中分离出两个单分量成分。

图5　仿真信号的FAEMD分解结果Fig.5 The decomposition results by applying FAEMD to the simulated signal in Fig. 3

3.3模式混叠

在进行信号的筛分过程中，有时会出现某个IMF分量中包含不同特征时间尺度或相近时间尺度的成分分步在不同的IMF中，这种情况称为模式混叠。

图6　含有间歇性分量的仿真信号Fig.6 The silumated signal with intermittency compoents

模式混叠归因于信号的间歇现象，它的产生和极值点的选择有关。图6所示的仿真信号，由正弦波s1、同周期发生的高频三角脉冲s2和趋势项s3组成，这是一种简单却典型的可产生模式混叠的信号。由于信号通常含有间歇型分量、脉冲干扰或噪声，使得基础成分(频率较低)的某个局部出现高频振荡。这种局部振荡和基础成分各自的包络线会很不相同，因此在绘制整体包络时，特别是在高低频成分的交界处会造成混乱现象，也就是模式混叠。它不仅能扭曲信号的时频分布，还能掩盖各IMF分量蕴含的真实物理意义。

为了解决模式混叠问题，Wu等提出了一种噪声辅助的数据分析方法——集合经验模态分解(Ensemble Empirical Mode Decomposition,EEMD)。 EEMD的思路是通过对信号添加高斯白噪声使得间歇的信号连续化，再进行EMD分解，最后对多次EMD分解结果做集合平均。由于白噪声随机产生，因此根据其统计特性，可以通过足够多次的试验抵消。

EEMD在使用时有两个参数需要设置，即高斯白噪声的幅值k和算法执行EMD的总次数M：

(1)k越小，越有利于分解精度的提高，但当小到一定程度时，可能不足以引起信号局部极值点的变化，就不能解决模式混叠问题，k太大则会淹没原始的信号特征，使得分解失去意义。建议k取原信号幅值标准差的0.01倍~0.5倍。

(2)M增大，显然噪声影响会减少，但M过大，会导致运算时间过长。当M取100~200次时，残留噪声引起的误差一般不足1%。

图7为EMD算法的分解结果，得到IMF分量C1-C4及趋势项。只有C4和趋势项具有真实的物理意义，分别代表了原始信号中的正弦分量和趋势项。C1出现严重的失真，多出C2和C3这两个无法解读其含义的成分。C1和C2都含有高频和低频成分，可以认为它们都是s1和s2的模式混叠。由于传统EMD无法将原信号中的间歇成分s2正确的分解出来，导致C1夹杂其他模态成分产生严重的失真，且剩余的间歇组分需要在后续的IMF中不断的被部分分解出来，使得信号中多出了C2和C3这种无意义的成分。由于C4幅值比前3个IMF大得多，所以受到的干扰较小。

图7　运用传统EMD算法分解结果Fig.7 The decomposition results by applying traditional EMD to the simulated signal in Fig. 6.

图8为FAEMD算法分解结果，得到IMF分量C1-C4及趋势项。res依旧保持了原有的趋势，C4较好的表现出S1的特征。前三项均反映出了S2的特征，有效的分解出了间歇信号。从能量的角度看，C1、C2和C3幅值依次减小，C2、C3可以分别认为是高频信号在上次分解中产生的能量泄漏。它们都代表了C1中损失的能量，这解释了C1中幅值小于3的原因。

图8　运用FAEMD算法分解结果Fig.8 The decomposition results by applying FAEMD to the simulated signal in Fig. 6.

FAEMD方法确实能高质量地分离出信号中的各个成分，有效抑制模态混叠的发生，客观展现了信号的真实物理内涵。之所以可以有效地避免间歇信号与低频信号混淆，是因为当高频分量存在时，窗口宽度较小，此时滤波法可以较准确的刻画低频信号从而对其保留，这样低频的信号可以留在下次分解

使用EEMD方法，虽然也可以消除模式混叠现象，得到类似于FAEMD方法的结果，但相比比下，EEMD 不足之处明显：

(1)EEMD运算时间较长，由于EEMD中包含了多个EMD循环，运算时间加长，处理数据效率低下。

(2)每个循环中引入的噪声不同，得到IMF数量可能不同，造成最终分析结果可能包含更多的IMF项，分散了真实分量应有的能量。

(3)由于EEMD所得的每个IMF是由M个EMD循环产生的IMF求平均所得，导致分析结果可能不能严格满足原本IMF定义。FAEMD不会有噪声引入的误差问题，也不会因为噪声幅值过大导致信号中微小的特征被淹没。

3.4运算效率

为解决EMD方法中存在的诸多弊端，学者们先后提出了很多技术。虽然这些弊端均得以解决，但却派生出另一个问题：运算时间大幅上涨。①EMD方法的筛分思想，需要不断地进行迭代和插值，这需要耗费大量的计算资源；②为了避免端点效应而采用的边界延拓技术增加了需要处理的数据量；③为了消除模态混叠效应而采用的集合经验模态分解法，需要数百次的EMD循环，因而增加数百倍运算时间。

对于任何一个信号都可以将其分解成为若干个固有模态分解和一个残余分量，其中每个固有模态分量包含了信号的各种不同频段成分，且它们在频域近似相互正交。为了评估EMD分解结果的优劣，文献[1]引入了“正交因子(Orthogonality Index，OI)”。在工程实际中的分解，OI值为0是几乎不可能的，只能通过改良EMD算法来获取尽可能低的OI值，OI值越低，代表分解的效果越来。一般来说，低于0.1的OI值是可以接受的。

我们以实验中测得的旋转机械松动故障实验信号作为测试数据，分别用EMD，EEMD和FAEMD三个MATLAB程序进行处理，得计算时间和正交因子。从表1中，可以清楚看到，FAEMD方法不仅具有最少的计算时间，而且得到了最低的正交因子，也就是最高质量的分解结果。

表1　传统EMDs与FAEMD计算时间及正交因子

4结论

FAEMD算法不仅能够显著提高原始EMD的计算速度，同时能很好的避免EMD的四个主要问题：终止准则、边界效应、模态混叠和运算效率。由于FAEMD自适应的IMF分解过程，可以有效避免原始EMD中IMF筛分过度所导致的幅值无意义问题，并无须考虑筛分的终止条件，实现了真正的自适应。由于FAEMD采用的全新包络均值求解方法，信号在边界处不会出现过冲或欠冲的拟合失真问题，避免了严重的端点效应。另外FAEMD在处理容易出现模态混叠现象的信号时，分解结果明显优于现有最好的EEMD算法。

参 考 文 献

[1] Huang N E, Shen Z, Long S R, et al. The empirical mode decomposition and the Hilbert spectrum for nonlinear and non-stationary time series analysis[J]. Proceedings of the Royal Society of London. Series A: Mathematical, Physical and Engineering Sciences, 1998,454(1971): 903-995.

[2] 陈林，宋海斌. 基于经验模态分解的地震瞬时属性提取 [J]. 地球物理学进展，2008，23(4): 1179-1185.

CHEN Lin, SONG Hai-bin. Seismic instantaneous attribute extraction based on empirical mode decomposition [J]. Progress in Geophysics, 2008, 23(4): 1179-1185.

[3] 李永勤，王清，邓亲恺. EMD及其在生物医学信号处理中的应用研究[J]. 生物医学工程学杂志, 2005, 22(5): 1058-1062.

LI Yong-qin, WANG Qing, DENG Qin-kai. Research on EMD and its application in biomedical signal processing [J]. Journal of Biomedical Engineering, 2005, 22(5): 1058-1062.

[4] 时培明, 丁雪娟, 李庚,等. 一种EMD改进方法及其在旋转机械故障诊断中的应用[J]. 振动与冲击, 2013, 32(4): 185-190.

SHI Pei-ming, DING Xue-jun, LI Geng, et al. An improved method of EMD and its applications in rotating machinery fault diagnosis [J]. Journal of Vibration and Shock, 2013, 32(4): 185-190.

[5] 苏文胜, 王奉涛, 张志新, 等. EMD 降噪和谱峭度法在滚动轴承早期故障诊断中的应用 [J]. 振动与冲击, 2010, 29(3): 18-21.

SU Wen-sheng, WANG Feng-tao, ZHANG Zhi-xin, et al. Application of EMD denoising and spectral kurtosis in early fault diagnosis of rolling element bearings [J]. Journal of Vibration and Shock, 2010, 29(3): 18-21.

[6] Rilling G, Flandrin P, One or two frequencies? the empirical mode decomposition answers[J]. IEEE Transactions on Signal Processing, 2008, 56(1): 85-95.

[7] 郭明威, 倪世宏, 朱家海, 等. 振动信号中HHT/EMD端点延拓方法研究[J]. 振动与冲击, 2012, 31(8): 62-66.

GUO Ming-wei, NI Shi-hong, ZHU Jia-hai, et al. HHT/EMD end extension method in vibration signal analysis [J]. Journal of Vibration and Shock, 2012, 31(8): 62-66.

[8] 蔡艳平, 李艾华, 石林锁，等. EMD端点效应的改进型混沌延拓方法及其在机械故障诊断中的应用[J]. 振动与冲击, 2011, 30(11): 46-52.

CAI Yan-ping, LI Ai-hua, SHI Lin-suo, et al. Processing method for end effects of EMD based on improved chaos forecasting model and its application in machinery fault diagnosis [J]. Journal of Vibration and Shock, 2011, 30(11): 46-52.

[9] Wu Z, Huang N E. Ensemble empirical mode decomposition: a noise-assisted data analysis method[J]. Advances in Adaptive Data Analysis, 2009,1(1): 1-41.

[10] Rilling G, Flandrin P, Gonalves P, et al. Bivariate empirical mode decomposition[J]. Signal Processing Letters, IEEE, 2007, 14(12): 936-939.

[11] Damerval C, Meignen S, Perrier V. A fast algorithm for bidimensional EMD. Signal Processing Letters[J]. IEEE, 2005,12(10): 701-704.

[12] Sharif B, Adhami R R, Khan J F. Fast and adaptive bidimensional empirical mode decomposition using order-statistics filter based envelope estimation[J]. EURASIP Journal on Advances in Signal Processing,2008,2008:728356.

Basic principle of a fast and adaptive empirical mode decomposition and its performance evaluation

ZHOUYi,LIHong-guang

(State Key Laboratory of Mechanical System and Vibraiton, Shanghai Jiaotong University, Shanghai 200240, China)

Abstract:Empirical mode decomposition is an effective signal decomposition method, especially for non-stationary and non-linear signals. But it is found that there exist a lot of drawbacks along with the deepening of studies. Therefore, a one-dimensional signal processing method, called fast and adaptive empirical mode decomposition, was proposed according to the Bhuiyan’s study. It has been proved by numerous numerical simulations that the method can not only overcome the drawbacks of traditional methods to obtain high quality decomposition results but also enhance the computing efficiency.

Key words:empirical mode decomposition; fast and adaptive empirical mode decomposition; numerical simulation

中图分类号:Tp12；Tp13.3

文献标志码：A

DOI:10.13465/j.cnki.jvs.2016.03.003

通信作者李鸿光 男，博士，教授，博士生导师，1972年生

收稿日期：2013-12-19修改稿收到日期：2014-01-28

基金项目：国家自然科学基金资助项目(11372176)

第一作者 周义 男，博士生，1985年生