基于塊空間中低維流形恢復的圖像去噪

張宗良，顏永豐

(西北農(nóng)林科技大學信息工程學院，陜西楊凌712100)

0 引言

圖像去噪具有極高的理論和實用價值，近年來數(shù)字圖像的急劇增多更加突顯了圖像去噪的重要意義。傳統(tǒng)的去噪方法通常先假設圖像具有某種良好的性質，如:具有較小的總變分［1-2］，或者在小波基中稀疏，然后通過最小化與該性質對應的目標函數(shù)來去噪。

最近幾年利用圖像塊之間的相似性來去噪取得了極大的成功。文獻［3］最早研究了基于塊處理的圖像去噪，并提出了非局部均值 (NLM)算法。該算法通過加權累加非局部塊來估計當前所求塊。權值由非局部塊與當前所求塊的相似度決定。在非局部塊方法的基礎上，文獻［4］提出了一種全新的算法BM3D。該算法先用塊匹配算法找到那些與當前塊相似的塊，再將這些相似塊組成一個柱體。由于該柱體具有極大的冗余性，對該柱體進行三維離散余弦變換、閾值截斷、余弦逆變換后便可去噪。文獻［5］的算法與BM3D算法的不同之處是它對柱體進行的是PCA變換。

上述基于塊相似的去噪簡單有效，但是只有少量的文章［6-9］給出其背后的數(shù)學解釋?？偨Y這些文章我們知道無噪圖像中的相似塊會在圖像塊空間中形成低維流形，而在帶噪圖像中低維流形受到噪聲污染而變形。復原這些低維流形便可去噪。本文據(jù)此提出了一個多步迭代的去噪算法。與大多數(shù)類BM3D算法的兩步迭代相比，本文算法的多步迭代更具有一般性。

1 圖像塊空間中的低維流形

本文只考慮灰度數(shù)字圖像。一幅圖像是從平面區(qū)域到像素域的映射u:D→［0，255］ R，其中D Z2是圖像所在的平面區(qū)域。ux，y∈［0，255］是位置 (x，y)∈D上的像素值，［0，255］是像素的取值范圍即“像素域”。記位置(x，y)的鄰域為Nx，y。本文涉及到的都是長方形或方形的鄰域。有兩種常見的對(x，y)取鄰域的方法，一種是(x，y)在鄰域的左上角，另一種是(x，y)在鄰域的中心。例如若采用第一種取法，則(x，y)的大小為I×J的鄰域Nx，y={(i，j)}，i∈［x，x+I－1］∩Z ，j∈［y，y+J－1］∩Z 。圖像塊是尺寸遠小于原圖像的子圖像，即原圖像上的一個小塊。位置(x，y)上的圖像塊是指原圖像在鄰域Nx，y上的一幅子圖像px，y=uNx，y。這里有個細節(jié)問題，當 (x，y)靠近圖像的邊界時，鄰域Nx，y有可能超出圖像所在的平面區(qū)域即Nx，yD。此時可對圖像進行擴展，本文采用的擴展方法是使擴展部分與原圖像沿邊界對稱。接下來我們進一步的記由圖像塊px，y的像素按行列順序排列而成的向量為塊向量px，y，稱塊向量所在的空間為塊空間。

下面我們考察圖像塊向量在塊空間中的分布是什么樣的。真實自然的圖像具有很多良好的性質。例如圖像中含有大量相同或相似的塊，這些相同或相似的塊向量將分布在塊空間中一個半徑很小的球中。如圖1中的標為1、2的塊分別與標為R1、R2的塊相似。還有些圖像塊由于光照強度的漸變等原因而有規(guī)律的變化，這些圖像塊向量也將按一定的規(guī)律分布在塊空間中。由此可見，特定區(qū)域內的所有圖像塊向量將在塊空間中形成若干個流形［6-8］。這里的“流形”可以理解為是由相似塊或有規(guī)律的塊組成的一個有規(guī)律可循的組織或結構。圖2是一個二維塊分布示意圖，在該圖中，塊空間是二維的即每個塊只含二個像素。從左圖可看到若干條直線 (近似)，每條直線可看成一個低維的流形。右圖是對應的受到噪聲污染后的流形，可以看到噪聲的影響直觀上是使直線變粗。可見，只要我們能從受噪聲污染的流形 (圖1b和圖2b)中恢復出原本低維的流形(圖1a和圖2a)，再將恢復后流形中的點即圖像塊放回到該圖像塊原來所在的位置便能達到去噪效果。

接下來看看真實圖像中的情況。同樣在塊大小等于2(即一個塊只含左右兩個像素)的情況下畫出塊分布圖。圖4(a)、(b)、(c)分別是圖3中實線條框出的甲、乙、丙三個區(qū)域及其中的圖像塊在塊空間中的分布圖。左邊的是無噪的情況，右邊是加入方差為5的高斯白噪聲的情況。從圖4可以看出塊向量在塊空間中的分布確實有一定的規(guī)律，或聚成小圓，或近似排成直線。從圖4還能明顯看出噪聲對分布圖的影響，或使小圓變大，或使直線變粗。由此驗證了低維流形假設的正確性。以上是圖像塊大小等于2時的情況，我們相信圖像塊更大時在塊空間也會呈現(xiàn)出一定的規(guī)律。

2 算法描述

輸入:含噪圖像u，塊尺寸參數(shù)n，搜索范圍尺寸參數(shù)M，距離閾值T，高斯函數(shù)平滑度參數(shù)h1、h2，間隔參數(shù)step，初始化偏移圖像ur=u。

第一步:塊匹配。對偏移圖像中的每個圖像塊 (記為當前塊)，在一定的搜索范圍內找到與之相似的塊。即視那些與當前塊距離小于一定閾值的塊為相似的塊。形式化的描述塊匹配如下。簡潔起見，以下用k表示位置(x，y)。令分別表示帶噪圖像u和偏移圖像ur中所有塊的集合，其中分別是將帶噪圖像和偏移圖像中第k塊 (圖像塊以k為中心，大小為(2n+1)×(2n+1))按列堆疊起來的向量，K是塊的總個數(shù)。特別的我們稱為位置k上的參考塊，簡稱參考塊。接下來對每個位置k，搜索所有與相似的塊q。理想情況下，搜索范圍應該是整幅帶噪圖像即Ω。但是這樣的搜索代價太大，本文將搜索范圍縮小到以位置k為中心的(2M+1)×(2M+1)大小的區(qū)域，并記這個區(qū)域內塊的集合為Ωk。易知Ωk為Ω的子集。從而所有與相似的塊的集合，即所在的流形可定義為:

圖4 圖3中甲乙丙區(qū)域及其中的塊分布

其中T是給定的距離閾值，距離小于該閾值的兩塊被認為是相似的。是塊p與塊q的距離。還有一點需要說明的是，在實際中不必對每個塊都進行塊匹配，也就是可以隔幾行隔幾列取塊進行塊匹配。本文中行間隔和列間隔取相等的值，記為step。

第三步:加權累加?，F(xiàn)在需要將復原后低維流形中的圖像塊放回到它原來的位置。由于搜索范圍的重疊和圖像塊的重疊，這一步將需要兩次的加權累加。在前面一步中我們得到了去噪后的低維流形對應的低秩矩陣Ak，現(xiàn)在將Ak中的塊加上偏移量放回到原來的位置。放回的時候還有講究。由于在第一步中搜索范圍Ωk有重疊，因此將會出現(xiàn)同一位置上的塊被多次處理的情況，即該位置上會有多個處理結果。我們需要對這些結果加權累加以得到該位置上圖像塊的去噪結果。令給定位置k上的塊被處理過J次，第j次的結果為則第j次對應的權為

其中參數(shù)h1用來控制高斯函數(shù)的平滑度。得到權值后，進行累加，則位置k上的塊的去噪結果k為

類似的，參數(shù)h2用來控制高斯函數(shù)的平滑度表示取絕對值。最終去噪圖像在位置k處的像素值為

3 實驗結果

本文采用兩種客觀標準來評價去噪效果。第一種是峰值信噪比

式中:u0——原始無噪圖像——無噪圖像的像素總數(shù)。——去噪后的圖像，·F——Frobenius范數(shù)。PSNR越大去噪效果越好。另一種是SSIM指數(shù)［13］，同樣越大越好。為簡單起見，在實驗中固定參數(shù)n=1，M=10，h1=1，h2=25，step=3，T=2σ，σ為所加噪聲的方差。本文的噪聲都是與原始無噪圖像不相關的加性高斯白噪聲。

圖5展示了對噪聲方差為20的House圖像的去噪過程，圖像長寬都為256，進行了4次迭代。從該圖可以看出經(jīng)過迭代后去噪效果逐漸提高。

圖6是對噪聲方差為20的House和Lena含噪圖像的去噪對比圖，圖像長寬都為256。其中圖6(c)、(f)是本文算法迭代5次時的結果。由圖可以看出本文所提方法確實能達到很好的去噪效果，視覺效果已與當前去噪效果最好的算法 BM3D［14］差別不大。

表1和表2分別是對加入不同方差噪聲后的House和Lena圖像的去噪結果，圖像長寬都為128，迭代次數(shù)不一。表中數(shù)據(jù)上一行是PSNR值，下一行是SSIM值。初步結果表明本文所提算法去噪效果良好，與BM3D算法差距較小。值得一提的是前面所用的實驗參數(shù)還不是最優(yōu)的，如何選擇好的參數(shù)是我們進一步研究的內容。

表1 圖像“House”去噪

表2 圖像"Lena"去噪

4 討論

前面提到過，傳統(tǒng)方法通過最優(yōu)化含約束項的目標函數(shù)來實現(xiàn)去噪。本文的方法也有類似的形式，約束是無噪的圖像塊空間中具有低維流形。令分別是圖像塊空間中無噪流形和被噪聲污染流形的集合，其中M是流形的個數(shù)。則要得到降噪的圖像塊可通過解下面的最優(yōu)化問題

其中Dm和Am分別是重新確定坐標原點后對應的矩陣，R(Am)是Am的秩，p指定范數(shù)的類型，λm是調節(jié)參數(shù)。進一步的，我認為NLM和BM3D算法的成功可以看成是低維流形假設的特殊情況:復原塊空間中的0維流形，即他們試圖將受噪聲污染變成的小球復原成一個點。

5 結束語

本文進一步分析了基于塊處理的圖像去噪。前人的低維流形假設指出:真實無噪圖像在塊空間中普遍存在低維流形，而噪聲的污染使這些低維流形變形。本文用實際例子給出了低維流形的幾何解釋。并根據(jù)該假設給出了通過復原塊空間中的低維流形來去噪的一個迭代算法。在算法上給出了偏移圖像的概念。初步實驗表明本算法能達到良好的去噪效果，與當今最好的去噪算法有可比性。本算法還有很大改進的空間，未來的工作主要是尋找更好的塊匹配和低維流形恢復方法，考慮在其他領域的應用如圖像的超分辨率重建等。

［1］Beck A，Teboulle M.Fast gradient-based algorithms for constrained total variation image denoising and deblurring problems［J］.IEEE Transactions on Image Processing，2009，18(11):2419-2434.

［2］Louchet C，Moisan L.Total variation as a local filter［J］.SIAM Journal on Imaging Sciences，2011，4(2):651-694.

［3］Buades A，Coll B，Morel J M.A review of image denoising algorithms，with a new one ［J］.SIAM Journal on Multiscale Modeling and Simulation，2005，4(2):490-530.

［4］Dabov K，F(xiàn)oi A，Katkovnik V，et al.Image denoising by sparse3-D transform-domain collaborative filtering［J］.IEEE Transactions on Image Processing，2007，16(8):2080-2095.

［5］ZHANG L，DONG W，ZHANG D，et al.Two-stage image denoising by principal component analysis with local pixel grouping ［J］.Pattern Recognition，2010，43(4):1531-1549.

［6］Peyre G.Manifold models for signals and images［J］.Computer Vision and Image Understandin，2009，113(2):1077-3142.

［7］Niyogi P，Smale S，Weinberger S.A topological view of unsupervised learning from noisy data ［J］.SIAM Journal on Computing，2011，40(3):646-663.

［8］NI J，Turaga P，Patel V，et al.Example-driven manifold priors for image deconvolution ［J］.IEEE Transactions on Image Processing，2011，20(11):3086-3096.

［9］Wetzler A，Kimmel R.Efficient beltrami flow in patch-space［J］.Scale Space and Variational Methods in Computer Vision，2012(6667):134-143.

［10］Wright J，Ganesh A，Rao S，et al.Robust principal component analysis:Exact recovery of corrupted low-rank matrices via convex optimization［C］//Vancouver:Twenty-Fourth Annual Conference on Neural Information Processing Systems，2009:1-44.

［11］LIN Z，CHEN M，MA Y.The augmented lagrange multiplier method for exact recovery of corrupted low-rank matrices［R］.Urbana-Champaign:University of Illinois at Urbana-Champaign Technical Report，2010:1-23.

［12］LIN Z，CHEN M，MA Y.Inexact augmented lagrange multiplier(ALM)method［CP/OL］.［2012-05-04］.http://perception.csl.uiuc.edu/matrix-rank/sample_code.html.

［13］WANG Z，Bovik A C，Sheikh H R，et al.Image quality assessment:From error visibility to structural similarity［J］.IEEE Transcactions on Image Processing，2004，13(4):600-612.

［14］Dabov K.Block-matching and 3D filtering(BM3D)algorithm and its extensions［CP/OL］.［2012-05-04］.http://www.cs.tut.fi/～ foi/GCF-BM3D/.