基于重加权l1范数的边缘保持图像平滑算法

宋昱 孙文赟

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

边缘保持图像平滑的目的是在保持图像中主要边缘的前提下对其中的纹理和细节进行平滑。边缘保持图像平滑是很多计算机视觉与图形学应用的基础和关键步骤。自然图像中通常包含很多纹理和细节,可能对很多计算机视觉和图形学算法(如边缘检测[1]、图像分割[2]、图像分类、 基于内容的图像编辑等)造成困难。使用边缘保持图像平滑算法对图像进行预处理通常是有利的,可以降低边缘检测、图像分割、图像分类与基于内容的图像编辑等任务的难度;另外,边缘保持图像平滑也可以作为图像增强算法的一个中间步骤[3],在平滑了原始图像后,通过从原始图像中减去平滑图像得到一幅细节图像,然后对细节图像进行放大,再加到平滑图像上,可以得到一幅增强后的图像。

现有的边缘保持图像平滑算法可以大致分为2类。第1类算法是基于像素点的局部邻域统计特征对像素点进行滤波,这些算法被称为基于局部统计特征的算法；基于局部统计特征的算法根据中心像素点和该像素点邻域中其他像素点的关系对中心像素点进行处理。双边滤波器是一种基于局部统计特征的算法[4- 5],该算法通过对邻域中的像素加权平均得到滤波后的像素,像素灰度距离越大或者像素点空间距离越大,权重越小。引导滤波采用和双边滤波器类似的方法对图像滤波[6- 7]。局部拉普拉斯滤波器是一种基于局部统计特征的算法[8],该算法是在对中心像素进行点操作之后,构建该像素点的局部拉普拉斯金字塔,通过组合局部拉普拉斯金字塔得到像素点的滤波结果。文献[9]中的算法是根据区域协方差平滑像素点,该算法和非局部均值图像去噪算法类似[10]。基于偏微分方程(PDE,Partial Differential Equation)的算法(例如各向异性扩散或者一致增强扩散[11- 13]),也是基于局部统计特征的算法,这类算法根据像素灰度距离以迭代的方式对像素点进行滤波,那些与邻域中像素灰度更接近的像素被以较大力度进行平滑,那些与邻域中像素灰度差异较大的像素被以较小的力度进行平滑。基于边缘回避小波的平滑算法也是基于局部统计特征的算法[14],该算法首先对原始图像进行边缘回避小波变换,然后对各个细节子带进行抑制,最后通过使用边缘回避小波逆变换得到平滑图像。最近提出的树滤波算法使用中心像素点与邻域中像素点的树距离和双边距离对图像进行平滑[15],通过使用该算法,可以有效地滤除图像中的强纹理。基于局部统计特征的算法有一个主要的优势——因为算法中涉及的主要计算都限制在各个像素点的局部邻域中,基于局部统计特征的算法的计算复杂度较低。而这类算法的一个主要缺点是每一像素点的滤波结果仅依赖该像素点局部邻域的统计特征,在滤波过程中没有充分考虑整幅图像提供的信息。

第2类算法是对整幅图像同时进行滤波,这类算法被称为基于全局统计特征的算法,通常需要构建一个目标函数,该目标函数以平滑图像与原始图像作为变量,通过最小化目标函数得到平滑后的图像；目标函数通常包含2种类型的能量项,其中一项是数据忠诚项,该项用以衡量平滑图像与原始图像之间的差值,另外一项是正则项,该项用以描述平滑后图像的特征。该类算法中最有代表性的算法是加权最小二乘(WLS,Weighted Least Square)算法[16],该算法在目标函数中有2项能量项,数据忠诚项用以衡量平滑图像与原始图像差值的l2范数的平方,正则项用以衡量平滑图像的加权梯度,可以通过求解一个线性方程组得到平滑图像；算法的平滑结果具有如下特征,那些高对比度的边缘被保留了,那些低对比度的纹理被平滑了。为了进一步提升原始WLS滤波器的运行速度,文献[17] 提出了一种快速加权最小二乘算法。文献[18]的算法基于局部极值对图像进行平滑,该算法首先基于局部邻域定位极值,然后根据全局统计特征构建最小和最大包络,最后通过对最小和最大包络求均值得到平滑图像。l0梯度最小化也属于这类算法[19],该算法的目标函数中有2项能量项,数据忠诚项用以衡量平滑图像与原始图像之差的l2范数的平方,正则项用以衡量平滑图像的l0梯度,通过迭代求解2个子最小化问题可以得到平滑图像。最近,研究者又提出了一种基于l1范数的图像平滑算法[20],该算法的目标函数中有3项能量项,数据忠诚项用以衡量平滑图像与原始图像之差的l2范数的平方,第一项正则项用以衡量平滑图像中彼此在对方邻域中的像素对差值的l1范数,第二项正则项用以衡量平滑图像中代表性像素对差值的l1范数,这些代表性像素从一个超像素区域中选取得到；该算法使用了稀疏表示中的概念进行边缘保持图像平滑[21- 23],算法目标函数是一个凸函数,通过最小化目标函数直至收敛得到平滑图像。基于全局统计特征的算法与基于局部统计特征的算法的特征是相反的；因为在最小化过程中,需要求解一个大型稀疏线性方程组,所以这类算法的计算复杂度较高；但是因为在滤波过程中考虑了整幅图像的全局信息,因此基于全局统计特征算法的平滑效果较好。除前述的相关算法外,最近研究者还提出了一些基于学习的边缘保持图像平滑算法[24- 28],这些算法利用深度学习方法有监督或者无监督的训练一个神经网络,利用神经网络平滑图像。

因为在算法的目标函数中考虑了中心像素点和邻域中其他像素点差值的l1范数,基于l1范数的算法可以很好地平滑图像;WLS算法和l0梯度最小化算法只考虑了中心像素点水平与垂直方向的梯度[16,19],与WLS算法和l0梯度最小化算法相比,基于l1范数的算法考虑的邻域信息更加完整。事实上,衡量差值的最优的范数是l0范数,但是,使用l0范数的目标函数是非凸的,最小化该目标函数非常困难甚至是不可能的。l1范数是l0范数的一个凸近似,基于l1范数的算法的滤波结果仍然包含很多纹理；受Candes等[29]提出的重加权l1范数的启发,本研究使用重加权l1范数替代l1范数,提出一种基于重加权l1范数的图像平滑算法,并通过将该算法与其他很多最新的边缘保持图像平滑算法进行比较,来验证该方法的性能。

1 基于重加权l1范数的图像平滑算法

1.1 原始的基于l1范数的图像平滑算法

基于l1范数的图像平滑算法的目标函数是3项能量项的加权之和[20],可表示为

E=αE1+βEg+θEa

(1)

其中,E1、Eg和Ea分别表示局部稀疏正则项、全局稀疏正则项和数据忠诚项,α、β和θ是这3项的权重参数。

(1)局部稀疏正则项

局部稀疏正则项可以表示为

(2)

其中：xi、xj为平滑图像中像素点i、j的RGB色彩空间中的3维色彩向量；Nh(i)是一个局部邻域窗口,窗口的中心点像素是i,窗口大小为h×h；wij表示像素点i和像素点j的相似性,可以定义为

(3)

其中：fi=[κ·liaibi]T,fj=[κ·ljajbj]T,fi、fj是像素点i处的Lab颜色空间中加权后的颜色向量;κ和σ是常数。

原始图像像素点i的CIELab色彩空间中的3维色彩向量可以表示为[liaibi]T。每一通道都归一化至[0,1]之间。当κ<1时,该能量项对光照变化的敏感性会降低。给出这一设置的目的是为了抑制光照变化。在该能量项中计算了平滑图像中彼此在对方邻域中的像素对的差值的l1范数。因为权重函数是原始图像中像素对差值的单调递减函数,那些在原始图像中有着较小颜色差异的像素对被赋以较大的权重,这会导致平滑图像中相应像素点l1差值的减小,从而使平滑中有着较大颜色差值的像素对的个数会比较稀疏。

为了给出能够最小化目标函数的算法,必须使用矩阵-向量的方法对式(2)中的能量项进行改写。假设图像中共有n个像素,式(2)中共有ml个相邻像素对。平滑图像所有像素的R通道值可以表示为向量zr；类似的,平滑图像所有像素的G通道值和B通道值可以分别表示为zg和zb。令zr、zg和zb是一个矩阵z的第1、第2和第3列,也即,z=[zrzgzb]是一个n×3的矩阵。令L=Lij是一个mi×n的矩阵。如果像素点i和像素点j构成了第k对像素对,那么Lki=wij、Lkj=-wij。通过使用矩阵-向量乘法,可以将式(2)中的能量项重新写为

E1=‖Lz‖1

(4)

(2)全局稀疏正则项

因为一幅图像中通常仅有有限数量的不同光照值,所以可以利用超像素计算全局稀疏正则项。通过使用文献[30]中的算法产生一个固定数量(这一数量表示为ns)的超像素。从超像素区域中选出代表性的像素点,选取的准则是该像素点的颜色值与超像素区域的颜色均值最为接近。这些代表性像素点构成一个集合Sr,该集合是原始图像像素点集合的一个子集。全局稀疏正则项可以定义为

(5)

在该能量项中考虑了所有代表性像素对之间的相互作用。和式(2)类似,原始图像中具有较小颜色差值的像素对被赋以较大的权重,这样会更容易使这些像素点在平滑图像中具有相同的颜色,从而能够提升平滑图像中像素点颜色的全局稀疏性。

所有代表性像素对的总数可以表示为mg(mg=ns(ns-1))。ns表示代表性像素点的个数。G=Gij是一个mg×n的矩阵(此处的n和前文的n是相同的变量,都表示图像中像素点的个数),如果像素点i和像素点j构成了第k对像素对,那么Gki=wij,Gkj=-wij。再次使用矩阵-向量乘法,式(5)中的能量项可以重新写为

Eg=‖Gz‖1

(6)

(3)数据忠诚项

如果仅考虑上述能量项,那么显然存在一个平凡解,即所有的像素都取相同的颜色值。为了避免得到上述平凡解,需要增加一个限制条件,即平滑图像与原始图像的差异应该尽可能的小。数据忠诚项可以表示为

(7)

用矩阵-向量形式写成的目标函数表示为

(8)

可以采用交替方向乘子方法(ADMM,Alternating Direction Method of Multipliers)或者分裂Bregman方法最小化目标函数[31- 32]。

1.2 利用重加权l1范数改进原始的基于l1范数的图像平滑算法

1.2.1 重加权l1范数和重加权l1范数最小化

l0范数最小化问题(记为问题P0)可以写成

(9)

s.t.y=Φx

其中:‖x‖0=|i:xi≠0|;Φ是一个m×n的矩阵,其行数小于列数,也即m

问题P0是非凸的,从而很难求解。为了求解问题P0,通常的做法是凸化该问题,并且求解该凸问题以近似原问题P0的解。如果使用l1范数而不是l0范数,那么问题就变为一个凸问题。

对应于问题P0的l1范数最小化问题(记为问题P1)可以表示为

(10)

s.t.y=Φx

尽管l1范数是l0范数最紧的凸近似,问题P1的解中仍然包含很多非零元素。为了使解更加稀疏,可以考虑问题P1的一个重加权版本。考虑“加权”l1范数最小化问题(记为问题WP1)可以表示为

(11)

s.t.y=Φx

其中,w1,w2,…,wn是正的权重。在下面,可以方便的将目标函数表示为‖Wx‖1,其中W是一个对角矩阵,其对角线上的元素是w1,…,wn,其他元素是0。在文献[29]中,Candes等指出,如果权重与x0中的数值成反比,那么问题WP1和问题P0就是等价的,其中x0是问题P0的解。因为事先不可能得到x0,那么可以将权重设置为反比于解x0的一个近似。Candes等[29]提出一个简单的迭代算法,其在估计x0和重新计算权重之间交替进行,算法如下。

步骤2 求解如下的重加权l1范数最小化问题

xl=argmin‖Wlx‖1

(12)

s.t.y=Φx

步骤3 对于i=1,…,n,按照下式更新权重

(13)

步骤4 达到收敛条件或者l超过一个最大的迭代次数lmax时,停止算法;否则,对l增加1,并且回到步骤2。

在步骤3中,有一个参数ε>0,引入ε的原因是为了保证算法的稳定性,即保证在xl中如果有一个零值存在,那么在下次迭代中,可以在对应位置存在非零值。就像文献[29]中的实验给出的那样,求解重加权l1范数最小化问题,得到的解比求解问题P1得到的解更加稀疏。

1.2.2 将重加权l1范数最小化应用于基于l1范数的图像平滑算法

重新查看式(8)中的目标函数,可以发现其本质上是对基于l0范数的图像平滑算法的一种近似。其对应的l0范数图像平滑算法的目标函数可以写为

(14)

原来的l1范数被替换为l0范数。如果平滑图像是通过求解目标函数式(14)得到的,而不是通过求解目标函数式(8)得到的,那么可以预见的是在平滑图像中应该有更少的具有不同颜色值的像素对。换句话说,式(14)中的目标函数的解应比式(8)中目标函数的解更加“稀疏”。因为式(14)中的目标函数是非凸的,很难求解,所以可以使用文献[29]中提出的重加权l1范数最小化来求解。文中,使用重加权l1范数最小化来求解式(14)中目标函数的近似解。式(14)中的目标函数可以重新写为

(15)

每一个单品拿出来，都营造出以女孩为主的消费者争先追逐的热潮。用上这款面膜，前男友看了要后悔；涂上这些口红，你就是生活的女王。

(16)

变量z的解与前面的解是一样的。变量u1的解可以写为

(17)

这是一个重加权l1范数最小化问题,变量u1有闭式解,可以写为

(18)

加权软阈值算子Sεw(x)的一般定义为

(19)

(20)

在求解完变量u1和u2后,根据下式更新权重

(21)

(22)

权重和当前的解(u1与u2)的值成反比,下一次迭代的权重是通过当前迭代时解的值计算得到的。用于求解目标函数(式(15))的重加权l1范数最小化算法如算法1所示。

算法1:用于最小化目标函数(式(15))的优化算法

输入：原始图像zin,参数α、β和θ

while not converged do

1.固定其他变量,通过下式更新z

zk+1=(2θI+μkLTL+μkGTG)-1·

(2θzin+μkLT(u1,k+Y1,k/μk)+

μkGT(u2,k+Y2,k/μk))

2.固定其他变量,通过下式更新u1

3.固定其他变量,通过下式更新u2

4.对于i=1,…,3,ml,j=1,…,3，通过下式更新权重

5.对于i=1,…,3,mg,j=1,…,3通过下式更新权重

6.更新拉格朗日乘子

Y1,l+1=Y1,l+μl(u1,l-Lzl)

Y2,l+1=Y2,l+μl(u2,l-Gzl)

7.更新惩罚因子μ,μl+1=min(ρμl,μmax)

8.检查是否满足收敛条件

‖u1,l-Lzl‖∞<εADMM

‖u2,l-Lzl‖∞<εADMM

9.l=l+1

end while

输出：平滑图像z=zl+1

比较算法1和文献[20]中的求解算法,可以看出有两处不同。第1处不同是变量u1和u2的解不同；在文献[20]的算法中,使用了普通的软阈值算子；与之相对的,在算法1中,使用了一个加权软阈值算子。第2处不同是在算法1中,权重是迭代更新的,但是在文献[20]的算法中,没有权重。可以预见到,通过使用重加权l1范数最小化,解可以更加稀疏,也就意味着在平滑图像中具有不同颜色值的像素对的数量更少。

2 实验结果与分析

为了衡量所提算法的性能,将所提算法和很多最新的边缘保持图像平滑算法进行了比较,用于比较的算法包括原始的基于l1范数的图像平滑算法(l1)[20]、l0梯度最小化图像平滑算法(l0)[19]、加权最小二乘图像平滑算法(WLS)[16]、快速加权最小二乘图像平滑算法(Fast WLS)[17]、使用扩散图的加权最小二乘图像平滑算法(Diffusion WLS)[28]、局部拉普拉斯滤波器(LL)[8]、基于边缘回避小波的图像平滑算法(EAW)[14]、基于局部极值的图像平滑算法(LE)[18]、基于树滤波的图像平滑算法(TF)[15]。上述算法中,l1、l0、Fast WLS和Diffusion WLS是基于全局统计特征的图像平滑算法,LL、EAW、 LE和TF是基于局部统计特征的图像平滑算法。在大量图像上进行了实验,限于篇幅,文中选取了其中6幅有代表性的图像作为说明。这些图像均取自BSDS300数据集[34],分辨率分别是341×512、560×368、321×481、481×321、321×481以及321×481。实验环境是Intel(R)Core(TM)i5- 4460S CPU @2.9 GHz,8 GB RAM,Matlab R2015a。在给出实验结果之前,先给出了所提算法和比较算法的相关参数,如表1所示。

表1 不同算法的参数

图像平滑结果如图1-图12所示。为了更加清楚地比较不同图像平滑算法的性能,从每幅图像中选取了一个代表性的区域(图像中白色矩形框所示区域),将该区域进行放大并示于相应的图像平滑结果之后。

图1 图像1及不同边缘保持图像平滑算法的图像平滑结果

图2 从图1各图中选取的代表性区域的放大图

图3 图像2及不同边缘保持图像平滑算法的图像平滑结果

图4 从图3各图中选取的代表性区域的放大图

图5 图像3及不同边缘保持图像平滑算法的图像平滑结果

图6 从图5各图中选取的代表性区域的放大图

图7 图像4及不同边缘保持图像平滑算法的图像平滑结果

图8 从图7各图中选取的代表性区域的放大图

图9 图像5及不同边缘保持图像平滑算法的图像平滑结果

图10 从图9各图中选取的代表性区域的放大图

图11 图像6及不同边缘保持图像平滑算法的图像平滑结果

图12 从图11各图中选取的代表性区域的放大图

因为没有关于“真实平滑图像”的明确定义,所以客观评价指标(峰值信噪比(PSNR)、结构相似度(SSIM)等)不能使用。文献[35]中给出了一种得到“真实平滑图像”的方法,即从现有的多种图像平滑结果中主观选出一幅最满意的图像作为“真实平滑图像”,并以此为标准来定义“真实平滑图像”,该方法具有一定的局限;由于现有方法不一定能取得较好的结果,如果从现有方法的平滑结果中选择其中一种平滑结果作为真实图像,则该图像不能作为“真实平滑图像”，而只是“主观的真实平滑图像”。下面主要采用主观评价方法对算法平滑结果进行比较。边缘保持图像平滑算法应该具有如下的特征：首先,在平滑图像中,主要的边缘应该被保留；第二,在平滑图像中,低对比度和高对比度的纹理应该被平滑。可以通过上述两点判断平滑图像的质量。通过观察实验结果发现,算法WLS、 Fast WLS、Diffusion WLS的滤波效果并不好,这些算法通过最小化一个目标函数得到平滑图像,目标函数中包括一项数据忠诚项和一项正则项,倾向于平滑那些有着较小图像梯度的像素点,也会平滑低对比度的边缘,且不能有效平滑高对比度的纹理(可以很清楚的从图5与图6(c)-6(e)中看出);算法LL的图像平滑效果也不好,该算法不能很好地滤除图像中高对比度的纹理；算法LL使用一个全局阈值区分高对比度的边缘和低对比度的纹理,因为算法会将高对比度的纹理判断为边缘,故该算法不能平滑高对比度的纹理；算法EAW的图像平滑结果也不佳,在平滑图像中有很多高对比度和低对比度的纹理,该算法可以保留主要的边缘,但是算法处理纹理的能力不强；算法LE的图像平滑结果同样不好,该算法依赖局部极小值和极大值估计一个最小和最大包络,平滑图像是2个包络的平均,当图像中的纹理非常复杂时,采用平均的方法不能有效地滤除这些纹理；与前述算法相比,算法l0和TF的图像平滑效果较好,这2个算法可以平滑一些高对比度的纹理并保留主要的边缘,但是,当图像中的纹理具有很高的对比度时,算法l0的滤波效果就会变得不好(这可以从图9(b)和图10(b)中看出来)；算法l1的图像平滑能力一般(从各图的图(j)中可以看出),尽管在平滑图像中保留了主要的边缘,但是平滑图像中却存在很多纹理,因为算法l1使用l1范数描述平滑图像中像素对颜色差值且基于l1范数的算法并不能保证算法解的稀疏性,从而在平滑图像中有很多具有不同颜色值的像素对；文中所提算法的图像平滑结果较好,该算法可以有效地滤除纹理并保留主要边缘,因为该算法中使用了重加权l1范数,得到的解更加“稀疏”,在平滑图像中有着更少的具有不同颜色值的像素对,证明使用重加权l1范数可以明显改善算法l1的图像平滑结果。

采用5分制对10种平滑算法在6幅实验图像上的平滑结果进行了评分。5分制中分数和平滑效果的对应关系为：5,很好；4,好；3,一般；2,差；1,很差。评分结果如表2所示,同时在表2中(第2列)给出了图像的主要特征。从实验结果可以看出,对于含有强纹理和复杂纹理的图像,文中所提算法可以取得最优的结果,而对于含有弱边缘的图像,文中所提算法会滤除其中的一些弱边缘。

表2 10种图像平滑算法的得分

4 结论

文中提出了一种新的基于重加权l1范数的边缘保持图像平滑算法。原始的基于l1范数的图像平滑算法使用l1范数衡量平滑后图像中局部邻域内像素点差值的稀疏性,其平滑结果中仍然存在较多的纹理。为了减少平滑后图像中的纹理,提出的算法将原始的基于l1范数的图像平滑算法与重加权l1范数最小化进行了结合,通过使用重加权方法使得解更加稀疏,即使得平滑后图像中局部邻域内像素点差值更小。实验结果表明,和原始的基于l1范数的图像平滑算法以及其他的图像平滑算法相比,该算法可以更加有效地平滑图像。所提算法的平滑图像中纹理残留较少。