基于图像块l0梯度最小化的边缘保持平滑算法

宋 昱，孙文赟

1)深圳大学电子与信息工程学院，广东深圳 518060；2)深圳大学深圳市媒体信息内容安全重点实验室，广东深圳 518060；3)深圳大学广东省智能信号处理重点实验室，广东深圳 518060

边缘保持图像平滑指在滤除纹理和细节的同时保持图像的边缘结构信息，它是很多计算机视觉与图形学，如边缘检测[1]、图像分割[2]、图像分类和基于内容的图像编辑等的基础，被广泛用于图像增强[3]、色调映射、图像抽象和基于内容的纹理编辑等领域．

边缘保持图像平滑算法可分为传统算法和基于学习的算法．传统算法有2种分类方式：一是根据平滑过程中的处理范围分类，二是根据所用信息的类型分类．根据平滑过程中的处理范围，可将传统算法分为基于局部信息的算法和基于全局信息的算法．前者通过对每个像素的邻域信息进行处理得到输出像素；后者对整幅图像同时处理．根据所用信息的类型，传统算法可被分为基于像素灰度距离的算法和基于其他准则的算法．前者根据图像梯度或像素灰度距离决定被平滑内容；后者采用其他衡量准则决定被平滑内容．双边滤波器[4]通过对附近的像素进行加权平均对图像进行平滑．域变换滤波器[5]和引导滤波[6]是双边滤波器和相关算法的改进．各向异性扩散和一致增强扩散[7-8]根据偏微分方程(partial differential equation, PDE)对像素点进行迭代滤波．局部拉普拉斯(local Laplican, LL)滤波[9]是一种对双边滤波和基于PDE的滤波器的改进．FATTAL[10]使用图像变换来平滑图像．文献[4-10]中的算法可以滤除小尺度的细节并保留大尺度的边缘，但不能滤除强纹理并同时保留弱边缘．KASS等[11]利用一个邻域中的平滑局部直方图对图像进行平滑．BAO等[12]提出树滤波算法，使用中心像素点与邻域中像素点的树距离和双边距离进行滤波．文献[13]中的算法使用一组图像块的局部统计特性对像素点进行滤波．边窗口滤波算法[14]使用一个边窗口对每个像素点进行平滑．因为文献[11-14]中的算法不依赖于图像梯度，所以这些算法对强纹理和弱边缘的处理效果比基于像素灰度距离的算法更好．然而，基于局部信息的算法对复杂图像的平滑效果并不好．因为基于局部信息的算法中所用到的计算主要限制在每个像素点的邻域中，所以该算法通常运行速度较快，且计算复杂度比较低． 但是，此类算法中，像素点的滤波结果依赖邻域像素的统计特性，整幅图像的信息并未被充分考虑．

基于全局信息的图像平滑算法对整幅图像同时滤波．算法通常先构建一个包含一项数据忠诚项与一项正则项的目标函数，再通过最小化目标函数得到平滑后的图像．数据忠诚项衡量了平滑后图像与输入图像之差的l2范数的平方，而正则项描述了平滑后图像的特征．不同算法含有相同的数据忠诚项，但是正则项不同．加权最小二乘(weighted least square, WLS)滤波算法[15]衡量了平滑后图像的加权平方梯度．基于l1范数的图像滤波算法[16]的目标函数中有2项正则项．l0梯度最小化算法[17]衡量了平滑图像梯度的l0范数．文献[15-17]中的算法都可以滤除小尺度的纹理并保留大尺度的边缘，但不能同时滤除强纹理并保留弱边缘．基于局部极值的平滑算法[18]通过计算最小和最大包络的均值得到平滑图像．基于相对全变差的算法[19]中，目标函数的正则项衡量了平滑图像的相对全变差．基于迭代全局优化的算法[20]在迭代全局优化的框架中使用了选择性的梯度抑制和引导梯度修正对图像进行平滑．文献[18-20]中的算法对于强纹理的滤波效果优于基于全局信息和像素灰度距离的算法，但它们在复杂图像上的滤波效果依然不佳．基于全局信息的算法因为需要在优化过程中求解一个大型稀疏线性方程组，所以计算复杂度通常较高．它与基于局部信息的算法不同，在滤波过程中考虑了由整幅图像提供的信息．

近年来，深度学习因具有端到端学习的特点，被应用于边缘保持算法．基于深度学习的边缘保持算法可被分为2类：第1类算法通过使用其他传统算法(基于局部信息和全局信息的算法)的平滑结果作为监督信号有监督的训练一个深度神经网络[21]；第2类算法使用基于输入图像和平滑图像的目标函数无监督的训练一个深度神经网络．通过使用大量的自然图像，在训练过程中对目标函数进行最小化从而得到平滑图像[22]．神经网络可以在训练过程中学习到平滑准则．无监督算法的运行速度通常比传统算法更快．

本研究提出一种基于图像块的l0梯度最小化算法(image patch basedl0gradient minimization algorithm, 简称IP-l0算法)，不直接对整幅图像而是对输入图像中的图像块进行图像平滑处理．算法根据图像块的不同特性调整算法的平滑力度，先在每个图像块中实现l0梯度最小化，再进行融合，从而将全局算法转变为局部算法，得到了更好的滤波结果．

1 基于图像块的l0梯度最小化算法

1.1 l0梯度最小化算法的缺点

l0梯度最小化属于基于全局信息和像素灰度距离的算法．该算法根据像素灰度距离选择被平滑区域，像素灰度距离较小的区域被平滑，而像素灰度距离较大的区域被保留．换句话说，该算法含有一个假设，纹理和细节的像素灰度距离应比边缘的小．但在实际应用中，上述假设并非总是成立，在一些有着复杂背景的自然图像中，纹理和细节的像素灰度距离可能比边缘大很多．此时，若直接使用l0梯度最小化算法将会平滑边缘，但也保留了纹理和细节．图1是一个成功应用l0梯度最小化算法的例子．由图1(a)可见，原始图像不包含强纹理，且边缘有较高的清晰度．采用l0梯度最小化处理后，小尺度的纹理被平滑了，同时大尺度的边缘被保留了，见图1(b)．图2为一个失败的l0梯度最小化算法应用例子．由图2(a)可见，原始图像包含很多强纹理，且部分边缘很弱．图像平滑处理后，强纹理被保留了，而弱边缘被平滑了，见图2(b)．由图2(b)可见，若使用较小的权重参数，如λ=0.01， 则边缘和纹理会被同时保留．换句话说，不管如何调整l0梯度最小化算法的权重参数，都不可能在滤除强纹理的同时保留弱边缘．为更直观的观察平滑处理对效果，选取原图中任一行的像素点(图1(a)和图2(a)中黄色线所覆盖的像素)，分析这些像素在图像平滑前后的灰度值分布，结果如图1(c)和图2(c)．其中，c为像素点的列数；I(c)和S(c)分别为原始图像和平滑图像在此位置的灰度值．

图1 使用l0梯度最小化算法对弱纹理、强边缘图像平滑结果Fig.1 The image smoothing results after applying l0 gradient minimization to images with strong textures or weak edges

图2 使用l0梯度最小化算法对强纹理、弱边缘图像平滑结果Fig.1 The image smoothing results after applying l0gradient minimization to images with strong textures or weak edges

本研究提出的IP-l0算法将整个输入图像分割为不同的图像块，通过计算其相对全变差来描述图像块的特征，来调节对应的图像块的l0梯度最小化参数，将全局算法转变为局部算法，也获得了更好地滤波效果．

1.2 l0梯度最小化算法

l0梯度最小化算法的目标函数[17]为

(1)

其中，Ip为输入图像I中像素p的值，p=(x,y)是其坐标索引；Sp为平滑图像S的像素灰度值；λ是平衡数据忠诚项和正则项的权重参数，数据忠诚项衡量了平滑后图像与输入图像之差的l2范数的平方，而正则项描述了平滑后图像的特征；C(S)为平滑图像的非零梯度的个数，数学表达式[17]为

C(S)=#{p||∂xSp|+|∂ySp|≠0}

(2)

其中， ∂xSp和∂ySp为平滑图像在x方向和y方向上的梯度； #{·}为集合内元素的个数．

目标函数(1)的解无法通过直接计算得到，参照文献[17]通过引入辅助变量扩展原来的项，将原始目标函数分解为2个子问题，再交替优化，最终求出目标函数的一个近似问题的解．在l0梯度最小化的情形中，这2个子问题都有闭式解．求解目标函数(1)的方法简单介绍如下：参照文献[17]引入对应于平滑图像梯度∂xSp和∂ySp的2个辅助变量hp和vp， 并将目标函数写为

β[(∂xSp-hp)2+(∂ySp-vp)2]}

(3)

其中，C(H,S)=#{p||hp|+|vp|≠0}，H={hp},V={vp}；β是一个用于控制辅助变量(hp,vp)和对应变量(∂xSp, ∂ySp)相似性的参数．当β取较大值时，式(3)是式(1)的一个很好的近似．将式(3)分解为2个单独求解的子问题，子问题1求解S， 子问题2求解(hp,vp)， 具体求解过程见文献[17]．在每一次迭代中，子问题中有一组变量的值是固定的，其值来自上一次迭代．

1.3 图像块统计特征的选择

如果直接对输入图像使用l0梯度最小化，那么会出现小尺度的细节被平滑，而大尺度的特征被保留的现象．XU等[19]提出的相对全变差，可以很好地区分纹理和边缘．因此，本研究使用相对全变差描述图像块的特征．

相对全变差R包括窗口全变差D和窗口固有变差L． 其中，D是图像梯度绝对值的线性卷积，则像素p在x方向和y方向上的窗口全变差分量[19]分别为

(4)

其中，q为像素p的局部邻域Np的像素点．由式(4)可见，图像梯度的绝对值被叠加，因此该式可描述窗口Np中的变化．权重函数gp, q与双边滤波中的类似，按照空间距离的定义，其数学表达[19]为

(5)

其中，σ是反映窗口空间尺度的参数．如果图像中包含显著的纹理，那么在纹理区域和边缘附近的像素，其窗口全变差会比较大．这说明那些在显著区域的像素，会有比较大的窗口全变差．

窗口全变差D本身并不能区分纹理和边缘，因此文献[19]引入窗口固有变差L，则像素p在x方向和y方向上的窗口固有变差分量为

(6)

由式(6)可见，窗口固有变差L在运算时不对梯度取绝对值，因此可得到以下结果：如果窗口中的所有梯度符号相同，梯度作用就会累积；反之，如果窗口中的梯度符号相反，梯度作用会抵消．对于纹理区域中的像素点，因梯度符号不一致，其窗口固有变差值会比较小．对于边缘附近的像素点，因梯度符号相同，其窗口固有变差值会比较大．相对全变差可定义为窗口全变差和窗口固有变差的比值，因此像素p的相对全变差[19]为

(7)

其中，经验常数ε取一个较小的值，目的是防止除数为0．

因为纹理区域内像素点的窗口固有变差值较小且窗口全变差值较大，而在边缘附近像素点的这两个值都较大，所以纹理区域内像素点的相对全变差值较大，而在边缘附近像素点的相对全变差值较小．可见，相对全变差的性质可用来决定输入图像中图像块的平滑力度．

1.4 将图像块统计特征引入l0梯度最小化

定义一个图像块集合系统，其图像域Ω被图像块集合系统覆盖，数学表达式为

(8)

其中，Nx为以像素x为中心的图像块Ix中像素点的下标集合，所有Nx的集合可覆盖整个图像；Ω′是所有x的集合．

本研究对每一个图像块Ix使用l0梯度最小化算法，其权重参数根据图像块的相对全变差来调整，则每个图像块Ix的目标函数为

(9)

(10)

目标函数(9)可通过交替优化进行求解．算法步骤与全局l0梯度最小化类似．子问题1的解为

其中，hx和vx是辅助变量，用于近似对应的图像梯度；F()和F-1()分别表示傅里叶变换和傅里叶逆变换；*为共轭运算符．

子问题2的解为

(hx, p,vx, p)=

(12)

其中，hx, p和vx, p是在像素点p处的图像块的辅助变量．权重参数λx为与图像块Ix的相对全变差r(Ix)成正比，即

λx=kr(Ix)

(13)

其中，k为比例参数．将λx设置为正比于图像块的相对全变差的原因如下：当图像块中存在很多纹理时，平滑力度要大，从而令权重参数较大，此时图像块的相对全变差r(Ix)也取一个较大的数值；当图像块中仅存在边缘时，平滑力度要小，从而令λx较小，此时， 图像块的r(Ix)也取一个较小的数值．按照这种方式设置权重参数，那些有着强纹理的图像块被平滑，而那些仅有边缘的图像块会被保留．

对所有图像块处理后，需要将这些图像块融合才能得到最终的平滑图像．因为图像块之间有重叠，而这些重叠区域的像素点会有多个平滑结果，取平均后可得到最终输出图像的像素灰度值Sp． 取平均再进行融合的步骤可写为

(14)

2 实验结果与分析

2.1 算法性能比较

为衡量IP-l0算法性能，比较IP-l0算法、l0梯度最小化算法[17]、LL滤波器算法[9]、基于相对全变差(relative total variation, RTV)算法[19]和基于树滤波(tree filtering, TF)算法[12]的边缘保持图像平滑结果．上述算法中，LL滤波算法属于基于局部信息和像素灰度距离的算法；TF算法属于基于局部信息和其他度量准则的算法；l0梯度最小化算法属于基于全局信息和像素灰度距离的算法；RTV算法属于基于全局信息和其他度量准则的算法．

限于篇幅，本研究从BSDS300分割数据集(hppts://www2.eecs.berkeley.edu/Research/Projects/ CS/vision/bsds/)选取含强纹理和含弱边缘的图像各1幅，两幅图像分辨率均是481×321像素．实验环境为：Intel®Core(TM) i5-4460S CPU@2.9 GHz，8 Gbyte RAM，Matlab R2015a．表1中给出了5种算法的参数设置，所有参数均取使算法达到最优效果时的数值．因为没有“真实平滑图像”的明确定义，无法使用如峰值信噪比或结构相似度等客观指标进行评价，本研究仅对不同算法进行主观比较．5种算法对两类图像的平滑结果见图3和图4．为更清楚地比较不同平滑算法的结果，从每幅图像中选择一个代表性区域(红色矩形部分)进行放大和平滑处理．同时为更直观的观察到相对全变差变化情况，输入图像的相对全变差图像见图3(f)和图4(f)．

表1 不同算法中的参数设置

图3 不同算法对强纹理图像的平滑结果Fig.3 The image smoothing results by using different methods for image with strong textures

图4 不同算法对弱边缘图像的平滑结果Fig.4 The image smoothing results by using different methods for image with weak edges

由图3可见，对于强纹理图像， LL算法对纹理的滤除效果并不好，原因是该算法属于基于局部信息和像素灰度距离的算法，强纹理的像素灰度距离比边缘的更大，因此在进行图像平滑时这些纹理也同时被保留了．TF算法采用了树距离平滑纹理，纹理滤除效果优于LL算法，但边缘保持效果不佳．使用TF算法会令图像平滑后边缘被扭曲，平滑图像中存在“阶梯”效应．RTV算法虽然可以得到相对较好的纹理平滑结果，但平滑效果太强以至某些边缘被模糊了．IP-l0算法的图像平滑结果优于其他算法，不仅图像中的强纹理被滤除了，边缘也被很好的保留下来．从全变差图像可见，有强纹理的区域具有较大的全变差，而边缘区域的全变差较小，利用这些信息可有效调整滤波力度，对纹理进行平滑并保留弱边缘．

滤除图像纹理的难点是保留弱边缘．由图4可见， LL算法可以保留弱边缘，但不能有效滤除纹理；RTV算法的边缘保持效果不佳，绝大部分的弱边缘都被平滑了；TF算法的边缘保持能力较好，但保持效果不如IP-l0算法．对比l0梯度最小化算法和IP-l0算法的平滑结果发现， IP-l0算法可以更好地保留弱边缘，这是因为在弱边缘区域，相对全变差的数值较小，令算法的平滑力度较小，可有效保留弱边缘．

综上可知，IP-l0算法的图像平滑结果优于l0梯度最小化算法，且在平滑强纹理和保留弱边缘的任务中，IP-l0算法是5种算法中综合能力最优的．在对强纹理图像的平滑实验中，在保留边缘的情况下，IP-l0算法的纹理平滑结果最优．在对弱边缘图像的平滑实验中，在滤除纹理的情况下，IP-l0算法的边缘保持能力最优．

2.2 IP-l0算法与基于深度神经网络平滑算法的比较

近年来，随着深度学习算法的广泛研究，有学者提出基于审读神经网络的图像平滑算法．例如，文献[21]采用传统平滑算法的结果作为监督图像训练深度神经网络，文献[22]采用无监督的方式训练深度神经网络．以下将对比分析采用基于深度学习算法和IP-l0算法进行图像平滑的效果．为便于对比，原图采用文献[22]的一幅大小为500×374 像素和一幅大小为600×444像素的图像．IP-l0算法的参数设置不变，文献[22]算法参数见原文．3种算法的图像平滑结果见图5和图6．

图5 含有弱边缘图像的平滑结果比较Fig.5 The comparison of different image smoothing results for an image with weak edges

图6 含有强纹理图像的平滑结果比较Fig.6 The comparison of different image smoothing results for an image with strong textures

由图5和图6可见，文献[22]算法和IP-l0算法都可得到较好的图像平滑结果，但后者的边缘保持与纹理平滑结果更佳．观察图5(b)和(c)发现，两算法的纹理平滑效果类似，但文献[22] 算法的边缘保持能力仍然不如IP-l0算法．图5(b)中右上部分门上的边缘和左下部分地板上的边缘都被平滑了，而在图5(c)中这些边缘是清晰的．在图6(d)中，女人头巾附近的纹理被更好地平滑了，女人脸部左侧与背景之间的边缘被更好的保留了．但是，在图6(c)中，采用文献[22]算法进行图像平滑处理后，头巾部分存在更多未被平滑的纹理，而脸部左侧和背景之间的边缘相对模糊．在图5(a)中的瓷砖周围有一些弱边缘，文献[21]的滤波结果中没有保留这些边缘，滤波效果显然不如IP-l0算法．图6(a)中头巾附近存在较强的纹理，文献[21]算法没有很好的滤除这些强纹理．这是因为文献[21]算法需要采用传统平滑算法的结果作为监督图像训练神经网络，得到的网络具有对应平滑算法的性能．当采用l0梯度最小化算法的结果作为监督信号时，文献[21]算法的滤波输出具有类似l0梯度最小化算法的性能和缺点，不能很好的保留弱边缘信息，对于强纹理的滤除也显得不足．

结 语

提出基于图像块的l0梯度最小化边缘保持图像平滑算法，将l0梯度最小化应用于图像块，权重参数依据图像块的特征正比于图像块的相对全变差，从而令平滑力度也根据图像块的特征实现自适应调整．实验结果表明，与l0梯度最小化以及其他一些最新的边缘保持图像平滑算法相比，IP-l0算法有效结合了基于全局信息和基于局部信息两种算法的优势，平滑结果最优．

未来仍需深入研究基于全局信息的算法和基于局部信息的算法之间的关系．此外，本研究使用了相对全变差的概念，能较好地区分纹理与边缘．下一步可以探索其他能更好地区分纹理与边缘的度量准则，例如采用监督学习的方式，以期得到更好的图像平滑结果．