呂照美，劉新國

(中國海洋大學(xué)數(shù)學(xué)科學(xué)學(xué)院，山東 青島 266100)

1 研究內(nèi)容

本文研究三約束矩陣最小二乘問題

(1)

X∈εpsd∩P∩B

其中

εpsd={X∈Rn×n；XT=X且λmin(X)≥ε},

B={X∈Rn×n；L≤X≤U},

〈A,B〉表示A與B的內(nèi)積,(A,B)= tr(ATB)，tr(A)表示A的跡。記號A≤B表示Aij≤Bij，其中1≤i≤m，1≤j≤n。

對應(yīng)的，

L=L1+L2+…+Lt,

U=U1+U2+…+Ut,

(2)

L，U有上述分裂，其中Li和Ui分別是L和U在Si上的投影，并記

(3)

上述問題出現(xiàn)于統(tǒng)計和數(shù)量經(jīng)濟學(xué)等領(lǐng)域。統(tǒng)計學(xué)中的一個經(jīng)典例子為：求一個到樣本協(xié)方差矩陣最近的對稱正定并且有一定結(jié)構(gòu)的矩陣。這種有特殊結(jié)構(gòu)的協(xié)方差矩陣[1]經(jīng)常出現(xiàn)在物理、生物、心理學(xué)和社會科學(xué)模型。經(jīng)濟學(xué)中的一個例子[2]是獲得效用函數(shù)問題，在這種情形下，擬合矩陣必須有界且對稱，同時它的最小特征值必須大于某一特定的正數(shù)。

由于Dykstra交替投影算法本身收斂較慢且內(nèi)迭代中SPG方法非單調(diào)，使得整個算法收斂更慢。且文獻(xiàn)[12]中設(shè)置內(nèi)外迭代次數(shù)上限均為5 000次，數(shù)值實驗表明，SPG達(dá)到迭代次數(shù)上限時存在沒有收斂或沒有達(dá)到精度要求的情況。文獻(xiàn)[12]中要求A,BT列滿秩并假定可行集非空以保證解的存在性與唯一性，但并未給出判斷可行集是否為空的方法以及在可行集上得到的點是否為最優(yōu)解的條件。

在本文中,我們把Majorize-Minimize(以下簡稱MM)方法與Dykstra算法結(jié)合，提出求解式(1)的一種交替投影迭代方法。MM方法是一種單調(diào)收斂方法，Han[13]利用對偶性推導(dǎo)了Dykstra方法，并分析了其收斂性。

2 理論結(jié)果

2.1 問題的等價形式

令BP=B∩P，則

(4)

其中

(5)

這里，li和ui如式(3)定義，顯然式(4)等價于式(1)。

2.2 可解性

[tf(X)+(1-t)f(Y)]-f(tX+(1-t)Y)=

(6)

定理1令Ω=εpsd∩BP,設(shè)Ω≠φ，那么

(1)f(X)在Ω上有最小點。

證明 顯然，B為有界閉凸集，P為凸集，而由

λmin(tX+(1-t)Y)≤tλmin(X)+(1-t)λmin(Y)，

t∈[0,1]，XT=X，YT=Y，εpsd為閉凸集，如果Ω≠φ，則Ω為有界閉凸集，而f(X)為連續(xù)可微函數(shù)，因此f(X)在Ω上有最小點。

證明畢。

推論1如果Ω≠φ，而A和BT列滿秩，那么式(1)的解存在且唯一[12]。

下面分析Ω≠φ的條件。

記δi=λmin(Wi)，υi=λmax(Wi)則

{yTWiy:y∈Rn，‖y‖2=1}=[δi,υi]。

從而，如果

則

定理2如果

(7)

那么Ω≠φ。

2.3 最優(yōu)性條件

定理3X*∈Ω為最小點的充要條件為

tr{(X-Y)TAT(AX,B-C)BT}≥0,?X∈Ω。

(8)

證明 如果X*為最小值點，則對?t∈[0,1]有

從而(8)式成立，否則取充分小，則得矛盾。

假設(shè)(8)式成立，則?X∈Ω，

證明畢。

注2式(8)不易驗證，下面給出一個充分條件。

推論2如果X*∈Ω，且

tr{(X-X*)AT(AX*B-C)BT}≥0，?X∈Bp，

(9)

則X*是原問題的最小點。

li≤αi≤ui(i=1,…,t)。

(10)

其中

bi=〈AGiB,C〉。

又

一方面，如果式(10)成立，則

(11)

另一方面，如果式(11)成立，則

有式(10)成立。

因此，推論2的等價形式為

(12)

則X*為原問題的最優(yōu)解。

3 算法

MM算法是求解最優(yōu)化問題的重要方法。為了求解問題(1)，設(shè)當(dāng)前迭代點為X(k)，那么

2tr{BT(X-X(k))TAT(AX(k)B-C)}+

(13)

顯然，gk(X(k))=f(X(k))，因此有MM算法：

(14)

tr{BT(X-X(k))TAT(AX(k)B-C)}=

tr{(X-X(k))TAT(AX(k)B-C)BT}

且

2tr {(X-X(k))TAT(AX(k)B-C)BT} =

其中Ck與X無關(guān)，而

(15)

由此可知式(14)等價于

(16)

因此提出算法MM-Dykstra Method如下。

步驟2：計算

步驟3:計算

計算ck+1=‖X(k+1)-X(k)‖，如果ck+1

Boyle和Dykstra證明了Dykstra算法的收斂性[14]。

由式(14)，gk(X(k+1))≤gk(X(k))≡f(X(k))，且gk(X(k+1))≥f(X(k+1))，所以{f(X(k))}單調(diào)遞減。因此，在每次迭代僅需要求解兩個子問題。

其一是

(17)

其中Q是一個正交矩陣，

Λ(1)=diag(μ1,…,μs),

Λ(2)=diag(μs+1,..,μn),

μ1≥…≥μn≥ε>μs+1≥…≥μn,

則式(17)的解為

(18)

因此，

Pεpsd(X)=QΛ(i)QT。

(19)

另一個是

(20)

Z=Z1+…+Zt，

(21)

Zi與Gi有相同的結(jié)構(gòu)，從而有

(22)

(23)

其解為

(24)

則

注3有多種取法。

由

其中，

aij=〈AGiB,AGjB〉，

4 初始點

由于MM算法是一種投影梯度法，一般而言收斂較慢。Lopez和Raydan[16]指出，Dykstra算法有時收斂也很慢，可以使用文獻(xiàn)[10]中的加速技巧，這里我們考慮選取特殊的來實現(xiàn)對算法的加速。

因此可以如下構(gòu)造初始迭代點

(2)取

(Gi是Hadamard積)，?i滿足li≤?i≤ui且

|?i|=min{|αi|:li≤αi≤ui}。

5 數(shù)值實驗

在例1和例2中，L為零矩陣，U為所有元素均為3的矩陣，TOLin=1×10-7，TOLout=1×10-8，MM算法的最大迭代次數(shù)為5 000，P≡Rn×n，且在本文中

即在兩種算法比較的過程中所用的停止準(zhǔn)則一致。

例1 對于文獻(xiàn)[12]中的實驗例子：ε=0，

B=I，

利用文獻(xiàn)[11]中算法4.1得到的結(jié)果與本文算法的結(jié)果對比見表1。

表1 MM_Dykstra算法與SPG_Dykstra算法的時間、誤差、殘差比較

其中在兩種算法中X*表達(dá)相同，均為

例2 取C=50×rand(m,p)，ε=0.1，此處以10階為例，計算結(jié)果見表2。

表2 MM_Dykstra算法與SPG_Dykstra算法10階矩陣的誤差與時間比較

從例1和例2兩個例子的對比中，我們能夠看到MM_Dykstra算法無論是在殘差還是在時間上都有了一定改善。

6 結(jié)語

本文提出了MM_Dykstra算法來求解多約束矩陣方程最小二乘問題，該方法將三約束問題轉(zhuǎn)換為雙約束問題，有效地減少了交替迭代次數(shù)。數(shù)值實驗結(jié)果表明，與SPG_Dykstra方法相比，本文算法具有更快的求解效果，結(jié)果精度也相對更高。