基于块空间中低维流形恢复的图像去噪

张宗良，颜永丰

(西北农林科技大学信息工程学院，陕西杨凌712100)

0 引言

图像去噪具有极高的理论和实用价值，近年来数字图像的急剧增多更加突显了图像去噪的重要意义。传统的去噪方法通常先假设图像具有某种良好的性质，如:具有较小的总变分［1-2］，或者在小波基中稀疏，然后通过最小化与该性质对应的目标函数来去噪。

最近几年利用图像块之间的相似性来去噪取得了极大的成功。文献［3］最早研究了基于块处理的图像去噪，并提出了非局部均值 (NLM)算法。该算法通过加权累加非局部块来估计当前所求块。权值由非局部块与当前所求块的相似度决定。在非局部块方法的基础上，文献［4］提出了一种全新的算法BM3D。该算法先用块匹配算法找到那些与当前块相似的块，再将这些相似块组成一个柱体。由于该柱体具有极大的冗余性，对该柱体进行三维离散余弦变换、阈值截断、余弦逆变换后便可去噪。文献［5］的算法与BM3D算法的不同之处是它对柱体进行的是PCA变换。

上述基于块相似的去噪简单有效，但是只有少量的文章［6-9］给出其背后的数学解释。总结这些文章我们知道无噪图像中的相似块会在图像块空间中形成低维流形，而在带噪图像中低维流形受到噪声污染而变形。复原这些低维流形便可去噪。本文据此提出了一个多步迭代的去噪算法。与大多数类BM3D算法的两步迭代相比，本文算法的多步迭代更具有一般性。

1 图像块空间中的低维流形

本文只考虑灰度数字图像。一幅图像是从平面区域到像素域的映射u:D→［0，255］ R，其中D Z2是图像所在的平面区域。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。这里有个细节问题，当 (x，y)靠近图像的边界时，邻域Nx，y有可能超出图像所在的平面区域即Nx，yD。此时可对图像进行扩展，本文采用的扩展方法是使扩展部分与原图像沿边界对称。接下来我们进一步的记由图像块px，y的像素按行列顺序排列而成的向量为块向量px，y，称块向量所在的空间为块空间。

下面我们考察图像块向量在块空间中的分布是什么样的。真实自然的图像具有很多良好的性质。例如图像中含有大量相同或相似的块，这些相同或相似的块向量将分布在块空间中一个半径很小的球中。如图1中的标为1、2的块分别与标为R1、R2的块相似。还有些图像块由于光照强度的渐变等原因而有规律的变化，这些图像块向量也将按一定的规律分布在块空间中。由此可见，特定区域内的所有图像块向量将在块空间中形成若干个流形［6-8］。这里的“流形”可以理解为是由相似块或有规律的块组成的一个有规律可循的组织或结构。图2是一个二维块分布示意图，在该图中，块空间是二维的即每个块只含二个像素。从左图可看到若干条直线 (近似)，每条直线可看成一个低维的流形。右图是对应的受到噪声污染后的流形，可以看到噪声的影响直观上是使直线变粗。可见，只要我们能从受噪声污染的流形 (图1b和图2b)中恢复出原本低维的流形(图1a和图2a)，再将恢复后流形中的点即图像块放回到该图像块原来所在的位置便能达到去噪效果。

接下来看看真实图像中的情况。同样在块大小等于2(即一个块只含左右两个像素)的情况下画出块分布图。图4(a)、(b)、(c)分别是图3中实线条框出的甲、乙、丙三个区域及其中的图像块在块空间中的分布图。左边的是无噪的情况，右边是加入方差为5的高斯白噪声的情况。从图4可以看出块向量在块空间中的分布确实有一定的规律，或聚成小圆，或近似排成直线。从图4还能明显看出噪声对分布图的影响，或使小圆变大，或使直线变粗。由此验证了低维流形假设的正确性。以上是图像块大小等于2时的情况，我们相信图像块更大时在块空间也会呈现出一定的规律。

2 算法描述

输入:含噪图像u，块尺寸参数n，搜索范围尺寸参数M，距离阈值T，高斯函数平滑度参数h1、h2，间隔参数step，初始化偏移图像ur=u。

第一步:块匹配。对偏移图像中的每个图像块 (记为当前块)，在一定的搜索范围内找到与之相似的块。即视那些与当前块距离小于一定阈值的块为相似的块。形式化的描述块匹配如下。简洁起见，以下用k表示位置(x，y)。令分别表示带噪图像u和偏移图像ur中所有块的集合，其中分别是将带噪图像和偏移图像中第k块 (图像块以k为中心，大小为(2n+1)×(2n+1))按列堆叠起来的向量，K是块的总个数。特别的我们称为位置k上的参考块，简称参考块。接下来对每个位置k，搜索所有与相似的块q。理想情况下，搜索范围应该是整幅带噪图像即Ω。但是这样的搜索代价太大，本文将搜索范围缩小到以位置k为中心的(2M+1)×(2M+1)大小的区域，并记这个区域内块的集合为Ωk。易知Ωk为Ω的子集。从而所有与相似的块的集合，即所在的流形可定义为:

图4 图3中甲乙丙区域及其中的块分布

其中T是给定的距离阈值，距离小于该阈值的两块被认为是相似的。是块p与块q的距离。还有一点需要说明的是，在实际中不必对每个块都进行块匹配，也就是可以隔几行隔几列取块进行块匹配。本文中行间隔和列间隔取相等的值，记为step。

第三步:加权累加。现在需要将复原后低维流形中的图像块放回到它原来的位置。由于搜索范围的重叠和图像块的重叠，这一步将需要两次的加权累加。在前面一步中我们得到了去噪后的低维流形对应的低秩矩阵Ak，现在将Ak中的块加上偏移量放回到原来的位置。放回的时候还有讲究。由于在第一步中搜索范围Ωk有重叠，因此将会出现同一位置上的块被多次处理的情况，即该位置上会有多个处理结果。我们需要对这些结果加权累加以得到该位置上图像块的去噪结果。令给定位置k上的块被处理过J次，第j次的结果为则第j次对应的权为

其中参数h1用来控制高斯函数的平滑度。得到权值后，进行累加，则位置k上的块的去噪结果k为

类似的，参数h2用来控制高斯函数的平滑度表示取绝对值。最终去噪图像在位置k处的像素值为

3 实验结果

本文采用两种客观标准来评价去噪效果。第一种是峰值信噪比

式中:u0——原始无噪图像——无噪图像的像素总数。——去噪后的图像，·F——Frobenius范数。PSNR越大去噪效果越好。另一种是SSIM指数［13］，同样越大越好。为简单起见，在实验中固定参数n=1，M=10，h1=1，h2=25，step=3，T=2σ，σ为所加噪声的方差。本文的噪声都是与原始无噪图像不相关的加性高斯白噪声。

图5展示了对噪声方差为20的House图像的去噪过程，图像长宽都为256，进行了4次迭代。从该图可以看出经过迭代后去噪效果逐渐提高。

图6是对噪声方差为20的House和Lena含噪图像的去噪对比图，图像长宽都为256。其中图6(c)、(f)是本文算法迭代5次时的结果。由图可以看出本文所提方法确实能达到很好的去噪效果，视觉效果已与当前去噪效果最好的算法 BM3D［14］差别不大。

表1和表2分别是对加入不同方差噪声后的House和Lena图像的去噪结果，图像长宽都为128，迭代次数不一。表中数据上一行是PSNR值，下一行是SSIM值。初步结果表明本文所提算法去噪效果良好，与BM3D算法差距较小。值得一提的是前面所用的实验参数还不是最优的，如何选择好的参数是我们进一步研究的内容。

表1 图像“House”去噪

表2 图像"Lena"去噪

4 讨论

前面提到过，传统方法通过最优化含约束项的目标函数来实现去噪。本文的方法也有类似的形式，约束是无噪的图像块空间中具有低维流形。令分别是图像块空间中无噪流形和被噪声污染流形的集合，其中M是流形的个数。则要得到降噪的图像块可通过解下面的最优化问题

其中Dm和Am分别是重新确定坐标原点后对应的矩阵，R(Am)是Am的秩，p指定范数的类型，λm是调节参数。进一步的，我认为NLM和BM3D算法的成功可以看成是低维流形假设的特殊情况:复原块空间中的0维流形，即他们试图将受噪声污染变成的小球复原成一个点。

5 结束语

本文进一步分析了基于块处理的图像去噪。前人的低维流形假设指出:真实无噪图像在块空间中普遍存在低维流形，而噪声的污染使这些低维流形变形。本文用实际例子给出了低维流形的几何解释。并根据该假设给出了通过复原块空间中的低维流形来去噪的一个迭代算法。在算法上给出了偏移图像的概念。初步实验表明本算法能达到良好的去噪效果，与当今最好的去噪算法有可比性。本算法还有很大改进的空间，未来的工作主要是寻找更好的块匹配和低维流形恢复方法，考虑在其他领域的应用如图像的超分辨率重建等。

［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，Foi 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/.