一种改善的HHT端点效应抑制方法

　　摘要：
　　经验模态分解（EMD）是Hilbert-Huang 变换（HHT）的核心算法。在EMD 分解中，利用三次样条函数拟合信号上下包络在数据两端不可避免地会出现发散现象，即端点效应。端点效应严重影响EMD 分解质量。在解决该问题已有算法基础上，提出了镜像延拓与多项式拟合结合的算法，使端点处极值取值有据可依，镜像延拓更加完善。实验结果表明提出的方法能够有效抑制端点效应，提高EMD 分解质量。
　　
　　关键词： 数字信号处理；EMD 端点效应； 多项式拟合； 镜像延拓
　　
　　中图分类号：TN911.72
　　
　　0 引言
　　
　　Hilbert-Huang 变换是一种适合处理非平稳、非线性信号的新型信号处理方法[1]。HHT采用经验模态分解方法（EMD）对信号进行平稳化处理，产生一系列具有不同特征时间尺度的固有模态函数（IMF），每个IMF 表示原始信号的一个固有振动模态，能很好地体现非平稳信号的局部特性。基于这些IMF 分量进行Hilbert 变换，可以得到每一个IMF 随时间变化的瞬时频率和振幅，最后求得的 Hilbert 谱能准确反映出该物理过程能量在时间空间上的分布规律。HHT 方法的质量很大程度上取决于EMD 的分解质量，而EMD 分解采用三次样条插值来获得信号的瞬时平均，使得该方法存在特殊的端点效应，在运用EMD 对非平稳信号进行分解时，数据的两端会产生发散现象，并且这种发散结果会逐渐向内“污染”整个序列，使得分解结果严重失真[5]。为了改善HHT 的端点效应，提出镜像延拓和多项式拟合相结合的方法，避免了单独采用多项式拟合存在端点误差较大及单独采用镜像延拓处理端点不是极值点的序列效果不佳的问题。
　　
　　1 EMD 分解方法
　　
　　在 HHT 中，为了计算瞬时频率，定义了IFM，每个IMF 是单分量的幅度或频率调制信号。它必须满足两个条件：（1）在整个数据段内，极值点和过零点数目必须相等或至多相差一个；(2)在任意时刻, 由局部极大值点和局部极小值点构成的上下包络平均值为零。对于不满足 IMF 的信号，用EMD 方法对其分解。分解过程如下：
　　
　　（1）设时间序列为x(t) ,确定x(t)的局部极大值集max x 和极小值集min x 。
　　
　　（2）分别根据max x 和 min x 作三次样条插值，拟合上下包络为u(t)和v(t)。则上下包络的平均曲线为m(t)：
　　
　　m(t) = [u(t) + v(t)] / 2 (1)(3)求得原始信号x(t)与局部极值差值，记为11 h ：
　　
　　11 h = x(t) m(t) （2）(4)若11 h 满足IMF 的要求，则令其为x(t)的第一个IMF 成分。否则以11 h 代替x(t)重复以上三步，直到1k h 满足IMF 要求。求得与原信号差值并代替原信号，作：
　　
　　1 c = 1k h ，1 1 r (t) = x(t) c ， 1 x(t) = r (t) (3)（5）重复以上四步，直到n r 或n c 比预定值小或( ) n r t 成为单调函数时，x(t)的EMD分解结束。信号变为各IMF 分量和一个趋势项之和：
　　
　　各 IMF 按瞬时频率由高到低分解得到，每个IMF 随信号的本身变化而变化，因而EMD方法是自适应的信号处理方法。但是由于分解过程中采用三次样条插值，使得EMD 分解产生端点效应。
　　
　　2 端点效应的抑制
　　
　　2.1 EMD 方法的端点效应
　　
　　应用 EMD 方法每一次筛选过程中，上、下包络是由信号的局部极大值和极小值通过三次条插值拟合给出的，但由于所分析信号的长度有限，信号的两端不能确定就是极值，则在三次样条插值时端点效应不可避免。在信号高频分量，由于时间尺度和极值点距离小，端点效应仅局限于信号两端很小部分；但对于低频分量，由于时间尺度和极值间距离大，对于多次分解的EMD 会使端点效应放大，使得IMF 没有实际意义。
　　
　　仿真信号为：y = cos(0.04π t) + cos(0.08π t) + cos(0.01π t)，不经过端点处理直接利用三次样条插值求得曲线上、下包络线如下图所示：可见包络曲线在端点处出现严重的失真，甚至有上下包络相互交织的趋势，这必然会对IMF 分解产生影响。
　　
　　2.2 镜像延拓与多项式拟合结合对端点效应的抑制
　　
　　2.2.1 镜像延拓
　　
　　镜像延拓就是在信号的左右两侧具有对称性的极值点上各放一面镜子，以此极值点为对称点，对镜面内的信号进行延拓，得到长度为两倍镜内信号长度的周期性信号[7]。对称镜像延拓的数学表达式为：
　　
　　实际中没必要对信号进行无限制延拓，只需要延拓后的端点效应不影响考察区间[0,T]即可。
　　
　　图 2 表示经过镜像法端点抑制处理后的包络，没有出现上下包络交织的趋势，且包络更加平滑。但由于未对端点值本身采取其它处理手段，拟合仍缺乏一定准确性。　　
　　
　　2.2.2 多项式拟合
　　
　　对于离散点，可以拟合一个函数表示这些数据的变化趋势，称为数据拟合。当拟合函数为多项式时，称为多项式拟合。求解多项式拟合的一般过程为[8]：
　　
　　(1)根据具体问题，确定拟合多项式的次数。
　　
　　2.2.3 镜像延拓与多项式拟合结合抑制端点效应
　　
　　由于当端点非极大极小值时，就无法确定它的终点值，而直接利用镜像延拓如果去掉端点或者直接作为极值处理，又会造成数据丢失或产生误差。利用镜像延拓与多项式拟合结合的方法，可以避免这一问题的发生。算法的执行过程为：（1）取出极值点序列在最左端的三个极值，拟合多项式函数，计算出函数对应于最左端序列的函数值，把此函数值作为极值点序列在此端点处的函数值；同理计算出极值点序列在最右端的函数值作为极值点序列在右端点处的函数值。（2）把新的端点近似值加入极值点序列。（3）对调整端点值后的数据进行镜像延拓处理。（4）EMD 分解，求解IMF。
　　
　　该方法在对数据镜像延拓前，先对数据进行多项式拟合预处理，使得端点处的取值有据可依，虽然只是近似处理，但正如赵进平等人所认为的[3]“延拓是为了提供一种条件，使得包络完全由端点内的数据确定”，并为后续的镜像延拓提供保证。经过两种方法结合的包络曲线如图3 所示：
　　
　　图 2 与图3 对比，可见上包络在端点处有明显的变化，有向上延伸的趋势。图4 为同一信号多个周期的包络曲线，可以看到在图3 端点处的下个周期内，将出现极大值点，上包络应该有向上递增趋势，可见图3 中包络预测的正确性。
　　
　　3 信号仿真及分析
　　
　　在信号仿真实验中，取频率分别为25HZ 和45HZ 的两个频率比较相近的信号叠加作为仿真信号： y = 2sin(50π t) + 5sin(90π t) ，分别用镜像延拓和本文提出的方法抑制端点效应，EMD 分解结果如下两图所示：
　　
　　从EMD 分解结果图上可以看出，两种方法均能在IMF1,IMF2 两个分量上基本分解出原始信号信息，即频率分别为25HZ 和45HZ 的两个正弦信号。但只经过镜像延拓后的EMD分解出6 个IMF 分量，而后者只分解出代表原信号信息的2 个分量。说明端点效应可能在前者EMD 分解时产生了虚假分量，尽管IMF3~IMF6 频率变化很小，但依然会一定程度上影响分解准确性。为了进一步分析，继续对二者得到的IMF 分量做边际谱，得到结果如下：
　　
　　从边际谱中可以看到两种方法都能基本正确的分解原始信号，得到频率分别为25HZ 和45HZ 的分量信息。但从边际谱图上可以看到，后一种方法分解质量较前一种有明显提高。
　　
　　如图7 所示，在小于25HZ 的部分伴随有虚假分量，而图8 则要好得多。这也就解释了经过后者端点延拓处理的EMD 分解出来的IMF 会少4 个分量，是因为前者端点效应产生了虚假分量。
　　
　　4 结论：
　　
　　介绍了EMD 分解中端点效应对IMF 产生的不良影响，在镜像延拓法基础上提出了镜像延拓和多项式拟合结合抑制端点效应的方法。仿真结果表明该方法能够正确分解IMF，提高EMD 分解质量，更加有效的抑制端点效应。