CN111553863A - 一种基于非凸全变分型正则化的图像增强方法 - Google Patents

Abstract

本发明提出了一种基于非凸全变分型正则化的图像增强方法，用以解决现有图像增强方法很少关注噪声存在，未考虑图像的稀疏性，无法得到高质量复原结果的问题。本发明采用基于指数变换的数据拟合项和用其中一个正则项来描述噪声先验信息，然后利用非凸全变分

(拟)范数的正则项来惩罚分段常数的稀疏性。根据图像结构，指数p在（0,1）中选取。最后基于Retinex理论，计算光照部分，并去除光照得到反射部分，即图像的本质特征，以达到偏移场校正的目的。本发明不但考虑了噪声对图像的干扰，而且通过正则项惩罚了图像的稀疏性，同时，用指数变换代替对数变换，提高了对比度，从而使得模型更具有鲁棒性。

Description

一种基于非凸全变分型正则化的图像增强方法

技术领域

本发明涉及图像处理的技术领域，尤其涉及一种基于非凸全变分型(TV)正则化的图像增强方法，主要涉及图像的偏移场校正，可用于医学图像的光照不均匀和彩色图像的阴影消除等。

背景技术

给定一幅图像，在不同的光照条件下，人类视觉***趋向于看到相同的颜色，这一现象被称为颜色恒常性。换句话说，它确保了在光照变化的情况下，物体的颜色和灰度等信息保持不变。最典型的一个例子是Adelson棋盘阴影错觉，如图2(a)所示，从视觉上看区域A比区域B更暗。然而，如果将该图像直接以绘图软件打开，用取色功能来确认，就会发现其实这两个区域的颜色完全一样。也就是说，在数字图像中，这两个区域的灰度值是相同的，如图2(b)所示。此外，在图2(c)中用同样颜色的、贯穿区域A和区域B的长条来对比，可以发现区域A和区域B的灰度是相同的。这一现象是由不同的光照条件所造成的，称之为Retinex效应。事实上，这种现象也会经常出现在医学图像中，即灰度不均匀。如在磁共振成像过程中，由于射频线圈产生不均匀的磁场等因素，造成磁共振图像具有灰度不均匀的特点，在彩色图像中的彩色阴影效应等。

Retinex理论认为，图像是由光照部分和反射部分的乘积组成，光照部分是由光源引起的对图像的干扰，反射部分代表的是图像的本质信息。Retinex理论通过计算光照部分，并将光照部分从原图中去除来得到图像的反射部分，即图像的本质特征。

关于图像增强问题，大多数现有的变分模型主要还是将原始图像进行对数变换，从而将乘积形式转化为加法形式。这种变换可以减少计算的复杂度，同时模拟人眼对光照强度的感知。但是，现存的方法很少关注噪声存在的情况，当观测到的图像被噪声干扰时，这些方法就会缺少鲁棒性，无法得到高质量的复原结果。此外，现存的方法也很少基于

(拟)范数来描述图像的稀疏性。

发明内容

针对现有图像增强方法很少关注噪声存在，未考虑图像的稀疏性，无法得到高质量复原结果的技术问题，本发明提出一种基于非凸全变分型正则化的图像增强方法，考虑了噪声对图像的干扰，通过正则项惩罚图像的稀疏性，提高图像的复原结果。

为了达到上述目的，本发明的技术方案是这样实现的：一种基于非凸全变分型正则化的图像增强方法，其步骤如下：

步骤一：输入包含噪声和灰度不均匀的待复原的原始图像，对原始图像进行对数变换；

步骤二：计算对数变换后图像的梯度映射以消除噪声；

步骤三：利用Retinex理论估计经过步骤二处理的复原图像的光照部分并去除光照部分；

步骤四：对步骤三中处理后的图像进行指数变换，获得最终的增强图像。

所述步骤一中原始图像包含噪声和偏移场的数学模型描述为：

S(x,y)＝R(x,y)×L(x,y)+n(x,y)；

其中，S(x,y)是从数字图像设备获得的原始图像在坐标位置(x,y)处的像素值；L(x,y)是原始图像的光照部分在位置(x,y)处的像素值；R(x,y)是原始图像的反射部分在位置(x,y)处的像素值，n(x,y)表示高斯白噪声；

运用对数变换将乘法运算变成加法运算：

v＝r+l；

其中，反射对数部分r＝log(R)，光照对数部分l＝log(L)和光照反射对数部分v＝log(R×L)；R为像素值R(x,y)对应的反射部分，L为像素值L(x,y)对应的光照部分。

所述步骤二中计算对数变换后图像的梯度映射的方法为利用数据拟合项

衡量指数变换e^v与原始图像S的拟合程度，利用正则项

描述了噪声的先验信息，其中，

对应指数变换e^v的梯度。

所述步骤三中利用高阶全变分

拟范数的正则项来估计光照部分，其具体形式为：

其中，

表示Hessian算子，F表示Frobenius范数以及q∈(0,1)，l表示光照部分L的对数变换。

所述所述步骤三中基于Retinex理论提出的模型为：

其中，α，β，γ是正的正则化参数，μ是罚参数；

是数据拟合项，

是正则项，

是惩罚项，

是适定项，参数θ取很小的值；S是原始图像S(x,y)的简写，v为光照反射对数部分，r为反射对数部分，l为光照对数部分，

表示Hessian算子，指数p∈(0,1)。

所述基于Retinex理论提出的模型是多变量优化问题，变量v,r,l是耦合在一起的，利用交替极小化算法进行求解，即将同时求解多个变量的问题转化为依次求解单个变量的过程，在求解其中一个变量的过程中固定其它两个变量。

所述交替极小化算法的求解步骤为：

(1)初值选取：根据图像特征选取合适的参数α,β,γ,μ,θ；并令v⁰＝r⁰＝log(S)；

(2)迭代：

(3)满足终止条件迭代终止：

max(R(rⁿ,lⁿ))≤10^-5，

其中，R(rⁿ,lⁿ)表示相对误差，具体形式为：

其中，rⁿ、lⁿ分别表示对数部分r,l第n次迭代的结果，vⁿ⁺¹、rⁿ⁺¹和lⁿ⁺¹分别是关于对数部分v,r,l的最优解，p,q∈(0,1)。

所述步骤(2)中的每个子问题vⁿ⁺¹、rⁿ⁺¹、lⁿ⁺¹利用交替方向乘子法和迭代加权

最小化算法进行求解。

所述采用交替方向乘子法求解子问题

的方法为：引入辅助变量u＝e^v和

其增广拉格朗日函数为：

其中，Λ₁，Λ₂是拉格朗日乘子，r₁,r₂是惩罚参数。

对于子问题

引入辅助变量

基于交替方向乘子法的迭代形式为：

利用迭代加权

最小化算法求解子问题w的形式为：

η^k+1＝ρη^k；

其中，T(w^k,η^k,p)＝p(||w^k||_2,1+η^k)^p-1，ρ和p∈(0，1)；r₃是惩罚参数，Λ₃是拉格朗日乘子，shrink()表示软阈值算子，η的作用是为了防止分母为零，ρ表示步长；

对复原结果r进行指数变换R＝e^r，获得最终的增强图像R。

与现有技术相比，本发明的有益效果：

1、本发明基于Retinex理论，数据拟合项用指数变换代替已有技术的对数变换，然后再做指数变换，通过指数变换，可以解决对数变换带来的对比度降低的问题，提高图像复原的质量。

2、本发明采用基于非凸全变分

(拟)范数的正则项来惩罚分段常数，以及基于非凸高阶全变分

(拟)范数来惩罚空间光滑区域，在梯度意义下能有效地保持本原图像的稀疏性和偏移场或反射率图像的光滑性。

3、在本发明中采用交替极小化算法，将原问题转变成一些容易求解的子问题。本发明采用交替方向乘子法(ADMM)来求解子问题，即通过引入辅助变量，将原问题转变成更小的子问题，使得每个子问题都更容易求解。

4、本发明考虑同时包含噪声和灰度不均匀的图像，对于不同水平的噪声，都可以得到高质量的复原结果，展示了本发明的鲁棒性。

5、本发明可用于不同类型的图片，例如合成图像、仿真图像以及医学MR图像等。本发明能广泛应用于图像分割、图像配准等领域。

附图说明

为了更清楚地说明本发明实施例或现有技术中的技术方案，下面将对实施例或现有技术描述中所需要使用的附图作简单地介绍，显而易见地，下面描述中的附图仅仅是本发明的一些实施例，对于本领域普通技术人员来讲，在不付出创造性劳动的前提下，还可以根据这些附图获得其他的附图。

图1为本发明的流程图。

图2为Adelson棋盘阴影错觉图，其中，(a)为原图像，(b)为(a)的数字图像，(c)为(a)的无错觉图像。

图3为包含噪声和偏移场的MR图像。

图4为本发明和现有Retinex变分方法对图3进行复原的结果对比图，其中，(a)为本发明对图3进行去噪后的复原结果，(b)为本发明利用Retinex理论对图3去噪后的复原结果分解出的光照部分，(c)为本发明利用Retinex理论对图3去噪后的复原结果分解出的反射部分，即最终的复原结果，(d)为TVH1方法对图3的复原结果，(e)为HOTVL1方法对图3的复原结果。

图5为本发明和现有Retinex变分方法对“Brain”图像复原后的结果对比图，其中，(a)为原始Brain图像，(b)为TVH1方法复原的结果图，(c)为HOTVL1方法的复原结果图，(d)为本发明复原结果图，(e)为图(b)的矩形局部区域图，(f)为图(c)的矩形局部区域图，(g)为(d)的矩形局部区域图。

图6为本发明和现有Retinex变分方法对“Head”图像复原后的结果对比图，其中，(a)为Head原始图，(b)为TVH1方法复原的结果图，(c)为HOTVL1方法的复原结果图，(d)为本发明复原结果图。

具体实施方式

下面将结合本发明实施例中的附图，对本发明实施例中的技术方案进行清楚、完整地描述，显然，所描述的实施例仅仅是本发明一部分实施例，而不是全部的实施例。基于本发明中的实施例，本领域普通技术人员在没有付出创造性劳动前提下所获得的所有其他实施例，都属于本发明保护的范围。

如图1所示，一种基于非凸全变分型正则化的图像增强方法，基于图像Retinex理论，估计光照分量并去除光照部分，然后，根据噪声的先验信息，计算图像的梯度映射来消除噪声，最后得到高质量的复原结果，其具体步骤如下：

步骤一：输入包含噪声和灰度不均匀的待复原的原始图像，对原始图像进行对数变换后再进行指数变换。

人类视觉***感知到的颜色信息是由物体表面的反射性质决定的，而与光照无关。但由于光照的影响，会使得获得的图像存在灰度不均匀的现象，即图像存在偏移场。本发明从数字成像设备获得一幅如图3所示的待复原的原始图像S(x,y)，原始图像S(x,y)存在灰度不均匀现象且被噪声干扰。

包含噪声和偏移场的原始图像的数学模型可以描述为：

S(x,y)＝R(x,y)×L(x,y)+n(x,y) (1)

其中，S(x,y)是从数字图像设备获得的原始图像在坐标位置(x,y)处的像素值；L(x,y)是原始图像的光照部分在位置(x,y)处的像素值，即物体接受的光照强度；R(x,y)是原始图像的反射部分在位置(x,y)处的像素值，即物体的反射率，它是物体的本质信息，以及n(x,y)表示高斯白噪声。

由于数学模型(1)中的乘法运算会导致数值的不稳定，所以本发明运用对数变换将乘法运算变成下面的加法运算，即：

v＝r+l；

其中，反射对数部分r＝log(R)，光照对数部分l＝log(L)和光照反射对数部分v＝log(R×L)。R为像素值R(x,y)对应的反射部分，L为像素值L(x,y)对应的光照部分，n表示高斯白噪声n(x,y)的简写，S是像素值S(x,y)对应的原始图像。

然后再对对数变换后的模型做v＝r+l指数变换，于是，数学模型(1)就等价于：

S＝e^v+n； (2)

这里采用指数变换e^v的目的是通过指数变换以解决对数变换带来的对比度降低的问题，提高图像复原的质量。

步骤二：计算指数变换后图像的梯度映射以消除噪声。

数据拟合项

衡量了指数变换e^v与原始图像S的拟合程度，确保去噪后的复原图像与原始图像不要偏差太远。正则项

描述了噪声的先验信息，使得去噪后的图像具有一定的平滑性,其中，

对应e^v的梯度。图3经过步骤二处理得到的复原结果如图4(a)所示，可以直观地发现去噪后的复原图像中不含有噪声，达到了很好的去噪效果。

步骤三：利用Retinex理论估计经过步骤二处理的复原图像的光照部分并去除光照部分。

Retinex理论认为，光照引起的颜色变化是缓慢的，一般表现为平滑的光照梯度。相反地，由反射引起的颜色变化是迅速。通过考虑这两种变化不同的形式，可以将图像信息分为光照信息和反射信息，然后利用分解的光照信息来获得反射信息，既物体的本质信息。

基于Retinex理论，本发明假设光照部分L是空间光滑的，因此利用高阶全变分

(拟)范数的正则项来估计光照部分，其具体形式为：

其中，

表示Hessian算子，F表示Frobenius范数以及q∈(0,1)，l表示光照部分L的对数变换。经过步骤三，利用Retinex理论对图3进行分解，得到光照部分和反射部分，如图4(b)和(c)所示，图4(c)为图像的本质信息。

步骤四：对步骤三中处理后的图像进行指数变换，获得最终的增强图像。

从步骤三中得到复原结果r，经过指数变换R＝e^r，获得最终的增强结果R。

基于Retinex理论，本发明提出下面的模型：

其中，α，β，γ是正的正则化参数，μ是罚参数。

是数据拟合项，

是正则项，

是惩罚项，最后一项

是为了保证模型(3)是适定的，即解的唯一性，参数θ可以取很小的值。

本发明所提模型(3)是多变量优化问题，其中变量v,r,l是耦合在一起的，如果直接求解，则需要求解大规模的问题，从而导致数值求解困难。为了克服这一问题的，本发明交替极小化算法进行求解，即将同时求解多个变量的问题转化为依次求解单个变量的过程。

由于模型(3)是多变量优化问题，在求解其中一个变量的过程中固定其它两个变量，即使用经典的交替极小化算法求解，其具体步骤为：

1.初值选取：根据图像特征选取合适的参数α,β,γ,μ,θ；并令v⁰＝r⁰＝log(S)；

2.迭代：

3.满足终止条件迭代终止：

max(R(rⁿ,lⁿ))≤10^-5，

其中，R(rⁿ,lⁿ)表示相对误差，具体形式为：

rⁿ、lⁿ分别表示对数部分r,l第n次迭代的结果，vⁿ⁺¹、rⁿ⁺¹和lⁿ⁺¹分别表示关于对数部分变量v,r,l的最优解，p,q∈(0,1)。

针对每个子问题vⁿ⁺¹、rⁿ⁺¹、lⁿ⁺¹，利用交替方向乘子法(ADMM)和迭代加权

最小化算法进行求解。

由于子问题v中含有非线性函数e^v,

范数的存在，使得这个子问题是非凸和非光滑的，因此直接求解非常困难。本发明采用交替方向乘子法(ADMM)来求解该子问题，即通过引入辅助变量，将原问题转变成更小的子问题，使得每个子问题都更容易求解。即引入辅助变量u＝e^v和

然后写出增广拉格朗日函数，最后利用ADMM进行求解，具体的形式为：

其中，Λ₁，Λ₂是拉格朗日乘子，r₁,r₂是惩罚参数。

通过步骤二去除噪声，获得初始复原结果v，然后再利用Retinex理论将初始复原结果v分解成光照对数部分l和反射对数部分r。

对于子问题r，由于

(拟)范数的存在，使得该问题是非凸、非光滑和非Lipschitz的。为了有效地求解该问题，本发明采用迭代加权

最小化算法来求解这两个子问题。首先引入辅助变量

基于ADMM，可以写成下面的迭代形式：

利用迭代加权

最小化算法，关于w子问题的具体形式为：

η^k+1＝ρη^k；

其中，T(w^k,η^k,p)＝p(||w^k||_2,1+η^k)^p-1，ρ和p∈(0，1)。

本发明的实验平台为：因特尔酷睿处理器5500U，[email protected]赫兹,8G内存个人笔记本,使用Matlab2018进行程序设计。

具体实验一：

通过本发明和现有Retinex变分方法对图3进行复原并对比，复原结果如图4，其中，图4(a)为本发明对图3进行去噪后的复原结果，图4(b)为本发明利用Retinex理论对图3去噪后的复原结果分解出的光照部分，图4(c)为本发明利用Retinex理论对图3去噪后的复原结果分解出的反射部分，即图像的本质信息，本发明最终的复原结果，图4(d)为TVH1方法对图3的复原结果，图4(e)为HOTVL1方法对图3的复原结果。

由图3可以直观地发现，图像中间位置较左右两边偏亮，但实际上它们的亮度应该是一样的。从图4(c)可以直接看出，本发明方法得到的复原结果的中间位置和左右两边的亮度一样，有效地校正了偏移场，同时，复原结果中不含有噪声，这说明了本发明拥有更高质量的复原结果。从图4(d)、图4(e)展示的复原结果，TVH1方法的复原结果中间位置还是偏亮，没有很好地校正偏移场。HOTVL1方法虽然校正了偏移场，但是降低了原始图像的对比度，此外，复原图像中还含有噪声。

具体实验二：

分别用本发明和现有Retinex变分方法对一幅仿真图像和一幅真实的医学MR图像进行复原并对比，结果如图5和图6，其中：

图5是用现有Retinex变分方法和本发明方法对“Brain”图像复原后的结果对比图，其中：图5(a)是“Brain”图像，图5(b)是用TVH1方法对图5(a)复原后的结果图像，图5(c)是用HOTVL1方法对图5(a)复原后的结果图像，图5(d)是用本发明方法对图5(a)复原后的结果图像。图5(e)为图5(b)的矩形局部区域图，图5(f)为图5(c)的矩形局部区域图，图5(g)为图5(d)的矩形局部区域图。

图6是用现有Retinex变分方法和本发明方法对“Head”图像复原后的结果对比图，其中：图6(a)是“Head”原始图像，图6(b)是用TVH1方法对图6(a)复原后的结果图像，图6(c)是用HOTVL1方法对图6(a)复原后的结果图像，图6(d)是用本发明方法对图6(a)复原后的结果图像。

从图5(b)和图5(c)可以看出，TVH1方法和HOTVL1方法这两种方法虽然能够校正偏移场，但是图像中还含有噪声。相反，通过图5(d)，可以发现本发明不仅可以校正偏移场，而且图像中不含有噪声，去噪效果很好。

为了进一步展示本发明的有效性，放大图4(a)的局部区域对应的三种方法的复原结果，即图4(e)-图4(g)。通过局部细节图，可以很好的展示本发明的方法更具有鲁棒性。

通过观察图6(a)，可以发现由于偏移场和噪声的影响，使得图像的底部较暗，细节信息很难辨认，这些都会影响后续图像处理的有效性，比如图像分割、图像配准等。因此，校正偏移场和去除噪声是非常重要的。从图6(b)、图6(c)和图6(d)可以看出，TVH1方法虽然可以有效的去除偏移场，使图像底部的细节信息能够被清晰辨认，但是背景图会变得偏白，破坏了原始图像的背景，HOTVL1方法同样地可以校正偏移场，但是图像的上部(即脑部)的很多细节都缺失了，变得过度光滑。相反地，本发明方法的复原结果中，不仅复原了图像底部的细节信息，而且保持了脑部的细节，体现了本发明方法的有效性和鲁棒性。

以上所述仅为本发明的较佳实施例而已，并不用以限制本发明，凡在本发明的精神和原则之内，所作的任何修改、等同替换、改进等，均应包含在本发明的保护范围之内。

Claims (10)

1.一种基于非凸全变分型正则化的图像增强方法，其特征在于，其步骤如下：

步骤一：输入包含噪声和灰度不均匀的待复原的原始图像，对原始图像进行对数变换；

步骤二：计算对数变换后图像的梯度映射以消除噪声；

步骤三：利用Retinex理论估计经过步骤二处理的复原图像的光照部分并去除光照部分；

步骤四：对步骤三中处理后的图像进行指数变换，获得最终的增强图像。

2.根据权利要求1所述的基于非凸全变分型正则化的图像增强方法，其特征在于，所述步骤一中原始图像包含噪声和偏移场的数学模型描述为：

S(x,y)＝R(x,y)×L(x,y)+n(x,y)；

其中，S(x,y)是从数字图像设备获得的原始图像在坐标位置(x,y)处的像素值；L(x,y)是原始图像的光照部分在位置(x,y)处的像素值；R(x,y)是原始图像的反射部分在位置(x,y)处的像素值，n(x,y)表示高斯白噪声；

运用对数变换将乘法运算变成加法运算：

v＝r+l；

其中，反射对数部分r＝log(R)，光照对数部分l＝log(L)和光照反射对数部分v＝log(R×L)；R为像素值R(x,y)对应的反射部分，L为像素值L(x,y)对应的光照部分。

3.根据权利要求2所述的基于非凸全变分型正则化的图像增强方法，其特征在于，所述步骤二中计算对数变换后图像的梯度映射的方法为利用数据拟合项

衡量指数变换e^v与原始图像S的拟合程度，利用正则项||▽e^v||_2，1描述了噪声的先验信息，其中，▽e^v对应指数变换e^v的梯度。

4.根据权利要求1-3中任意一项所述的基于非凸全变分型正则化的图像增强方法，其特征在于，所述步骤三中利用高阶全变分l^q拟范数的正则项来估计光照部分，其具体形式为：

其中，▽²表示Hessian算子，F表示Frobenius范数以及q∈(0,1)，l表示光照部分L的对数变换。

5.根据权利要求4所述的基于非凸全变分型正则化的图像增强方法，其特征在于，所述所述步骤三中基于Retinex理论提出的模型为：

其中，α，β，γ是正的正则化参数，μ是罚参数；

是数据拟合项，

是正则项，

是惩罚项，

是适定项，参数θ取很小的值；S是原始图像S(x,y)的简写，v为光照反射对数部分，r为反射对数部分，l为光照对数部分，▽²表示Hessian算子，指数p∈(0,1)。

6.根据权利要求5所述的基于非凸全变分型正则化的图像增强方法，其特征在于，所述基于Retinex理论提出的模型是多变量优化问题，变量v,r,l是耦合在一起的，利用交替极小化算法进行求解，即将同时求解多个变量的问题转化为依次求解单个变量的过程，在求解其中一个变量的过程中固定其它两个变量。

7.根据权利要求6所述的基于非凸全变分型正则化的图像增强方法，其特征在于，所述交替极小化算法的求解步骤为：

(1)初值选取：根据图像特征选取合适的参数α,β,γ,μ,θ；并令v⁰＝r⁰＝log(S)；

(2)迭代：

(3)满足终止条件迭代终止：

max(R(rⁿ,lⁿ))≤10^-5，

其中，R(rⁿ,lⁿ)表示相对误差，具体形式为：

其中，rⁿ、lⁿ分别表示对数部分r,l第n次迭代的结果，vⁿ⁺¹、rⁿ⁺¹和lⁿ⁺¹分别是关于对数部分v,r,l的最优解，p,q∈(0,1)。

8.根据权利要求7所述的基于非凸全变分型正则化的图像增强方法，其特征在于，所述步骤(2)中的每个子问题vⁿ⁺¹、rⁿ⁺¹、lⁿ⁺¹利用交替方向乘子法和迭代加权l¹最小化算法进行求解。

9.根据权利要求8所述的基于非凸全变分型正则化的图像增强方法，其特征在于，所述采用交替方向乘子法求解子问题

的方法为：引入辅助变量u＝e^v和b＝▽u，其增广拉格朗日函数为：

其中，Λ₁，Λ₂是拉格朗日乘子，r₁,r₂是惩罚参数。

10.根据权利要求9所述的基于非凸全变分型正则化的图像增强方法，其特征在于，对于子问题

引入辅助变量w＝▽r，基于交替方向乘子法的迭代形式为：

利用迭代加权l¹最小化算法求解子问题w的形式为：

η^k+1＝ρη^k；

其中，T(w^k,η^k,p)＝p(||w^k||_2,1+η^k)^p-1，ρ和p∈(0,1)；r₃是惩罚参数，Λ₃是拉格朗日乘子，shrink()表示软阈值算子，η的作用是为了防止分母为零，ρ表示步长。

Images