劉勇進(jìn),湯婉紅

(福州大學(xué)數(shù)學(xué)與統(tǒng)計(jì)學(xué)院,福建 福州 350108)

0 引言

其中,xi∈R,i=1,…,n表示向量x∈Rn中的第i項(xiàng)分量.保序回歸問題在統(tǒng)計(jì)、數(shù)學(xué)規(guī)劃、生產(chǎn)計(jì)劃和庫存控制等領(lǐng)域具有重要的實(shí)際應(yīng)用.具體地,保序回歸問題應(yīng)用于分類器、廣告排序和質(zhì)量控制等等實(shí)際問題.由于保序回歸問題沒有考慮到估計(jì)量所在空間的已知信息.在估計(jì)量中加入先驗(yàn)信息可以帶來更好的估計(jì)性能,一般的方法是在該估計(jì)量上添加一個(gè)線性約束[2],常見的例子包括非負(fù)最小二乘回歸[3],保序回歸[1]等.基于以上背景,考慮在保序回歸問題上添加一個(gè)線性約束,即對(duì)于給定向量b∈Rn和τ>0,考慮以下優(yōu)化問題,即

(1)

其中,en表示分量全為1的n維列向量.可以觀察到問題(1)的最優(yōu)解是b在超平面交單調(diào)錐上的投影.

目前國內(nèi)外的學(xué)者設(shè)計(jì)了各種有效的算法求解保序回歸問題,主要包括: Van Eeden算法[4]、最小下集算法(minimum lower sets,MLS)[5]、池相鄰違反算法(pool adjacent violators,PAV)[4]等.但是求解超平面交單調(diào)錐上的投影問題(1)的算法比較少,所以設(shè)計(jì)求解問題(1)的有效算法是必要的.在文獻(xiàn)[6]中,求解有序加權(quán)l(xiāng)1范數(shù)投影問題[7]可以轉(zhuǎn)變?yōu)榍蠼馊缦聝?yōu)化問題,即

(2)

其中:λ∈Rn,λ中的分量不全為0且滿足λ1≥λ2≥…≥λn≥0和集合C={x∈Rn|x1≥x2≥…≥xn≥0,x≠0}.

Li等[7]通過隨機(jī)數(shù)據(jù)集上的數(shù)值實(shí)驗(yàn)表明了半光滑牛頓算法可以有效地求解問題(2)的對(duì)偶問題.由于問題(1)和問題(2)的形式相似,所以考慮用半光滑牛頓法求解問題(1)的KKT系統(tǒng).此外,在文獻(xiàn)[8]中,PAV算法能夠高效地求解保序回歸問題,所以也考慮用PAV算法求解問題(1)[9].

主要研究在超平面交單調(diào)錐上的投影問題.首先,介紹一些相關(guān)的預(yù)備知識(shí);接著,給出求解該問題的PAV法和半光滑牛頓法并分析算法的有效性;最后,根據(jù)兩種算法的數(shù)值結(jié)果驗(yàn)證: 在求解隨機(jī)數(shù)據(jù)集上的投影問題時(shí),PAV算法比半光滑牛頓算法更高效.

1 PAV算法

1.1 問題描述

Av(S)=(W(Li)Av(Li)+W(Ui)Av(Ui))/W(S) (i∈S)

(3)

其中:σ∈R是問題(1)中等式約束對(duì)應(yīng)的拉格朗日乘子,v∈Rn-1是問題(1)中不等式約束對(duì)應(yīng)的拉格朗日乘子.問題(1)的KKT條件為

有效集J的KKT系統(tǒng)的解可以通過求解有效集J的每個(gè)分塊相應(yīng)的子系統(tǒng)來獲得.定義對(duì)應(yīng)于有效集J的KKT系統(tǒng)的解為擬穩(wěn)定點(diǎn).下面的定理提供了KKT子系統(tǒng)解的顯式公式.

定理1方程組(4)的解為:xi=Av(B),i∈B和vi可以從以下3個(gè)等價(jià)表達(dá)式中的任意一個(gè)獲得,即

下面用數(shù)學(xué)歸納法證明定理的第二個(gè)結(jié)果.當(dāng)i=p時(shí),由(4a)得:vp=W(Lp)(Av(B)-Av(Lp));假設(shè)p≤i-1

vi=xi-bi+σ*+vi-1=(Av(B)-(bi-σ*))+W(Li-1)(Av(B)-Av(Li-1))=W(Li)(Av(B)-Av(Li))

綜上可知,(5a)第一個(gè)等式成立.另外,(5a)第二個(gè)等式和(5b)可由(5a)第一個(gè)等式和式(3)得到.證畢.

由定理1以及原始和對(duì)偶可行性的定義可得下述推論.

推論1令J是一個(gè)有效集,B是有效集J的一個(gè)分塊.擬穩(wěn)定點(diǎn)是原始可行的當(dāng)且僅當(dāng)對(duì)任意滿足B-≠?的分塊B,Av(B-)≥Av(B)成立.擬穩(wěn)定點(diǎn)是對(duì)偶可行的當(dāng)且僅當(dāng)對(duì)?i∈B以下3個(gè)等價(jià)條件每一個(gè)條件都成立:

Av(B)≥Av(Li),Av(Ui)≥Av(B),Av(Ui)≥Av(Li)

(6)

1.2 PAV法的算法框架和有效性

下面給出PAV法求解問題(1)的算法迭代框架.

算法1: PAV算法輸入: 令J={{1}, {2}, …, {n}}, B0={1}.1: if B+=?, then終止算法, 輸出xi=Av(B),i∈B,B?J.2: else3: if Av(B0)≥Av(B+), B0=B+; 轉(zhuǎn)步驟1;4: else5: J=J{B0, B+}∪{B0∪B+}, B0=B0∪B+. 6: while B-≠?&&Av(B-)

下面給出用于建立PAV算法有效性的引理.

接下來介紹PAV算法求解問題(1)的有效性定理.

定理2PAV算法在有限步內(nèi)得到問題(1)的最優(yōu)解,且算法復(fù)雜度為O(n).

證明: 由PAV算法的描述可知,算法終止時(shí)獲得的有效集對(duì)應(yīng)的擬穩(wěn)定是原始可行的.下面用數(shù)學(xué)歸納法證明PAV算法每一次產(chǎn)生的擬穩(wěn)定點(diǎn)均是對(duì)偶可行的.PAV算法第一次執(zhí)行時(shí),擬穩(wěn)定點(diǎn)滿足對(duì)偶可行.假設(shè)第k次執(zhí)行時(shí),擬穩(wěn)定點(diǎn)是對(duì)偶可行的.則第k+1次執(zhí)行時(shí),如果Av(B0)≥Av(B+),根據(jù)假設(shè),此時(shí)擬穩(wěn)定點(diǎn)仍是對(duì)偶可行的;如果Av(B0)

2 半光滑牛頓法

半光滑牛頓算法是目前用于求解中等規(guī)模或大規(guī)模的優(yōu)化問題最有效的求解方法之一,其收斂速度是超線性收斂甚至可達(dá)到二次收斂.本節(jié)主要介紹半光滑牛頓算法求解問題(1)的詳細(xì)過程.

2.1 問題重構(gòu)

令矩陣G∈R(n-1)×n滿足

Gx=[x1-x2,x2-x3,…,xn-1-xn]T(?x∈Rn)

那么,問題(1)等價(jià)于

(7)

問題(7)的拉格朗日函數(shù)為

(8)

由方程組(8)可以得到最優(yōu)拉格朗日乘子σ*為

(9)

求解問題(7)轉(zhuǎn)變?yōu)榍蠼庀率龇蔷€性方程組:

(10)

F(y)=0

(11)

2.2 半光滑牛頓方向的求解

J(y)=I-W(I-GGT)

下面給出關(guān)于雅可比矩陣J(y)的一些性質(zhì).

證明: 令[n-1]={1,2,…,n-1}.由于當(dāng)i∈C0且Wi,i取值為0或1時(shí),可以把指標(biāo)i視作其屬于指標(biāo)集C-或C+的情形考慮,因此當(dāng)i∈C0時(shí),以下證明只考慮Wi,i取值為(0,1)區(qū)間任一數(shù)值的情形.

情形1 當(dāng)C+=[n-1]時(shí),J(y)=I-W(I-GGT)=GGT.由于det(J(y))=det(GGT)=n≠0,所以,雅可比矩陣J(y)是非奇異的.

情形2 當(dāng)C-=[n-1]時(shí),J(y)=I.由于det(J(y))=det(I)=1≠0,所以,雅可比矩陣J(y)是非奇異的.

情形3 當(dāng)C0=[n-1]時(shí),J(y)=I-W(I-GGT)為嚴(yán)格對(duì)角占優(yōu)矩陣,所以,雅可比矩陣J(y)是非奇異的.

情形4 當(dāng)C+=?,C0和C-都不為[n-1]時(shí),J(y)為嚴(yán)格對(duì)角占優(yōu)矩陣,所以,雅可比矩陣J(y)是非奇異的.

情形5 當(dāng)C0=?,C+和C-都不為[n-1]時(shí),令s∈R(n-1),構(gòu)造齊次線性方程組為

(I-W(I-GGT))s=0

(12)

要證明雅可比矩陣J(y)是非奇異的,只需證明方程組(12)的解s為零向量.對(duì)式(12)進(jìn)行初等變換可得

?匌

(13)

情形6 當(dāng)C-=?,C+和C0都不為[n-1]時(shí),雅可比矩陣J(y)為三對(duì)角陣.通過一系列初等行列式變化使得矩陣J(y)上下兩條次對(duì)角線都為-1,即

從最后一列展開得:Dn-1=an-1Dn-2-Dn-3.而D1=a1>0;D2=a1a2-1>0;由于D2-D1=(a2-1)a1-1≥a1-1>0,所以D3=a3D2-D1≥2D2-D1>0;由于D3-D2=(a3D2-D1)-D2=(a3-1)D2-D1≥D2-D1>0,所以D4=a4D3-D2≥2D3-D2>0;由于D4-D3=(a4D3-D2)-D3≥D3-D2>0;所以D5=a5D4-D3≥2D4-D3>0;…;Dn-1=an-1Dn-2-Dn-3≥2Dn-2-Dn-3>0.由上述可知,det(J(y))≠0.所以,雅可比矩陣J(y)是非奇異的.

情形7 當(dāng)C-,C+和C0都不為[n-1]時(shí),由情形6可知,當(dāng)C-=?,C+和C0都不為[n-1]時(shí),雅可比矩陣J(y)的行向量組線性無關(guān).再結(jié)合情形5的證明思路,即可證明命題.證畢.

下面給出當(dāng)i∈C0,Wi,i=0時(shí),雅克比矩陣J(y)的LU分解.

證明: 令L,U∈R(n-1)×(n-1),L是主對(duì)角元全為1的下三角矩陣,U是上三角矩陣,使得

J(y)=I-W(I-GGT)=LU

成立.由于J(y)是三對(duì)角矩陣,通過計(jì)算可知: 矩陣U滿足Ui,j=0,j-i>1和矩陣L滿足Li,j=0,i-j>1.下面給出矩陣L和矩陣U其余元素值.具體地,令p∈Rn-1,若i∈C+,pi=1,否則pi=0.p可劃分成N個(gè)子向量,這些子向量要么是全1向量要么是零向量,并且任意兩個(gè)相鄰子向量是不同的,即

最后,給出矩陣U主對(duì)角線上一行元素的值: 若i∈C+,Ui,i+1=-1,否則Ui,i+1=0,i=1,…,n-2.

由于雅可比矩陣J(y)是非奇異且存在LU分解,從而得出雅可比矩陣J(y)的LU分解唯一的.證畢.

由命題2知,雅可比矩陣J(y)可以被分解成矩陣L和矩陣U的乘積,因此可以利用這兩個(gè)稀疏矩陣通過下面步驟獲得牛頓方向,即:Uz=-F,Ud=z.其中,z∈Rn-1是一個(gè)臨時(shí)變量.經(jīng)過運(yùn)算可得:Uz=-F,Ud=z.其中:z∈Rn-1是一個(gè)臨時(shí)變量.通過計(jì)算可得

z1=-F1,zi=-Fi-Li,i-1zi-1(i=2,…,n-1)

(14)

因此,牛頓迭代方向?yàn)?/p>

(15)

2.3 半光滑牛頓法的算法框架和收斂性

下面給出半光滑牛頓法求解問題(7)的算法框架.

算法2: 半光滑牛頓算法(SSN)輸入: b∈Rn, y0∈Rn-1+, ω∈(0, 0.5), μ∈(0, 1), ε>0, t=0.1: whileF(yt)>ε do2: 通過式(14)～(15)計(jì)算dt;3: 令ξt=ωmt, 其中mt為滿足下面條件的最小非負(fù)整數(shù)m:F(yt+ωmdt)2≤F(yt)2-2μωmF(yt)2.4: yt+1=yt+ξtdt; t=t+1;5: end while;6: 通過式(9)計(jì)算σ?;7: 輸出x=b+σ?en+GTy.

由文獻(xiàn)[12]中定理3.2和文獻(xiàn)[13]中定理2.4可得下述定理成立.

定理3若{yk}是由算法2生成的序列,則{yk}收斂到問題(7)的最優(yōu)解y*,且收斂速度是二階的.

3 數(shù)值實(shí)驗(yàn)與結(jié)果分析

表1 當(dāng)δ=10-3時(shí),PAV,SSN 算法在隨機(jī)數(shù)據(jù)集上的實(shí)驗(yàn)結(jié)果

表2 當(dāng)δ=1時(shí),PAV,SSN 算法在隨機(jī)數(shù)據(jù)集上的實(shí)驗(yàn)結(jié)果

4 結(jié)語

針對(duì)問題(1),充分利用 KKT 最優(yōu)性理論條件,給出求解該問題的PAV法和半光滑牛頓法并對(duì)算法進(jìn)行了有效性分析,最后通過大量的數(shù)值實(shí)驗(yàn)證明了在求解隨機(jī)數(shù)據(jù)集上的投影問題(1)時(shí),PAV算法比半光滑牛頓算法更高效.然而,PAV法不能推廣應(yīng)用于求解問題(2),而半光滑牛頓算法仍可用于求解問題(2).