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

(東莞理工學(xué)院 計(jì)算機(jī)科學(xué)與技術(shù)學(xué)院，廣東 東莞 523808)

0 引言

本文將考慮以下模型：

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

(1)

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

(2)

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

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

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

1 預(yù)備知識(shí)

τ>1時(shí)，問題(1)為 MCP 問題；當(dāng)

τ>2時(shí)，問題(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個(gè)n維向量對(duì)形成 Hessian 逆的估計(jì)，并且不需要貯存矩陣，被廣泛應(yīng)用于求解無約束和約束優(yōu)化問題[20-21].本文中采用文獻(xiàn)[20]的方法得到問題(1)的 Hessian 逆的估計(jì)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 新算法

本節(jié)中，將利用文獻(xiàn)[1]的積極集識(shí)別技術(shù)和 L-BFGS 算法，提出一種求解大規(guī)模l1、SCAD 和 MCP 問題的 L-BFGS 算法.設(shè)式(1)中的目標(biāo)函數(shù)滿足以下假設(shè).

假設(shè)2.1

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

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

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

成立.

2.1 搜索方向的計(jì)算及性質(zhì)

(4)

2.2 搜索方向的計(jì)算及性質(zhì)

對(duì)任意k，存在 0

(5)

成立.

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

(6)

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

成立.

引理2.2設(shè)假設(shè) 2.1和假設(shè) 2.2成立，d(k)由 式(4)定義，則有d(k)=0 當(dāng)且僅當(dāng)x(k)為式 (1)的一個(gè)穩(wěn)定點(diǎn).

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

和

成立.式中，σ=0 對(duì)應(yīng)于l1問題和 SCAD 問題，σ=1 對(duì)應(yīng)于 MCP 問題.

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

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

2.3 算法及其收斂性

下面給出新算法的詳細(xì)步驟.

算法2.1 L-BFGS 擬牛頓法

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

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

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

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

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

(7)

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

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

(8)

6.設(shè)k：=k+1.轉(zhuǎn)到步驟1.

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

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

(9)

成立.

證明因?yàn)?0<γ<1，由式(8)和Q(0)=1 得到

(10)

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

(11)

因此{(lán)C(k)} 是單調(diào)遞減的.進(jìn)一步由線性搜索式(7)得到 {x(k)}?Θ.由假設(shè)2.1的1)得到 {C(k)} 有下界，因此 {C(k)} 是收斂的.對(duì)式(10)兩邊求和，利用 {C(k)} 的收斂性得到結(jié)論.

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

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

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

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

考慮以下兩種情況.

(12)

成立.由中值定理可知，當(dāng)k∈K足夠大時(shí)，有下式成立，

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

3 數(shù)值實(shí)驗(yàn)

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

將算法應(yīng)用于壓縮感知式 (2).其目標(biāo)是從m個(gè)觀測數(shù)據(jù)中，重構(gòu)出一個(gè)長度為n的稀疏信號(hào)(m

其中ε=‖e‖.

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

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

4 總結(jié)

本文基于文獻(xiàn)[1]的積極集識(shí)別技術(shù)和 L-BFGS 算法，提出了一種用于求大規(guī)模l1、SCAD 和 MCP 問題的新方法.在適當(dāng)條件下，證明了新方法在采用非單調(diào)線性搜索時(shí)的全局收斂性.數(shù)值實(shí)驗(yàn)證明了新方法對(duì)于求解大規(guī)模l1、SCAD 和 MCP 問題具有優(yōu)秀的性能.