2.2 平滑去噪
图像噪声是指存在于图像数据中不必要的或多余的干扰信息,产生于图像的采集、量化或传输过程,对图像的后处理、分析均会产生极大的影响,因此一种好的去噪方法在去除噪声的同时,还需要保持图像的边界和细节。早期的去噪方法多为空间滤波(Spatial Filter)。Donoho D等人于1995年[10]、Simoncelli E P等人于1996年[11]提出的小波变换(Wavelet Transform)在该领域实现了巨大的突破。2005年非局部(Non-local)方法被首次提出,将图像去噪带入了新的阶段。近年来,随着深度学习的不断发展,基于神经网络的方法也在不断涌现。本节将分别介绍这四大类方法。
2.2.1 空间滤波
空间滤波由一个邻域和对该邻域内像素执行的预定义操作组成,滤波器中心遍历输入图像的每个像素点之后就得到了处理后的图像。每经过一个像素点,邻域中心坐标的像素值就替换为预定义操作的计算结果。若在图像像素上执行的是线性操作,则该滤波器称为线性空间滤波器,反之则称为非线性空间滤波器。
1.线性空间滤波器
常见的线性空间滤波器有平滑线性滤波和高斯滤波器[12]。平滑线性滤波器的输出是包含在滤波器模板邻域内像素的简单平均值,因此也称为均值滤波器。通常用下式来表示m×n 的滤波器对 M×N 的图像进行平滑线性滤波操作:
其中,w表示滤波器,f表示输入图像,简单来说,即进行加权均值的计算。若滤波器模板内所有系数相同,则称为盒装滤波器。
高斯滤波器的模板系数随着与模板中心距离的增大而减小。如图2-14所示,对于一个(2k+1)×(2k+1)大小的高斯滤波器模板,取模板中心为原点后得到各位置的坐标,将坐标分别带入高斯公式即可得到系数,公式如下:
图2-14 3×3滤波器模板示意图
通过公式我们不难发现,对于高斯滤波模板最重要的参数就是 σ, σ 决定了数据的离散程度,σ 越大,分布越分散,即模板内各系数差别不大,效果类似于平滑线性空间滤波;σ越小,分布越密集,即模板中心的系数高于其他部分。因此,相比于平滑线性空间滤波器,高斯滤波器的平滑效果更柔和,边缘保留也更好。
2.非线性空间滤波器
中值滤波是一种非线性数字滤波器。与上述滤波不同,其首先是对邻域内的值进行排序,中值作为输出。该滤波器消除椒盐噪声的效果比平滑线性空间滤波器更好。
双边滤波(Bilateral Filter)[12]不仅如高斯滤波一样考虑到像素间的欧氏距离,还关注到像素范围域中的辐射差异(如像素与中心间的相似程度、颜色强度和深度距离等),这两个权重域分别称为:空间域(Spatial Domain S)和像素范围域(Range Domain R)。权重的计算公式如下:
其中,σd和 σr是平滑参数,I(i, j) 和 I(k, l) 分别是像素 (i, j) 和 (k, l) 的强度函数。计算权重之后,对结果做归一化,公式如下:
其中,ID(i, j) 是去噪后的像素强度。
3.代码
Python的OpenCV包给出了上述四种滤波器的实现。
1. import numpy as np 2. import cv2 3. ######## 四种不同的滤波器 ######### 4. img = cv2.imread('test.jpg') 5. # 平滑线性滤波 6. img_mean = cv2.blur(img, (5,5)) 7. # 高斯滤波 8. img_Guassian = cv2.GaussianBlur(img, (5,5),0) 9. # 中值滤波 10. img_median = cv2.medianBlur(img, 5) 11. # 双边滤波 12. img_bilater = cv2.bilateralFilter(img,9,75,75)
2.2.2 小波阈值去噪
对于二维图像信号,可分别在水平和垂直的方向使用滤波器,从而实现二维小波多分辨率分解。图2-15为经过小波三层多分辨率分解之后子图像的划分。
图2-15 小波三层多分辨率分解示意图
其中,LL子带是在两个方向利用低通小波滤波器产生的小波系数,表示图像的近似;HL子带是在行方向利用低通小波滤波器后,再用高通小波滤波器在列方向卷积产生的系数,表示图像水平方向的奇异特性;LH子带是在行方向利用高通小波滤波器后,再用低通小波滤波器在列方向卷积产生的系数,表示图像垂直方向的奇异特性;HH子带是在两个方向利用高通小波滤波器产生的小波系数,表示图像的对角边缘特性。通常情况下,噪声部分包含在HL、LH和HH中,只要对这三部分作相应的小波系数处理,最后进行重构即可达到消噪的目的。
因此小波去噪的基本思路如下。1)二维图像的小波分解。选择一个小波和小波分解的层次N,计算信号s到第N层的分解。2)对高频系数进行阈值量化。对1~N的每一层选择阈值,并对该层的高频系数进行软阈值量化处理。3)二维小波重构。根据小波分解第N层的低频系数和经过修改的第一到第N层的各层高频系数,计算二维信号的小波重构。
阈值的选择是离散小波去噪中最关键的一步,阈值处理函数分为硬阈值和软阈值,对于硬阈值,如式(2.4);对于软阈值。
其中,w为小波系数;λ为给定阈值;sign(w)表征取值符号。当w>0时,sign(w)=1;当w<0时,sign(w)=-1;当w=0时,sign(w)=0。
2.2.3 非局部方法
1. NL-means
A. Buades等人于2005年首次将Non-Local理论应用于图像去噪[13],提出了基于图像所有像素的非局部均值算法(NL-means),他们对去噪方法Dh做了如下定义:
v=Dhv+n(Dh, v)
其中,v是带噪声的图像,n 是滤波参数,通常取决于噪声的标准偏差。理想情况下,Dhv 比 v 更加平滑,n(Dh, v) 为白噪声。给定离散噪声图像 v={v(i)|i∈I},对于像素 i,估计值 NL[v](i) 为图像中所有像素的加权平均值:
其中,权重 w(i, j) 取决于像素 i 、j 间的相似性,0≤w(i, j)≤1且。而所谓的相似性,则进一步取决于强度灰度级向量v(N)i 和v(N )j 的相似性。其中,Nk表示以像素 k 为中心固定大小的正方形邻域,该相似度被定义为加权欧式距离的递减函数:, a 是高斯核的标准差。与v(Ni) 相近的灰阶邻域像素在平均值上有较大的权重:
其中,Z(i) 是归一化常数,参数h表示影响过滤程度。通过控制指数函数的衰减,可以降低权重。
该方法不仅比较了单个点的灰度,而且在整个邻域中进行了几何构型的对比。这使得它具有比邻域筛选器更强的稳定性。从图2-16可以看出,像素 q3和像素p具有相同的灰度值,但因为邻域差别很大,所以最终权重 w(p, q3) 近似为零。
图2-16 NL-means算法思想
对于图2-16,高斯滤波算法与NL-means算法处理结果对比如图2-17所示。
图2-17 去噪算法结果对比
2. BM3D
BM3D(Block-Matching and 3D Filtering)[14]是当前效果最好的算法之一。它与NL-means算法具有相似之处在于同样运用到了非局部块匹配的思想,但其复杂度远高于NL-means。如图2-18所示,BM3D共包含两大步骤,每一步里又包括相似块分组、协同滤波和聚合。找到相似块后,不同于NL-means的均值处理,BM3D会对其做域转换操作,利用协同过滤以降低相似块自身的噪声,并在聚合阶段对相似块进行加权处理。
图2-18 BM3D算法流程图
步骤一:基础估计(Basic Estimate)
1)相似块分组。首先在噪声图像中选择参照块,在参照块周围适当大小区域内进行搜索,寻找若干个差异度最小的块,与参照块一起整合成一个三维的矩阵。距离函数的公式如下:
其中,ZXR表示以像素 R 为中心的参照块,Zx为搜索块,为步骤一的块大小,为 L2距离,Y′ 为硬阈值操作,为归一化后的二维线性操作。根据距离找到的相似块集合可用如下公式表示:
其中,为判断是否相似的超参数。
2)协同滤波。先使用归一化的3D线性变换降低相似块中的噪声,再使用对应的反变换得到处理后的相似块:
其中,为归一化的3D线性变换,为对应的反变换,Y为硬阈值操作,为处理后的相似块。这一步可以有效避免NL-means均值操作导致的相似块信息冗余和引入参照块自身包含的噪声这两大问题。
3)聚合。对上一步处理后得到的相似块进行加权平均值操作得到目标块的像素值:
其中,为权重,为式(2.6)硬阈值操作后非0系数的个数。
步骤二:最终估计(Final Estimate)
1)相似块分组。对步骤一初步处理后的图像,重新计算 L2距离,得到相似块集合:
其中,为判断是否为相似块的超参数,为相似块集合。
2)协同过滤。对第一步得到的相似块集合做域变换后,使用维纳收缩系数加权,再经过反变换得到新的集合:
其中,为维纳收缩系数。
3)聚合。与步骤一的聚合操作相似,最终结果与式(2.7)一致:
官方提供了MATLAB实现(http://www.cs.tut.fi/~foi/GCF-BM3D/)和C++ 实现(http://www.ipol.im/pub/art/2012/l-bm3d/),同时https://github.com/ericmjonas/pybm3d 提供了Python的接口pybm3d。
2.2.4 基于神经网络的方法
1. MLP
Harold Christopher Burger等人于2012年提出了使用多层感知机(Multi Layer Perceptron, MLP)实现图像去噪[15]的思想,从而探讨深度学习在文字识别任务上是否有必要。整体的网络结构非常简单,含有两层隐藏层的MLP可用如下公式表示:
f(x)=b3+W3 tanh(b2+W2tanh(b1+W1x))
训练时,以三步长将图像分割为图像块,随机选择一个原始图像块 x,加入高斯白噪声(White Gaussian Noise),生成对应的噪声图像 y,网络将根据 f(x) 和 y 误差(如:(f(x)-y)2)的反向传播进行参数更新。
2. LLNet
Kin Gwn Lore等人于2017年提出了The Low-Light Net(LLNet)[16]的概念。LLNet通过引入序贯相似性检测算法(Stacked Sparse Denoising Autoencoder, SSDA)的思想实现低噪度图片的自适应增强(增亮、去噪)。
首先介绍一下自动编码机(Auto-Encoder, AE)。AE属于非监督学习,由三层网络组成,其中输入层神经元数量与输出层神经元数量相等,中间层神经元数量少于输入层和输出层。对于每个训练样本,AE学习的目的是使输出信号与输入信号尽可能相同。流程可以用下式表示:
其中,s 是非线性函数(如sigmoid), W 和 W'分别为输入层到中间层和中间层到输出层的链接权值,b 和 b'分别为中间层和输出层的偏置。y 可视为 x 的有损压缩形式,信号 y经过解码层解码传到输出层,变为信号 z。计算x 和 z 的误差时可使用典型的平方误差,若x是位向量,则可以使用交叉熵。
使用AE对神经网络进行预训练,可以确定编码器的初始权重 W,然而受模型复杂度、训练集数据量以及数据噪声等问题的影响,通过AE得到的初始模型往往存在过拟合的风险。
降噪自动编码器(Denoising Auto-Encoder, DAE)就是在AE的基础上,为防止过拟合,对输入的数据加入噪声,从而提高模型的稳定性。DAE的示意图如图2-19所示。
图2-19 降噪自动编码器示意图
其中,x 为输入信号,DAE以特定概率将部分输入层节点的值置为0,进而得到含有噪声的输入信号 ,通过 计算 y 和 z,并对 z 与原始输入信号 x 做误差迭代,减小了测试样本与训练样本因分布不同而对模型造成的影响。不仅如此,DAE还可以多层堆叠形成栈式自编码器(Stacked Denoising Auto-encoders, SDA),即前一层的中间层输出是后一层的输入,通过逐层贪婪训练,对每层自编码都进行非监督训练,用训练好的第 K 层的输出作为第 K+1层自编码训练的输入。
前面提到过,中间层节点的个数一般小于输入输出层,因此中间层倾向于学习信号的内部规律。但如果设定的中间层节点数多于输入输出层,同时又希望自编码器能够学习到内部规律,就需要利用稀疏自编码器(Sparse Denoising Auto-encoder),即对中间层的激活函数增加稀疏性限制,保证同一时间内只有部分节点是活跃的。假设中间层的激活函数为sigmoid,那么节点输出1表示活跃,输出0表示不活跃,可以使用 KL 散度来衡量节点实际活跃程度与人为设定稀疏度 ρ 之间的相似性:
其中,m 为训练样本个数,为中间层第 j 个节点对第 i 个训练样本的响应输出,为中间层第 j 个节点的平均激活输出。一般设定 ρ 为0 .5或0 .01, KL散度越小意味着和 ρ 的相似度越高,因此可以将 KL 散度作为正则化项加入损失函数从而对网络的稀疏性进行约束。
LLNet正是利用了它的降噪能力和深度网络的复杂建模能力来学习低噪度图片的特征并生成含有最少噪声和更高对比度的图片。LLNet的结构如图2-20所示。
图2-20 LLNet框架结构
图2-20a是图2-20b和图2-20c中灰色模块的具体结构,自编码器含有多个隐藏层,编码部分采用非监督学习,解码部分的权重由编码部分转置得到,并根据误差的反向传播进行微调。图2-20b表示LLNet中的模块同时进行对比度增强和去噪。图2-20c中则是将两个模块分离开按顺序进行对比度增强和去噪,称为Staged LLNet(S-LLNet)。
LLNet的训练流程如图2-21所示,其中,训练图片被人工调低亮度并添加噪声。LLNet通过重构图像与原始图像误差的反向传播调整参数。
图2-21 LLNet训练流程