一种新的基于NIG模型的四元数小波图像去噪方法

朱 芳，刘 卫



一种新的基于NIG模型的四元数小波图像去噪方法

朱 芳1，刘 卫2

（1. 安徽新华学院通识教育学院，安徽 合肥 230088；2. 中国科学院合肥智能机械研究所，安徽 合肥 230031）

基于正态逆高斯模型和快速双边滤波，本文在四元数小波变换域提出一种新的图像去噪算法。首先将噪声图像进行四元数小波变换分解，其次应用快速双边滤波算法处理低频子带，高频子带采用正态逆高斯模型对其进行建模，并结合最大后验概率估计准则推导出相应的阈值函数，最后结合最优线性插值得到的阈值算法实现图像去噪。对提出的算法进行实验仿真，通过与已有优秀去噪算法比较，结果显示本文方法取得了不错的视觉效果，且在峰值信噪比和平均结构相似性上都得到一定的提高。

四元数小波变换；双边滤波；正态逆高斯模型；最大后验概率估计；最优线性插值；图像去噪

0 引言

当今社会已经进入数字化的信息时代，超大容量、实时地、高效率地获取各种有用信息已成为当代的典型特点。大量数据显示，人类有70%以上的信息都是来源于图像，自然图像在获取和传输过程中都会受到不同程度噪声的污染，如何快速有效地去除噪声已成为当今学者研究的重点。经典的去噪方法主要有空间域滤波和变换域滤波，都在一定程度上实现了噪声去除，但在降噪和保留细节折中方面难以令人满意。小波变换作为新型的时频分析方法，具有多尺度、多分辨率、低熵性及去相关性等特点，其中Donoho基于小波变换提出一系列的阈值去噪方法[1-2]。由于小波变换不具有平移不变性，使得去噪后的图像出现伪吉布斯现象，为此，1998年Kingsbury等提出双树复小波变换（DT-CWT）[3]并广泛应用于图像去噪中，DT-CWT达到近似平移不变性，同时方向特性也增强。Chang[4]等人利用广义高斯模型建模原始图像的小波系数，提出BayesShrink去噪算法。Portilla[5]等人则借助一种高斯尺度混合模型实现图像去噪。2002年，Sendur[6]等基于小波系数尺度间关系提出Bishrink算法。这些方法一定程度上达到了图像去噪的效果，但仍旧存在一些局限性。Cho D[7]在2005年基于Gaussian、GGD、Non-Gaussian提出了多元广义高斯模型，但去噪过程中涉及到较为复杂的参数估计。2006年，Cunha[8]等人提出了非下采样Contourlet变换，但NSCT计算复杂度较高，所需的运行时间较长，同时Guo[9]等人提出了Shearlet变换，但该方法不具备平移不变性，在图像融合时容易在奇异点附近产生伪吉布斯现象。由于传统的实小波变换缺乏平移不变性和方向选择性，同时双树复小波变换（DT-CWT）具有近似平移不变性，但相位信息不够丰富，在图像信号局部特征表征时容易产生相位歧义，2008年，Chan[10]等提出了四元数小波变换（QWT），是基于二维Hilbert变换给出的一种新的多尺度分析工具，具有近似平移不变性和3个相位信息，克服了实小波与复小波的不足，且广泛应用于图像分割、边缘检测和图像重构等。

常见的基于变换域的图像去噪方法主要集中在处理变换后的系数上面，处理的原则一般都是选择合适的统计模型来逼近分解后的系数，其中基于最大后验概率（MAP）准则的Bayesian估计是经典的方法之一。一般情况下，变换域中分解系数的先验概率密度函数中含有未知参数，其中Laplacian密度函数中只含有一个未知参数，是描述分解系数的标准模型，广义Gaussian[11]密度函数中涉及两个参数等，然而这些方法在方差不同的条件下不能准确建模分解后系数，因为它们忽略了系数层间和层内之间的关系。正态逆高斯模型（NIG）的密度函数中含有4个参数，其模型可以用来描述任意形状的曲线，可以准确建模不同程度拖尾的图像分解系数，因而用其作为QWT分解系数的先验模型是合适的。2011年，贾建[12]等人运用正态逆高斯模型（NIG）[13]建模NSCT变换后的系数，在去噪效果上取得了一定成效。

本文提出一种基于快速双边滤波[14]和NIG模型的QWT变换域图像去噪算法，充分考虑图像经过QWT变换后低频和高频部分，因为图像经过多尺度几何变换后，噪声主要集中在高频子带，但低频子带中仍旧含有一些噪声。新算法将QWT变换与NIG模型结合，对变换后的低频子带采用快速双边滤波算法，高频子带的系数用NIG模型进性建模，并利用Bayesian最大后验概率估计推导出阈值函数，最后结合最优线性插值阈值算法（OLI-shrink），实现图像去噪的目的。实验结果表明本文所提方法可以取得较好地去噪效果。

1 基于NIG模型的四元数小波变换（QWT）

1.1 四元数小波变换

W. R. Hamilton[15]在1843年基于复数理论提出四元数概念，也可作为特殊的Clifford代数。设四元数＝r＋i＋j＋kk且r、i、j、kÎ，r为四元数的实部，i、j、k为四元数的虚部，i、j、k为正交虚数单位且满足下式：

四元数的共轭记为：*＝r－qi－jj－kk，其范数记为：

四元数小波变换的思想是将图像进行多分辨率分解，得到不同空间、频率的子图像，再对系数选择合适的方法进行处理。设二维信号(,)的四元数解析信号[11]为：

设一维小波函数的尺度函数和小波函数分别为h()，h()：

记(,)＝h()h()＋g()h()＋h()g()＋g()g()为2(2;)空间的四元数小波尺度函数，同时记：

为2(2;)空间的四元数小波函数。则对于"(,)Î2(2;)，令c,k,m＝((,),,k,m(,))，d,k,m＝((,),,k,m(,))（＝1, 2, 3,,,Î），则称d,k,m（＝1,2, 3）为(,)的离散四元数小波变换。在文献[16]中详细讲解了双数复小波及四元数小波的分解和合成算法。

1.2 QWT域图像去噪算法

多分辨率分析是图像去噪中常用的经典方法，在以往的图像处理算法中，很多学者只单纯处理小波变换后的高频子带部分，而忽略低频子带系数认为其不含任何噪声，使得处理后的图像视觉效果不是很理想。因此，本文对图像进行QWT分解后，兼顾低频和高频子带系数特征，首先对低频子带采用快速双边滤波[17-18]算法；然后在高频子带基于MAP估计准则用NIG模型对分解系数进行建模，得到相应的阈值公式；最后结合OLI-Shrink算法实现图像去噪。

1.2.1 NIG模型

NIG模型[12]主要是由逆高斯和不同均值高斯分布混合构建而成，理论上NIG模型可以描述任何形状的曲线，则可以作为QWT变换系数的先验分布，且其概率密度函数记为：

由于图像通常由一些平滑区域和细节丰富区域组成，对于平滑区域，其四元数小波变换系数趋近于零值，对于纹理丰富区域，其变换系数会具有较大的幅值，因此非高斯分布特性是图像变换系数的一个重要统计特征，它们呈现出零均值，重拖尾的特点。NIG模型可以灵活选择模型参数用来准确描述具有不同拖尾形状的数据属性，因此NIG模型可以作为QWT高频系数的先验分布并对其进行建模。为了检验使用NIG模型作为四元数变换后系数分布的先验模型的合理性，我们分别使用纹理结构较为丰富的图像（Barbara）及纹理较少、对比度较低的红外图像（Plane）作为拟合示例，图1是结果图。图中显示的是Barbara图像及Plane图像经3层QWT分解的第3层水平方向的实部系数直方图分布结果以及分别用NIG和GGD拟合这个概率密度函数结果。从图1中我们可以看出，NIG模型可以准确地对这一高频系数分布进行拟合，特别是重拖尾部分。NIG模型分布在零值附近呈现出尖锐的峰值，在峰值两端表现出很长的拖尾特点，这与高频系数的分布特点非常接近，这种高度拟合情形也同样出现在其他的高频子带中。在其他测试图像中也出现了相似的系数分布直方图拟合结果，这再次证明了我们假定NIG为图像的系数先验分布的合理性。

图1 不同先验模型拟合结果

1.2.2 贝叶斯最大先验概率估计

设原始图像(,)，含噪图像(,)，(,)是均值为零，方差为2的高斯白噪声，且满足＝＋。经过QWT变换后，得到系数：

＝＋(6)

由假设知噪声服从高斯分布，则其概率密度函数为：

通过求导并结合文献[15]得出的估计解为：

利用(9)式进行图像去噪时需要在每个子带估计参数和，可以自适应的在所有子带中去除噪声，同时(9)式类似于经典的软阈值函数，阈值函数2的值依据QWT变换后系数。

1.2.3 参数估计

用(9)式进行图像去噪时需要估计NIG模型的参数、及噪声方差2，其中2可以用第一层小波分解对角子带估计得到：

式中：HH1为第一层小波分解对应的高频子带，参数、根据不同的分解系数进行估计。

1.3 阈值算法

基于小波变换的阈值去噪主要通过设定一个特定的值，将大于阈值的系数去除以达到去噪的效果，可以应用到很多变换域中，比如四元数小波变换域。经典的阈值函数有硬阈值函数，软阈值函数等，使用硬阈值公式时令分解后系数小于特定阈值的部分为零，硬阈值准则为：

软阈值准则为：

软阈值函数比硬阈值函数更加实用有效，但对于大的系数很难选择合适的最优阈值，Fathi和Nilchi[19]应用最优线性插值给出一种新的阈值算法（OLI-Shrink）：

1.4 算法的主要步骤

Step 1：对含噪图像进行3层QWT变换，得到高频子带和一系列的低频子带；

Step 2：低频子带采用快速双边滤波算法进行处理；

Step 3：采用(10)式估计噪声标准，正态逆高斯模型的参数、采用(11)式估计；

Step 4：利用(9)式估计阈值函数中的变量，再通过(14)式估计系数；

Step5：最后进行逆QWT变换，得到去噪图像。

2 实验结果与分析

为了验证本文提出的去噪算法在去除加性高斯白噪声时的可行性和有效性，我们分别对自然图像Lena、Barbara、红外图像Plane及遥感图像（Pentagon）进行测试，加入零均值、方差为2的高斯白噪声。在仿真实验中，将本文所提的去噪算法与Mihcak等人[20]提出的LAWML算法，非下采样Contourlet变换[21]（NSCT）算法、基于双树复小波（DT-CWT）的双变量收缩（BiShrink）[22]算法、四元数小波变换域的Bayesian阈值算法（BQWT）[23]以及SURE-LET[24]进行比较；然后使用峰值信噪比（PSNR）、平均结构相似性[25]（MSSIM）以及去噪后的视觉效果结合起来评价本文方法。其中峰值信噪比定义为：

式中：、分别为原始图像和降噪图像；2为图像大小，下面给出实验结果。

图2～图4为含高斯噪声标准差为20的自然图像（Barbara）、红外图像（Plane）及遥感图像（Pentagon）采用不同方法得到的去噪效果图。从图中可以看出，相比较于其他去噪方法，本文方法表现出更加优秀的抑制噪声的性能，得到的图像视觉效果更好。LAWML方法和SURE-LET方法去除噪声不够彻底，图像中仍存在一些噪声；NSCT方法可以有效地去除图像噪声，但是纹理信息丢失严重，边界信息不能很好地得到保留；BiShrink方法和BQWT方法得到的去噪图像出现了严重的块效应，在边界处存在较为严重的伪Gibbs现象；本文所提算法可以有效地去除噪声，图像的边缘及纹理信息可以有效地保留，伪Gibbs现象得到很好的抑制。

图2 噪声方差为20的图像（Barbara）去噪后局部放大图

图3 噪声方差为20的红外图像（Plane）去噪后局部放大图

图4 噪声方差为20的遥感图像（Pentagon）去噪后局部放大图

从图5可以看出，随着噪声方差的增加，PSNR值存在一个下降的趋势。对于不同的含噪图像，在所有方法中本文方法取得了最大的PSNR值。与LAWML及SURE-LET方法相比，本文方法获得的PSNR值有较大的提高，特别是对红外图像Plane和遥感图像Pentagon，本文方法获得的PSNR值的优势较为明显；相比于其他几种去噪方法，本文方法去噪后的PSNR也有不同幅度的提升，这主要是因为所提算法处理高频噪声时的阈值合适，此外对低频噪声也有相关对策，这些处理都为抑制噪声作了贡献。

平均结构相似性可以用来评估去噪算法保留边界信息的能力，因此我们在实验中也使用了平均结构相似性作为客观度量。不同方法处理不同噪声等级图像得到的客观数据如表1所示。在表1中，本文去噪后的平均结构相似性是最高的，这表明了本文方法去噪后的图像保持原图像的结构能力是最高的。

3 结束语

本文在QWT小波域基于快速双边滤波和NIG模型提出一种新的图像去噪算法，充分考虑小波变换后低频和高频子带中的噪声情况。低频子带采用快速双边滤波算法进行处理，去除低频子带中的噪声；高频子带采用最优线性插值阈值算法（OLI-shrink）进行处理，其中该算法中的阈值函数是通过对子带系数应用NIG模型建模并结合MAP估计准则推导给出的，由于NIG模型中参数的自适应性，使得OLI-shrink算法中阈值也具有相应的自适应性，提高了去噪后图像质量。实验结果显示，本文算法在主观视觉效果和客观评价（PSNR，MSSIM）上都取得明显的改进，达到一定的去噪效果并充分保留原始图像的细节信息，说明本文算法的有效性。

图5 不同图像去噪性能对比

表1 不同方法对不同图像去噪后的平均结构相似性（MSSIM）

[1] Donoho D L, Johnstone I M. Ideal spatial adaptation via wavelet shrinkage[J]., 1994, 81(3): 425-455.

[2] Donoho D L, Johnstone I M. Adapting to unknown smoothness via wavelet shrinkage[J]., 1995, 90(432): 1200-1224.

[3] Selesnick I W, Baraniuk R G, Kingsbury N G. The dual-tree complex wavelet transform[J]., 2005, 22(6): 123-151.

[4] CHANG S G,YU B, Vetterli M. Adaptive wavelet thresholding for image denoising and compression[J]., 2000, 9(9): 1532-1546.

[5] Portilla J, Strela V, Wainwright M, et al. Image denoising using scale mixtures of gaussians in the wavelet domain[J]., 2003, 12(11): 1338-1351.

[6] Sendur L, Selesnick I W. Bivariate shrinkage functions for wavelet-based denoising exploiting interscale dependency[J]., 2002, 50(11): 2744-2756.

[7] Cho D, Bui T D. Multivariate statistical for image denoising using wavelet transforms[J]., 2005, 20(1): 77-89.

[8] Cunha A L, ZHOU J P, DO M N. The nonsubsampled contourlet transform: theory design and application[J]., 2006, 15(10): 3089-3101.

[9] GUO K, Labate D. Optimally sparse multidimensional representation using shearlets[J]., 2007, 39(1): 298-318.

[10] Chan W L, Choi H, Baraniuk R G. Coherent multiscale image processing using dual-tree quaternion wavelets[J]., 2008, 17(7): 1069-1082.

[11] Moulin P, LIU J. Analysis of multiresolution image denoising schemes using generalized Gaussian and complexity priors[J]., 1999, 45(3): 909-919.

[12] 贾建, 陈莉. 基于正态逆高斯模型的非下采样Contourlet变换图像去噪[J]. 电子学报, 2011, 39(7): 1563-1568.

JIA Jian, CHEN Li. Using normal inverse Gaussian model for image denoising in NSCT domain[J]., 2011, 39(7): 1563-1568.

[13] Barndorff-nielesn O E. Normal inverse Gaussian distribution and stochastic volatility modeling[J]., 1997, 24(1): 1-13.

[14] Chaudhury K N. Acceleration of the shiftable O(1) algorithm for bilateral filtering and non-local means[J]., 2013, 22(4): 1291-1300.

[15] Corrochano E B. Multi-resolution image analysis using the quaternion wavelet transform[J]., 2005, 39(1): 35-55.

[16] Selesnick I W, Kingsbury N G, Baraniuk R B. The dual-tree complex wavelets transform[J]., 2005(11):123-151.

[17] Bulow T. Hypercomplex spectral signal representations for the processing and analysis of images[D]. Kiel, Germany: Christian Albrechts Univ., 1999.

[18] Bhuiyan M I H, Ahmad M O, Swamy M N S. Wavelet based despeckling of medical ultrasound images with the symmetric normal inverse Gaussian prior[C]//,, 2007: I-721-I-724.

[19] Fathi A, Naghsh-Nilchi A R. Efficient image denoising method based on a new adaptive wavelet packet thresholding function[J]., 2012, 21(9): 3981-3990.

[20] Mihcak M K, Kozintsev M, Ramchandran I K. Low-complexity image denoising based on statistical modeling of wavelet coefficients[J]., 1999, 6(12): 300-303.

[21] 刘炳良. 一种小波域改进双边滤波的水果图像去噪算法[J]. 红外技术, 2014, 36(3): 196-199.

LIU Bingliang. A new improved bilateral filtering algorithm for the fruit image based on wavelet transform domain[J]., 2014, 36(3): 196-199.

[22] Crouse M S, Nowak R D, Baraniuk R G. Wavelet-based statistical signal processing using hidden Markov models[J]., 1988, 46(4): 886-902.

[23] YIN Ming, LIU Wei, SHUI Jun, et al. Quaternion wavelet analysis and application in image denoising[J]., 2012, 18(6): 587-612.

[24] Luisier F, Blu T, Unser M.A new SURE approach to image denoisng: interscale orthonormal wavelet thresholding[J]., 2007, 16(3): 593-606.

[25] ZHOU W, Bovik A C, Sheikh H R, et al. Image quality assessment: from error visibility to structural similarity[J]., 2004, 13(4): 600-612.

New Image Denoising Algorithm Based on NormalInverse Gaussian Model in Quaternion Wavelet Domain

ZHU Fang1，LIU Wei2

(1.,,230088,; 2.,,230031,)

This paper proposes a novel image denoising algorithm based on the normal inverse Gaussian model and fast bilateral filtering in the quaternion wavelet transform domain. The quaternion wavelet transform is utilized to decompose a noised image. The fast bilateral filtering algorithm is used to deal with the low frequency sub-band coefficients. The normal inverse Gaussian model for the prior model is used to describe the distributions of the image’s high frequency coefficients. Its corresponding threshold function is then derived using Bayesian maximum a-posteriori probability estimation theory. Finally, an optimal linear interpolation thresholding algorithm is employed to guarantee a gentler thresholding effect. Experimental results show that the proposed method outperforms other existing state-of-the-art denoising methods in terms of peak signal-to-noise and mean structural similarity.

quaternion wavelet transform，bilateral filtering，normal inverse Gaussian model，Bayesian maximum posterior estimation，optimal linear interpolation，image denoising

TP391

A

1001-8891(2017)10-0928-08

2017-03-08；

2017-04-07.

朱芳（1987-），女，硕士，讲师，研究方向为数字图像处理。

安徽省高等学校自然科学研究重点项目（KJ2016A310）；安徽新华学院《概率论与数理统计A》教改课程项目（2015jgkcx11）；安徽科技学院校级项目（ZRC2016499）。