一種基于二維GARCH模型的圖像去噪方法

李俊澤，袁小芳，張振軍，王耀南，王國(guó)鋒

（1.湖南大學(xué) 電氣與信息工程學(xué)院，湖南 長(zhǎng)沙 410082； 2.中國(guó)公路工程咨詢集團(tuán)有限公司，北京 100097）

?

一種基于二維GARCH模型的圖像去噪方法

李俊澤1，袁小芳1，張振軍1，王耀南1，王國(guó)鋒2

（1.湖南大學(xué) 電氣與信息工程學(xué)院，湖南 長(zhǎng)沙 410082； 2.中國(guó)公路工程咨詢集團(tuán)有限公司，北京 100097）

提出了一種基于小波系數(shù)統(tǒng)計(jì)模型的圖像去噪方法。該方法利用二維廣義自回歸條件異方差（2D-GARCH）模型對(duì)小波系數(shù)進(jìn)行建模，這種小波系數(shù)模型能夠更好地利用小波系數(shù)“尖峰厚尾”的分布特性和系數(shù)間的相關(guān)性等重要特性。利用基于果蠅優(yōu)化算法的極大似然估計(jì)（ML Estimation based on FOA）代替?zhèn)鹘y(tǒng)的線性規(guī)劃方法求解模型參數(shù)，提高了建模的準(zhǔn)確性。在此基礎(chǔ)上再采用最小均方誤差估計(jì)（MMSE Estimation）對(duì)未受噪聲污染的原始圖像的小波系數(shù)進(jìn)行估計(jì)。實(shí)驗(yàn)結(jié)果表明，與當(dāng)前主流的去噪方法相比，該算法能更有效地去除圖像中的噪聲，獲得更高的峰值信噪比（PSNR）和較好的視覺效果。

小波變換；統(tǒng)計(jì)建模；二維GARCH模型；果蠅優(yōu)化算法；圖像去噪

在智能交通領(lǐng)域，交通場(chǎng)景圖像質(zhì)量的高低直接影響對(duì)場(chǎng)景中車輛、行人、障礙物等目標(biāo)的檢測(cè)和識(shí)別。然而，圖像在獲取和傳輸?shù)倪^(guò)程中往往容易被噪聲污染，從而影響圖像的質(zhì)量。有效的去噪方法不僅能夠抑制圖像中的隨機(jī)噪聲，同時(shí)還能較好地保持重要的圖像特征和細(xì)節(jié)。如果能夠了解圖像的某個(gè)先驗(yàn)統(tǒng)計(jì)模型，便能很大程度上地改善圖像處理的效果[1]。在圖像去噪領(lǐng)域，基于小波變換的方法已經(jīng)成為了主流[2]。例如Donoho等提出閾值去噪算法[3]以及Chang等在此基礎(chǔ)上改進(jìn)得到的閾值萎縮算法[4-6]。這些去噪算法主要包括3個(gè)步驟：首先，對(duì)含噪圖像進(jìn)行小波變換；然后，對(duì)各細(xì)節(jié)子帶上小波系數(shù)進(jìn)行閾值萎縮處理；最后，通過(guò)小波逆變換重構(gòu)去噪后圖像。然而閾值萎縮處理主要存在2個(gè)不足：1） 如何選擇一個(gè)最合適的閾值一直是個(gè)難以解決的問題；2） 在不同的尺度上，小波系數(shù)的先驗(yàn)統(tǒng)計(jì)模型并不能和信號(hào)很好地匹配[7]。研究發(fā)現(xiàn)，相較于目前廣泛用于小波系數(shù)建模的廣義高斯模型（GGM），二維廣義自回歸條件異方差（2D-GARCH）模型能夠更準(zhǔn)確地描述圖像小波系數(shù)的統(tǒng)計(jì)特性，是一種更加優(yōu)秀的小波系數(shù)先驗(yàn)?zāi)Ｐ停时疚氖褂?D-GARCH模型對(duì)小波系數(shù)進(jìn)行建模。1D-GARCH模型[8]在金融時(shí)間序列建模中已經(jīng)得到了廣泛的應(yīng)用，Noiboar等最早將1D-GARCH模型擴(kuò)展為二維模型并用于圖像的畸點(diǎn)探測(cè)[9]，Amirmazlaghani等利用2D-GARCH模型抑制合成孔徑雷達(dá)圖像上的斑點(diǎn)噪聲[10-11]。本文利用基于果蠅優(yōu)化算法的極大似然估計(jì)準(zhǔn)確求解2D-GARCH模型的參數(shù)，獲得了各細(xì)節(jié)子帶的先驗(yàn)統(tǒng)計(jì)模型，在此基礎(chǔ)上再通過(guò)最小均方誤差估計(jì)對(duì)原始圖像的小波系數(shù)進(jìn)行估計(jì)，從而達(dá)到去噪目的。

1 2D-GARCH模型的小波系數(shù)建模

1.1 2D-GARCH模型

假設(shè)yij是二維的隨機(jī)變量且均值為零即E（yij）=0，若yij有如下形式則稱其服從純2D-GARCH（p1,p2,q1,q2）模型[9]：

（1）

（2）

式中：

,（kl）≠（0,0）}

Λ2={kl|0≤k≤p1,0≤l≤p2,（kl）≠（0,0）}

hij是yij的條件方差，εij是獨(dú)立同分布的二維標(biāo)準(zhǔn)正態(tài)分布，即εij～Ν（0,1），（p1,p2,q1,q2）為模型的階數(shù)。從式（2）可以看出，每一個(gè)確定的空間位置上的變量值yij的條件方差hij是由與之鄰近的變量值yi-k,j-l和其條件方差hi-k,j-l決定的。現(xiàn)定義隨機(jī)變量yij鄰近位置上的變量值和條件方差的集合為

根據(jù)式（1）、（2）可知yij的條件分布為

yij|ψij～Ν（0,hij）

（3）

由于條件方差hij非負(fù)，那么模型參數(shù)需滿足式（4）：

α0>0

αkl≥0,kl∈Λ1

βkl≥0,kl∈Λ2

（4）

同時(shí)為了確保2D-GARCH（p1,p2,q1,q2）過(guò)程具有廣義平穩(wěn)性[9]，需保證

（5）

式（4）、（5）是2D-GARCH（p1,p2,q1,q2）模型有意義且具有廣義平穩(wěn)性的充要條件，下文的所有推導(dǎo)都在此約束條件下進(jìn)行。

然而統(tǒng)計(jì)發(fā)現(xiàn)各子帶上的小波系數(shù)的均值E（yij）≠0，那么首先對(duì)yij進(jìn)行相應(yīng)的去均值處理。令zij表示二維線性回歸的新息：

zij=yij-rijTb

（6）

式中：rij=[的均值服從式（6）所描述的模型，方差服從式（1）、（2）描述的模型，其條件分布為

yij|ψij～Ν（rijTb,hij）

（7）

根據(jù)式（7）可知yij的條件概率密度函數(shù)的表達(dá)式為

（8）

式中：

（9）

求解上述模型，即需要估計(jì)出模型中的未知參數(shù)Γ={{α0,α01,…,αq1q2,β01,…,βp1p2},b}。本文采用極大似然估計(jì)求解模型參數(shù)，現(xiàn)定義樣本空間：

M×N:Φ={ij|1≤i≤M,1≤j≤N}

根據(jù)式（8）、（9）可知關(guān)于參數(shù)Γ的對(duì)數(shù)極大似然函數(shù)為

（10）

為了完成對(duì)各細(xì)節(jié)子帶上小波系數(shù)的建模，必須求出使得上述對(duì)數(shù)極大似然函數(shù)取得最大值時(shí)的參數(shù)Γ。

1.2 2D-GARCH模型參數(shù)估計(jì)

由于較高的模型階數(shù)會(huì)導(dǎo)致計(jì)算量過(guò)大且更容易受到舍入誤差的不利影響，因此本文采用2D-GARCH（1,1,1,1）模型對(duì)各細(xì)節(jié)子帶建模，即未知參數(shù)為：Γ={{α0,α01,α10,α11,β01,β10,β11},b}。利用傳統(tǒng)的線性規(guī)劃方法求解式（10）的極大值問題時(shí)往往容易陷入局部最優(yōu)解，對(duì)子帶的準(zhǔn)確建模產(chǎn)生較為不利的影響。為了盡可能地消除這種影響，本文采用果蠅優(yōu)化算法[12]求解式（10）。果蠅優(yōu)化算法屬于進(jìn)化算法的一種，對(duì)高維目標(biāo)函數(shù)具有較好的尋優(yōu)能力和較快收斂速度。

基于果蠅優(yōu)化算法的2D-GARCH（1,1,1,1）模型參數(shù)估計(jì)的原理如下：

1） 在滿足約束條件的前提下初始果蠅群體位置InitΓ_axis.其中Γ={{α0,α01,α10,α11,β01,β10,β11},b}。

2） 賦予果蠅利用嗅覺搜索食物的隨機(jī)方向與距離Γi=Γ_axis+RandomValue。

3） 將Γi代入式（10）求出果蠅個(gè)體的味道濃度

Smelli=LF（Γi）

4） 找出果蠅群體中味道濃度最高的果蠅

5） 保留最佳味道濃度與各參數(shù)的坐標(biāo)，此時(shí)果蠅群體利用視覺往該方向飛去：

Smellbest=bestSmell

6） 進(jìn)入迭代尋優(yōu)，重復(fù)執(zhí)行步驟2）～4），并判斷味道濃度是否優(yōu)于前一迭代味道濃度，若是則執(zhí)行步驟5）。

為了檢驗(yàn)?zāi)Ｐ蛥?shù)估計(jì)的準(zhǔn)確性，采用圖像小波分解后第2層的水平子帶上的小波系數(shù)進(jìn)行參數(shù)估計(jì)實(shí)驗(yàn)。表1給出了該子帶模型各參數(shù)的估計(jì)結(jié)果，圖1是歸一化的小波系數(shù)直方圖和參數(shù)值為表1中各估計(jì)值的2D-GARCH（1,1,1,1）模型直方圖。

表1 模型參數(shù)估計(jì)值Table 1 Estimated value of the model’s parameters

圖1 第2層水平子帶的小波系數(shù)直方圖

2 圖像去噪算法

本文去噪算法的具體步驟如下，流程圖如圖2。

圖2 本文去噪算法流程圖

1） 含噪圖像小波分解：對(duì)被噪聲污染的圖像進(jìn)行N層離散小波分解，得到不同尺度和方向上的子帶。

2） 2D-GARCH（1,1,1,1）模型參數(shù)求解：利用果蠅優(yōu)化算法求解式（10），得到各個(gè)尺度上的細(xì)節(jié)子帶的模型參數(shù)。

3） 原始圖像小波系數(shù)估計(jì)：在上一步估計(jì)出各尺度上的細(xì)節(jié)子帶的模型參數(shù)的基礎(chǔ)上，利用最小均方誤差估計(jì)原始圖像的小波系數(shù)。

4） 圖像重構(gòu)：通過(guò)離散小波逆變換重構(gòu)去噪后的圖像。令x和n分別表示原始圖像和加性噪聲，y=x+n表示觀測(cè)到的含噪圖像。首先，對(duì)含噪圖像進(jìn)行任意尺度的離散小波變換，從而得到了不同尺度和方向上的子帶。令Xij，Nij和Yij分別代表x、n和y的任意子帶上的小波系數(shù)，那么有

Yij=Xij+Nij

（11）

zij=yij-rijTb

（12）

（13）

（14）

（15）

（16）

本文采用最小均方誤差估計(jì)，則式（16）可轉(zhuǎn)化為

（17）

即

（18）

（19）

根據(jù)上文推導(dǎo)可知：

（20）

（21）

（22）

將式（20）～（22）代入式（19）中可得

（23）

最后，通過(guò)離散小波逆變換重建去噪后的圖像。

3 實(shí)驗(yàn)結(jié)果分析

實(shí)驗(yàn)采用圖3（a）所示含有不同強(qiáng)度高斯白噪聲的交通場(chǎng)景灰度圖像（Cars 256×256），并使用Daubechies（Db4）小波基[10-11]對(duì)圖像進(jìn)行二層小波分解。為了檢驗(yàn)本文算法的去噪效果，文中與維納濾波法[14]、BayesShrink方法[15]、改進(jìn)的局部窗口自適應(yīng)估計(jì)（LAWML）方法[16]進(jìn)行對(duì)比，并根據(jù)峰值信噪比（PSNR）來(lái)評(píng)價(jià)各算法的去噪效果。其中峰值信噪比定義為

表2 4種去噪算法的去噪性能比較

Table 2 The comparison of four different denoising methods PSNR

dB

圖3 σ=20時(shí)Cars圖像4種算法的去噪比較

從表2中可以看出，本文方法的峰值信噪比明顯高于其他方法。當(dāng)噪聲強(qiáng)度增大時(shí)，去噪效果下降的比較明顯，這是因?yàn)楫?dāng)噪聲強(qiáng)度過(guò)大時(shí)可能會(huì)對(duì)模型的適用范圍產(chǎn)生一定的影響。

圖3是當(dāng)噪聲標(biāo)準(zhǔn)差σ=20的時(shí)候各種算法的去噪效果對(duì)比圖。圖3（c）、（d）、（e）分別是利用Wiener濾波法、BayesShrink方法和改進(jìn)的LAWML方法對(duì)圖像去噪后的效果，出現(xiàn)了一定程度圖像模糊、部分邊緣輪廓辨認(rèn)困難的現(xiàn)象。圖3（f）是采用本文方法的去噪效果，可以看出用本文方法對(duì)圖像去噪后較好地保持了圖像的邊緣輪廓特征，且具有較好的視覺效果。

在算法復(fù)雜度方面，若圖像的大小為N×N，那么Wiener濾波法、BayesShrink方法和改進(jìn)的LAWML方法的時(shí)間復(fù)雜度均為O（N2）。本文算法的步驟1）、3）和4）的時(shí)間復(fù)雜度均為O（N2），步驟2）用到果蠅優(yōu)化算法求解各細(xì)節(jié)子帶的模型參數(shù)，令T0代表果蠅個(gè)體一次迭代需要的時(shí)間，M為種群數(shù)量，K為迭代次數(shù)，那么步驟2）的時(shí)間復(fù)雜度為T=K·M·T0。

在算法耗時(shí)方面，由于需要估計(jì)出每個(gè)細(xì)節(jié)子帶的2D-GARCH（1,1,1,1）模型參數(shù)，所以本文算法在計(jì)算量上相較于幾個(gè)對(duì)比算法要大一些，但是本系統(tǒng)面向交通的圖像處理硬件平臺(tái)采用的是多核分布式處理器網(wǎng)絡(luò)以及流水線式的任務(wù)分配，所以本文算法并不影響整個(gè)系統(tǒng)的實(shí)時(shí)性。

4 結(jié)束語(yǔ)

二維GARCH模型是一種新的小波系數(shù)統(tǒng)計(jì)模型，能夠更好地利用小波系數(shù)“尖峰厚尾”的分布特性和層內(nèi)系數(shù)間的相關(guān)性等重要特性。本文采用果蠅優(yōu)化算法的極大似然估計(jì)來(lái)求解2D-GARCH模型參數(shù)，完成了對(duì)小波系數(shù)的準(zhǔn)確建模，在此基礎(chǔ)上再采用最小均方誤差估計(jì)原始圖像的小波系數(shù)，從而得到了一種新的2D-GARCH模型的圖像去噪方法。該方法在峰值信噪比和視覺效果上要好于目前一些主流的去噪方法。然而該方法也存在一些需要改進(jìn)的地方，比如在求解各細(xì)節(jié)子帶的模型參數(shù)時(shí)存在計(jì)算量較大的缺點(diǎn)，今后將在此基礎(chǔ)上進(jìn)一步完善和改進(jìn)。

[1]SIMONCELLI E P. Modeling the joint statistics of image in the wavelet domain[C] //SPIE’s International Symposium on Optical Science, Engineering, and Instrumentation. International Society for Optics and Photonics.Denver， USA, 1999: 188-195.

[2]ANTONIADIS A, BIGOT J, SAPATINAS T. Wavelet estimators in nonparametric regretssion: A comparative simulation study[J]. Journal of Statistical Software, 2001, 6（6）: 1-83.

[3]DOHONO D L. Denoisingby soft-thresholding[J]. IEEE Transactions on Information Theory, 1995, 41（3）: 613-627.

[4]CHANG S, YU B, VATTERELI M. Wavelet thresholding for multiple noisy image[J]. IEEE Transactions on Image Processing, 2000, 9（9）: 1631-1635.

[5]CHANG S, YU B, VATTERELI M. Spatially adaptive wavelet thresholding with context modeling for imaged noising[J]. IEEE Transactions on Image Processing, 2000, 9（9）: 1522-1531.

[6]CHANG S, YU B, VATTERELI M. Adaptive wavelet thresholding for image denoising and compression[J]. IEEE Transactions on Image Processing, 2000, 9（9）: 1532-1546.

[7]ACHIM A, BEZERIANOS A, TSAKALIES P. SAR image denoising via Bayesian wavelet shrinkage based on heavy tailed modeling[J]. IEEE Transactions on Geoscience and Remote Sensing, 2003, 41（8）: 1773-1784.

[8]BOLLERSLEV T. Generalized autoregressive conditional heteroscedasticity[J]. Journal of Econometrics, 1986, 31（3）: 307-327.

[9]NOIBOAR A, COHEN I. Two-dimensional GARCH model with application to anomaly detection[C] //13th European Signal Processing Conference. Istanbul, Turkey， 2005: 1594-1597.

[10]AMIRMAZLAGHANI M, AMIRNDAVAR H. Speckle suppression in SAR image using the 2D GARCH model[J]. IEEE Transactions on Image Processing, 2009, 18（2）: 250-259.

[11]AMIRMAZLAGHANI M, AMIRNDAVAR H. Two novel Bayesian multiscale approaches for speckle suppression in SAR images[J]. IEEE Transactions on Geoscience and Remote Sensing, 2010, 40 （7）: 2980-2993.

[12]PAN W T. A new fruit fly optimization algorithm: taking the financial distress model as an example[J]. Knowledge Based Systems, 2012, 26（2）: 69-74.

[13]DOHONO D L, JOHNSTONE I M. Ideal spatial adaptati on via wavelet shrinkage[J]. Biometrika, 1994, 81（3）:425-455.

[14]LEE J S. Digital image enhancement and noise filtering by use of local statistics[J]. IEEE Transactions on Pattern Analysis and Machine Intelligence, 1980, 2（2）: 165-168.

[15]胡海平，莫玉龍.基于貝葉斯估計(jì)的小波閾值圖像降噪方法[J].紅外與毫米波學(xué)報(bào)， 2002， 21（1）: 74-76. HU Haiping, MO Yulong. Method of wavelet threshold denoising based on bayesian esitimation[J].Journal of Infrared and Millimeter Waves, 2002, 42（9）: 74-76.

[16]謝杰成，張大力，徐文立.一種小波去噪方法的幾點(diǎn)改進(jìn)[J].清華大學(xué)學(xué)報(bào)：自然科學(xué)版，2002，42（9）：1269-1272. XIE Jiecheng, ZHANG Dali, XU Wenli. Several improvements for a wavelet denoising method[J]. Journal of Tsinghua University:Science and Technology, 2002, 42（9）: 1269-1272.

李俊澤，男，1988年生，碩士研究生，主要研究方向?yàn)閳D像處理。

袁小芳，男，1979年生，副教授，主要研究方向?yàn)橹悄芸刂评碚撆c應(yīng)用、電動(dòng)汽車控制、新能源發(fā)電，發(fā)表學(xué)術(shù)論文30余篇。

張振軍，男，1981年生，講師，碩士生導(dǎo)師，主要研究方向?yàn)闄C(jī)器視覺與智能交通、大規(guī)模機(jī)器學(xué)習(xí)與海量數(shù)據(jù)分析。

A method of image denoising based on two-dimensional GARCH model

LI Junze1， YUAN Xiaofang1， ZHANG Zhenjun1， WANG Yaonan1， WANG Guofeng2

（1. School of Electrical and Information Engineering, Hunan University, Changsha 410082, China; 2. China Highway Engineering Consulting Corporation, Beijing 100097, China）

An image denoising method based on the statistical model for wavelet coefficients is proposed. It uses a two-dimensional Generalized Autoregressive Conditional Heteroscedasticity （2D-GARCH） model for modeling the wavelet coefficients. A novel wavelet coefficients model is also used to make better use of the important characteristics of wavelet coefficients such as "sharp peak and heavy tailed" marginal distribution and the dependencies between the coefficients. It utilizes maximum likelihood estimation based on fruit fly optimization algorithm （ML Estimation based on FOA） to estimate the model parameters instead of using traditional linear programming in order to improve the accuracy of the modeling. The minimum mean square error estimation （MMSE Estimation） is applied to estimating the parameters of the wavelet coefficients of the original image that is not affected by noise. Experimental results showed that compared to the present widely-used denoising methods the proposed method is more effective in image denoising, and it may achieve higher peak signal-to-noise ratio （PSNR） and good visuality.

wavelet transform; statistical modeling; two-dimensional GARCH model; FOA; image denoising

2014-03-25.

日期：2015-01-13.

國(guó)家“863”計(jì)劃資助項(xiàng)目（2012AA112312）.

張振軍. E-mail：zhenjun@hnu.edu.cn.

10.3969/j.issn.1673-4785.201403066

http://www.cnki.net/kcms/doi/10.3969/j.issn.1673-4785.201403066.html

TP751.1

A

1673-4785（2015）01-0062-06

李俊澤，袁小芳，張振軍，等. 一種基于二維GARCH模型的圖像去噪方法[J]. 智能系統(tǒng)學(xué)報(bào)， 2015, 10（1）: 62-67.

英文引用格式：LI Junze，YUAN Xiaofang，ZHANG Zhenjun，et al. A method of image denoising based on two-dimensional GARCH model[J]. CAAI Transactions on Intelligent Systems, 2015, 10（1）: 62-67.