一種慣性交替極小化算法及其應用

吳雙雙，唐玉超

(南昌大學數學系，江西 南昌 330031)

1 引言與研究背景

在本文中，我們主要考慮以下形式的凸極小化問題：

(1)

其中H1和H2是實Hilbert空間，f:H1→(-∞,+∞]和g:H2→(-∞,+∞]是正則下半連續(xù)凸函數，G是實Hilbert空間，A:H1→G和B:H2→G是非零有界線性算子，以及b∈G。我們假設f是σ-強凸的，其中σ>0。下面，我們介紹幾個信號和圖像處理中的具體問題，可以看作是問題(1)的特例。

例1(全變分圖像去噪問題)全變分圖像去噪問題最早由Rudin、Osher和Fatemi[1]提出，其形式如下：

(2)

其中u∈Rm×n是觀測到的噪聲圖像，‖x‖TV表示全變分，λ>0是正則參數。(2)通常被稱為ROF模型。因為模型中的全變分‖x‖TV可以表示為凸函數φ和一階差分算子L:Rm×n→Rm×m×Rm×n的組合，即‖x‖TV=φ(Lx)[2-3]。因此，(2)可以表示為下面的形式：

(3)

如果令y=Lx，則(3) 可以寫成如下形式：

s.t.Lx-y=0

(4)

例2(具有約束的全變分圖像去噪問題) 考慮到數字圖像中像素值往往分布在[0,1]或[0,255]之間。設C是一個非空閉凸集，考慮如下具有約束的ROF模型[4-6]：

s.t.x∈C

(5)

通過使用示性函數δC(x)，以及‖x‖TV=φ(Lx)， (5)等價于：

(6)

s.t.Lx-y=0

(7)

這也是(1)的一個特例。除此之外，還有許多問題是(1)的特例，例如閉凸集交集上的投影問題[7-8]，強凸魯棒主成分分析(RPCA)問題[9-10]和資源分配問題[11]等。

交替方向乘子法(ADMM)是求解(1)的一種常用方法，最早可以追溯到Gabay和Mercier[12]、Glowinski和Marroco[13]以及Gabay[14]等人的工作。眾所周知，ADMM等價于Douglas-Rachford分裂算法[15]應用于(1)的對偶問題。為了充分利用(1)中f(x)的強凸性，基于鄰近梯度算法(也稱為向前向后分裂算法)，Tseng[16]提出了一種交替極小化算法(AMA)求解問題(1)的對偶問題，具體算法如下：

(8)

為了加速AMA(8)，Goldstein等[17]基于快速迭代收縮閾值算法 (FISTA)[18]提出了一種求解(1)的對偶問題的算法，即Fast AMA：

(9)

由于Fast AMA (9) 是由FISTA推導建立的，其步長參數范圍比原始AMA (8)小，而文獻[17]并沒有證明序列{xk}和{yk}的收斂性。因此，本文的主要目的是提出一種用于求解(1)的慣性交替極小化算法，該算法源于慣性鄰近梯度算法。我們不僅分析所提出算法的收斂性，而且證明所提出的算法收斂到原始問題(1)和對偶問題的最優(yōu)解。作為應用，我們將所得算法用于求解一類復合凸極小化問題。為驗證所提出算法的有效性和優(yōu)越性，我們考慮具有約束的全變分圖像去噪模型，該模型包含非空閉凸集約束和全變分正則項。

本文的主要貢獻：(ⅰ)提出一種慣性的交替極小化算法求解凸極小化問題(1)，其中AMA(8)可以看作本文所提算法的特例。(ⅱ)證明所提出算法同時收斂到問題(1)和其對偶問題的最優(yōu)解。(ⅲ)應用所提算法求解具有約束的全變分圖像去噪模型，數值結果表明了所提算法的優(yōu)越性。

本文具體安排如下：在第2節(jié)，我們回顧凸分析中的一些基本的定義和重要結論。在第3節(jié)，我們介紹所提算法并證明其收斂性。在第4節(jié)，我們利用提出的算法來解決一類復合凸優(yōu)化問題。在第 5節(jié)，通過數值實驗說明所提出算法的性能。最后，我們給出全文小結。

2 預備知識

C上的投影用PC表示，定義為PC(x)=argminy∈C‖x-y‖。

定義1[20]設f:H→(-∞,+∞]，

(ⅰ)如果

f(λx+(1-λ)y)≤λf(x)+(1-λ)f(y),?x,y∈domf,λ∈[0,1]

那么f是凸的。

f在x處的次微分定義為:

?f(x)≡{v∈H|f(y)≥f(x)+〈v,y-x〉,?y∈H}

注意f是凸的就意味著?f是單調的，即

〈x-y,u-v〉≥0,?x,y∈H,?u∈?f(x),?v∈?f(y)

此外，如果f是σ-強凸的，則?f是σ-強單調的，即?f-σI是單調的。

定義2[20]設f∈Γ0(H)，則f的Fenchel共軛定義為

由f∈Γ0(H)，可以得到(?f)-1=?f*。

定義3[20](鄰近算子)設λ>0，f∈Γ0(H)，λf的鄰近算子定義為：

如果f(x)=δC(x)，那么proxλδC(x)=PC(x)。

Moreau[21]首次提出了鄰近算子的概念，它在設計和解決凸優(yōu)化問題的鄰近算法中起著非常重要的作用。下面的引理展示了如何利用一個正則下半連續(xù)凸函數f的凸共軛f*或其逆來計算它本身的鄰近算子。

引理1[20](Moreau's decomposition)設λ>0，f∈Γ0(H)，則

定義4[20](Moreau envelope)設f∈Γ0(H)，λ>0，則Moreau envelope可以表示為：

下面的引理表明Moreau envelope在H上是可微的。

2.1 慣性鄰近梯度算法

考慮下面凸極小化問題：

(10)

其中f:H→R是具有L-Lipschitz連續(xù)梯度的凸可微函數，L>0，g:H→(-∞,+∞]是一個正則的、下半連續(xù)凸函數。鄰近梯度算法 (PGA) 是求解(10)的一種常用的方法。PGA 可以看作是在凸極小化問題中向前-向后分裂算法[22-24]的一種特殊情況。此外，慣性鄰近梯度算法是一種有效的加速鄰近梯度算法被廣泛應用于求解問題(10)。下面，我們給出關于慣性鄰近梯度算法的收斂性定理，具體證明可以參考文獻[25]。

定理1設g∈Γ0(H)，f:H→R是具有L-Lipschitz連續(xù)梯度的凸可微函數，其中L>0。假設 argmin(f+g)≠?，x0∈H，定義

(11)

如果{αk}和{γk}滿足下列條件之一：

則以下結論成立：

(b){xk}弱收斂到argmin(f+g)的一個極小點。

2.2 最優(yōu)性條件

對于前面我們所考慮的凸極小化問題(1)，它的拉格朗日函數定義為

L(x,y,p)=f(x)+g(y)+〈p,b-Ax-By〉

(12)

我們通過極小化關于x和y的拉格朗日函數(12)來獲得對偶問題，如下所示：

(13)

其中f*和g*分別是f和g的Fenchel共軛。如果

L(x*,y*,p)≤L(x*,y*,p*)≤L(x,y,p*),?(x,y,p)∈H×G×K

則稱(x*,y*,p*)是拉格朗日函數L的鞍點。眾所周知，當且僅當 (x*,y*)是(1)的最優(yōu)解，p*是(13)的最優(yōu)解時，(x*,y*,p*)是拉格朗日函數L的鞍點，并且原始問題(1) 的最優(yōu)值和對偶問題(13)的最優(yōu)值是相同的。(x*,y*,p*) 通常稱為(1)和(13) 的原始-對偶最優(yōu)解。

引理3[11](x*,y*,p*)是拉格朗日函數L的鞍點，當且僅當

A*p*∈?f(x*)

B*p*∈?g(y*)

A*x*+B*y*=b

(14)

最后，在本文中我們將充分利用以下余弦等式：

2〈a-b,c-a〉=‖b-c‖2-‖a-b‖2-‖a-c‖2,?a,b,c∈H

(15)

3 慣性交替極小化算法收斂性分析

在本節(jié)中，我們提出一種用于求解(1)的慣性交替極小化算法，它可以看作是原始AMA算法(8)的推廣。下面我們給出慣性交替極小化算法的具體格式。

算法1慣性的交替極小化算法

當αk=0時,所提出的IAMA退化成原始的AMA(8)。下面，我們證明所提的慣性交替極小化算法的收斂性。

定理2考慮問題(1)，假設拉格朗日函數(12)的鞍點集是非空的，{xk},{yk}和{pk} 是由算法1生成的序列，如果{αk}和{γk}滿足下列條件之一：

(16)

(17)

則存在拉格朗日函數L的一組鞍點(x*,y*,p*)，使得：

(a)pk?p*；

(b)xk→x*；

(c)Byk→By*；此外，如果存在β>0，使得B*BβI，則yk→y*；

證明(a)由引理4可知，f*(A*p) 是連續(xù)可微的，并且具有Lipschitz連續(xù)梯度。在對偶問題(13)中，極小化的是一個光滑函數和一個凸函數的和，因此我們可以應用慣性鄰近梯度算法來求解(13)，即：

(18)

(19)

即

(20)

將h(p)代入(18)的第二個式子，有

(21)

由(21)可得到其一階優(yōu)化條件為

(22)

設yk+1∈?g*(B*pk+1),我們有

0∈?g(yk+1)-B*pk+1

(23)

以及

(24)

這是我們所提出算法的最后一步。由(23)和(24)，我們可以得到

(25)

這是慣性交替極小化算法的第三步。根據定理 1，我們可以得到pk?p*。

(b),(c)由一階優(yōu)化條件和算法可得到

(26)

和

(27)

根據?f強單調和?g單調，則有

(28)

和

(29)

將(28)和(29)相加，可以得到

(30)

根據pk+1的定義，以及Ax*+By*=b，可以得到

(31)

將(31)代入(30)，可以得到

(32)

對(32)利用余弦等式，我們有

(33)

接著在上述不等式兩邊乘以2γk，經過移項和整理可以得到

(34)

(35)

設φk=‖pk-p*‖2，上式可以重新改寫為

(36)

將(36)重新整理，可以得到

(37)

將(34)再次進行移項，可以得到

(38)

對(38)兩邊同時進行求和，則

(39)

在(39)兩邊取極限，即當k→+∞時，我們能夠得到

即

xk→x*,Byk→By*=b-Ax*

如果β>0，使得B*BβI，那么yk→y*。

(d)下面我們證明目標函數值f(xk+1)+g(yk+1)收斂于(1)的最優(yōu)解。由(12)所定義的拉格朗日函數L我們可以得到：

L(xk+1,yk+1,p*)≥L(x*,y*,p*)

(40)

把它展開來，即：

f(xk+1)+g(yk+1)≥f(x*)+g(y*)+〈p*,Axk+1+Byk+1-b〉

(41)

(42)

由算法1的一階優(yōu)化條件我們可以得到

(43)

由次微分的定義，可以得到：

(44)

把上面兩個式子相加得到f(x*)+g(y*)≥f(xk+1)+g(yk+1)+〈A*(pk+αk(pk-pk-1)),x*-xk+1〉+〈pk+1,By*-Byk+1〉

(45)

(46)

最后，結合(42)和(45)，可以得出

(47)

注1在定理 2中，我們證明了慣性交替極小化算法在兩種條件下的收斂性。

(ⅱ)在第二種情況下，為了方便，我們固定步長參數γk，記為γ，根據文獻[25]中的定理 2， (17)可以寫為

因為序列{αk}是非遞減的，所以0≤αk≤α(γ)，其中

α(γ)=1+

(48)

在本文中，我們取ε=10-6。

4 應用

在本節(jié)中，我們應用所提出的算法1來求解一類復合凸極小化問題，此類問題在信號和圖像處理中具有廣泛的應用。具體來說，我們考慮如下形式的凸極小化問題。

問題4.1假設g∈Γ0(H)、h∈Γ0(G)、b∈H，并且L:H→G是一個非零有界線性算子。考慮如下凸極小化問題：

(49)

我們對所考慮的凸極小化問題 4.1做以下假設：(i)，最優(yōu)解存在；(ii)，滿足一定的正則性條件，例如0∈sqri(L(domg)-domh)。 容易看出問題 4.1 是問題(1) 的一個特例。下面，我們給出具體的收斂性定理。

定理3在問題4.1的設定中，設p0∈G，定義

(50)

其中{γk}和{αk}滿足定理2的條件。則下列說法成立:

(ⅰ){pk}弱收斂到(51)的對偶問題的一個最優(yōu)解p*，其對偶問題為：

(51)

(ⅱ)x*=proxg(u+L*p*)是(51)的唯一最優(yōu)解，并且{xk}強收斂到x*。

s.t.Lx-y=0

(52)

易見(52)是(1)的一種特殊情況。因此，我們可以應用算法1來求解(52)。從而我們得到如下迭代格式的算法：

(53)

由xk+1的一階優(yōu)化條件，可以得到

(54)

再由(53)中yk+1的定義和引理1，可以得到

yk+1=

(55)

將(55)代入到(53)的pk+1中，可以得到

(56)

即有(50)。

(ⅰ)根據定理2，我們知道{xk}強收斂到x*,{(yk,pk)} 弱收斂到 (y*,p*)，其中(x*,y*,p*) 是下面拉格朗日函數的鞍點：

L(x,y,p)=f(x)+h(y)+〈p,-Lx+y〉

(57)

由引理3，我們得到

(58)

根據(13)我們可以得到(52)的對偶問題是：

(59)

(ⅱ)由(58) 知，x*=proxg(u+L*p*)是(49)的唯一解，且{xk}強收斂到x*。

5 數值實驗

在本節(jié)中，我們應用所提算法求解具有約束的全變分圖像去噪模型(5)，并與其他算法進行比較。所有實驗均在Windows 7和MATLAB R2016a下進行，該筆記本電腦的配置為 Intel Core i7-6700 CPU @3.40GHZ和8G RAM。在下文中，我們將具有約束的全變分圖像去噪模型(5)稱為C-L2-TV。通常存在兩種類型的全變分：一種是各向同性全變分 (ITV)，另一種是各向異性全變分 (ATV)。它們都可以用凸函數φ和一階差分算子L的組合來表示，即‖x‖TV=φ(Lx)。本文所提出的算法可應用于 ITV和 ATV。在下面的實驗中，我們只使用ATV作為(5)中的正則函數，類似地可以處理ITV的情況。我們定義約束集合C={x|0≤x≤255}。

我們用峰值信噪比 (PSNR) 和結構相似性指數 (SSIM) [26]來衡量恢復圖像的質量，其定義如下：

和

實驗測試圖像包括“Cameraman” 、 “Building”、“Boat”和“Goldhill”，如圖1所示。

(a) 256×256 “Cameraman”圖像(b) 493×517 “Building”圖像 (c) 512×512 “Boat”圖像 (d) 512×512 “Goldhill”圖像

表1 不同圖像在不同噪聲水平下的最優(yōu)正則參數選擇Table 1 Optimal regularization parameter selection for different images under different noise levels

我們同時考慮了兩種不同慣性參數條件下的IAMA，分別稱為IAMA_1和IAMA_2。對原始圖像分別添加均值為0和標準方差分別為σ=15,25和50的高斯噪聲。為了獲得最佳的去噪圖像，我們調整了正則化參數，所得數值結果見表2。

在第一個實驗中，我們比較本文提出的慣性AMA算法與原始的AMA算法在不同步長參數下的恢復效果。對于模型(5)，σ=1，‖L‖2=8，由注1，我們有

為了直觀比較，我們給出α(γ)的圖像，見圖2。

圖2 由γ決定的慣性參數αk的上界Figure 2 The upper bound of the inertial parameter αk determined by γ

接下來我們設置了四種不同的步長參數來進行試驗，γk分別選取：γk=0.1,γk=0.15,γk=0.2,γk=0.225。在IAMA_2中，根據注1我們取：αk=0.26,αk=0.207,αk=0.132,αk=0.078。這里的測試圖像為“Building”，我們對原始圖像添加均值為0 和標準方差為σ=50的高斯噪聲。

所得數值結果如表2所示，由表中數據我們可以看出來，在不同的步長參數下，iAMA_1和原始的AMA在峰值信噪比 (PSNR)、結構相似性指數 (SSIM) 以及迭代步數的數值結果都是一致的。在給定停止條件下，隨著γk在取值的增大，3種算法所需迭代次數逐漸減少。同時，IAMA_2算法比原始的AMA算法收斂稍快。

表2 慣性AMA和原始AMA在不同參數關于PSNR、SSIM 和 Iter的數值結果Table 2 Numerical results of the proposed inertial AMA and the AMA in terms of the PSNR,SSIM and Iter with different parameters

在第2個實驗中，我們測試了慣性AMA和原始AMA在對于不同圖像且在不同噪聲下的恢復效果。我們統(tǒng)一設置γk=0.1。在IAMA_1中，我們按照注1對慣性參數來進行更新。在IAMA_2中，我們設定慣性參數αk=0.26。

表3展示了4幅圖像在不同噪聲水平下去噪效果的數值結果。由表中數據可以看出來，IAMA_1和原始的AMA所選取的圖像在不同的噪聲水平下的峰值信噪比(PSNR)、 結構相似性指數(SSIM)以及迭代步數 (Iter) 幾乎都是一樣的。對于所有的測試圖像，在不同噪聲水平下，IAMA_2的迭代步數都是要比原始的AMA迭代步數少。在4幅測試圖像中，噪聲水平比較大的時候，即σ=50時，IAMA_2的恢復效果都要比原始的AMA恢復效果稍好。在“Cameraman”圖像中，對于設定的3種噪聲水平，IAMA_2 的恢復效果都要比原始的AMA恢復效果稍好。總體來說，我們所提出的算法比原始的AMA算法具有一定的優(yōu)越性。接下來，我們展示“Cameraman”圖像在不同噪聲水平下，本文所提算法與原始AMA算法的恢復效果，如圖3所示。

表3 所提算法在兩種不同的慣性參數設置下與AMA的數值結果Table 3 Numerical results of the AMA and the proposed algorithm under two different inertial parameter

(a)σ=15 (b)σ=25 (c)σ=50

(d)AMA (e)AMA (f)AMA

(g)IAMA_1 (h)IAMA_1 (i)IAMA_1

(j)IAMA_2 (k)IAMA_2 (l)IAMA_2

6 結論

為解決具有線性等式約束的兩塊可分離凸極小化問題(1)，我們提出了一種慣性交替極小化算法，所提出的算法推廣了 Tseng[16]提出的經典交替極小化算法。在實Hilbert空間中，我們證明了所提出的算法在兩種不同慣性參數條件下的收斂性。進一步，我們應用所提出的算法求解一類復合凸極小化問題，該問題在圖像處理，特別是圖像去噪中具有廣泛的應用。為了驗證所提出算法的有效性和優(yōu)越性，我們應用所提出的算法求解具有約束的全變分圖像去噪模型。數值實驗表明慣性交替極小化算法比原始交替極小化算法的迭代步數少。