董卫东，彭宏京

(南京工业大学 计算机科学与技术学院，南京 211816)

1 引 言

图像恢复是图像处理和分析领域的重要分支，其目的在于增强在图像获取和传输过程中受到破损的图像质量，模型通常建立在被系统和噪声模糊的图像基础上，寻找合适的数学模型并设计优化算法，求解得到恢复图像，其过程可描述为线性模型：

f=Au+η

(1)

其中f表示观测图像，A为线性算子，例如图像修复中的映射算子或图像去卷积过程中的卷积算子等，η代表均值为0，方差为δ2的高斯白噪声.

在图像恢复领域，总变分(TV)最小化模型应用最早且最广泛[1]，自然图像通常被建模为分片光滑函数，变分模型通过建立合适的能量泛函形式，并转换成对应的偏微分方程设计算法求解，但变分模型在修复过程中扩散系数不能根据图像结构特征的变化而自适应改变，因此会随着求解偏微分方程迭代次数的增加，出现明显的“阶梯效应”.小波及其多分辨分析理论的兴起，极大地推动了图像表示方式的发展，通过选择合适的小波基对图像进行稀疏表示，并设计分裂算法对能量泛函进行求解，不仅降低了图像数据处理的规模，使得算法效率大大提升，而且提高了图像恢复效果[2].自然图像结构特征信息复杂，很难用单一的小波基变换对图像进行完备表述.随着国内外对图像数据处理的重视，最近几年紧框架理论的快速发展，其在稀疏表示信号和描述图像特征方面的优势，使其在图像处理分析领域的应用得到重视[3].多分辨率分析得到的紧框架系统具有冗余性，即图像序列和稀疏表达系数之间不是一对一的映射关系，从而观测图像在紧框架系统下的稀疏表达形式并不唯一，实际应用中可根据具体的图像特征选择合适的表达方式，不仅提高紧框架变换对于噪声和误差的鲁棒性，使图像处理系统更为稳健，而且冗余表示增强了信号逼近的灵活性，提高了对图像复杂结构的稀疏表示能力.紧框架理论中主要包含三种表达方式：基于分析的紧框架表达形式[4]，基于合成的紧框架表达形式[5]和基于分析与合成平衡的表达形式[6].本文主要研究基于分析的紧框架表达形式，其可看作TV模型的推广，将TV模型中的梯度约束项做在紧框架下的离散化表示.同TV能量泛函正则化函数类似，基于分析的紧框架表达模型建立以惩罚图像稀疏性的正则化项，一般以紧框架变换系数的L0范数为约束项[7]，通过求解最小L0范数来精确重构图像，然而求解L0范数问题是NP-hard难题，并且是非凸优化问题，无法建立较为理想的优化算法.Rchan等提出基于L1范数的松弛优化算法[8]，通过建立框架变换系数L1范数的正则化项来近似表达L0范数，在一定条件下，最小化L1范数和最小化L0范数问题是等价的，从而证明其近似稀疏表达的合理性.紧框架域下建立的能量泛函模型，通过近几年发展较快的split bregman方法[9]进行设计求解，其不仅具有较快的收敛速度，而且能够大大改善图像恢复质量.

近年来压缩感知理论在图像处理领域的应用获得快速发展[10]，信号建立变换域系数L1范数的惩罚项，即是对信号在变换域中稀疏性的恢复.信号稀疏表示后，其能量大多集中在少数较大的系数中，而其他大部分系数都为零或接近于零[11].同样地，对于图像经过框架变换后，原始图像信息主要集中在少量较大的框架系数中，而剩余大部分系数都为零或接近于零.建立框架系数L1范数的约束项对图像进行稀疏性恢复时，少量较大的框架系数与大量为零或接近于零的系数对于稀疏性恢复的影响程度是不一样的，E.J.Candes等人证明了信号重构过程中，稀疏表示后不同的变换系数添加不同程度的惩罚因子可大大增强信号的稀疏性[12].传统框架变换正则化模型中不同框架系数对稀疏性恢复的影响因子是一致的，其忽略了变换后框架系数与原始图像结构信息之间的联系，从而不利于框架变换后系数稀疏性的较好恢复.

本文提出一种新的紧框架域下，对稀疏性约束项进行重加权的正则化模型：首先通过正交匹配追踪算法(OMP)[13]确定图像多层紧框架变换后多尺度不同方向上具有较大系数的集合，即稀疏支撑集，由于支撑集内的框架系数包含了原始图像主要结构信息，故对其加上较小的惩罚因子，而对于其他为零或接近于零的框架系数加上较大的惩罚因子，实现不同系数间的非均匀惩罚约束.模型的算法设计上，结合split bregman方法的思想，提出一种多步交替优化算法，将每一步交替迭代得到的恢复结果，应用在下一步的交替迭代过程中，充分利用图像的先验信息，大大提高了图像恢复的质量.在算法的验证过程中，通过真实的图像进行恢复实验，将本文紧框架域重加权正则化模型与基于TV的小波恢复模型[2]和传统紧框架域正则化模型[4]进行对比分析.

2 紧框架域重加权正则化模型

2.1 基于多分辨率分析的图像紧框架变换

对于可数集合X⊂L2(R)，被称为L2(R)上的紧框架，当且仅当：

(2)

<·>表示内积，可以得到标准正交基是一组特殊的紧框架.一般情况下，小波框架系统表示对确定集合ψ={φ1，φ2，…，φn}中元素的伸缩和平移，即{2l/2φm(2l·(-k))}⊂L2(R)，l和k为伸缩和平移因子，ψ称为生成集合.特别的是，若X为紧框架时，φ称为小框架.

(3)

若尺度函数对称，则生成的框架系统也是对称的，且由对称或反对称的尺度函数和小框架函数构造的紧支撑框架系统应具有线性相位，该特性十分重要.在图像数据分析中，当尺度函数和小波函数作为滤波函数，且滤波器具有线性相位时，能有效避免图像信号在分解与重构中的失真，从而保证图像恢复效果的改善[15].实际应用中，考虑到B样条函数的特性，本文在B样条函数基础上构造紧框架系统.如图1基于分段三次B样条构造的紧框架系统的尺度函数和小框架函数.

图1 分段三次B样条框架系统Fig.1 Framelet system of piecewise cubic

当对离散图像序列进行处理时，W∈m×n，m≥n代表分解算子，WT代表重构算子，根据UEP原则，则有WTW=I，且WWT≠I，则对于完整的图像重构过程描述为：

u=WTWu

(4)

为获取图像多方向上的细节特征，进一步对图像序列u作L层框架分解：

Wu=(…，Wl，ju，…)0≤l

(5)

其中Ω代表每一层分解中所包含的框架系数子带集，例如对于分段三次B样条构造的紧框架变换中|Ω|=25，小波变换模型中|Ω|=4；Wl，ju代表第l层分解的第j个子带图像的小波框架系数.下页图2(a)-图2(d)为分段三次B样条函数构造的紧框架系统对图像进行4层框架分解.

2.2 紧框架域L1范数正则化图像恢复模型

对系统和噪声模糊的图像，考虑其在紧框架变换下的稀疏性，建立紧框架变换框架系数L1范数的正则化项，得到图像恢复正则化模型表达式：

(6)

其中λ为保真项与正则化项间的调整参数，q=1或q=2分别对应于各向异性L1范数和各向同性L1范数，J.F Cai等证明了在选择正则化项时各向同性的L1范数明显优于各向异性L1范数[16]，因此本文优先选择q=2.

令α=Wu，将式(6)转换为条件约束优化问题，即：

(7)

图2 4层框架分解子带图像Fig.2 Four-level framelets decomposition of u in bands

2.3 紧框架域重加权正则化图像恢复模型及其分裂算法

2.3.1 紧框架域重加权正则化图像恢复模型

传统紧框架变换L1范数模型建立的稀疏正则化能量泛函本质上是紧框架变换后框架系数的稀疏性恢复，但模型算法实现过程中稀疏性的惩罚系数是一致的，并没有充分利用紧框架变换后框架系数的特性.在图像稀疏表达的惩罚项中，对少数包含大量原始图像信息的框架系数与其他系数应赋予不同程度的影响因子，以平衡所有系数对图像信号稀疏性恢复效果的影响.本文在传统的紧框架变换L1范数正则化模型基础上，通过引入框架分解系数模相关的权值函数κ，提出了一种新的紧框架变换L1范数的重加权正则化模型，即：

(8)

其中0≤l≤L-1，i={1，2，…，n}，n表示像素位置，j∈Ω.对于运算符⊗，若a=(a1，a2，…，ak)，b=(b1，b2，…，bk)，c=(c1，c2，…，ck)，则：a⊗b=(a1×b1，a2×b2，…，ak×bk)，a⊗b⊗c=(a1×b1×c1，a1×b1×c1，…，ak×bk×ck).λl，j为第l层第j个子带图像保真项与正则项的调整参数.

传统紧框架变换L1范数模型中，框架系数对稀疏性恢复的影响因子是一致的，即式(8)中κl，i，j始终为1.而在本文的重加权正则化模型中，对少量包含原始图像结构信息的框架系数和其它非零或接近于零的系数实行非均匀惩罚，令

(9)

其中|Wul，i，j|代表框架系数的模，ε为调和参数.

传统紧框架L1正则化，忽略了紧框架变换后框架系数与原始图像之间的联系，丢失了很多重要的结构信息，不利于框架变换后稀疏性的恢复，而本文模型将框架变换后稀疏支撑集内的多尺度不同方向的框架系数引入非均匀惩罚因子中，充分利用原始图像的多层结构信息，大大增强正则化模型稀疏性的恢复，从而更好地保护恢复图像边缘及一些微小的细节信息，改善图像恢复的质量.

2.3.2 重加权正则化图像恢复模型的分裂算法(RMSBA)

紧框架变换正则化模型算法由于权值因子始终为1，故只采用单步交替迭代优化得到最终恢复结果.本文模型在split bregman算法的思想基础上，结合恢复过程中权值因子的更新，提出一种重加权多步交替迭代优化算法(RMSBA)，将原始复杂的正则化模型(8)解耦为两个较简单的子问题进行重复交替优化：

1)固定权值函数κl，i，j(第一步κl，i，j=1)，利用紧框架变换正则化模型算法求解关于u的凸优化问题；

2)利用步骤1)得到的恢复图像u，自适应更新权值函数κl，i，j，转步骤1).

首先将模型(8)转化为约束优化问题，即：

利用增广拉格朗日乘子技术，引入辅助变量ν，根据增广拉格朗日原理，将有条件约束最优化问题(10)转化为易于处理的无条件约束最优化问题：

(11)

关于u的最小化过程，采用最小二乘算法：

(ATA+μI)uk+1=ATf+μWT(αk-νk)

(12)

(13)

可转化为等价形式：

(14)

由上述分析过程可得到本文模型的分裂算法(RMSBA)：

1)给定观测图像f，线性算子A，交替优化的步数s=1，令u为第s-1步交替优化后的输出图像，并且初始化α0=v0=Wus-1，κl，i，j=1，迭代终止误差ζ，外部和内部迭代最大次数分别为Smax，MaxIt

2)u的更新：uk+1=(ATA+μI)-1(ATf+μWT(αk-νk))

4)ν的更新：νk+1=vk+Wuk+1-αk+1，k=k+1

5)若err>ζ或者迭代次数k

6)寻找u在紧框架变换下的框架系数的稀疏支撑集，记：θ=supp(Wu)={i：Wul，i，j≠0} ，0≤l≤L-1

8)s=s+1

9)若s≤Smax，转步骤2)

10)输出最终恢复图像u

3 实验结果及分析

对于不同的图像恢复实验，采用不同方法生成线性算子A：

图3 卷积算子Fig.3 Convolution operator

(II)图像去卷积：fspecial(′gaussian′,[7,7],5)，如图3.

图4 原始图像Fig.4 Original image

3.1 实验1——图像修复

实验中采用线性算子(I)对原始图像进行处理.图5第一列为图像破损率从30%到60%时对应的破损图像，第二列至第四列为不同算法恢复结果的定性对比，下页表1为不同算法所获得恢复结果的定量对比，下页图6为不同PR值时，不同算法恢复结果PSNR的变化曲线，图7为不同PR值时，不同算法恢复结果SSIM的变化曲线.

图5 δ=5，不同PR值(30%-60%)时，不同算法的恢复性能比较Fig.5 δ=5，the PR value from first row to four row rangesfrom 30% to 60%，respectively，compariation of restoration by different algorithms

表1 不同算法恢复图像的SSIM和PSNR指标对比Table 1 SSIM and PSNR comparison of restoration image by different algorithms

表2 多步交替优化过程中图像恢复的SSIM和PSNR指标对比Table 2 SSIM and PSNR comparison of restoration image in multi-stage alternating optimization

从图5获得的恢复结果中看出，TV小波正则化模型在一些平滑区域很容易引起“阶梯效应”，不能很好的修复破损图像；文献[4]算法较文献[2]算法有较好的恢复效果，但是在一些细小的边缘处易产生模糊，例如相机支架之间细小长条的边缘处并不清晰，而本文算法能够很好地解决边缘模糊的现象，图像的细小边缘恢复效果较为理想.

图6 δ=5，PR值(30%-60%)时，不同模型PSNR值的变化曲线Fig.6 δ=5，changing curve of different models PSNR values

从表1可知，本文算法和文献[4]算法比文献[2]算法有更高的SSIM，PSNR值，多层紧框架分解稀疏表示之后，可以更加精确的描述图像的特征，从而可以获得更好的图像恢复效果.而本文算法比文献[4]算法具有更高的SSIM，PSNR值，提高了2-4db，表明本文算法较文献[4]更加有效.

分析图6和图7可以对不同算法稳定性进行对比.不同PR值时，对本文5幅原始图像的修复效果来看，本文算法较其他两种算法一直具有更高的SSIM，PSNR值，从曲线的近似下降速率可以看出，紧框架域的正则化模型较TV正则化更加的稳定.而本文的算法稳定性比文献[4]算法也有一定的提升.

图7 δ=5，PR值(30%-60%)时，不同模型SSIM值的变化曲线Fig.7 δ=5，changing curve of different models SSIM values

图8 PR=0.5，δ=5时，多步交替优化的恢复效果Fig.8 Restoration performance comparison of multistage alternating optimization

图8为PR=0.5，δ=5时，本文所采用的多步交替迭代优化过程.图8(a)为第一步交替迭代过程，图9(a)是其局部放大图像，可以看出相机支架区域出现明显的边缘修复不完全的现象.图8(b)以第一步交替迭代的修复结果作为先验信息，更新权值函数，进行交替迭代，得到的修复结果由图9(b)的局部放大图像中可以看出有了明显的改善.再进行图8(c)第三步后，相机支架部分从图9(c)中看出边缘部分模糊现象已基本消失.从表2中可知，多步交替优化过程可以提高修复图像的SSIM和PSNR值，从图像修复效果直观分析，对于相机支架中间较粗的支架修复效果越来越清晰，边缘模糊现象也基本消失，大大改善图片修复效果.

图9 PR=0.5，δ=5时，多步交替优化恢复图像的局部放大Fig.9 Local areas of multistage alternating optimization restoration

3.2 实验2——图像去卷积

实验中采用卷积算子(II)对原始图像进行处理.图10(a)，(e)为待恢复的模拟图像，(b)(c)(d)和(f)(g)(h)为不同算法恢复结果的定性对比，表3为不同算法所获得恢复结果的定量评价指标对比.

图10 不同算法恢复性能的比较Fig.10 Restoration performance comparison of different algorithm

从图10(a)和图10(e)可知，图像视觉效果不理想，文献[2]算法恢复的图10(b)对于中间及其相连的两根短的支架并没有恢复，且两侧黑色支架恢复边缘也十分模糊，在平滑区域不满足“连通性准则”，不能很好地恢复图像.图10(f)对于船头上的文字部分并不能够恢复出来，左上角处船杆的恢复也较为模糊.文献[4]算法恢复的图10(c)虽然改善了黑色支架边缘模糊的现象，但是对于相机支架之间相连的细小长条并不能清晰的恢复出来，图10(g)中也不能够清晰地恢复出船头上的文字.而本文算法的修复图10(d)支架部分恢复十分清晰，边缘也得到较好保护，图10(h)基本恢复出船头的文字，几乎接近于原始图像，充分说明本文算法恢复性能比较理想.

从表3可知，文献[2]算法获得的SSIM和PSNR数值最小，而本文算法获得的SSIM和PSNR数值最大，进一步验证了本文对框架系数进行非均匀惩罚的思想，对于图像紧框架变换后稀疏性增强的有效性.

4 结 论

结合压缩感知理论，在传统紧框架变换L1范数正则化模型的基础上，对紧框架变换后的稀疏系数实行非均匀惩罚，提出一种新的紧框架域L1范数的重加权正则化模型.将多层框架分解系数的模作为先验信息引入权值因子，充分利用了紧框架变换后框架系数所包含的原图像结构信息.由多分辨率结构得到的紧框架模型保证了能量泛函求解快速算法的存在.本文在经典split bregman的算法思想基础上，提出一种多步重复交替优化算法(RMSBA).仿真实验表明，该模型相比于传统的TV小波模型和紧框架变换L1范数正则化模型，能够有效改善图像恢复过程边缘细节模糊的不足，大大提高图像恢复的质量，且模型具有较好的稳定性，能够取得较高的峰值信噪比和结构相似性测度.另外本文所采用的紧框架系统对于图像结构成分的稀疏表示十分有效，但对于图像纹理部分的稀疏表示结果并不理想，因此研究能同时稀疏表示图像结构和纹理成分的框架系统是今后的重要方向.

[1] Chambolle A，Lions P L.Image recovery via total variation minimization and related problems[J].Numerische Mathematik，1997，76(2)：167-188.

[2] Chan R，Wen Y W，Yip A.A fast optimization transfer algorithm for image inpainting in wavelet domains[J].IEEE Transactions on Image Processing，2009，18(7)：1467-1476.

[3] Chan R，Riemenschneider S，Shen L，et al.Tight frame：an efficient way for high-resolution image reconstruction[J].Applied and Computational Harmonic Analysis，2004，17(1)：91-115.

[4] Cai J F，Osher S，Shen Z.Split bregman methods and frame based image restoration[J].SIAM Journal on Multiscale Modeling and Simulation，2009，8(2)：337-369.

[5] MAT Figuriredo，R D Nowark.A bound optimization approach to wavelet-based image deconvolution[C].IEEE International Conference on Image Processing，2005，2：782-785.

[6] Cai J F，Chan R，Shen Z.Simutaneouscartoon and texture inpainting[J].Inverse Problems Images，2010，4(3)：379-395.

[7]Zhang Y，Dong B，Lu Z.L0 minimization of wavelet frame based image restoration[J].Mathematics of Computation，2013，82(282)：995-1015.

[8] Chan R，Chan T，Shen L，et al.Gorithms for high-resolution image reconstruction[J].SIAM Journal on Scientific Computing，2003，24(4)：1408-1432.

[9] Goldstein T,Osher S.The split bregman algorithm for L1 regularized problems[J].SIAM Journal on Imaging Science，2009，2(2)：323-343.

[10] Zhu J，Chen C W.A compressive sensing image reconstruction algorithm based on low-rank and total variation regularization[J].Journal of Jinling Institute of Technology，2015，31(4)：23-26.

[11] Dornoosh Zonoobi，Ashraf A Kassim.On the reconstruction of sequences of sparse signals-the weighted-cs[J].Journal Visual Communication and Image Representation，2013，24(2)：196-202.

[12] Candes E J，Wakin M B，Boyd S P.Enhancing sparsity by reweighted L1 minimization[J].Journal of Fourier Analysis and Applications，2008，14(5-6)：877-905.

[13] Sahoo S K，Makur A.Signal recovery from random measurements via extended orthogonal matching pursuit[J].IEEE Transactions on Signal Processing，2015，63(10)：2572-2581.

[14] Daubechies I，Han B，Ron A，et al.Framelets：MRA-based construction of wavelet frames [J].Applied and Computational Harmonic Analysis，2003(1)，14：1-46.

[15] Ron A，Shen Z.Affine system in L2(Rd)：the analysis of the analysis operator[J].Journal of Functional Analysis，1997，148(2)：408-447.

[16] Cai J F，Dong B，Osher S，et al.Imagerestorations：totalvariation，wavelet frames and beyond[J].Journal of the American Mathematical Society，2012，25(4)：1033-1089.