中央民族大学理学院，北京

收稿日期：2023年1月26日；录用日期：2023年2月21日；发布日期：2023年2月28日

摘要

针对分数阶全变差图像修补模型的梯度下降算法计算速度缓慢的问题，本文受到增广拉格朗日方法在图像去噪领域成功应用的启发，对基于分数阶全变分的图像修补模型也采用增广拉格朗日方法来进行求解。首先将原问题转化为等价的对应于增广拉格朗日泛函的鞍点问题，然后使用交替方向法将鞍点问题分解成u子问题和p子问题的序列来进行求解。接着对u子问题应用快速傅里叶变换(FFT)来求解，对p子问题应用Shrinkage算子来求解。参数的选取对数值仿真有很大的影响，本文对参数的选择进行了大量的实验并选取了较好的参数，实验表明数值算法是收敛的，并对含有文字遮挡及有人工划痕的实验图像有较好的修补效果。

关键词

图像修补，分数阶全变差，增广拉格朗日方法，参数选取

A Fast Algorithm Based on the Fractional Order Total Variation Image Inpainting Model

Qiang Liu, Zhaoxia Liu^*

College of Science, Minzu University of China, Beijing

Received: Jan. 26^th, 2023; accepted: Feb. 21^st, 2023; published: Feb. 28^th, 2023

ABSTRACT

In view of the slow calculation speed of gradient descent algorithm of fractional order total variation image inpainting model, inspired by the successful application of augmented Lagrangian method in the field of image denoising, this paper also uses the augmented Lagrangian method to solve the image inpainting model based on the fractional order total variation. Firstly, the original problem is transformed into an equivalent saddle point problem corresponding to the augmented Lagrangian functional, and then the saddle point problem is decomposed into a sequence of u subproblem and p subproblem by the alternate direction method. Then the fast Fourier transform (FFT) is applied to solve the u subproblem and the shrinkage operator is applied to solve the p subproblem. The selection of parameters has a great influence on the numerical simulation. This paper has carried out a lot of experiments on the selection of parameters and selected better parameters. The experiments show that the numerical algorithm is convergent, and has a good inpainting effect on the experimental images with text occlusion and artificial scratches.

Keywords:Image Inpainting, Fractional Order Total Variation, Augmented Lagrangian Method, Parameter Selection

Copyright © 2023 by author(s) and Hans Publishers Inc.

This work is licensed under the Creative Commons Attribution International License (CC BY 4.0).

http://creativecommons.org/licenses/by/4.0/

1. 引言

数字图像修补是利用图像的未缺损已知信息对图像中的待修复区域进行填充，以便得到可信的修复效果 [1] ，它本质上是一个插值问题 [2] 。数字图像修补作为低层图像处理技术之一，在文物的修复、图片上文本的去除、老旧视频和照片的修复以及人脸识别中眼镜的自动移除等诸多领域有广泛应用。

2000年Bertalmio、Sapiro、Caselles和Ballester [3] 提出的一种基于三阶PDE的图像修复模型(BSCB模型)；后来，Chan和Shen等人 [4] 通过借鉴ROF去噪模型，提出了TV修补模型及曲率驱动扩散(Curvature Driven Diffusion，CDD)模型，克服了TV修复模型所不能解决的连通性问题 [5] 。

2011年，张军等人提出了一类分数阶多尺度的图像去噪变分模型 $B{V}_{\alpha }-M{H}^{-s}$ [6] ：

$\underset{u\in B{V}_{\alpha }}{\min }\left\{E\left(u\right)=J_{\alpha }\left(u\right)+\frac{\lambda }{2}{\displaystyle \underset{j=0}{\overset{L}{\sum }}{2}^{-2j{s}_j}{\Vert {\left(f-u\right)}_j\Vert }_X^2},1\le \alpha \le 2,0\le s_j\le 1,j=0,1,\cdots ,L\right\}$ (1)

其中 ${\left(f-u\right)}_j,j=0,1,\cdots ,L$ 是 $f-u$ 在不同尺度下的分量，该文献详细介绍了分数阶梯度及其共轭算子的理论与数值离散方法。

2012年，张益等人受文献 [7] [8] [9] 的启发，提出了一类可处理噪声图像的分数阶p-Laplace TV图像修复模型 [10] ：

$\min J_{\lambda }\left[u\right]=\frac{1}{p}{\displaystyle {\int }_{E\cup \Omega }{\left|{\nabla }^{\alpha }u\right|}^p\text{d}x\text{d}y}+\frac{{\lambda }_E}{2}{\displaystyle {\int }_{E\cup \Omega }{\left|u-u_0\right|}^2\text{d}x\text{d}y},\alpha \in R^{+},p\in \left[1,2\right]$ (2)

其中， ${\lambda }_E=\lambda \cdot 1_E\left(x,y\right)=\{\begin{array}{l}\lambda ,\left(x,y\right)\in E\\ 0,\text{otherwise}\end{array}$ ，u在 $\partial \Omega$ 上满足边界条件 $\frac{\partial u}{\partial n}=0$ 。此模型是采用梯度下降算法进行计算，速度较慢。2018年，张俊等人针对分数阶TV泊松去噪模型 [11] ：

$\underset{u>0}{\min }{\displaystyle {\int }_{\Omega }\left|{\nabla }^{\alpha }u\right|\text{d}x}+\alpha {\displaystyle {\int }_{\Omega }\left(u-f\log u\right)\text{d}x}$ (3)

应用增广拉格朗日快速算法等优化算法进行求解。2019年，张龙等人对去除脉冲噪声的混合变分 $L^1$ 保真项模型采用增广拉格朗日方法进行实现 [12] ，取得了良好的视觉效果。2020年，刘洪琛等人使用增广拉格朗日方法成功地求解了抑制混合噪声的 $L^2$ 和Kullback-Leibler保真项混合模型 [13] 。台雪成等人首次使用增广拉格朗日方法求解欧拉弹性模型并应用于图像修补领域 [14] 。

由于已有分数阶变分修补模型采用梯度下降流方法进行数值离散，计算速度慢，需要计算大量的系数，复杂性高。受文献 [11] 的启发，本文对分数阶全变差修补模型提出增广拉格朗日快速算法，分为u子问题和p子问题，u子问题应用快速傅里叶变换(FFT)进行求解，p子问题应用Shrinkage算子进行求解。由于参数的选取对数值实验有很大影响，我们对参数进行了优化。最后的数值实验表明，本文算法数值上是收敛的，达到了一定的修补效果。

本文第二节提出分数阶全变差图像修补模型的增广拉格朗日算法及子问题的具体求解方法；第三节介绍对参数的优化及数值实验，从客观评价指标峰值信噪比(PSNR)和结构相似性指数(SSIM)及主观视觉效果可以看出本文提出的算法有较好的修补效果。

2. 模型及快速算法

本文要求解的分数阶全变差图像修补模型( $FOTV-L_2$ )的形式为：

$\underset{u>0}{\min }{\displaystyle {\int }_{E\cup \Omega }\left|{\nabla }^{a}u\right|\text{d}x}+{\alpha }_E\left(x\right){\displaystyle {\int }_{E\cup \Omega }{\left(u-f\right)}^2\text{d}x}$ (4)

其中， $a>0$ ， ${\alpha }_E\left(x\right)=\{\begin{array}{l}\alpha ,x\in E\\ 0,\text{otherwise}\end{array}$ ， $\Omega$ 是修复(开)区域， $\partial \Omega$ 是它的边界，E是在 $\partial \Omega$ 周围的拓展区域，使得 $\partial \Omega$ 位于 $E\cup \Omega$ 的内部。将图像u看作是一个大小为 $N\times N$ 的实矩阵，如果设空间 $X=R^{N\times N}$ ， $Y=X\times X$ ，则有 $u\in X,{\nabla }^{a}u\in Y$ 。

到目前为止，已有的求解分数阶图像修补模型的数值算法都是采用梯度下降流的方法，速度缓慢。为了克服这个缺点，我们对问题(4)，引入辅助变量 $p$ ，可转化为一个带约束的优化问题：

$\{\begin{array}{l}\underset{u>0,p}{\min }{\Vert p\Vert }_1+{\alpha }_E{\displaystyle {\int }_{E\cup \Omega }{\left(u-f\right)}^2\text{d}x}\\ \text{s}\text{.t}\text{.}p={\nabla }^{a}u\end{array}$ (5)

对问题(5)我们采用增广拉格朗日算法进行数值求解。首先写出(5)式的增广拉格朗日泛函如下：

$L\left(u,p;\lambda \right)={\displaystyle {\int }_{E\cup \Omega }\left|p\right|\text{d}x}+\langle \lambda ,p-{\nabla }^{a}u\rangle +\frac{\mu }{2}{\Vert p-{\nabla }^{a}u\Vert }_2^2+{\alpha }_E{\displaystyle {\int }_{E\cup \Omega }{\left(u-f\right)}^2\text{d}x}$ (6)

于是，对原问题(4)的求解转化为求解鞍点问题：

寻找 $\left(u^{\ast },p^{\ast };{\lambda }^{\ast }\right)\in X\times Y\times Y$ ，使得 $L\left(u^{\ast },p^{\ast };\lambda \right)\le L\left(u^{\ast },p^{\ast };{\lambda }^{\ast }\right)\le L\left(u,p;{\lambda }^{\ast }\right)$

$\forall \left(u,p;\lambda \right)\in X\times Y\times Y$ (7)

对问题(7)的求解，我们给出如下定理：

定理2.1 形如(7)式的鞍点问题在空间 $X\times Y\times Y$ 中至少存在一个解。

这个结论的证明与文献 [15] 中的方法完全类似，故不再重复。

于是我们可以导出求解鞍点问题(7)的一般算法如下：

第一步：输入 $u^0=0,p^0=0,{\lambda }^0=0$ ；

第二步：对 $k=0,1,2,\cdots$ 固定拉格朗日乘子 ${\lambda }^k$ ，近似计算：

$\left(u^{k+1},p^{k+1}\right)\approx \underset{u,p}{\arg \min }L\left(u,p;{\lambda }^k\right)$ (8)

第三步：更新拉格朗日乘子：

${\lambda }^{k+1}={\lambda }^k+\mu \left(p^{k+1}-{\nabla }^a{u}^{k+1}\right)$ (9)

不过，在实际求解中我们会发现，问题(8)中变量u和 $p$ 的耦合导致直接求解非常困难，解决这个问题的方法是应用交替方向法。于是问题(8)转化为两个子问题，鞍点问题(7)转化为如下问题：

$\{\begin{array}{l}{u}^{k+1}=\underset{u}{\arg \min }L\left(u,p^k;{\lambda }^k\right)\\ p^{k+1}=\underset{p}{\arg \min }L\left(u^{k+1},p;{\lambda }^k\right)\\ {\lambda }^{k+1}={\lambda }^k+\mu \left(p^{k+1}-{\nabla }^a{u}^{k+1}\right)\end{array}$ (10)

先介绍u-子问题的求解，写出具体形式：

$\begin{array}{c}{u}^{k+1}=\underset{u}{\arg \min }L\left(u,p^k;{\lambda }^k\right)\\ =\underset{u}{\arg \min }\langle {\lambda }^k,p^k-{\nabla }^{a}u\rangle +\frac{\mu }{2}{\Vert p^k-{\nabla }^{a}u\Vert }_2^2+{\alpha }_E{\displaystyle {\int }_{E\cup \Omega }{\left(u-f\right)}^2\text{d}x}\end{array}$ (11)

(11)式对应的欧拉–拉格朗日方程为：

$-{\left({\nabla }^a\right)}^{\ast }{\lambda }^k-\mu {\left({\nabla }^a\right)}^{\ast }\left(p^k-{\nabla }^{a}u\right)+2{\alpha }_E\left(u-f\right)=0$ (12)

我们用快速傅里叶变换(FFT)从(12)式中求解出 $u^{k+1}$ ：

$u^{k+1}=F^{-1}\left(\frac{F\left[\left({\alpha }_{E}f+{\left({\nabla }^a\right)}^{\ast }{\lambda }^k\cdot \frac{1}{2}\right)/\mu +\frac{1}{2}{\left({\nabla }^a\right)}^{\ast }{p}^k\right]}{\frac{{\alpha }_E}{\mu }+\frac{1}{2}F\left({\left({\nabla }^a\right)}^{\ast }\left({\nabla }^a\right)\right)}\right)$ (13)

再求解p-子问题：

$\begin{array}{c}{p}^{k+1}=\underset{p}{\arg \min }L\left(u^{k+1},p;{\lambda }^k\right)\\ =\underset{p}{\arg \min }{\Vert p\Vert }_1+\langle {\lambda }^k,p-{\nabla }^a,u^{k+1}\rangle +\frac{\mu }{2}{\Vert p-{\nabla }^a{u}^{k+1}\Vert }_2^2\\ =\arg \min {\Vert p\Vert }_1+\frac{\mu }{2}{\Vert p-\left({\nabla }^a{u}^{k+1}-\frac{{\lambda }^k}{\mu }\right)\Vert }_2^2\end{array}$ (14)

观察(14)式的最后一行，可推出其具有闭形式解：

$p^{k+1}=Shrinkage\left({\nabla }^a{u}^{k+1}-\frac{{\lambda }^k}{\mu },\frac{1}{\mu }\right)$ (15)

将式(13) (15)的结果代入(10)式就可以直接计算 $\lambda$ 。

为了加快算法收敛的速度，我们对 ${\lambda }^{k+1}$ 的计算公式中的惩罚因子 $\mu$ 进行一次更新迭代：

${\lambda }^{k+1}={\lambda }^k+{\mu }^k\left(p^{k+1}-{\nabla }^a{u}^{k+1}\right)$ (16)

${\mu }^{k+1}=\theta \cdot {\mu }^k$ (17)

其中参数 $\theta$ 是引入的加速因子，需要根据特定的待修补图像选取合适的值。

综上所述，我们得出求解 $FOTV-L_2$ 图像修补模型的增广拉格朗日算法如下：

Step 0：为参数赋值： $\mu ,K,a,\alpha ,\theta$ (K为计算分数阶梯度时所用到的邻域数)

Step 1：输入：待修补图像f， $p=p_0,\lambda ={\lambda }_0$ ，最大迭代步数N

Step 2：迭代：

For i = 1：N = 1000

I.利用公式(13)计算 $u^{k+1}$ ；

II.利用公式(15)计算 $p^{k+1}$ ；

III.利用公式(16)更新 ${\lambda }^{k+1}$ ；

IV.利用公式(17)更新惩罚参数 ${\mu }^{k+1}$ ；

If

$\frac{{\Vert u^{k+1}-u^k\Vert }_2}{{\Vert u^k\Vert }_2}<1.0\times {10}^{-4}$ (18)

Break

End for

Step 3：输出 $u=u^{k+1}$

3. 实验结果与分析

本文实验均在Intel(R) Core(TM) i7-6500U CPU @ 2.50 GHz 2.59 GHz处理器上运行，操作系统为64位Windows 10，编译程序所使用的软件为MATLAB R2019b。本文所使用的三张测试图像是大小均为256 × 256的添加了文字遮挡的灰度图像Cameraman (如图4(b)所示)、添加了人工涂痕的灰度图像Barbara (如图6(b)所示)以及添加了人工涂痕的灰度图像Baboon (如图8(b)所示)。

3.1. 参数的优化方法

由于算法1中的五个参数(惩罚参数 $\mu$ 、分数阶梯度计算中的邻域数K、分数阶阶次a、数据保真项参数值 $\alpha$ 以及加速因子 $\theta$ )取值不同会极大地影响图像修补的效果甚至影响算法收敛性，故使用算法1对上述三张测试图像进行数值实验的最重要的环节便是为五个参数选取合适的数值。

本文对参数的优化采取“控制变量、逐步搜索”的方法，判断参数是否达到最优除了观察图像修补的视觉效果以外，更客观的判据便是经典的图像质量评价指标：峰值信噪比(PSNR)和结构相似性指数(SSIM)。具体的优化方法是：先优化第一个参数 $\mu$ ，此时按经验为其他四个参数选取合适的初值，之后为 $\mu$ 选取第一个范围(通常是10³~10⁴，步长为1000)来循环地执行算法1，每执行一次算法1便会在算法停止步数(由停止准则决定)内确定一个最大PSNR值(或最大SSIM值)，然后画出参数 $\mu$ 关于最大PSNR值的曲线，再求出第一个参数范围内的最大参数值和对应的PSNR值(SSIM值)，通过画出的曲线来确定下一步的参数取值范围(第二个参数范围的步长通常为下一数量级100，但要避免选取局部极值附近的范围)然后再循环执行算法1，以此类推，直到第一个参数 $\mu$ 的数值精确到小数点后第二位(因为此时再去对 $\mu$ 调优不会使评价指标值增加超过0.0001)。此时已经优化了第一个参数 $\mu ={\mu }^{\ast }$ ，接下来固定 $\mu ={\mu }^{\ast }$ 以及 $a,\alpha ,\theta$ 的初值，对K取值从1到256，步长为1来循环地执行算法1，求出最优的K值 $K^{\ast }$ 以及对应的PSNR值(SSIM值)。其余三个参数的优化方法与 $\mu$ 的优化完全类似，只需注意：1) 参数a的最大取值一般不超过200，否则程序无法运行；2) a通常合理的调参范围在0~4附近，否则修补结果会越来越暗直到变为全黑图像；3) 参数 $\theta$ 只有在 $\theta >1$ 时会起到加速收敛的作用，在 $0\le \theta \le 1$ 时会使得计算速度缓慢并得到失真图像(全黑色以及“雪花”图像)。最后我们得到了一组最优参数以及最大的PSNR值(SSIM值)。综上所述，我们把上述参数优化方法按观察的评价指标分为了按PSNR最大调参和按SSIM最大调参两大类。

经过实验，我们画出了三张测试图像在两种参数调优实验下的评价指标随着迭代步数变化的曲线图，具体的最优参数值将在3.2节介绍，需要指出的是每张图像在两种调参方法下会得到两组不同的最优参数值。为了图片标题的简洁性，用“方法1”表示“按PSNR最大调参”，用“方法2”表示“按SSIM最大调参”。

图1(a)、图1(b)分别展示了Cameraman图像上应用算法1按PSNR最大调参方法下所得出的PSNR与SSIM随迭代步数的变化曲线，算法运行11步停止(最大迭代步数为1000，满足停止准则(18)时迭代停止)，同时达到最大PSNR值为23.3032 (大于未修补图像的PSNR：19.4686)，此时SSIM值也逐步增加直至平稳为0.6134 (小于未修补图像的SSIM：0.8912)。图1(c)、图1(d)是摄影师图像在另一种调参方法下画出的曲线，实验运行了12步停止，但是图1(c)告诉我们按SSIM调参无法使PSNR值得到提升，因为第1步是最大的PSNR值等于23.0023 (大于未修补图像的PSNR：19.4686)，之后逐步下降到很低的水平。图1(d)表明在第12步达到最大SSIM为0.7437 (小于未修补图像的SSIM：0.8912)。这个实验证明无法在同一套最优参数下使得Cameraman图像同时取得最佳的PSNR与SSIM。

图2绘制了Barbara图像应用算法1在两种调参方法下的评价指标曲线，纵观四张曲线图可以发现按一种评价指标最大调参并没有影响到另一种评价指标值的提升。图2(a)、图2(b)展示了按PSNR调参的PSNR与SSIM曲线，实验运行到11步停止，在第3步达到最大的PSNR值为27.4416 (大于未修补图像PSNR：24.4007)，11步达到的SSIM值为0.8670 (小于未修补图像SSIM：0.9562)。图2的(c)、(d)绘制了Barbara图像在另一种调参方法下的PSNR与SSIM曲线，可以看出与Cameraman图像不同，按SSIM调参没有影响PSNR值的提升。实验运行14步停止并达到最大SSIM为0.9255 (仍小于未修补图像SSIM：0.9562)，14步的PSNR值为25.6002 (大于未修补图像PSNR：24.4007)。本实验证明对Barbara图像而言，无法在同一套最优参数下同时达到最佳的PSNR与SSIM。

图3展示了Baboon (狒狒)图像应用算法1在两种调参方法下的评价指标随迭代步数的变化曲线。图3(a)、图3(b)分别是在按PSNR调参下的PSNR曲线与SSIM曲线。实验运行13步停止，达到最大PSNR值为22.6092 (大于未修补图像PSNR：18.7732)，达到SSIM值为0.4832 (小于未修补图像SSIM：0.8525)。图3(c)、图3(d)展示了按SSIM调参下的PSNR与SSIM曲线，可以看出对狒狒图像而言按SSIM调参不能让PSNR值稳步提升，而会先快速上升再逐步下降。实验运行到14步停止，在14步达到最大SSIM为0.8089 (仍小于未修补图像SSIM：0.8525)，在第2步达到PSNR为20.5345 (大于未修补图像PSNR：18.7732)。本实验证明对Baboon图像来说无法在同一套最优参数下同时达到最佳的PSNR和SSIM。

综上所述，对本文所选取的三张测试图像来说，对模型(4)应用算法1求解时使用按PSNR调参与按SSIM调参两种方法来优化参数无法得出使PSNR和SSIM同时达到最大值的参数组。此外，在通常情况下，按PSNR调参时得到的PSNR值大于按SSIM调参得到的PSNR峰值，按SSIM调参时得到的SSIM值也大于按PSNR调参时得出的SSIM峰值。最后需要指出，按SSIM调参得到的所谓最优SSIM值仍小于未修补图像的SSIM，而不论按哪个评价指标来调参，所得的PSNR值都使未修补

图1. 摄影师图像的评价指标曲线

图像的PSNR得到了提升。三张测试图像在两种调参方法下的最大PSNR与SSIM值的对比汇总在了表1和表2中。

图2. 芭芭拉图像的评价指标曲线

图3. 狒狒图像的评价指标曲线

3.2. 数值实验

本节介绍三张测试图像(图4(b)、图6(b)和图8(b))上的数值实验及最优参数的具体结果，与此同时我们绘制了相对误差及各种能量泛函值随着迭代步数的变化曲线图，旨在说明算法1在我们选取的参数下是数值收敛的。

待修补的Cameraman图像(图4(b))是在原始干净未缺损的图4(a)上添加文字遮挡形成的，图像掩膜如图4(c)所示。按照PSNR调参的最优参数组为： $\mu =1.40,K=2,a=0.89,\alpha =0.41,\theta =1.77$ ，对应的最佳修补结果如图4(d)所示，相对误差曲线如图4(f)所示，可以看出算法1是收敛的。按照SSIM调参的最优参数组为： $\mu =0.26,K=2,a=0.84,\alpha =0.38,\theta =1.65$ ，对应的最佳修补结果如图4(e)所示，其相对误差曲线如图4(g)所示，验证了算法1的收敛性。图5(a)~(d)展示了按PSNR调参的各级能量曲线，从曲线的趋势可以看出u子问题能量逐渐下降直到稳定，但是p子问题能量值的变化趋势反而增大了，但是图5(c)、(d)符合问题(7)和原问题(4)能量求极小的要求，算法依然是收敛的。图5(e)~(h)展示的是另一种调参方法下的能量曲线，我们发现除了p子问题能量先下降后上升以外，其他能量依然符合求极小的要求，图5(g)、(h)印证算法1是收敛的。

待修补的Barbara图像(图6(b))是在原图(图6(a))的基础上添加了人工涂痕形成的，相应的图像掩膜如图6(c)所示。按PSNR调参得到最优参数组： $\mu =1.02,K=3,a=0.80,\alpha =0.42,\theta =1.76$ ，对应的最佳修补结果如图6(d)所示，相对误差曲线如图6(f)所示，表明算法收敛。按SSIM调参的最优参数组为： $\mu =0.33,K=3,a=0.27,\alpha =0.32,\theta =1.67$ ，相应修补结果如图6(e)所示，相对误差曲线如图6(g)所示，表明相对误差满足了停止准则(18)。图7展示了两种调参方法下的各种能量泛函的曲线图，与前述分析类似，不论是图7的第一行图像(按PSNR调参)还是第二行图像(按SSIM调参)，p子问题的能量趋势是上升直至平稳，但通过观察图7的(a)(c)(d)(e)(g)(h)的六张子图我们会发现原问题及其等价的鞍点问题依然满足求极小的要求，曲线先下降后趋于平稳，算法依然是收敛的。

图8展示的是Baboon(狒狒)图像修补实验及其相对误差曲线。在原图(图8(a))的基础上添加人工涂痕后得到待修补图像(图8(b))，图像掩膜如图8(c)所示。按PSNR调参的最优参数组为： $\mu =5.05,K=13,$ $a=0.60,\alpha =0.39,\theta =1.66$，对应的最优修补结果如图(d)所示，相对误差曲线如图8(f)所示，可以看出误差逐渐减小直到满足停止准则。按SSIM调参的最优参数结果为： $\mu =0.60,K=4,a=0.50,\alpha =0.37,$ $\theta =1.59$，对应的最佳修补结果如图8(e)所示，相对误差曲线如图8(g)所示，依然满足停止准则。图9的8张曲线图中前四张展示了按PSNR调参的u子问题能量、p子问题能量、增广拉格朗日能量及原问题能量随迭代步数的变化曲线，后四张展示的是按SSIM调参的相应能量曲线。通过观察可以看出p子问题能量先上升后平稳的趋势没有影响到原问题中极小值点的求解，图9(c)(d)(g)(h)曲线的变化趋势表明了数值算法的收敛性。

综观图4至图9我们发现：按照PSNR调参相比按SSIM调参可以更好地修补图像，但是会小程度地模糊部分像素值和纹理，按照SSIM调参虽然修补效果不佳但是会更大程度地保护已知区域的像素值及纹理，且我们可以看出无法在同一套最优参数下使得PSNR和SSIM值同时达到最大。

为了说明本文算法的有效性，与经典的TV图像修补模型的梯度下降算法 [16] 进行了对比。需要指出的是，在对Cameraman (摄影师)进行实验的时候，方法1与方法2会得到同一套参数(即 $\lambda =0,a=1.0\times {10}^3$ )，修补结果如图10(a)所示，而对另外两张测试图像(Barbara和Baboon图像)按照方法1与方法2调参会得到两组不同的参数：Barbara图像按照方法1调参得到参数 $\lambda =0,a=1.0\times {10}^4$ ，按方法2调参得到参数 $\lambda =0,a=16.19$ ；Baboon图像按方法1调参得到参数 $\lambda =0,a=147.13$ ，按照方法2调参得到参数 $\lambda =0,a=44.11$ 。Barbara与Baboon图像的修补结果见图10(b)~(e)。三张图像的TV修补评价指标值见表3，通过和表1、表2的对比我们发现，本文所提的分数阶图像修补模型的增广拉格朗日算法在PSNR上优胜于传统的TV算法，但是在SSIM指标上逊色于TV算法。

图4. 摄影师图像的实验结果与相对误差曲线

图5. 摄影师图像修补的各种能量曲线

图6. 芭芭拉图像的实验结果与相对误差曲线

图7. 芭芭拉图像修补的各种能量曲线

图8. 狒狒图像的实验结果与相对误差曲线

图9. 狒狒图像修补的各种能量曲线

图10. 三张测试图像的TV修补结果

表1. 三张测试图像使用方法1的评价指标对比

表2. 三张测试图像使用方法2的评价指标对比

表3. 三张测试图像的TV修补评价指标

4. 总结

本文针对目前分数阶全变差图像修补模型都是采用梯度下降流算法进行数值求解，速度慢的缺点，我们提出新的增广拉格朗日快速算法，将问题转化为拉格朗日泛函的鞍点问题，应用交替方向法，将问题分为u子问题和p子问题进行求解，u子问题采用FFT进行求解，p子问题采用Shrinkage算子进行求解。另外由于参数的选取对数值实验有很大影响，我们对参数进行了优化。最后进行了数值仿真以及与TV修补进行了对比，实验表明算法是收敛的且三张测试图像的PSNR值相比于TV修补得到了提高，但修补的直观效果仍有待进一步提高。

基金项目

国家自然科学基金面上项目，项目号：11871457。

文章引用

刘 强,刘朝霞. 基于分数阶全变差图像修补模型的快速算法
A Fast Algorithm Based on the Fractional Order Total Variation Image Inpainting Model[J]. 应用数学进展, 2023, 12(02): 752-763. https://doi.org/10.12677/AAM.2023.122077

参考文献