数字图像处理（十六）非局部均值去噪

一、前言

在之前我们已经介绍过许多图像去噪的方法，比如均值滤波、中值滤波、高斯滤波等。今天我们要介绍一种新的去噪方法——非局部均值去噪（the non local means, NL-means）。那么NL-means相较于之前的方法有什么不同的地方呢？
我们先思考之前的图像去噪方法，用均值滤波进行举例。均值滤波在对噪声图像中的某个位置的像素进行处理时，需要以该像素为中心，选取3x3或者5x5邻域大小的像素求和后取平均值作为该像素的像素值。可以看到，使用均值滤波进行去噪时，我们只使用了该像素的邻域区域（也就是局部区域），中值滤波类似。
NL-means来自于这篇论文“A non-local algorithm for image denoising”（发表于2005年）。该方法充分利用了图像中的冗余信息，在去噪的同时能最大程度得保持图像的细节特征。

二、NL-means

先观察下面的这张图像：

可以看到$p、q_1、q_2$三个区域比较相似，那么$q_1、q_2$相较于$p$可以说是重复区域，也就是说图像中存在着冗余信息。那么在对$p$区域进行去噪时，能否利用$q_1、q_2$这两个区域呢？NL-means就是基于这样的想法提出来的。

1.两个邻域块的相似度

如何衡量$A$与$B$两个邻域块的相似度呢？一种常见的衡量方法是MSE(均方误差)，当然还有其他的方法，例如PSNR，余弦相似度等等。

MSE的计算公式如下：
$$MSE(A,B)=\frac{1}{mn}\sum _{i=0}^{m-1}\sum _{j=0}^{n-1}{\left[A(i,j)-B(i,j)\right]}^2$$
其中$m$和$n$分别为两个邻域块的宽和高。

2.NL-means原理

使用NL-means对图像进行去噪，假设A点为当前待滤波点，我们需要为A设置两个窗口：搜索窗口和邻域窗口。搜索窗口的大小为$D\times D,D=2\times D_s+1$，邻域窗口的大小为$d\times d,d=2\times d_s+1$。如下所示：

那么对A点进行NL-means滤波后的像素值为：
$$\text{NL-means}=\sum _{B}w(A,B)\times I(B)$$
其中$I(B)$为B点的像素值；$w(A,B)$为$A、B$两个邻域块的相似度，其计算公式如下：
$$\begin{gathered} w(A,B)=\frac{1}{sum}{e}^{-\frac{MSE(A,B)}{h^2}} \\ sum=\sum _B{e}^{-\frac{MSE(A,B)}{h^2}} \end{gathered}$$
可以看到，$w(A,B)$的计算公式与$A,B$两点所在邻域窗口的相似度有关。其中h也是一个重要的参数，h越大去噪效果越好，但是图像越模糊，反之h越小去噪效果越差，但是去噪后的失真度越小。

3.数学理论推导

假设无噪声像素块的值为$f(x,y)$，噪声为$n(x,y)$，这里我们认为噪声是不相关的。那么带噪的像素块的值为$g(x,y)=f(x,y)+n(x,y)$。将$k$个相似的带噪像素块进行叠加后取均值得到
$$\bar{g}(x,y)=\frac{1}{k}\sum _{i=1}^k{g}_i(x,y)=\frac{1}{k}\sum _{i=1}^k{f}_i(x,y)+n_i(x,y)$$

对$\bar{g}(x,y)$取期望如下：
$$E\left[\bar{g}(x,y)\right]=E\left[\frac{1}{k}\sum _{i=1}^k{f}_i(x,y)+n_i(x,y)\right]=\frac{1}{k}\sum _{i=1}^{k}E\left[f_i(x,y)\right]+E\left[n_i(x,y)\right]$$
其中假设噪声的期望为0，$f_i(x,y)$为相似块，认为其期望相等，所以
$$E\left[\bar{g}(x,y)\right]=f(x,y)$$
对$\bar{g}(x,y)$取方差如下：
$$D\left[\bar{g}(x,y)\right]=D\left[\frac{1}{k}\sum _{i=1}^k{f}_i(x,y)+n_i(x,y)\right]=\frac{1}{k^2}D\left[\sum _{i=1}^k{f}_i(x,y)+\sum _{i=1}^k{n}_i(x,y)\right]=\frac{1}{k^2}D\left[\sum _{i=1}^k{f}_i(x,y)\right]+\frac{1}{k^2}D\left[\sum _{i=1}^k{n}_i(x,y)\right]$$
因为$D\left[{\sum }_{i=1}^k{f}_i(x,y)\right]$=$D\left[kf(x,y)\right]=k^{2}D\left[f(x,y)\right]=k^2{\sigma }_f^2$，$D\left[{\sum }_{i=1}^k{n}_i(x,y)\right]={\sigma }_{\eta 1}^2+{\sigma }_{\eta 2}^2+\cdots +{\sigma }_{\eta k}^2$，所以
$$D\left[\bar{g}(x,y)\right]={\sigma }_f^2+\frac{1}{k^2}({\sigma }_{\eta 1}^2+{\sigma }_{\eta 2}^2+\cdots +{\sigma }_{\eta k}^2)$$
又因为相似无噪声像素块的方差为0，所以最后的推导出的方差公式为：
$${\sigma }_{\bar{g}}^2=\frac{k}{k^2}{\sigma }_{\eta }^2=\frac{1}{k}{\sigma }_{\eta }^2$$
可以看到$\bar{g}(x,y)$的期望为$f(x,y)$，方差为$\frac{1}{k}{\sigma }_{\eta }^2$。所以说，通过非局部均值滤波的方法可以进行图像去噪。

4.代码链接

github

参考链接

非均质滤波的原理
https://www.cnblogs.com/luo-peng/p/4785922.html
https://zhuanlan.zhihu.com/p/355267754