基于圖結(jié)構(gòu)濾波的圖像去噪①

張哲浩,葛華勇,孫家慧

(東華大學(xué) 信息科學(xué)與技術(shù)學(xué)院,上海 201620)

1 引言

受各種記錄設(shè)備所固有的物理限制,大多數(shù)真實(shí)世界中獲取到的圖像都包含了一定程度的噪聲和失真[1],下面的模型可以描述期望圖像退化至觀測圖像的過程:

其中,y∈RN是由觀測圖像的像素組成的N維列向量,z∈RN是對應(yīng)的期望圖像,e∈RN假定是均值為0,方差為 σ2的加性高斯白噪聲,噪聲主要發(fā)生在圖像采集和傳輸過程中的模擬電路中,是實(shí)際圖像中最常見的噪聲. 另外,A∈RN×N是由模糊核產(chǎn)生的降級矩陣,在我們所關(guān)注的去噪問題中,A=I.

一類常見的去噪算法基于某些特定變換域,例如傅里葉變換、快速傅里葉變換、離散余弦變換、小波變換,這類算法主要利用稀疏性,使用較少的非零系數(shù)來表示信號. 確定變換域后還需要確定閾值,閾值函數(shù)一般有硬閾值和軟閾值兩種選擇[2]. 硬閾值方法將小波系數(shù)全部切除,其余的全部保留,軟閾值方法為了避免頻域中的階躍,將保留的系數(shù)都減去了閾值.

另一類方法基于空域,直接在原圖像上進(jìn)行處理,這類方法還可分為局部和非局部兩種. 常見的局部方法包括高斯濾波器、最小均方濾波器、雙邊濾波器、維納濾波器,局部算法的基本思想是利用像素之間的相關(guān)性,而像素的相關(guān)性會由于噪聲水平變高而嚴(yán)重受損,因此局部方法在高噪聲條件下表現(xiàn)不佳. 非局部方法中最常用的是非局部均值濾波(NLM),它利用了圖像的非局部相似性,是一種典型的逐點(diǎn)去噪算法[3].

圖信號處理的方法為圖像去噪提供了一種新思路.非局部圖變換算法通過圖信號處理在圖像塊中構(gòu)造出數(shù)據(jù)自適應(yīng)的變換域用于去噪,并在深度圖像去噪問題中表現(xiàn)優(yōu)秀[4]. Meyer 等人利用圖拉普拉斯矩陣的低維特征向量進(jìn)行去噪[5]. Yan 等人提出了GNNLG,在非局部方法的基礎(chǔ)上,結(jié)合圖拉普拉斯先驗(yàn)項(xiàng)進(jìn)行去噪[6].

目前基于圖信號的去噪的研究中,主要方法都是通過圖像塊匹配后,使用圖信號的方法對每一個(gè)圖像塊分別進(jìn)行去噪,關(guān)于圖像塊的去噪比較類似于傳統(tǒng)方法中基于變換域的方法,本質(zhì)上是利用圖傅里葉變換得到一個(gè)新的變換域. 而本文所提出的算法雖然同樣基于圖信號處理,但更類似于基于空域方法,直接針對圖像本身設(shè)計(jì)濾波器得到去噪圖像.

在本文的方法中,首先選擇適當(dāng)?shù)南嗨菩远攘?計(jì)算出圖像對應(yīng)的權(quán)重矩陣以及拉普拉斯矩陣,根據(jù)所得矩陣的濾波效果,提出了一個(gè)基于圖信號的目標(biāo)函數(shù),這個(gè)目標(biāo)函數(shù)使用權(quán)重矩陣來擬合數(shù)據(jù)保真度項(xiàng),并同時(shí)加入拉普拉斯先驗(yàn)項(xiàng)以確?；謴?fù)信號的平滑性.我們隨后說明了所提出的目標(biāo)函數(shù)的合理性,并最終使用迭代的方式對優(yōu)化問題進(jìn)行快速求解.

本文的其余安排如下: 首先在第2 節(jié)中介紹準(zhǔn)備工作,然后在第3 節(jié)中提出并詳細(xì)描述了我們的算法,并分別在第4 節(jié)和第5 節(jié)中介紹實(shí)驗(yàn)和結(jié)論.

2 準(zhǔn)備工作

2.1 圖信號處理的相關(guān)概念

根據(jù)圖信號處理的相關(guān)理論,可以將圖像建模為無向加權(quán)圖G(V,E,W)上的信號[7,8]. 通過將圖像中的每個(gè)像素都看作一個(gè)頂點(diǎn),可以得到頂點(diǎn)集合V; 邊緣集合E中的每一條邊連接了一對頂點(diǎn); 每條邊(i,j)∈E對應(yīng)一個(gè)非負(fù)邊緣權(quán)重wi,j,權(quán)重越大代表邊緣兩端的頂點(diǎn)i,j∈V的相關(guān)性越強(qiáng). 由于G為無向加權(quán)圖,因此有wi,j=wj,i,顯然W是一個(gè)對稱矩陣. 一個(gè)3×3圖像對應(yīng)的圖結(jié)構(gòu)與相應(yīng)的權(quán)重矩陣如圖1所示. 圖像的像素值被認(rèn)為是圖G頂點(diǎn)上的圖信號,y=[y1,y2,···,yN]T.于是,對于每一個(gè)頂點(diǎn)i∈V,本文定義了一個(gè)向量來描述位于該頂點(diǎn)的位置和像素強(qiáng)度信息,即:

圖1 圖像的圖結(jié)構(gòu)與其對應(yīng)的權(quán)重矩陣

其中,xi,1,xi,2分別代表頂點(diǎn)i在二維平面的空間坐標(biāo),yi代表頂點(diǎn)i的像素強(qiáng)度.

確定圖權(quán)重的一種常用方法是使用距離的指數(shù)衰減函數(shù),即:

其中,dist(ri,rj)代表頂點(diǎn)ri和rj之間的某種距離度量,h2是給定的參數(shù),這樣的選擇在直覺上也是合理的,因?yàn)閞i與rj距離越接近,權(quán)重越接近1; 反之,距離越疏遠(yuǎn),權(quán)重也就越小. 本文認(rèn)為較大權(quán)重wi,j所連接的頂點(diǎn)之間應(yīng)該具有相似的圖信號值.

關(guān)于常用的距離度量在表1 中給出,所列舉距離度量的直觀圖形解釋如圖2 所示. 由于經(jīng)典回歸濾波器對感興趣信號的底層結(jié)構(gòu)缺乏更強(qiáng)的適應(yīng)性,雙邊濾波器在低信噪比情況下無法提供有效性能,因此在后文中均選擇非局部均值濾波器作為權(quán)重矩陣中的距離度量dist(ri,rj). 在非局部均值濾波器中,和分別指 以頂點(diǎn)ri和rj為中心的圖像塊.

表1 常用的距離度量

圖2 常用距離度量的圖形解釋

3 基于圖濾波的圖像去噪算法

3.1 圖信號的譜分析

圖信號處理中的一個(gè)重要算子是圖拉普拉斯算子L=D?W,其中,D是一個(gè)對角矩陣,它的第i個(gè)對角元素等于W第i行的所有元素之和.

對拉普拉斯矩陣進(jìn)行特征分解,有Lu=λu,λ 為特征值,u為相應(yīng)的特征向量. 由于L的每一行與常向量的內(nèi)積為0,假設(shè)u為一個(gè)常向量,則Lu=0u=→0,另外對于無向圖結(jié)構(gòu),L是半正定的實(shí)對稱矩陣,于是可以寫作L=UΛUT,其中Λ=diag(λ1,λ2,···,λN)是一個(gè)由L的特征值所組成的對角矩陣,U由特征向量構(gòu)成,每一列代表一個(gè)特征向量,并且滿足U?1=UT. 根據(jù)特征向量的定義,拉普拉斯矩陣L的每一個(gè)特征向量uk都滿足Luk=λkuk.

接下來說明,對于圖信號z,zTLz表示信號z在L所延拓的空間內(nèi)的平滑度.

由Luk=λkuk,左乘ukT,得ukTLuk=ukTλkuk,特征向量可以進(jìn)行常數(shù)倍的放縮,因此經(jīng)過適當(dāng)?shù)募s定,令ukTuk=1,可以得到:

由于拉普拉斯矩陣L的特征向量u1,u2,···,uN可以構(gòu)成一個(gè)完備正交集,因此可以將圖信號z分解為不同特征向量的系數(shù)和. 分解后z對應(yīng)的低頻分量越多,zTLz的值越小,對應(yīng)的高頻分量越多,zTLz的值也就越大. 在極端情況下,如果z=u1=[1,1,···,1],那么此時(shí)zTLz=0,因此如果希望優(yōu)化后的信號更加平滑,可以在目標(biāo)函數(shù)中加入zTLz作為先驗(yàn)項(xiàng),這樣傾向于得到平滑度更高的z作為最優(yōu)解.

換一個(gè)角度,如果將L看作濾波矩陣,那么它作用在圖信號z時(shí)得到Lz,L的效果相當(dāng)于一個(gè)自適應(yīng)的高通濾波器,假如z為直流向量,那么濾波后的結(jié)果為.基于這樣的性質(zhì),通過計(jì)算圖像所對應(yīng)的拉普拉斯矩陣,就可以過濾圖像所對應(yīng)的圖信號,形成自適應(yīng)濾波.

同理,由于L=D?W,W和L的特征向量是相同的,對應(yīng)的特征值逆序,因此W也可以認(rèn)為是一個(gè)自適應(yīng)低通濾波器.

3.2 問題建模與優(yōu)化

本文提出下面的優(yōu)化問題:

其中,第一項(xiàng)為數(shù)據(jù)保真度項(xiàng),用來衡量觀測輸入y和期望輸出z之間的差異; 第二項(xiàng)是數(shù)據(jù)自適應(yīng)差分項(xiàng),有助于獲得更平滑的輸出z.

對于目標(biāo)函數(shù)E(z)=(y?z)TW(y?z)+ηzTLz,關(guān)于E(z)的優(yōu)化可以通過計(jì)算函數(shù)梯度直接得到其最優(yōu)解的數(shù)學(xué)表達(dá)式.

首先計(jì)算E(z)的梯度:

考慮到權(quán)重矩陣W∈RN×N,矩陣求逆的時(shí)間復(fù)雜度過高,不便于計(jì)算,因此采用梯度下降法對問題進(jìn)行求解.

由式(6),迭代更新表達(dá)式為:

其中,μk代表梯度下降法每一輪更新的步長.

如果每一輪都要求出最優(yōu)的 μk,就會轉(zhuǎn)換為對μk的優(yōu)化問題,顯然會使計(jì)算量會變大,于是可以采用A rmijo 線搜索確定每一輪的步長.

3.3 算法流程

總結(jié)歸納后得到算法1.

算法1. 基于圖結(jié)構(gòu)濾波的圖像去噪算法1)非局部均值核計(jì)算輸入圖像的權(quán)重矩陣 ;W L y W 2)根據(jù)權(quán)重矩陣 計(jì)算得到拉普拉斯矩陣 ;3)循環(huán)執(zhí)行下列迭代直到收斂:zk=zk?1?μk?E(zk?1)︿z=(W+ηL)?1Wy 4)得到輸出圖像近似等于.

算法1 中,步長 μk通過Armijo 線搜索確定.

4 實(shí)驗(yàn)分析

為了驗(yàn)證本文所提出的方法的性能,使用圖3 中的6 張灰度圖像作為測試圖像,如果要推廣至彩色圖像去噪,將算法分別應(yīng)用于各個(gè)通道即可. 在這些圖像上添加加性高斯白噪聲,均值為0,噪聲方差從0.2 到400,從而得到噪聲圖像. 對照算法分別選用NLM[3],BM3D[9],WNNM[10],這3 種算法都是去噪效果較好且比較常用的算法. NLM 的基本思想是基于圖像塊構(gòu)建圖像的逐點(diǎn)估計(jì). WNNM 基于加權(quán)核范數(shù)進(jìn)行去噪.BM3D 基于圖像塊進(jìn)行三維濾波. 所有實(shí)驗(yàn)都是基于Matlab,運(yùn)行在LENOVO 82DN 環(huán)境下,所使用的處理器是Intel Core i5-10210U.

圖3 用于評估算法性能的一組圖像

其中,n為每像素的比特?cái)?shù),對于灰度圖像,n=8,MSE表示當(dāng)前圖像和參考圖像之間的均方誤差.PSNR的單位是 dB,該值越大,代表失真越小.

對于圖像恢復(fù)效果的衡量,使用峰值信噪比(PSNR)作為性能指標(biāo)[11]. 定義為:

表2 中給出了不同方法的PSNR指標(biāo). 噪聲方差處于0.2 和400 之間. 可以看到,與其他的去噪算法相比,本文的方法在所有圖像中的PSNR指標(biāo)都是最佳的.BM3D 和NLM 在各個(gè)條件下呈現(xiàn)出相似的去噪水平,而我們的方法在噪聲方差為0.2 時(shí)比BM3D 的表現(xiàn)平均提升了26.8 dB. 可以觀察到本文的方法在低噪聲水平下具有更明顯的優(yōu)勢,這是由于在本文的目標(biāo)函數(shù)中的數(shù)據(jù)保真度項(xiàng)通過權(quán)重矩陣的擬合,可以更加自適應(yīng)地還原出原圖像的細(xì)節(jié)信息.

表2 基于圖濾波去噪算法,NLM,BM3D,WNNM 的PSNR 性能 (dB)

視覺上的主觀比較如圖4,以Goldhill 疊加方差為400 的噪聲為例. 使用本文的方法所生成的圖像紋理更加清晰. 相比WNNM,本文的方法沒有造成過度模糊,相比NLM,本文的方法顯然能更好地去除噪聲,相比BM3D,本文的方法消除了其中的振鈴效應(yīng). 因此,使用本文的方法可以獲得更好的細(xì)節(jié),極大抑制了噪聲且避免了振鈴,同時(shí)紋理信息清晰可見.

圖4 Goldhill 圖像的去噪示例

5 結(jié)論與展望

本文提出了一種基于圖信號的圖像去噪通用方法.這個(gè)方法基于以下兩點(diǎn)分析: (1)將基于NLM 濾波計(jì)算得到的圖結(jié)構(gòu)用于描述圖像的內(nèi)在框架. (2)基于圖像計(jì)算的權(quán)重矩陣和拉普拉斯矩陣分別可以看作自適應(yīng)低通和高通濾波器,可用于表征圖信號的特征. 我們創(chuàng)新性地提出了一種基于圖信號的目標(biāo)函數(shù)用于去噪,對于所提出的目標(biāo)函數(shù),給出了理論的最優(yōu)解和迭代逼近算法. 實(shí)驗(yàn)表明,本文所提出的方法在數(shù)值指標(biāo)和視覺效果上都優(yōu)于BM3D 和WNNM 等前沿去噪算法.