随着信息科技水平的发展, 生活学习研究中处理的数据量越来越大, 大数据因此受到了广泛关注.而大数据需要专业技术, 以有效地从大量含有意义的数据中提取有用信息.趋势噪声分解技术就是其中的关键技术.趋势噪声分解是从原始信号里去除噪声找出趋势的过程.趋势表示的是事物或局势的发展动向, 在金融市场当中表示股票、债券、期货、外汇等的运动方向.噪声表示的是不同频率和不同强度无规则地组合在一起的声音或信号.当我们关注于事物或局势的趋势时, 噪声将影响我们对趋势的把握和判断.在日常生活及实际数据中, 噪声对趋势的污染是不可避免的, 严重影响了我们的学习、生活、工作, 甚至是科学研究.因此我们需要从原始被污染的信号中去除干扰的噪声, 找出相关的趋势.而趋势噪声分解的最终目标是将信号分解为趋势和噪声两部分.
对于原始信号的趋势噪声分解, 处理信号的不同, 人们采用的方法也不同.比较常见的有Fourier变换、短时Fourier变换、Wagner-Ville分布、小波变换等分析方法.但是这些方法都不能很好的处理非线性和非平稳的数据.而在实际应用当中, 我们碰到的数据往往具有非线性非平稳的特征, 因此经验模态分解(Empirical Mode Decomposition, EMD)在趋势提取与噪声去除中得到了广泛的应用. 1998年, 美籍华人Huang Norden等人提出了Hilbert-Huang变换(Hilbert-Huang Transform, HHT)的新方法[1]:在创造性的提出本性模态函数(Intrinsic Mode Function, IMF)的基础上, 提出了经验模态分解, 将任意信号分解为有限若干本性模态函数之和, 最终通过Hilbert变换得到Hilbert谱, 从而为复杂的非线性非平稳信号的分析提供了行之有效的工具[2]. 2004年, Wu[3]通过大量的实验, 证实了EMD方法的滤波性质, 为EMD用于去噪提供了理论依据.进一步Flandrin[4]基于Wu的理论基础构建滤波器, 选择产生相应的IMF并重构, 实现了基于EMD的去噪.
EMD算法在用于去噪与趋势提取时, 其一般方法是:在EMD分解出多个IMF后, 由于EMD的滤波性质, 产生的IMF序数越小(越早产生), 频率越高, 序数小的IMF代表高频成分, 序数大的IMF代表低频成分, 而噪声是高频的.那么当信号被噪声污染时, 序数小的高频IMF中一般包含的是噪声; 这样序数较小的高频IMF就作为噪声去除, 而序数较大的就保留作为信号.当然EMD用于去噪也存在着弊端, 有待完善, 首先高频信号也可能包含趋势的能量, 而低频的也有可能包含噪声能量, 按这个方法有可能造成相互的能量损失; 其次到目前为止没有特定的准则来判定舍弃的IMF的数目.
在对EMD去噪的研究中, Boudraa[5]提出了连续均方误差的方法来判别舍弃的IMF数量, 而孙伟峰等[6]在2008年改进了此方法, 较好的解决了这一问题; 在2009年, 文献[7]对比小波分析, 较为详细的研究了EMD算法, 并将其应用于电力系统的趋势提取与噪声去除.文献[8]利用EMD对国际原油价格数据进行趋势分解, 对原油价格波动原因进行分析, 并进行预测.
在利用EMD算法去噪时, 作为判断舍弃IMF数量的准则至关重要, 无论对EMD去噪的算法怎么改进, 最后大多都需要一个精确的准则来判断舍弃IMF数量, Boudraa[5]提出了连续均方误差的方法能够较好的解决这一问题, 但是仍然存在缺陷, 因此本文从噪声的平稳性质出发, 基于非平稳性度量, 构造一个全新的判别准则[9].
非平稳性度量[10-13]方法以遍历论和信息论为思想基础, 研究非平稳性的度量问题, 提出了其度量指标NS, 给出了该指标的近似算法. NS取值于区间[0, 1], 数据流越平稳, 该值越接近于0, 反之越接近于1, 因此NS能够反映出数据的非平稳性程度.
噪声是平稳的, 而趋势则是非平稳的, 如果一组数据是非平稳的, 那么有理由认为这组数据含有趋势, 同样的, 如果一组数据是平稳的, 那么没有理由认为这组数据是趋势.因此可以通过判断数据的平稳程度, 来判断数据是否噪声.而度量一组数据是否平稳以及其平稳的程度, 通过非平稳性度量, 能够很好的进行度量.
本文主要结合非平稳性度量(Non-stationary measure, NS), 研究利用经验模态分解(Empirical Mode Decomposition, EMD)算法进行趋势噪声分解, 提出基于非平稳性度量的准则来判定舍弃IMF的数目.通过数值模拟证明了该准则克服了连续均方误差准则的缺陷, 在不同噪声强度和复杂趋势下, 都能够达到很好的去噪效果.
EMD算法是基于一个简单的假设:任何信号都是由有限若干IMF组成的.
本性模态函数(IMF)的定义有如下两条:
1) 所有极值点的数目与零值点的数目相等或至多相差一个;
2) 局部极大值形成的包络与局部极小值形成的包络均值为0.
依据IMF的定义, EMD通过以下的方法将每一个IMF从数据中分离出来:
首先, 寻找原始信号$X(t)$的极值点, 分别利用三次样条函数拟合$X(t)$极大值及极小值, 形成上包络$e_1(t)$, 及下包络$e_2(t)$, 并计算上下包络的均值
利用原始信号$X(t)$减去上下包络的均值得到$h_1(t)$, 即$h_1(t)=X(t)-m_1$.此时$h_1(t)$, 一般不满足IMF的两个条件, 然后将$h_1(t)$当成数据继续进行上述过程, 直至其变为IMF, 将其记为$h^{(1)}(t)$.
接着原始信号$X(t)$减去$h^{(1)}(t)$, 得到了剩余信号$r_1(t)$, 即$r_1(t)=X(t)-h^{(1)}(t)$.对$r_1(t)$进行与$X(t)$一样的操作, 在$r_1(t)$中继续提取其他的IMF.直至使得某个剩余信号$r_n(t)$为单调函数或常值, 则提取IMF完毕.
此时信号可表示为
在此方法的过程中, IMF的第二个条件上下包络均值为0, 很难到达, 很可能使产生IMF的过程出现死循环, 因此一般引入一个近似条件
当指标${SD}_k$小于某一个预设值时, 认为满足了第二个条件, 一般这个值设置为0.2到0.3.
综合上述过程, 可以归纳算法如下:
$\bullet$STEP1:找原始信号$X(t)$的所有极值点;
$\bullet$STEP2:分别利用三次样条函数拟合$X(t)$极大值及极小值; 形成上包络$e_1(t)$, 及下包络$e_2(t)$, 并计算上下包络的均值$m_1$;
$\bullet$STEP3:利用原始信号$X(t)$减去上下包络的均值得到$h_1(t)$, 即$h_1(t)=X(t)-m_1$;
$\bullet$STEP4:判断$h_1(t)$极值点数目与零点数目是否至多相差一个, ${SD}_k$是否小于0.2, 若都满足, 则为IMF, 若有一个不满足, 继续1-3步, 直到满足为止, 记产生的IMF为$h^{(1)}(t)$;
$\bullet$STEP5:原始信号$X(t)$减去$h^{(1)}(t)$, 得到了剩余信号$r_1(t)$;
$\bullet$STEP6:判断$r_1(t)$是否为单调函数或常数时, 若是算法停止, 若不是, 重复1-5步, 直至使得某个剩余信号$r_n(t)$为单调函数或常值, 算法停止;
$\bullet$OUT:最后信号被分解为有限个IMF与一个剩余项的和, 即$X(t) = \sum\limits_{i = 1}^n {{h^{(i)}}(t)} + {r_n}(t)$.
在提出EMD算法后, Wu及Huang[3]通过大量实验, 证实了EMD方法的滤波性质.由于EMD算法的自适应滤波性质, 分解形成的IMF的频率会随着产生的先后顺序依次减小. Flandrin[4]基于EMD的这个特性, 提出了可以利用产生的IMF构造滤波器组.这种不需要任何基函数的自适应滤波器, 为EMD算法处理非平稳非线性的信号提供了基础.设EMD产生的IMF总数为$n$, ${IMF}_{i}(t)$为产生的第i个IMF, 此处为表述方便, 最后的残差看成产生的最后一个IMF.则其滤波器组简单表述如下:
1) 低通滤波器:去掉EMD最开始产生的$k$个频率较高的IMF, 由剩下的IMF重构信号, 即
2) 高通滤波器:去掉EMD最后产生的$h$个频率较低的IMF, 由剩下的IMF重构信号, 即
3) 带通滤波器:去掉EMD最先产生的$k$个频率较高的IMF和最后产生的$h$个频率较低的IMF, 由剩下的IMF重构信号, 即
4) 带阻滤波器:保留EMD最先产生的$k$个频率较高的IMF和最后产生的$h$个频率较低的IMF, 由这些IMF重构信号, 即
噪声是高频的.那么当信号被噪声污染时, 信号经EMD分解出来的IMF中, 序数小的高频IMF中一般包含的是噪声; 这样序数较小的高频IMF就作为噪声去除, 而序数较大的就保留作为信号.结合EMD的滤波特性, 通过其构造的低通滤波器, 就可以将信号中的噪声去除, 达到去噪的效果.
利用EMD去噪也存在着弊端, 主要的问题有以下两个方面:
1) 高频信号可能包含趋势的能量, 而低频的也有可能包含噪声能量, 按这个方法有可能造成相互的能量损失.
2) 通过EMD低通滤波器进行去噪, 到目前为止没有特定的准则来判定舍弃的IMF的数目, 即低通滤波器中的$k$是多少.
针对EMD去噪弊端的第二个问题, 经典的判定舍弃IMF数目的准则为连续均方误差(Continuous Mean Square Error, CMSE)准则.
首先Wu及Huang[3]在对白噪声进行大量的研究实验知道, 白噪声分解产生的每一个IMF的能量密度$E_k$与其平均周期$\overline{T}_k$乘积是一个常数, 且产生的IMF的平均周期, 也是前一个IMF的两倍.即
其中IMF的能量密度定义为$E_k=\frac{1}{L}\sum\limits_{ i= 1}^L [{IMF}_k(t_i)]^2$, $L$为数据量.
可以看出如果是噪声的分解出来的IMF, 那么随着产生序数的增大, 其能量密度是依次降低的, 根据这个特征, Boudraa[5]提出了连续均方误差的准则, 来判断EMD去噪过程中要舍弃的IMF的数量, 该准则旨在度量两个连续重构信号($\widehat{X}_k, \widehat{X}_{k+1}$)之间的均方误差, 定义为
找到使CMSE值最小的第$j_s$个IMF, 那么就由这个IMF以后的IMF重构信号, 即
其中的分界点$j_s$定义为
也就是寻找能量密度的全局最小值, 作为分界点, 分离噪声IMF及趋势IMF.但是当趋势能量较小时, 包含趋势的IMF有可能比包含噪声的IMF能量密度低, 这样选取全局最小值作为分界点就不太合适, 因此文献[6]改进了这一准则, 选取第一个局部最小值作为分界点, 较好的解决了这个问题, 即
但是进行的改进仍然没有完全解决这一问题, 在趋势的能量较低, 噪声又太强时, 包含趋势的首个IMF有可能比包含噪声的IMF能量低, 因此改进的连续均方误差准则仍然会造成误判.对于连续均方误差准则的缺陷, 本文提出一种全新的判别准则, 来解决这个问题.
噪声是平稳的, 而趋势则是非平稳的, 如果一组数据是非平稳的, 那么认为这组数据含有趋势; 同样的, 如果一组数据是平稳的, 那么不能认为这组数据是趋势.度量一组数据是否平稳, 即数据的非平稳性度量.
非平稳性度量$NS$是一个统计量, 其抽样分布取决于数据序列的分布, 以及初始划分块数$N$, 参数$\lambda$, ${P_0}$ (即$p^*$) (取决于稳定集合假设检验的显著性水平).
噪声的判别实质上等价于下面这样一个假设检验.
原假设${H_0}$:待判定的随机序列$X$为噪声, 而对立假设${H_1}$:随机序列$X$不为噪声, 即带有趋势; 在原假设成立, 并给定参数$N$, $\lambda$, ${P_0}$的前提下, 非平稳性度量值$NS(X)$服从某一个抽样分布, 在这个抽样分布下构造小概率的区间, 如果$X$为噪声, $NS(X)$集中在较小值, 这样就可以构造$NS(X)$大于某个值的小概率区间(单边检验), 若这个小概率为$\alpha$, 由分位数的定义, 这个值为$NS(X)$的$\alpha $分位数, 记为$NS_\alpha$, 即$P(NS(X) > NS_\alpha ) = \alpha $.
如果$NS(X)\leqslant NS_\alpha$, 说明正常噪声情况, 不能认为$X$带有趋势, 判定为噪声;
如果$NS(X)>NS_\alpha $, 小概率事件发生了, 否定原假设, $X$认为带有趋势;
这样$NS_\alpha$就作为判定数据序列是否噪声的一个阈值, $\alpha$是犯第一类错误的概率, 也就是$X$本为噪声, 却认为带有趋势的概率.
在原假设成立, 并给定参数$N$, $\lambda$, ${P_0}$的前提下, $NS(X)$的抽样分布很难得出的, 对于噪声的判定的阈值也就很难得出, 但是可以通过多次抽样得到$NS(X)$的经验分布, 从而确定阈值.
趋势是非平稳的, 如果一组数据是非平稳的, 那么认为这组数据含有趋势; 再进一步, 如果一组平稳的数据加上一组数据是平稳的, 不能认为加上的这组是趋势; 如果一组平稳的数据加上一组数据变得不平稳了, 那么认为加上的这组数据为趋势.
EMD算法对一组带噪声的趋势利用EMD算法分解后, 产生的IMF及剩余项$r_n$, 频率是随着产生的次序依次从高到低递减, 我们知道噪声是高频的又是平稳的, 那么产生次序靠前的高频平稳的IMF认为是噪声的组成部分, 而剩下的IMF及$r_n$就是趋势组成部分, 这样就可以按产生的次序依次将IMF引入组合成噪声, 并观察组合成的噪声的非平稳性(非平稳性由非平稳性度量值$NS$衡量), 在没有趋势引入时这个组合成的噪声是平稳的, 如果加入某个IMF使得组合噪声的变得不平稳了(通过数值模拟得到当$NS$值大于0.1时[9]), 则认为加入的这个IMF为趋势组成部分, 并在此处断开, 以前的IMF组合成噪声, 以后的(包括此IMF)组合成趋势, 此IMF的序号为分界点$j_s$, 即
在去噪中会出现的两类问题, 通过辅助判别都可以较好的解决.
第一, EMD算法对带噪趋势进行分解的时候, 噪声会导致产生虚假的趋势.解决方法:噪声产生的虚假趋势由于在数量级上和噪声一个水平, 所以可以被依次组合的噪声吸收, 即这个虚假的趋势IMF组合在噪声中后, 仍然平稳, 非平稳值很小, 这样就不会把它当成趋势组成部分.利用原方法可以很好的解决.
第二, 噪声掩盖了较弱的趋势, 导致非平稳性度量找不出来, 即在加入某个趋势组成部分的IMF后, 噪声的非平稳值仍然很小, 这样就无法分离出这个趋势.解决方案:噪声掩盖了较弱的趋势, 但不能掩盖其本身的非平稳性, 这个时候我们采用从分界点处回溯的方法, 通过验证IMF本身的非平稳性度量值NS, 找回丢失的趋势.即修正准则为
解决这两个问题对应的方法是矛盾的, 比如虚假趋势通过回溯的方法, 很可能被当成真实趋势, 噪声掩盖的是真实的趋势还是本身产生的虚假趋势, 不好判断.通过多次尝试及实验, 得到了这两个方案应用的经验准则:
1) 在趋势较为简单或噪声较弱时, 容易产生第一类问题, 适合运用第一种解决方案;
2) 在趋势较为复杂或噪声较强时, 容易产生第二类问题, 适合运用第二种解决方案.
进行实验的数据首先前提要是非平稳的, 即非平稳性度量值$NS$较大, 非平稳的数据才会带有趋势.
进行去噪的数据序列的数量, 从减少抽样成本的考虑, 应该是尽可能的少, 但是较少的数据量不能刻画太复杂的趋势, 数据量与趋势的复杂程度有关, 从而对$NS$度量值造成影响; 因此数据量不能过少, 也不能过多, 按照趋势的复杂度与数据量相适应的原则取值; 我们取数据量为700, 并对与其相适应的趋势进行研究.
非平稳性度量参数中初始化分块数$N$, 根据拟合优度检验的经验公式
给出, $L$为数据量, 在数据量为700的时候, $N = 25$.
参数$\lambda $, ${P_0}$, 决定非平稳性度量中稳定集合的假设检验, 在此假设检验给定显著性水平为0.05的时候, 设定$\lambda = 1.96$, ${P_0} = 0.9258$.
信噪比是用来衡量信号与噪声的强度的比值, 本文中分别定义了输入信噪比($SNR_{in}$)及输出信噪比($SNR_{out}$)来衡量模拟去噪的好坏, 为不带噪声趋势(信号)与噪声的标准差之比, 即
其中$x_0$为输入的真实趋势, $\varepsilon$为输入的真实噪声; $\widehat x_{0}$为输出的拟合趋势, $\widehat \varepsilon$为输出的去除噪声.由定义可以看出, 输入信噪比($SNR_{in}$)和输出信噪比($SNR_{out}$)越接近去噪效果越好.
拟合评价值是用来衡量拟合好坏的值, 定义为
$x(t)$为原始数据,
$\widehat x_0(t)$为去噪后的趋势, 这个定义来源于残差平方和, 由于残差平方和在小于1及大于1时会放大及缩小误差, 因此采用绝对值, 并为便于比较, 除以原数据绝对值和, 达到归一化目的.这个值不能太大也不能太小, 过小就过度拟合了, 使得去噪没有意义, 过大就拟合得不好, 使趋势造成损失; 在实验时通常与真实的$ERR$值的数量级进行比较.
数据量及参数设定值, 在不加说明的情况下, 与此设定值一样.
本节考察不同噪声强度下, 用两种准则选择IMF, 从而进行比较.
例1 设置带噪趋势函数
其中${\varepsilon_t} \sim N(0, \sigma^2)$, 且相互独立, 分别加入$SNR_{in}$为2, 1.5, 1.
利用两种准则对组成趋势的IMF进行选取得到去噪结果见表 1.表中分界点指此序号以后的IMF (不包含分界点IMF)组合成趋势, 即若分界点为$j_s$, 则去噪后的趋势: $\widehat X_0 (t) = \sum\limits_{i = j_s + 1}^n {IMF_i(t)} + r_n(t)$, 其中$n$为EMD分解产生的IMF数量.
在三种信噪比下, 得到的两种准则去噪效果对比图为图 1-3, 计算的对比数据为表 2.通过例1, 对两种准则在输入信噪比为2, 1.5, 1 (噪声由弱到强)的情况下, 通过图像与数据进行去噪效果的比较, 非平稳性度量准则在不同强度的噪声下, 都能选取到最佳的IMF, 组合成较好的趋势; 与连续均方误差准则选取一致, 都能达到很好的去噪效果.
本节对前一节的趋势函数做出了一些变换, 使得趋势函数的变化频率加快, 振幅减小, 得到新的趋势函数, 并利用两种准则在不同信噪比下进行计算对比.
例2 设置带噪趋势函数
其中${\varepsilon_t} \sim N(0, \sigma^2)$, 且相互独立.分别加入$SNR_{in}$为2, 1.5, 1时, 利用两种准则对组成趋势的IMF进行选取得到去噪结果见表 3.得到的两种准则去噪效果对比图为图 4-6.计算的对比数据见表 4.
通过例2, 可以看出在对趋势加快变化频率, 减小增幅后, 在信噪比较低($SNR_{in}=1, 1.5$), 噪声强度较大时, 基于连续均方误差准则的判别就会产生错判, 不能达到理想的去噪效果, 而基于非平稳性度量准则的判别, 在较低信噪比的情况下, 依然能够很好的组合IMF, 形成趋势.
通过分析知道, 由于对趋势加快频率, 减少振幅后, 趋势产生的IMF能量减弱, 在较高能量的噪声下(低信噪比), 包含趋势的首个IMF能量(CMSE的值)就会低于包含噪声的最后一个IMF的能量, 利用连续均方误差就会产生错判, 不能很好的通过IMF组合趋势, 而非平稳性度量基于噪声的平稳性质, 因此不受能量的影响, 依然能够正确判断.
本文介绍了一种基于经验模态分解和非平稳性度量的新的趋势噪声分解方法, 并给出了具体的算法.同时进行了相关的模拟实验, 对比分析了改进的基于非平稳性度量准则的趋势噪声分解方法与基于连续均方误差准则下的趋势噪声分解方法.通过模拟分析表明, 当趋势中包含较多小趋势, 并且趋势振幅较小时, 趋势产生的IMF能量将减弱, 在较高能量的噪声下(低信噪比), 包含趋势的首个IMF能量(CMSE的值)就会低于包含噪声的最后一个IMF的能量, 利用连续均方误差准则就会产生错判, 不能很好的通过IMF组合趋势, 而非平稳性度量是基于噪声的平稳性质, 因此不受能量的影响, 可以取得比连续均方误差准则下更好的趋势噪声分解效果.