求解非凸正則化問題的 L-BFGS 算法

陳鴻升，葉建豪，胡子健，程萬友

(東莞理工學院 計算機科學與技術學院，廣東 東莞 523808)

0 引言

本文將考慮以下模型：

minφ(x)：=f(x)+r(x)，

(1)

式中，f：Rn→R為連續可微的函數.問題(1)的一個重要應用形式是壓縮感知問題：

(2)

求解l0問題的一般做法是將問題轉化為凸l1問題.常見的用來求解l1問題的算法有：迭代收縮閾值算法(ISTA)[2-3]、內點法[4]、投影梯度法[5]、交替方向乘子法和可分離增廣拉格朗日收縮算法(split augmented Lagrangian shrinkage algorithm，SALSA)[6]、約化空間算法[7]、坐標下降方法[8]、積極集方法[9]等.

為了克服l1罰函數的弊端，人們提出了許多非凸非光滑的罰函數.如 log-sum[10]、cappedl1[11]、lp[12-13]、SCAD[14]、MCP[15]等.這些罰函數被證明比l1罰函數具有更好的稀疏性，且能降低大系數的偏差估計[14，16].人們提出了一些求解這些罰函數的方法，第一類是運用光滑技術的方法，如：光滑序列二次規劃(SQP)和光滑信賴域牛頓方法(TRN)[13]用來求解包含lp、SCAD 和 MCP 問題在內的一類非凸非光滑問題；第二類方法是迭代重加權最小化算法[10，14-15]，這類算法通過求解一系列的l1子問題或l2子問題來實現問題的求解，但這類方法所生成序列的收斂性通常是不明確的；第三類方法基于 CD規劃思想.文獻[17]中提出了一種通用的 ISTA 算法，并在 SCAD、MCP 和 cappedl1罰函數上進行了說明，證明了該算法能收斂到穩定點.更多方法，請參考文獻[18-19].

有限內存的擬牛頓方法(L-BFGS)收斂速度快，而且僅需要貯存最近的m個n維向量對形成Hessian 逆的估計，并且不需要貯存矩陣，被認為是求解大規模優化問題[20-21]的最有效方法之一.本文利用文獻[1]的積極集識別技術和 L-BFGS 算法提出一種求解大規模l1、SCAD 和 MCP 問題的 L-BFGS 算法.在積極集集合上算法的搜索方向與文獻[1]的方向相同，自由空間集合上使用了 L-BFGS 的搜索方向.在適當的條件下，證明了使用非單調技術[22]的算法是全局收斂的.數值實驗證明提出的算法是有效的.

1 預備知識

τ>1時，問題(1)為 MCP 問題；當

τ>2時，問題(1)為 SCAD 問題.

為定義新算法的搜索方向，將B(x) 劃分為4部分：

B11(x)={i：0λ}，B12(x)={i：0≤xi≤λ，xi-gi(x)<-λ}，

B13(x)={i：-λ≤xi<0，xi-gi(x)<-λ}，B14(x)={i：-λ≤xi≤0，xi-gi(x)>λ}.

特別的，

L-BFGS 算法收斂速度快，而且僅需要貯存最近的m個n維向量對形成 Hessian 逆的估計，并且不需要貯存矩陣，被廣泛應用于求解無約束和約束優化問題[20-21].本文中采用文獻[20]的方法得到問題(1)的 Hessian 逆的估計H(k)，即：

(3)

式中：γ(k)>0；s(k)=x(k+1)-x(k)；y(k)=g(x(k+1))-g(x(k))；S(k)=[s(k-m)，…，s(k-1)]；

D(k)=diag[(s(k-m))Ty(k-m)，…，(s(k-1))Ty(k-1)].

2 新算法

本節中，將利用文獻[1]的積極集識別技術和 L-BFGS 算法，提出一種求解大規模l1、SCAD 和 MCP 問題的 L-BFGS 算法.設式(1)中的目標函數滿足以下假設.

假設2.1

1)水平集 Θ：={x∈Rn：φ(x)≤φ(x(0))} 有界.

2)存在Θ 的一個鄰域N，使得f是二階可微的函數，且f的一階導數滿足 Lipschitz 連續性.即存在常數L>0，使得

‖g(x)-g(y)‖≤L‖x-y‖，?x，y∈N

成立.

2.1 搜索方向的計算及性質

(4)

2.2 搜索方向的計算及性質

對任意k，存在 0

(5)

成立.

利用柯西不等式結合式(5)得到

(6)

引理 2.1設假設 2.1和假設 2.2成立，d(k)由式(4)定義，則有

成立.

引理2.2設假設 2.1和假設 2.2成立，d(k)由 式(4)定義，則有d(k)=0 當且僅當x(k)為式 (1)的一個穩定點.

引理 2.3設假設 2.1和假設 2.2成立，d(k)由 式(4)定義，β∈(0，1]，則有

和

成立.式中，σ=0 對應于l1問題和 SCAD 問題，σ=1 對應于 MCP 問題.

下面的定理表明，若x(k)不是式 (1)的一個穩定點，則d(k)為φ(x) 在x(k)處的下降方向.該定理保證了算法 2.1是正確定義的.由于定理的證明相似于文獻[1]的引理 4.1，此處省略.

定理 2.1設假設 2.1和假設 2.2成立，d(k)由 式 (4)定義.若d(k)≠0，則有

2.3 算法及其收斂性

下面給出新算法的詳細步驟.

算法2.1 L-BFGS 擬牛頓法

輸入：x(0)∈Rn、γ，η，δ∈(0，1)，C(0)=φ(x(0))，Q(0)=1，k：=0.

輸出：最終解xk、最終函數值φ(xk)、迭代次數k.

1.若滿足終止條件，則終止算法.

2.通過式(4)計算d(k).

3.計算β(k)=max{ηj：j=0，1，…} 滿足φ(x(k)+β(k)d(k))≤C(k)-δ(β(k)‖d(k)‖)2.

(7)

4.計算x(k+1)=x(k)+β(k)d(k).

5.計算Q(k+1)=γQ(k)+1，C(k+1)=(γQ(k)C(k)+φ(x(k+1)))/Q(k+1).

(8)

6.設k：=k+1.轉到步驟1.

在步驟 3 中，采用非單調線性搜索技術，與文獻[22]的公式 (7)相同，通過公式(8)對C(k)進行更新.與文獻[1]中的算法相比，在d(k)的計算上采用了 L-BFGS算法，這使得文中的算法能夠適用于大規模問題的求解.

引理2.4設假設 2.1和假設 2.2成立，{x(k)} 為算法 2.1產生的序列，則有

(9)

成立.

證明因為 0<γ<1，由式(8)和Q(0)=1 得到

(10)

由式(10)及式(7)和式(8)得到

(11)

因此{C(k)} 是單調遞減的.進一步由線性搜索式(7)得到 {x(k)}?Θ.由假設2.1的1)得到 {C(k)} 有下界，因此 {C(k)} 是收斂的.對式(10)兩邊求和，利用 {C(k)} 的收斂性得到結論.

接下來的引理 2.5 和定理 2.2表明，{x(k)} 的每一個聚點均為問題 (1)的穩定點.由于引理 2.5的證明相似于文獻[1]的引理 4.5，此處省略.

引理 2.5設假設 2.1和假設 2.2成立，d(k)由式(4)定義.若x(k)→x*，d(k)→0，則x*為問題(1)的一個穩定點.

定理2.2設假設 2.1和假設 2.2成立，{x(k)} 為算法 2.1產生的序列，則 {x(k)} 的任一聚點x*均為問題(1)的穩定點.

證明由于 {x(k)}?Θ，且水平集 Θ 有界，因此 {x(k)} 至少存在一個聚點.假設x*為 {x(k)} 的任意一個聚點，則存在一個無限集合K使得

考慮以下兩種情況.

(12)

成立.由中值定理可知，當k∈K足夠大時，有下式成立，

其中，τ(k)∈(0，1).將式(12)帶入最后一個不等式，可得

3 數值實驗

本節中用數值實驗來測試提出的新算法，并將新算法與 GIST 算法[17]和 FPC_AS算法[23]進行比較.GISTA 和FPC_AS的代碼可以找到.分別用 GISTA_SCAD 和 GISTA_MCP 表示代碼 GISTA 的 “SCAD” 和 “MCP” 形式.實驗的代碼通過 MATLAB R2018a 編寫，并在一臺搭載有Windows 10 系統的 PC 上運行 (CPU：i5-6300HQ，內存：16G).所有參數都使用默認值，初始點均為零向量.

將算法應用于壓縮感知式 (2).其目標是從m個觀測數據中，重構出一個長度為n的稀疏信號(m

其中ε=‖e‖.

本節中取n=215，m=round(0.1×n)，xs中非零元素個數為T=1，2，…，50.“CPU time” 表示以秒為單位的 CPU 時間，“err” 表示誤差(err=‖x*-x(k)‖∞).

為了更直觀地展示數值實驗的結果，文中采用文獻[24]中所描述的性能曲線圖，將數值結果中的 “CPU time”“err”分別用圖1～圖4進行展示.從圖1～圖4可以看出，GISTA_SCAD 和 GISTA_MCP 總是需要花費更多的 CPU 時間；GISTA_SCAD、GISTA_MCP 和 FPC_AS 具有相近的誤差.在大多數問題中，LBFGS_L1、LBFGS_SCAD 和 LBFGS_MCP 花費的 CPU 時間更少，且具有更小的誤差.

4 總結

本文基于文獻[1]的積極集識別技術和 L-BFGS 算法，提出了一種用于求大規模l1、SCAD 和 MCP 問題的新方法.在適當條件下，證明了新方法在采用非單調線性搜索時的全局收斂性.數值實驗證明了新方法對于求解大規模l1、SCAD 和 MCP 問題具有優秀的性能.