双小波非凸稀疏正则化去噪算法研究

马 敏， 王 涛

(中国民航大学 电子信息与自动化学院，天津 300300)

1 引 言

航空发动机工作环境复杂，滑油监测数据易受噪声干扰。噪声不仅会导致信号在外观形态上的变换，而且会影响信号后续处理的清晰度和准确度。因此，有效的数据去噪算法至关重要。

数据去噪是信号处理中的常见问题[1]。Donoho和Johnstone在小波变换的基础上提出小波阈值去噪算法，运算逻辑简单，去噪效果良好，成为主流的去噪算法[2～6,12～19]。对于基于小波变换的去噪算法，由特定的小波基函数生成的小波表示域确定小波系数的稀疏性，阈值函数用于获得噪声信号稀疏表示系数，影响重构信号的稀疏性。因此，小波基函数和阈值函数的选择都会影响数据的去噪效果[7]。针对阈值函数的选取，文献[8]提出了一种带有调节因子的阈值函数，阈值函数具有自适应性，但去噪效果不太理想；秦冬冬等提出一种基于多层阈值函数的去噪方法，阈值函数可根据采样的长度的变化进行调节，但未分析分解层数对去噪效果的影响[9]。针对小波基函数的选取，许罗鹏等牺牲部分对称性，放宽正交性为双正交性，选出了去噪效果更好的小波基函数[10]；郝泳涛等为了提高小波重构精度，在分解和构造时，分别构造不同的小波基函数，去噪效果得到改善[11]。

数据去噪的原则是在滤除噪声的同时，保留原始数据的细节信息。本文利用两个小波基函数同时去噪，重新构造双小波去噪(double wavelet denoising, DWAD)模型，利用交替方向乘子法(alternating direction method of multipliers, ADMM)和Majorized-Minmizetion(MM)方法交替求解目标函数。将该算法应用于滑油监测数据，取得了较好的去噪效果。

2 数据集构造

利用12电极的ECT系统采集航空发动机不同磨损状态下的滑油监测数据[16]，将归一化的数据按时间轴和空间轴依次排列，构造数据集。以某1磨损状态为例，数据集构造方法如图1所示(E1-E2表示激励E1传感器时，E1传感器和E2传感器之间的电容值，数据集时间轴长度节点为1 024，空间轴长度节点为66)。

图1 数据集构造方法

3 双小波非凸稀疏正则化去噪算法

3.1 问题描述

通常，含有噪声的真实信号为：

y=x+n

(1)

式中：y∈RN为观测信号；x∈RN为真实信号；n表示噪声信号，通常假设n服从高斯分布。基于小波变换的去噪模型可以表示为：

(2)

式中：Cy表示观测信号y的小波系数，通过软阈值函数去噪，利用小波逆变换得到去噪信号。

信号可看作不同尺度下小波基的线性组合，观测信号和原始信号都可用2个小波基来表示：

(3)

式中：C1、C2、C1y、C2y为原始信号和观测信号在2个小波域下的小波系数；W1、W2表示小波变换矩阵。

在实验室条件下获得原始ECT信号，通过加性高斯白噪声(additive white gaussian noise, AWGN)对数据破坏作为观测信号。以磨损状态1中E1～E2传感器之间的噪声数据为例，由图2可以看出，在Daubechies(db4)和Symlets(sym4)小波基函数下，噪声数据的小波系数具有不同的分布，这表明可在特定的小波域下滤波，在其它小波域下保留信号的细节。

图2 不同小波基下小波系数

3.2 DWAD算法

为尽可能地保留数据细节信息，提高去噪性能，提出一种双小波非凸稀疏正则化去噪算法，对信号进行非凸优化去噪。DWAD的约束优化问题如下：

(4)

式中：λ1、λ2表示正则化系数；P1、P2表示2个表示域下表示系数的罚函数，表示如下：

(5)

式中：φ表示罚函数；N表示系数数量；参数a控制罚函数的非凸性，且a∈[0,1/λ]。

式(4)是当罚函数的参数a∈[0,1/λ]时，具有约束的凸优化函数。可采用基于凸收敛理论的优化方法来求得式(4)的最优解。ADMM适合求解分布式凸优化问题，根据增广拉格朗日方法，用二次惩罚方法将约束问题转化为无约束优化问题，得到增广拉格朗日量：

(6)

式中：μ>0，式(4)可以通过迭代最小化式(6)中的C1和C2来解决。ADMM中的每次迭代分3次求解，交替求得式(4)的最优解。

(7)

(8)

(9)

(10)

式中：C1n是前一次迭代的最优解；替代函数g1n(C1,C1n)由MM方法设计，定义为：

(11)

式中：

(12)

式(11)中，K1不依赖于C1，参数α应满足：

(13)

小波变换矩阵满足单帧条件。因此，式(7)等同于以下问题：

(14)

最后，通过阈值函数，直接获得公式(7)的最优解：

(15)

式中：阈值函数Th( )由相应的罚函数定义；参数λ1为第一小波域下的正则化参数。

式(8)中最优变量与变换矩阵W2耦合，类似于式(7)，可用MM方法来解决，替代子目标函数为：

(16)

式中：

(17)

由于小波变换矩阵满足单帧条件，式(8)可以简化为：

(18)

使用由特定罚函数指定的阈值函数来生成最优解：

(19)

式中：λ2是第二个小波域下的正则化参数。

不同小波域下原始信号的表示系数具有不同的稀疏性。对于表示系数稀疏性强的小波域，设置小正则化参数。两个小波域的正则化参数设为：

λ1=γλ
λ2=(1-γ)λ

(20)

式中：γ表示比例系数，且每次迭代都会更新，定义如下：

(21)

通过阈值函数交替求解2个子问题，在2个小波域下获得滤波后的小波系数。原始信号的不同细节信息在不同的小波域下保留，采用2个去噪信号差对式(9)中的d值进行更新，2个目标子函数的解可推导出相同的去噪信号。用小波变换矩阵分别得到去噪信号：

(22)

DWAD算法中，小波基函数的选择有2个原则：第1是小波系数的稀疏性，选择L1范数较小的小波基函数作为观测信号；第2是选择两个相关系数较低的小波基函数，有利于细节信息的保留。

对于用于求解DWAD模型的ADMM，参数μ影响算法的收敛速度。为了提高收敛速度，每次迭代都会更新μ：

(23)

4 实验结果分析

4.1 DWAD算法对连续信号的去噪性能分析

首先，测试DWAD算法对于ECT连续信号的去噪性能。将噪声信号的信噪比(SNR)设置为16，利用sym4和db4两个小波基函数分别生成了2个5尺度小波变换。设置迭代次数Iter=20，迭代参数μ=0.8。

使用不同尺度下的通用阈值滤波。在db4和sym4小波域下，采用经典的软阈值方法和罚函数法(对数函数、有理函数和反正切函数)与DWAD算法进行比较。对于软阈值方法，阈值设置为2.5σ。对于罚函数法和DWAD方法，阈值和参数分别设置为λ=3.0，a=1/λ。

图3为软阈值法去噪和无罚函数DWAD算法去噪结果(横轴为数据个数，纵轴为原始数据、噪声数据和重构数据数值)。结果表明，该方法具有较好的去噪效果和平滑度。

图3 软阈值方法和无惩罚函数DWAD的去噪数据

使用信噪比(SNR)来评价去噪性能，均方根误差(RMSE)和结构相似度指数(SSIM)评价原始信号的细节保留能力，测量指标如表1所示。通过使用罚函数，提高了去噪和细节保留的性能。Atan罚函数获得最佳去噪效果，当使用相同的罚函数时，DWAD算法在各项指标上均优于传统的小波去噪方法。采用DWAD算法，信噪比提高约2.3 dB。

Atan罚函数法和带Atan罚函数的DWAD法去噪结果如图4所示，DWAD算法在很好地抑制噪声的同时，很好地保留了原始信号的细节。

4.2 DWAD算法对阶跃信号的去噪性能分析

在haar小波域内，噪声信号具有较小的L1范数，haar和sym4小波基函数的相关系数较低。因此，利用haar和sym4小波基函数分别生成了2个5尺度的非抽取小波变换矩阵。同样，设置参数μ=0.8，Iter=20。

表1 基于db4和sym4基函数的比较算法的测量指标

图4 Atan惩罚方法的去噪信号和Atan的DWAD信号

将硬阈值、软阈值、罚函数法与DWAD算法比较。原始信号与噪声信号(σ=4)分别如图5(a)和图6(a)，通过不同的方法在haar小波域下的去噪信号如图5(b)～图5(d)。为ym4小波域下的相应结果如图6(b)～图6(d)。采用DWAD算法去噪后的结果如图5(e)和图6(e)。DWAD算法的RMSE最低，去噪效果最好，且haar小波域下的结果比sym4具有更佳的RMSE。

图5 haar小波域下不同方法的去噪信号

利用具有不同标准偏差的AWGN损坏信号，测试不同方法对具有不同噪声水平的数据的去噪性能。表2显示了20次迭代实验的去噪信号的平均信噪比(SNR)，结果表明：DWAD算法对噪声水平较低的信号具有更好的性能，DWAD算法对不同噪声水平的数据的平均信噪比提高了至少4.2 dB。

图6 sym4小波域下不同方法的去噪信号

表2 不同表示域的几种算法的平均信噪比

5 结 论

本文提出了一种改进的双小波非凸稀疏正则化去噪算法，在滤除噪声的同时，保留了原始数据的细节信息。针对连续和阶跃的滑油监测数据，DWAD算法均具有良好的去噪性能。本文尚未讨论小波分解层数、参数μ和Iter的设置对该算法去噪效果的影响，这为以后的研究指明了方向，笔者将对其做进一步研究。