時閃閃,宇振盛
(上海理工大學 理學院,上海 200093)
本文考慮如下復合函數優化問題

式中:f:Rn→R 是連續可微函數;g:Rn→R 是凸的非光滑函數;可行集Ω={x∈Rn:aiTx=bi,i∈ME,aiTx≥bi,i∈MI},ME={1,···,m},MI={m+1,···,p},ai∈Rn,bi∈R,i=1,···,p。
從現實生活中的問題抽象出來的數學模型,可以分為光滑和非光滑兩大類。對于光滑問題,由于其梯度及Hessian 矩陣易于求出,可以使用經典的優化算法求解。而對非光滑問題,函數是不可微的,因此無法直接使用經典的優化算法求解。此外,非光滑復合函數的極小化問題在信息理論[1]、無線通信[2-4]、圖像恢復[5-7]、信號處理[8-9]、變量選擇[4,10]等領域都有廣泛的應用。
當可行集 Ω=Rn時,問題(1)就是無約束優化問題。目前,已經有很多學者提出了各種非常有效的算法[11-18]解決無約束復合函數問題。比較成熟的算法有迭代收縮閾值算法[12-13]、臨近牛頓法[14]、分裂算法[15-16]、鄰近點算法[17-18]等。
實際生活中的很多問題都要求滿足一些特定的條件,所以研究帶約束的復合函數問題十分有必要。例如,壓縮感知理論中的稀疏信號重構問題就是經典的復合函數優化問題。稀疏重構主要解決這兩類問題:在對非負信號處理時,有非負約束的限制;進行圖像恢復時,灰度強度取值在0~255 之間。在信號恢復、圖像恢復模型里加入約束,可以提高恢復的準確率。同樣地,對于實際問題,考慮約束限制,能更加精確地求出符合實際條件的解,從而節約成本。
當可行集 Ω={l≤x≤u}時,問題(1)就是界約束的復合函數問題。Bian 等[6]提出光滑二次正則化方法,用SQR 表示解決圖像恢復中出現的界約束非Lipschitz 連續優化問題。有效集方法是求解約束優化問題的常用方法,采用有效集策略,可以加快算法的收斂速度。張濤[19]提出了一種求解界約束優化問題的新積極集信賴域算法。
生產生活中的大多數問題,往往具有更加復雜的約束,因此研究復雜約束的復合函數優化問題更具有現實意義。當可行集 Ω為多面體約束時,問題(1)就是多面體約束復合函數優化問題。簡單約束復合函數優化問題的理論已經發展得十分成熟,但是關于多面體約束的復合函數優化問題,研究工作仍處于基礎階段,還有很多問題尚未解決。Liu 等[4]提出了光滑序列二次規劃框架(SSQP)求解一類多面體上復合函數Lq極小化問題。實際上,有效集方法已成功地應用于大規模的線性約束光滑優化問題的求解。William 等[20]提出了求解多面體約束優化問題的新的有效集算法(PASA)。Zhang 等[1]提出了求解線性約束非凸非Lipschitz 優化的光滑有效集算法(SASA),且用于求解簡單約束非光滑復合函數優化問題效果良好,即

也就是本文實際求解的模型:

它是Liu 等[4]在求解一類多面體上復合函數Lq極小化問題中求解的模型。
以上算法在求解大多數問題時表現良好,但它們不能很好地解決一般多面體約束復合函數問題,例如,SQR 算法[6]不能求解問題中f(x)的復合項Lq和一般多面體約束問題。SSQP 算法[4]的目標函數不具一般性,且由于投影收縮算法是不精確求解子問題,這會使得算法收斂速度慢。PASA 算法[20]需要求解投影梯度信息,SASA算法[1]對線性約束優化器有特殊條件限制,這兩種算法都是兩階段算法,程序實現難度較大。
序列二次規劃(sequential quadratic programming,SQP)是求解非線性約束優化問題的主要方法之一,它在求解中等規模和小規模的非線性問題[21-22]中表現不俗。有效集算法 (active set algorithm,ASA)在求解不等式約束優化問題時十分有效,通過有效集策略,每次只需求解僅含等式約束的優化問題,從而可以降低維度,達到提高算法收斂速度的目的。
本文提出了將ASA 和SQP 結合的新算法,該算法主要借鑒文獻[4]中提到的求解多面體上最小化復合函數的光滑近似技巧及SSQP 框架,將非光滑目標函數轉化成光滑函數,并利用有效集策略求解二次規劃子問題。
定義1[23]函數:Rn×R+→R 稱為非光滑函數g:Rn→R 的近似光滑函數,若 ?μ>0,是連續可微的,且對 ?μ∈Rn,有=g(x)。
非光滑優化問題的光滑近似已經得到了廣泛的應用[23-24]。本文對絕對值函數
θ(t)=|t|
構造光滑函數[25]:


則式(5)是問題(1)的光滑近似函數。
不難得到 θμ(t,μ)具有的特殊性質。顯然地,對任意固定的 μ>0,有
θμ(t,μ)=θ(t),?t≥μ
和

此外,θq(t,μ)是連續可微的。除在t=和t=外,θq(t,μ)在其他處都是二次連續可微的。θq(t,μ)關于t的一階導數和二階導數分別為

接下來,給出一個非常重要的引理。
引理1對 ?q∈(0,1),μ∈(0,+∞),

引理1 的詳細證明參見文獻[4]中的引理 4.1。
本文受文獻[4]的啟發,采用該文獻中的SSQP框架,提出了序列有效集算法(SQP-ASA),在每一迭代步中利用有效集策略求解一個凸二次規劃(QP)子問題。此QP 子問題的目標函數是問題(1)光滑近似后式(5)的目標函數的局部上界。光滑參數的更新依賴于問題(1)的殘差。
對固定的光滑參數 μ>0,定義非光滑函數g(·)的光滑近似函數在xk處的二次近似為

類似地,定義光滑函數f(x)在xk處的二次近似為

假設f(x)的梯度是Lipschitz 連續的,即
‖?f(x)-?f(y)‖2≤L‖x-y‖2,?x,y∈Ω
定義

由于Q1(x,xk,μ)和Q2(x,xk)分別是和f(x)在xk處的二次近似,因此,可以把Q1(x,xk,μ)看作是在xk處的二次近似。
下面的引理給出了凸二次規劃(6)的目標函數Q(x,xk,μ)是式(5)的目標函數之間的大小關系,即Q(x,xk,μ)是在xk處的一個上界,這保證了問題(6)的最小值是問題(5)最小值的上界。
引理2對任意的xk,x,使得

該引理的證明可參考文獻[4]中的引理 4.2。
基于引理2,提出序列有效集二次規劃算法。
本文研究的問題(2)光滑后可以表示為

那么,問題(7)在xk處的QP 子問題為

求解問題(7)的算法思想為:在每一迭代步中,首先使用有效集方法求出其二次規劃子問題(8)的最優解d*,并將它作為問題(7)下一步的搜索方向dk,再利用價值函數通過Armijo 線搜索求出步長 αk,進而得到問題(7)的下一步迭代點xk+1,重復這一過程,直到求得問題(7)的最優解x*。具體求解過程如算法1 所示。

令x-xk=d,則可將QP 子問題(8)轉化與d有關的子問題求解。通過積極集策略,會有越來越多的xk逐漸滿足=bi,i∈Sk。這樣問題就轉化為了等式約束優化問題,再利用拉格朗日乘子算法對其進行求解。因此,本文的算法可以在大大減少問題維數的同時提高計算的效率。
接下來,利用文獻[26]中的有效集方法求解QP 子問題,如算法2 所示。


此處的價值函數 ψ(x,μ,?)=。這是因為由有效集方法得到的方向dk一定是可行方向。
為證明算法1 的全局收斂,現給出如下基本假設:
假設1對 ?x∈Ω,線性獨立約束限定 (LICQ)成立,即梯度ai(i∈Sk)線性無關。
由Hk的構造可知,Hk一定是對稱正定的。本文提出的SQP-ASA 算法通過求解如下形式的二次規劃問題

得到解d*,作為原問題變量x在第k次迭代過程中的搜索方向dk。在求解過程中利用有效集策略,通過逐次求解僅含帶等式約束的相關二次規劃問題(9)來逼近其最優解。
下面證明SQP-ASA 算法的全局收斂性,即證該算法產生的點列 {xk}的任何聚點x*都是問題(7)的KT點。
引理3由算法2 中的有效集方法求得的QP子問題方向dk一定是原問題(7)的可行下降方向。


因此,dk一定是原問題(7)的可行方向。
其次,證明方向dk是下降的。對 ?μ∈(0,1],利用QP 子問題(9)的KT條件,有

由式(10)的第一個式子可得

引理4對于上述二次規劃子問題(9),設其KT點對為。如果d(x)=0,則x為問題(7)的KT點。
證明不妨設

如果d(x)=0,依據有效集方法求解方向d(x)的過程,結合問題(7)的KT條件可知,≥0,否則,d(x)不是二次規劃問題(9)的最優解,與KT點對這個條件相矛盾。因此

這說明x是問題(7)的一個KT點。
定理1由序列有效集算法生成的序列{xk}的任何聚點x*都是原問題的KT點。
證明對 ?μ∈(0,1],當該算法是有限步終止時,由引理4 可知,算法產生的點列{xk}的任何聚點x*都是問題(7)的KT點。
下設該算法產生無窮點列{xk},x*為其任一給定聚點。由于 Sk∈M為有限集,不妨設存在無窮子集K,使xk→x*,Sk≡S,k∈K。為不失一般性,不妨設由SQP-ASA 算法產生的迭代點xk+1=xk+βkdk(?k∈K),而由引理3 可知,單調下降,從而由xk→x*,k∈K可以得到,k→∞。
一方面,由于步長 β是由Armijo 線搜索得到的,所以 β∈(0,1]。另一方面,通過引理3 知dk是下降方向,不難得到

故可知dk→0,k∈K。從而-bi=0,i∈S,S?M(x*)。那么

這說明x*是問題(7)的一個KT點。
利用Matlab R2016a 實現算法SQP-ASA。其中試驗環境為:Intel(R)Core(TM)i5-5200U CPU@2.20GHz 4.00GB RAM。通過求解下列問題的稀疏解來檢驗算法的有效性。
例 1非負箱約束稀疏信號恢復的模型如下

尋求該模型的稀疏解,可以轉化為求解模型:

分別運用本文提出的序列有效集算法(SQPASA)和文獻[4]中的序列投影收縮算法(SQPPG),求出該模型的稀疏解。相關參數取法如下:μ0=0.1,σ2=0.000 1,A=randi([-10,,10],50,50),b=randi([-10,10],50,1),L0=,Lmax=3 000,Lmin=1 137.5,精度參數 ?=10-5,初始點d0=254×randi([0,1],50,1),最大迭代次數為1 000。目標函數與迭代次數的關系如圖1 所示。

圖1 迭代次數比較Fig.1 Comparison of iteration times
這兩種方法的迭代次數和迭代時間的比較見表1。從表1 可以看出,對于同樣的初始點,SQPASA 算法比SQP-PG 算法的迭代次數少,運行效率也更高。

表1 例1 的數值結果Tab.1 Numerical results of example 1
例 2尋求如下函數

的稀疏解。可以轉化為求解模型

本文通過求解以下兩個模型在不同約束條件下的稀疏解,來驗證算法的有效性。兩個模型均來自文獻[26]。


分別使用SQP-ASA 算法和SQP-PG 算法,求出該模型的稀疏解。這里的相關參數取法如下:μ0=0.01,ξ=0.5,σ=0.1,η=1.2,Lmax=20,Lmin=3,L0=3,精度參數 ?=10-5,初始d0=[0.5-1]T,最大迭代次數為1 000。迭代過程如圖2 所示。

圖2 迭代次數比較Fig.2 Comparison of iteration times
這兩種方法的迭代次數和迭代時間的比較見表2。從表2 可以看出,對于同樣的初始點,SQP-ASA 算法比SQP-PG 算法的迭代次數少,運行效率也更高。


表2 模型1 的數值結果Tab.2 Numerical results of model 1

分別使用SQP-ASA 算法和SQP-PG 算法,求出該模型的稀疏解。此時的相關參數取法如下:μ0=0.01,ξ=0.5,σ=0.1,η=1.2,Lmax=20,Lmin=2.181,L0=2.618,精度參數 ?=10-5,初始d0=[0-3]T,最大迭代次數為1 000。目標函數與迭代次數的關系如圖3 所示。

圖3 迭代次數比較Fig.3 Comparison of iteration times
這兩種方法的迭代次數和迭代時間的比較見表3。從表3 可以看出,對于同樣的初始點,SQPASA 算法比SQP-PG 算法的迭代次數少,運行效率也更高。

表3 模型2 的數值結果Tab.3 Numerical results of model 2
提出了一種新的求解多面體約束非光滑復合函數優化問題的序列有效集算法。首先,將非光滑目標函數轉化為光滑目標函數;然后,在每一迭代步中用有效集方法求解一個二次規劃子問題得到方向,利用價值函數通過線搜索取得步長,重復這一過程直到求得原問題的解;最后,證明了該算法的全局收斂性。此外,通過實例進行了數值驗證,結果表明,該方法較序列投影收縮方法具有一定優勢。