1、2024年3月Mar.,2024DOI:10.15960/ki.issn.1007-6093.2024.01.001求解一类线性等式约束凸优化问题的加速方法*孟辛晴1张文星1,t摘要具有线性约束的凸优化问题是数学规划中的一类经典问题。本文将借助对偶理论,研究求解一类具有线性等式约束的凸优化问题的加速算法。由于此类问题的对偶问题是一个具有两块可分离结构的凸优化问题,我们基于Goldstein等人在加速交替方向乘子法方面的重要工作,提出了一种在弱化条件下求解线性等式约束凸优化问题的加速方法。我们的方法与Goldstein等人的加速交替方向乘子法的不同之处为:1)目标函数仅要求具有凸性(而不必强凸)
2、;2)罚参数仅要求0(而不受目标函数的利普希茨常数、强单调系数的限制)。基于上述弱化的条件,我们证明了所提的加速交替方向乘子法依然具有收敛性和O(1/k)的收敛率。我们将条件弱化后的加速交替方向乘子法用于求解一个图像重建问题。数值实验结果表明,条件弱化后的加速交替方向乘子法依然具有较好的数值效果。关键词线性等式约束,对偶,可分离结构凸优化,交替方向乘子法,Nesterov加速技术中图分类号_0 2 2 1.22010数学分类号90 C47,90C26,90C30An accelerated method for solving a class oflinear equality constra
3、ined convex optimization*MENG XinqinglAbstract The linearly constrained convex optimization is a class of quintessentialproblem in optimization community.In this paper,we will devote to the algorithmicstudy of such problems on the perspective of duality.By introducing auxiliary variables,the dual of
4、 this problem can be reformulated into a separable structured optimization.Byvirtue of the recent advances in accelerated alternating direction method of multiplier,weanalyze the convergence and establish the O(1/k2)convergence rate of Goldstein et alsaccelerated algorithm under some mild conditions
5、.Specifically,1)the objective suficesto be convex instead of strongly convex;2)the penalty parameter merely requires 0instead of a bounded interval involving Lipschitz constant and strongly convex modulus.Some numerical experiments on image reconstruction problem were conducted to demon-strate the p
6、erformance of algorithm on solving linearly constrained convex optimization,which is computationally appealing.Keywords linear-equality constraint,duality,separable structured convex opti-mization,alternating direction method of multipliers,Nesterovs acceleration techniqueChinese Library Classificat
7、ion O221.22010 Mathematics Subject Classification 90C47,90C26,90C30这筹学学报(中英文)Operations Research TransactionsZAHNG Wenxingl.t第2 8 卷第1期Vol.28No.1收稿日期:2 0 2 1-0 7-2 5*基金项目:国家自然科学基金(No.11971003),中央高校基本业务费(No.ZYGX2019J090)1.电子科技大学数学科学学院,四川成都6 117 31;School of MathematicalSciences,U n i v e r s i t y o f
8、 El e c-tronic Science and Technology of China,Chengdu 611731,Sichuan,China+通信作者E-mail:2其中:Rn(-o 0,+o o)为一个具有利普希茨连续梯度的正常闭凸函数*;ARIxn 为给定的矩阵;bER为给定的向量。问题(1)是一类经典的具有线性等式约束的凸优化问题,这类问题被广泛地应用于信号压缩感知 1、计算机视觉中的图像修复和重建 2 、统计学习 3、系统辨识与实现 4 等各个科学领域中,并发挥着重要的作用。针对问题(1)的特殊情况,我们借助对偶理论,研究求解问题(1)的对偶问题的加速算法。具体地,我们假设目
9、标函数的临近点函数(见第1节中的定义4)有显式解或者能用已有的软件包进行高效的求解。借助共轭函数,问题(1)的Fenchel对偶形式可以表示为(2)其中sERn和tER为问题(1)的对偶变量;*为函数的共轭函数(见第1节中的定义3)。为了方便后面的讨论,不妨记F(s):=0*(s),G(t):=b T t。因此,我们可以得到问题(2)如下的等价形式(3)问题(3)是一个典型的具有可分离结构的线性等式约束优化问题,其对应的拉格朗日函数定义为L(s,t,A):=0*(s)+bTt-T(s+ATt),V(s,t,X)E2:=IR R R,其中入Rn为拉格朗日乘子。进一步地,问题(3)的增广拉格朗日函
10、数定义为Lp(s,t,):=0*(s)+b6Tt-AT(s+ATt)+号lls+ATt2,V(s,t,)E 2,其中0称为罚参数。交替方向乘子法(Alternating Direction Method of Multiplers,简称ADMM)最早是由 Glowinski与Marroccol5、G a b a y 与Mercier6于2 0 世纪7 0 年代提出的,适用于求解非线性的椭圆型方程问题。后来,交替方向乘子法的思想被广泛应用于求解具有可分离结构的凸优化问题 7-8 。借助增广拉格朗日函数,应用交替方向乘子法求解特殊问题(3)时的迭代格式为:给定任意的初值(to,入o)和固定的罚参数
11、0,通过依次计算Sk+1=argmin(F(s)-Ts+lls+ATsERnth+1=argmin/G(t)-入T(ATt)+tER!(入k+1=入k-(sk+1+ATt+1),*对于函数Q:R(-8 0,+o o l,如果它的有效定义域 dom(0):=E R|0(a)+o0)是非空集,即dom(0),则称函数是正常的。如果它的上镜图(或上图)epi(0):=【(,)RR/o(a)是RnR中的闭集,则称函数 是闭的(或下半连续的)。孟辛晴,张文星本文致力于求解如下具有线性约束的凸优化问题min(ac)s.t.Ac=b,min 0*(s)+bTts.t.$+ATt=O,min F(s)+G(t
12、)s.t.s+ATt=O。BllSk+1+AT28卷(1)(4)(5)(6)1期实现变量(s,t)和拉格朗日乘子入的更新。一方面,交替方向乘子法可以看作是对偶分解方法和增广拉格朗日乘子法的结合;另一方面,从算子分裂的角度,它也可以看作Douglas-Rachford算子分裂方法应用于求解可分离凸优化问题(3)的对偶问题。它是在对偶上升法的基础上进行了改进,在每次迭代过程中交替地求解若干个优化的子问题。因此,其迭代格式既有分解子问题(即将较复杂或者规模较大的问题分解为若干个简单易解的小规模子问题)的特点,又能保证算法在理论上的收敛性,这主要得益于增广拉格朗日乘子法和Douglas-Rachfor
13、d算子分裂法的一系列优点。大量的数值实验和应用领域的研究表明,交替方向乘子法的迭代速度相对较快,而且数值性能比较稳定。因此被广泛应用于求解工程领域中精度要求不是很高的大规模结构型凸优化问题。目前,许多国内外学者相继对交替方向乘子法进行了大量的研究,包括非精确方法、并行方法、序列方法、预估-校正方法、多块可分离结构的方法、动态罚参数选取方法、带有加速技术的方法、适用于分布式计算的方法、带有随机优化技术的方法等。该领域在近二十多年来发展迅速,产出了大量研究成果,限于篇幅,我们在此不一一列出,感兴趣的读者可以参考Glowinski9,Boyd10和Eckstein11 等学者在交替方向乘子法方面所做
14、出的开创性或综述性工作,或者参考国内南京大学何炳生教授主页上的关于交替方向乘子法的讲义t。尽管经典的交替方向乘子法(6)在理论方面已经比较成熟(如收敛性 5,12 、收敛率等都已被证明 13),并且在实际中得到了广泛的应用,但是随着Nesterov 的梯度加速技术的产生 14,促使学者们设计了一系列的加速方法,不仅提高了算法的数值性能,而且更加适用于求解大规模的实际应用问题。例如,当交替方向乘子法中的子问题为线性问题时,Zhang 等 15】提出了基于广义最小残差算法(Generalized Minimum Residual Method,简称GMRES)的加速算法来求解二次规划问题。Zhan
15、g等 16 在ADMM中使用安德森(A n d e r s o n)加速技术,即在交替方向乘子法的第二个原始变量和对偶变量的更新过程中使用安德森加速。Ouyang等 17 研究了线性化交替方向乘子法的加速迭代格式,针对不同的线性化目标函数或者正则项,作者将Nesterov的梯度加速技术 14 应用于线性化的交替方向乘子法中,提出了加速的线性化交替方向乘子法(Accelerated Linearized ADMM,简称A-LADMM)。特别地,他们的加速方法理论上可以应用于可行域无界的情形。Goldstein等 18 同样采用Nesterov的梯度加速技术,研究了当目标函数均为强凸函数时的加速交
16、替方向乘子法(Accelerated ADMM,简称A-ADMM)。特别地,当罚参数受控于一个跟强凸系数有关的区间时,该加速方法具有全局收敛性。Bai等 19 设计了非精确的加速增广拉格朗日乘子法的迭代格式,并用于求解矩阵优化问题。He和Yuan20的加速增广拉格朗日乘子法、Luo 和Yang21 的加速对称迭代格式的交替方向乘子法、Lei 和Xie22的加速交替极小化方法均是基于Nesterov的梯度加速技术,在求解原始优化问题(1)的算法中,对拉格朗日乘子入的更新过程进行加速处理。而本文则是从凸优化问题(1)的对偶形式出发,利用加速交替方向乘子法进行求解,并且利用对偶问题的特殊结构,进一步
17、地将加速交替方向乘子法中的若干收敛性条件进行了弱化。另一方面,上述三个加速方法是针对不同的算法的加速版本(具体地,分别是增广拉格朗日乘子法、对称格式的交替方向乘子法和交替极小化方法)。针对问题(2 的特殊结构,本文提出了相应的加速交替方向乘子方法,与文献 18 中的方法相比,我们弱化了算法收敛性分析中对目标函数的凸性要求(我们仅需要其中一个函数是强凸的),同时我们也弱化了文献 18 中的算法对罚参数的要求。在收敛性分析上,我们获得了理论上更好的收敛性结果(即对偶残差dk由Ildkll2=-1aii表示 R空间中标准内积。给定非空闭集合Rn,我们用o()表示定义在上的所有正常的闭凸函数。定义1设
18、函数f:R (-,+,其中 是一个非空闭凸集。对任意,yE,常数 E0,1,若函数f 满足f(a+(1-)y)f()+(1-)f(y),则称函数f为集合上的凸函数。特别地,若对于任意y,不等式(7)总是严格成立,此时称函数为集合上的严格凸函数;若存在常数 0,使得函数号2是凸函数,则称f为凸集上的强凸函数,其中c称为强凸系数。进一步地,函数的凸性与其(次)梯度的单调性密切相关。定义2 设T是定义在闭凸集R上的一个映射。(1)T在X上是单调的,如果满足:T(u),u-u)0 使得IIT(u)-T(u)Il LIlu-ull,Vu,uE X,特别地,当L=1时,称T为非扩张映射。引理1设连续可微函
19、数f:R n(-0,+)具有利普希茨连续梯度,Lf0为利普希茨常数,则有f(a)f(a)+(Vf(a),a-a)+Ia-l2,Va,a E R。孟辛晴,张文星1in28卷(7)(8)1期共轭函数与函数的临近点算子在优化算法的设计中具有广泛的应用 2 3,其定义为:定义3函数的共轭函数f*:R n R 定义为由上述定义可知,对任意a,yER,不等式f(a)+f*(y)(亦称Fenchel不等式)恒成立。特别地,共轭函数的(次)微分还满足以下的重要结论 2 3:Vo()V0*(g)。值得注意的是,共轭函数*始终是一个闭凸函数(即使函数f不是闭凸函数)。引理2 设函数f:R (-0,+l,如果*是系
20、数为c0的强凸函数,则函数f是可微的且f是1-利普希茨连续的。定义 4函数f:R (-0,十)是一个正常的闭凸函数,函数的临近点算子定义为proxf(a):=argmin f(y)+Illy-all2),Va E Rn。特别地,e1-范数的临近点算子为软阈值算子 2 4。下面将介绍一个临近点算子的重要性质。定理1设常数0,函数:R(-,+)是一个正常的闭凸函数。记*是函数的共轭函数。对任意的 R,有=prox(a)+prox()。此定理是著名的Moreau 恒等式,它揭示了函数与它的共轭函数*的临近点算子之间的关系。更多关于单调算子、临近点算子、共轭函数的性质,感兴趣的读者可以参阅文献 2 3
21、 等专著。2加速交替方向乘子法Goldstein 等在文献 18 中提出了一种基于目标函数均是强凸函数的加速交替方向乘子法,并且证明了当罚参数满足一定的条件时,才能保证算法具有收敛性。本文考虑到问题(2)中目标函数的特殊结构(不仅仅局限于目标函数都是强凸函数),在保证收敛的情形下,还进一步地弱化了算法对罚参数的要求,给出了求解问题(2)的加速交替方向乘子法(A-ADMM)。根据上述工作,我们将基于以下假设来求解可分离凸优化问题(2):假设1(1)问题(2)的目标函数是可微的,并且具有利普希茨连续梯度;矩阵A是行满秩的;(2)问题(2)的解集记为2,是非空的。即问题(2)存在最优解(s*,t*)
22、,以及存在入*E R,使得*=(s*,*,*)E 是问题(2)对应的鞍点问题的解。注1对于实际应用领域中的线性约束优化问题,上述的假设条件并不苛刻。例如,在信息科学、图像处理等领域的凸优化问题中 2 5-2 8 ,问题(2)的线性约束中的系数矩阵通常具有形如A=A1,A2,,An】的分块结构,其中AERmxpi表示A的第i块子矩求解一类线性等式约束凸优化问题的加速方法f*(c)=sup(c,y)-f(y)。yeRnyERn5(9)6阵。而且图像处理中的很多优化模型对应的矩阵A中一般会有一个子矩阵A;为单位矩阵。因此,自然满足矩阵A行满秩矩阵的条件。在基于光滑正则项的图像处理模型(例如Huber
23、正则化、Pseudo-Huber正则化、对数正则化等)中,目标函数通常是可微的,而且梯度是利普希茨连续的。2.1算法的提出Goldstein 等在文献 18 中所提出的加速交替方向乘子法是在每一次更新变量的过程中,将Nesterov的梯度加速技术应用到临近点算子的求解过程中,使迭代序列能够更快地逼近最优解。用于求解问题(3)具体迭代格式的伪代码见算法1所示。算法1:加速交替方向乘子法求解问题(3)1 初值:任取 t-1=to E R,入-1=oE R,0,1=1。2 for k=0,1,2,.,doSk=argmin F(s)-XTs+号lls+ATtil2,3SERnth=argmin G(
24、t)-XT(ATt)+lls+ATt2),4tERI5入i=入k-(sk+ATtk),6Qk+1=1Vi+4a,27X+1=入k+9-(入元 入k-1)。89 end针对加速交替方向乘子法(即算法1),我们可以用原始-对偶残差来量化算法的迭代点与最优解之间的逼近程度。一方面,问题(3)的最优性条件可描述为(10a)0=VF(s*)-*,(10b)(0=VG(t)-A*;(10c)另一方面,加速交替方向乘子法中s-子问题(即算法1中的第3行)满足如下形式0=VF(sk)-入h+(sk+ATth)。将算法1中入的更新方式(即算法1中的第5行)代入上式中,我们有O=VF(sk)-入k+AT(th-t
25、h)。同理,算法1中t-子问题满足如下形式0=VG(th)-AAk-(sk+ATth)。将算法1中入的更新方式(即算法1中的第5行)代入上式,我们有0=VG(tk)-A入k o孟辛晴,张文星$*+ATt*=0,28卷(11)(12)1期我们不妨记rk:=Sk+ATtk,d k:=AT(tk-t),对比(11)(12)式和问题(3)的最优性条件(10)可知:当rk0和dk0时,则有(sh,ts,入)满足最优性条件(10)。因此,在算法1的收敛性分析中,我们可以用原始-对偶残差(rk,ds)来衡量迭代点与最优解之间的逼近程度。2.2收敛性分析对于经典的交替方向法,只要给定定义域内的任意值(s,t,
26、入),我们都可以通过迭代格式得到相应新的迭代点(s,t,入),即s:=argmin F(s)SERnt:=argmin/G(t)-T(ATt)+lls+AteRl(X:=入-(s+ATt)。对于问题(3)而言,等价于求解如下的鞍点问题:maxs-minr F(s)+G(t)-(,s+ATt)max(D():=-F*(A)-G*(AA)。(14)AERIsERn,tER!为了证明算法1的收敛性,我们先给出一些引理:引理对于任意 tR,R,定义:=-(s+ATt),则有 s=F*()和ATt=ATVG*(AX)。证明根据(13a)式中s-子问题的最优性条件,我们有0 E VF(s)-+(s+ATt
27、)=VF(s)-X,因此,我们可以得到XEVF(s)s=VF*()。同理,根据(13b)式中t-子问题的最优性条件,再结合(13c)式我们有0 E VG(t)-A+A(s+ATt)=VG(t)-AX。因此,我们有AXEVG(t)-D(k)求解一类线性等式约束凸优化问题的加速方法-入Ts+21ERI入,入一入7(13a)(13b)(13c)口112(15)8证明根据(13b)式中-子问题的最优性条件,我们有0=VG(t)-A+A(s+ATt)=b-A(-(s+ATt)。首先,我们将(13c)式代入(16)式中,得到A=6,即引理3中的第二个结论。其次,根据(12)式和函数G的定义,我们有A入=b
28、。最后,根据(13c)式以及算法1中t的更新方式,我们有 A=A入+(s+ATts)=A=b。进一步地,根据引理3,我们可以得到II-XII2=2IATt-ATt2=2/AT(VG*(AA)-VG*(AX)I?=0,即有入=入。对于 F(s)=0*(s)而言,我们有F*(n)-F*(X)=0(n)-0(A)(V(X),h-X)(V0(),-X)-/-X2,(17)其中第一个不等式利用了函数的凸性,第二个不等式利用了=入的等式。同理,对于G(s)=6Tt而言,我们可以得到G*(AR)-G*(AX)(VG*(AX),AK-AX)=)-D(k)=F*(r)-F*()+G*(Ak)-G*(AX)孟辛晴
29、,张文星28卷(16)(18)=(s+ATt,K-X)-一(19)2B其中第三个等式利用了引理3,最后一个等式利用了(13c)式。为了便于下面的讨论,我们不妨定义序列(mk=k入k(a k 一1)入k-1入 1%=o,其中(a%)-0,(入)%=。为算法1产生的选代序列。根据算法1中入的更新方式,我们可以很容易地验证(20)进一步地,序列(mk%=。满足如下引理。引理5算法1产生的迭代序列【(h,入,入%)%=。满足Ilmk+1l2-mill2 2a D(A*)-D(入%)-2 a+1 D(A*)-D(入+1)。口mh+1=mk+k+1(入+1-h+1)。(21)1期证明根据序列【mk%。的定
30、义以及恒等式(2 0),我们可以得到I/mk+1l2-1/ml/2=2ak+1(mk,a+1-A+1)+x+1/a+1-Aa+112=2h+1(x入-(a%-1)入a-1-*,入a+1-#+1)+/入a+1-#+1l2=2k+1K(as-1)(入-入x-1)+h+1入a-Qh+1Ah+入-入*,入a+1-Xa+1)+元+1/入a+1-k+1l2=2h+1+1+(1-h+1)入-入*,入+1-%+1)+%+1l/入+1-%+1l2=2k+1+1-*,入a+1-x+1)+11入A+1-其中第二、四个等式利用了算法1中入,的更新公式。现在,我们令(15)式中的:=入,入:=入k+1,入:=入k+1,
31、则有D(入k+1)-D(入k)另一方面,若令(15)式中的:=入*,入:=入+1,:=入+1,则有D(入h+1)-D(*)将(2 3)(2 4)式代入(2 2)式中,我们可以得到求解一类线性等式约束凸优化问题的加速方法(#+1 入n,入#+1-#+1)+11入a+1-Ah+1l29(22)(23)(24)2k+1(k+1-1)D(入k+1)-D(入k)+2k+1D(入k+1)-D(*)=2%+1D(A%+1)-2alk+1(h+1-1)D(A%)-2h+1D(A*)=2ai+1D(入x+1)-2+1D(Ax)+2ah+1D(入)-2 ak+1D(A)=2aa+1D(Ak+1)+2(ah+1-%
32、+1)D()-28(a+1-x)D(A*)=2ah+1D(Ak+1)-2a%,D(Ak)-2(ai+1-Q%)D(A*)=2 D(A*)-D(入)+2 ai+1 D(Ak+1)-D(*),其中第三个等式利用了算法1中步长 x+1 的更新公式,即 x+1=V,价形式i=0+k面基于上面的引理,下面我们将对加速交替方向乘子法的收敛性进行分析。首先,我们可以得到如下定理。定理2 算法1产生的选代序列【入%=。满足D(*)-D(入k),或者其等2口21/入1-入*1I2(25)(k+1)210证明根据引理5,不等式(2 1)可以等价为Ilmh+1I2+2%+1 D(A*)-D(入k+1)Ilml2+2
33、azD(入*)-D(入%),Vk 1。(2 6)这意味着1ml12+2D(入*)D(入)%=。是一个递减数列。因此,我们可以得到Ilmal/2+2az D(*)-D(入x)I/m1 I2+2a D(A*)-D(入x)。根据(2 4)式,我们有(即令(2 4)式中的k=0)D(入1)-D(*)不等式(2 6)可以变形为2a+1 D(A)-D(A+1)Iml/-mk+1l2+2az D(A*)-D(x)孟辛晴,张文星工211 11+(1 *,1 1)(1-*/11-*1)。2 I/ml/2+2a%D(A*)-D(入z)28卷(27)I/m1/2+2i D(A*)-D(1)=II入1-入*12 +2
34、 i D(*)-D(入1)I入1-*1,其中第一个等式利用了数列ms的定义,第三个不等式利用了数列(mal/2+2%(D(入*)-D(入)%=。的单调递减性,最后一个不等式利用了(2 7)式以及初值1=1这一事实。由算法1中步长k+1的更新公式,我们有k+1+ak2进一步地,根据不等式(2 8)(2 9),我们最终可以得到D(*)-D(入k+1)定理3算法1产生的迭代序列的原始-对偶残差满足 Ilrel12O(1/k2)和Idel12=0,其中rk和dk在(12)式下面给出了定义。证明在(15)式中,我们令:=入*,=入:=入,可以得到D(*)D(入)1入一入*12。根据定理2 中关于D(入)
35、的收敛性证明,我们有4入1-入*I2II1入%-入*12 2 D(*)-D(入x)(k;+1)2,其中第二个不等号利用了(2 5)式,这表明11入一入*12 O(1/k2)。进一步地,我们还可以得到I1入+1-入 12 2 (1入h+1-入*12 +I1*-入/12)0(1/k2)。(28)k-1k+1+1,Vk1。221/1-入*I2 211入1-*220元+1(k:+2)2(29)122/入1-入*I2(k+1)2口1期另一方面,由于 D()=-F*()-G*(A),其中 VF*=的利普希茨常数为 Le,A 是行满秩矩阵,所以D()也是一个具有梯度利普希茨连续性质的函数。在这里,我们不妨记
36、VD的利普希茨常数为L。根据引理1,我们有D(Ax)D(A)+(VD(A*),x-X)-1A-*I2=D(A)-I/x-I1,Vh 0,其中第二个不等式利用了最优性条件D()*)=0。利用(30)式在k+1时的情形,我们可以得到D(入k+1)-D(入h+1)D(入k+1)-D(A*)+求解一类线性等式约束凸优化问题的加速方法1入;+1-入*1211(30)21入+1-入*I2L2一*+%二k+1L11入-*I2+(ak-1)2Q元+10(1/k2),其中第一个等式利用了算法1中入h+1的更新公式,第二个不等式利用了D(入+1)D(*)(即入*是极大值点)。进一步地,在(15)式中取:=入+1,
37、入:=入+1,入:=入+1,则有D(入+1)-D(Ak+1)其中第一个等式利用了算法1中入+1更新公式,第二个等式利用了r的定义。最后,由(31)式可以得到 IIrh+12 O(1/k2)。对于残差dk而言,利用(10 c)式和引理4,则有IIdiI2=IIBAT(t-s)=IIBATVG*(A)-VG*(A,)I=0。2211入+1-入1+1/21(s+1+AT ta+1)2=/r+11,(31)3数值实验及结果工程应用领域中的很多实际问题都可以建模为形如问题(1)的优化模型,如信号处理中的压缩感知问题、图像处理中的图像重建以及恢复问题等。在本节中,我们将经典交替方向乘子法(ADMM)、加速
38、的交替方向乘子法(A-ADMM)以及加速的线性化交替方向乘子法(A-LADMM)应用于求解图像重建问题,对比它们在不同参数条件下的数值效果。本节数值实验中的所有的程序代码用MatlabR2020a编写,并在机器配置为Intel(R)Core(TM)i5CPU1.60GHz的个人笔记本电脑上完成。在压缩感知和图像处理中,问题(1)的一个重要应用是基追踪问题(Basis Pursuit)29。我们将考虑如下基于小波的图像重建问题(32)min Ill,s.t.SHWa=b,12其中bER表示观测到的图像(可能存在信息缺失、光照不均等现象);;SERIxm表示图像的采样矩阵或者信息缺失矩阵,根据矩阵
39、S的特征来刻画不同的图像处理问题(例如,如果S是一个对角矩阵且对角元素为0 或1,则模型(32)对应的是图像填补问题;如果S是一个对角矩阵且主对角线元素在区间(0,1)中,则模型(32)对应的是图像增强问题;如果S是由单位矩阵的某些行构成的行满秩矩阵,则模型(32)对应的是图像放缩问题;如果S是单位矩阵,则模型(32)对应的是图像去模糊问题);HERmxm是卷积核的矩阵表示(亦称卷积矩阵或模糊算子)。在二维图像的离散化过程中,如果采用了特殊的边界条件(例如,周期、镜像边界条件等),那么H是一个具有特殊结构的分块矩阵,我们可以通过快速傅里叶变换(FFT)或离散余弦变换((DCT)将它进行对角化;
40、WE Rmxn是一个小波变换的字典矩阵,表示图像在小波变换W下的系数。一般地,一幅图像在小波变换W下的系数是稀疏的,因此我们用 el-范数或光滑化的el-范数(如Huber范数)来获取稀疏的小波系数a。在模型(32)中,我们考虑用Huber范数来对小波系数进行正则化,更多关于正则化函数的内容,感兴趣的读者可以参阅文献 31。具体地,给定常数0,向量 E Rn 的Huber 范数定义为(33)-1其中 h:IR R 是 Huber函数,定义为 32(t2h(t):=2M2”根据文献 33 可知,一方面 Huber 函数是可微函数,而且满足 t一h(t)和h(t),也就是说Huber 范数是 e1
41、-范数的光滑逼近而且逼近误差为/2。另一方面,由于Huber范数可以表示为一个强凸的优化问题孟辛晴,张文星nIl1,=Z h(ci),如果围0是保证算法理论上收敛的必要条件,在数值实验中,罚参数对问题求解速度的影响往往很大。图2 中展示了不同罚参数对加速交*设C表示离散余弦变换(DCT),则-1表示它的逆变换。矩阵H在镜像边界条件下可以借助DCT变换进行对角化,即H=C-1C,其中是一个对角矩阵,对角元素为矩阵H的特征值 30 ;同理,在周期边界条件下,矩阵H可以通过傅里叶变换进行对角化。1期求解一类线性等式约束凸优化问题的加速方法13真实Peppers图像像素缺失2 0%的图像像素缺失40%
42、的图像像素缺失6 0%的图像真实House图像10010-2像素缺失2 0%的图像像素缺失40%的图像像素缺失6 0%的图像图1真实图像以及观测图像(有信息缺失)100-=0.1-8-=1=10-*=100-0=0.1-8-=1=10-*=10010-0图2 加速交替方向乘子方法取不同罚参数时的原始-对偶残差(rk,d s)随迭代变化的曲线图替方向法的原始-对偶残差 I(rk,ds)l关于迭代次数变化的曲线图。以加速交替方向乘子法为例,为了测试不同的罚参数的值对算法数值性能的影响,我们通过试误的方式选取了不同的罚参数。最终,我们选取 A-ADMM 的罚参数=0.4。由于Goldstein等人
43、18 的算法是针对普通的可分凸优化问题,并没有考虑到本文中优化问题的特殊结构,所以Goldstein等人 18 的工作中的参数选取范围比较局限。倘若按照Goldstein等人 18 的参数范围进行选取,数值效果要比按本文中的参数选取时的加速方法的效果差。因此,在这里我们分别选取经典的ADMM和A-LADMM的罚参数为=10以及=0.15。在数值实验中的各个算法,我们均采用了相同的初值,即用观测到的图像作为各个算法的初始向量。特别地,数值算法中跟图像变量无关的初始向量(如初始拉格朗日乘子等)全部取为零向量。在信息科学中,信噪比(SNR)通常作为衡量图像/信号质量的一个重要指标,其定义为其中u表示
44、真实图像,表示修复/重建得到的图像。我们分别例举了两种情况,即图1中的Peppers和House两幅观测图像在2 0%、40%像素缺失的时候,它们初始信噪比的2040选代次数左:Peppers图像的曲线,右:House图像的曲线SNR:=20 1g60801000ullTi-ull2040选代次数608010014值分别为6.8 2 dB以及3.90 dB。在图3中,我们绘制了问题(32)的目标函数值和恢复图像的信噪比随迭代次数变化的曲线,图4和5分别展示了基于不同比例信息缺失的观测X1062.2F2.01.81.61.41.202.4X1062.22.01.21.00孟辛晴,张文星30-AD
45、MM-A-ADMM+A-LADMM50100送代次数ADMM-A-ADMMA-LADMM50100选代次数图3问题(32)的目标函数值和图像信噪比值随迭代次数的变化曲线上:Peppers图像的曲线,下:House图像的曲线像素缺失图像ADMM方法28卷25P/105150200150200-ADMM-A-ADMM+A-LADMM25003530P/1052500A-LADMM方法5050100送代次数-ADMM-A-ADMM+A-LADMM100150送代次数A-ADMM方法150200200250250图4求解模型(32)恢复出的Peppers图像1期求解一类线性等式约束凸优化问题的加速方法
46、像素缺失图像ADMM方法15A-LADMM方法A-ADMM方法图5求解模型(32)恢复出的House图像图像而重建出的Peppers和House图像。从曲线图中的数据可以看到,在保证图像质量的情况下,加速的交替方向乘子法对于解决特定的基于小波的图像处理问题(1)是有效的,它能较快地达到较高分贝的信噪比。4结语本文研究了求解一类具有线性等式约束的凸优化问题的加速算法,借助对偶理论,将此类问题转化为 Fenchel对偶形式。由于其对偶问题是一个具有可分离结构的凸优化问题(目标函数为共轭函数与线性函数之和),结合加速交替方向乘子法的最新进展,我们在Goldstein 等学者的加速交替方向乘子方法的基
47、础上,通过对目标函数以及罚参数条件的弱化,证明了加速交替方向乘子法在弱化的条件下仍能保证算法的收敛性。最后,我们将算法应用于一个图像重建问题,验证了算法的有效性。我们接下来将进一步研究加速交替方向乘子法在求解某些特殊的多块可分离结构凸优化问题时的收敛性。参考文献1 Bruckstein A,Donoho D,Elad M.From sparse solutions of systems of equations to sparsemodeling of signals and images J.SIAM Review,2009,51:34-81.2 Chan T,Shen J.Image Pr
48、ocessing and Analysis:Variational,PDE,Wavelet,and StochasticMethods M.Philadelphia:SIAM,2005.3 Hastie T,Tibshirani R,Friedman J.The Elements of Statistical Learning:Data Mining,Inference,and Prediction M.New York:Springer,2008.164 Wets R.Stochastic Systems:Modeling,Identification and Optimization M.
49、Amsterdam:TheMathematical Programming Study,1976.5 Glowinski R,Marrocco A.Analyse numerique du champ magnetique dun alternateur parelements finis et sur-relaxation ponctuelle non lineaire J.Computer Methods in AppliedMechanics and Engineering,1974,3(1):55-85.6 Gabay D,Mercier B.A dual algorithm for
50、the solution of nonlinear variational problems viafinite element approximation J.Computers&Mathematics with Applications,1976,2(1):17-40.7 He B,Liao L,Han D,et al.A new inexact alternating directions method for monotone varia-tional inequalities J.Mathematical Programming,2002,92(1):103-118.8 Eckste