一種求解矩陣填充問(wèn)題的交替共軛梯度最小化法

郭佳浩，閆喜紅

(太原師范學(xué)院 數(shù)學(xué)系，山西 晉中 030619)

0 引言

數(shù)據(jù)缺失下的低秩填充問(wèn)題是近幾年研究熱點(diǎn)之一，在機(jī)器學(xué)習(xí)[1]、模式識(shí)別[2]、模型縮減[3]等科學(xué)工程領(lǐng)域有著重要的作用.一個(gè)矩陣填充問(wèn)題如果沒(méi)有約束條件，可以有無(wú)窮多解，事實(shí)上需要填充的矩陣一般具有低秩或者近似低秩的結(jié)構(gòu).對(duì)于低秩矩陣，就可以通過(guò)有效的算法準(zhǔn)確合理地恢復(fù)出原矩陣.由于Fazel在近似低秩矩陣[4]和Candes與Recht在矩陣填充方面做的開(kāi)拓性工作[5]，使得這一問(wèn)題從理論和算法方面都得到了深入的研究[3，4，6-11].低秩矩陣填充最初是由Candes和Recht提出的，其目的是將僅部分元素已知的采樣矩陣精確合理地填充成一個(gè)低秩矩陣.矩陣填充和低秩逼近技術(shù)依賴于低秩結(jié)構(gòu)的元素之間的依賴關(guān)系.在已知矩陣部分元素的前提下，合理精確地尋找與已知項(xiàng)一致的最低秩矩陣，其數(shù)學(xué)模型可表示為:

(1)

其中Z0∈Rm×n是采樣矩陣，Ω是已知樣本元素下標(biāo)集，PΩ是集合Ω上的正交投影.問(wèn)題(1)是非凸的，直接求解非凸問(wèn)題非常困難，因此很多學(xué)者用它的凸松弛模型取代(1)，其凸松弛核范數(shù)模型如下:

(2)

其中，‖Z‖*是矩陣Z的核范數(shù).相關(guān)的理論已經(jīng)證明，只要奇異向量與正則基之間存在弱相關(guān)性，就可以通過(guò)求解上述凸松弛問(wèn)題(2)得到(1)的解[5].針對(duì)凸松弛模型(2)，有許多算法(如硬閾值算法、迭代閾值算法等)可以直接進(jìn)行求解.但是這些算法的實(shí)現(xiàn)需要在每次迭代中計(jì)算部分奇異值分解(SVD).眾所周知，計(jì)算奇異值分解(SVD)成本很高.為了避免求解奇異值的高計(jì)算成本，Haldar[12]和Wen[13]提出如下分解模型:

(3)

其中，Z=XY，X∈Rm×r，Y∈Rr×n，r?min{m，n}.在模型(3)中，Z表示成X∈Rm×r和Y∈Rr×n的乘積，使得確保滿足低秩的要求.兩個(gè)針對(duì)模型(3)的特點(diǎn)，一個(gè)自然的求解方法就是交替最小化方法，如Power Factorization[12]和LMaFit[13].在這兩種交替最小化方法中，變量X和Y交替更新，每一步要精確求解子問(wèn)題得到X和Y.精確求解子問(wèn)題的計(jì)算量較高.為此，文獻(xiàn)[14]設(shè)計(jì)了交替最速下降法(ASD)和加速交替最速下降方向法(ScaledASD)來(lái)求解(3).在此算法中變量X和Y交替更新，且在每一步迭代過(guò)程中利用一步最速下降法更新X和Y， 從而提高算法的效率.眾所周知，共軛梯度方向較最速下降方向效果好.因此，本文在文獻(xiàn)[14]算法的基礎(chǔ)上，每一步利用共軛梯度方法來(lái)更新X和Y，提出了求解模型(3)的新方法.

本文結(jié)構(gòu)如下，在第一節(jié)中介紹了幾種相關(guān)的交替最小化法;在第二節(jié)中，本文提出了4種交替共軛梯度法(ACG);在第三節(jié)中將ASD算法和共軛梯度算法進(jìn)行了數(shù)值實(shí)驗(yàn)，并進(jìn)行了比較，驗(yàn)證了交替共軛梯度最小化算法的有效性.

1 交替最小化法

交替最小化法，也稱高斯-塞德?tīng)柣驂K坐標(biāo)下降法.由于其簡(jiǎn)單性、低內(nèi)存要求和靈活性被廣泛用于求解優(yōu)化問(wèn)題.本文考慮的模型(3)的函數(shù)是具有兩塊變量X∈Rm×r和Y∈Rr×n的非凸函數(shù)，直接求解較難.注意到模型(3)中關(guān)于X和Y都是最小二乘子問(wèn)題，一個(gè)自然的想法就是利用交替最小化方法求解，如算法1(PF)是一種將交替最小化方法應(yīng)用于模型(3)的矩陣填充算法.

算法1冪因式分解算法(PF).

步一:取初始矩陣PΩ(Z0)，初始點(diǎn)X0∈Rm×r，Y0∈Rr×n.

步二:重復(fù)計(jì)算以下步驟:

直到達(dá)到停機(jī)標(biāo)準(zhǔn).

算法1中的每個(gè)子問(wèn)題都是標(biāo)準(zhǔn)的最小二乘問(wèn)題.作為一種典型的交替最小化方法，算法1的極限點(diǎn)必然是不動(dòng)點(diǎn)[12].算法1中的最小二乘子問(wèn)題的求解決定了每次迭代計(jì)算成本.為了降低每次迭代子問(wèn)題的計(jì)算成本，文獻(xiàn)[13]提出了如下模型:

(4)

其中，投影到Ω上的要求從(3)中的目標(biāo)移動(dòng)到線性約束.應(yīng)用于模型(4)的交替最小化方法給出了低秩矩陣擬合算法(LMaFit).

算法2低秩矩陣擬合(LMaFit).

步一:取初始矩陣PΩ(Z0)，初始點(diǎn)X0∈Rm×r，Y0∈Rr×n.

步二:重復(fù)計(jì)算以下步驟:

3)Zi+1=Xi+1Yi+1+PΩ(Z0-Xi+1Yi+1).

直到達(dá)到終止標(biāo)準(zhǔn).

算法2中的每一步迭代中先固定Z，通過(guò)求解最小二乘問(wèn)題，LMaFit不斷地更新X和Y.注意LMaFit中的最小二乘子問(wèn)題有顯示表達(dá)式，而可以精確求解，降低了PF算法中的每次迭代計(jì)算成本.但是子問(wèn)題求精確解，計(jì)算成本仍很高.為了改進(jìn)計(jì)算效率，文獻(xiàn)[14]采用單步沿梯度下降方向的直線搜索代替求解最小二乘法問(wèn)題，并給出如下交替最速下降算法(ASD算法)[14].

算法3交替最速下降算法(ASD).

步一:取初始矩陣PΩ(Z0)，初始點(diǎn)X0∈Rm×r，Y0∈Rr×n.

步二:重復(fù)計(jì)算以下步驟:

2)Xi+1=Xi-txifYi(Xi).

4)Yi+1=Yi-tyifXi+1(Yi).

2 交替共軛梯度最小化法

共軛梯度法是各種優(yōu)化算法中非常重要的一種，最早是由Hestenes和Stiefle[15]在解決大規(guī)模線性方程組問(wèn)題時(shí)提出來(lái)的，在他們的基礎(chǔ)上Fletcher和Reeves(1964)將其引入到了非線性最優(yōu)化問(wèn)題中，創(chuàng)立了非線性共軛梯度法.共軛梯度法是一個(gè)典型的共軛方向法，其中每一個(gè)搜索方向都是互相共軛的.該方法的基本思想是基于前一迭代點(diǎn)的搜索方向?qū)Ξ?dāng)前迭代點(diǎn)的負(fù)梯度方向進(jìn)行修正來(lái)產(chǎn)生新的搜索方向，換言之，搜素方向是當(dāng)前迭代點(diǎn)的負(fù)梯度方向與前一搜索方向的線性組合即:

dk=-gk+βk-1dk-1，

其中，dk是當(dāng)前點(diǎn)的搜索方向，gk是當(dāng)前點(diǎn)的梯度，dk-1是前一迭代點(diǎn)的搜索方向.共軛梯度法是介于最速下降法與牛頓法之間的一種方法，它僅利用了一階導(dǎo)數(shù)的信息，既彌補(bǔ)了最速下降法收斂慢的不足，又克服了牛頓法需要存儲(chǔ)和計(jì)算Hesse矩陣并求逆的缺點(diǎn).由于這種方法不需要矩陣存儲(chǔ)、穩(wěn)定性高、又有較快的收斂速度和二次終止性等優(yōu)點(diǎn).目前在實(shí)際問(wèn)題中已經(jīng)得到廣泛的應(yīng)用.根據(jù)參數(shù)βk-1的多種表述形式有如下四種常用的共軛梯度算法.

算法4共軛梯度法(ACG).

步一:取x0∈Rn和終止參數(shù)ε≥0.計(jì)算g0，d0=-g0.

步二:如果‖gk‖≤ε，算法終止；否則進(jìn)行下一步.

步三:重復(fù)計(jì)算如下步驟:

xk+1=xk+αkdk，dk+1=-gk+1+βkdk，其中αk是步長(zhǎng)，βk-1可由如下四種公式計(jì)算:

直到滿足停機(jī)準(zhǔn)則.

將此方法應(yīng)用到模型(3)的子問(wèn)題中，基于βi的四種計(jì)算方法，建立了交替共軛梯度法，其中，βi的四種計(jì)算方法如下:

下面詳細(xì)介紹四種交替共軛梯度法.

算法5交替共軛梯度法1(ACG-CW算法)，β是由Crowder-Wolfe公式計(jì)算.

步一:取PΩ(Z0)，X0∈Rm×r，Y0∈Rr×n和終止參數(shù)ε.

X1=X0-tx0fY0(X0)，

Y1=Y0-tY0fX1(Y0)，

令dx0=-fY0(X0)，dY0=-fX0(Y0).

i=2，3，…

Xi+1=Xi+tXidXi，

算法6交替共軛梯度法2(ACG-FR算法)，β是由Fletcher-Reeves公式計(jì)算

步一:取PΩ(Z0)，X0∈Rm×r，Y0∈Rr×n和終止參數(shù)ε.

X1=X0-tx0fY0(X0)，

Y1=Y0-tY0fX1(Y0)，

令dx0=-fY0(X0)，dY0=-fX0(Y0).

i=2，3,…

Xi+1=Xi+tXidXi，

算法7交替共軛梯度法3(ACG-PR算法)，β是由Polar-Ribiere公式計(jì)算.

步一:取PΩ(Z0)，X0∈Rm×r，Y0∈Rr×n和終止參數(shù)ε.

X1=X0-tx0fY0(X0)，

Y1=Y0-tY0fX1(Y0)，

令dx0=-fY0(X0)，dY0=-fX0(Y0).

i=2，3,…

Xi+1=Xi+tXidXi，

算法8交替共軛梯度法(ACG-D算法)，β是由Dixon公式計(jì)算.

步一:取PΩ(Z0)，X0∈Rm×r，Y0∈Rr×n和終止參數(shù)ε.

X1=X0-tx0fY0(X0)，

Y1=Y0-tY0fX1(Y0)，

令dx0=-fY0(X0)，dY0=-fX0(Y0).

i=2，3，…

Xi+1=Xi+tXidXi，

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

在本小節(jié)中，我們通過(guò)數(shù)值實(shí)驗(yàn)將文獻(xiàn)[14]中的ASD算法與我們的ACG-CW算法、ACG-FR算法、ACG-PR算法、ACG-D算法進(jìn)行比較.這里所有的數(shù)值實(shí)驗(yàn)均是在MATLAB(R2013a)中進(jìn)行的，機(jī)器精度為10-64，機(jī)器類(lèi)型為英特爾核2.4 GB.我們的矩陣也是在相同的實(shí)驗(yàn)環(huán)境中隨機(jī)產(chǎn)生的，矩陣的維數(shù)從500到1 000，秩取10和15，停機(jī)精度為10-4，采樣率為10%，運(yùn)行時(shí)間以秒記錄.最大迭代次數(shù)為1 000，數(shù)值結(jié)果如表1.

表1 ASD算法與ACG的四種算法的比較

通過(guò)表1中數(shù)值實(shí)驗(yàn)結(jié)果，四種共軛梯度算法中ACG-CW和ACG-PR比ASD的迭代次數(shù)少，計(jì)算時(shí)間短.ACG-FR算法和ACG-D算法的優(yōu)勢(shì)在于大大提高了精度.因此，矩陣填充的交替共軛梯度法從精度、迭代次數(shù)、運(yùn)行時(shí)間方面優(yōu)化了之前的交替最速下降法.

4 結(jié)論

本文針對(duì)矩陣填充這個(gè)熱點(diǎn)問(wèn)題，基于分解模型，建立一種求解矩陣填充問(wèn)題的共軛梯度交替方向最小化算法.由于現(xiàn)有的算法大都需要計(jì)算矩陣的奇異值分解，計(jì)算量較大，計(jì)算成本較高，為了避免高成本計(jì)算，本文利用分解模型提出交替共軛梯度算法，且把此算法應(yīng)用到隨機(jī)產(chǎn)生的低秩矩陣填充問(wèn)題中，數(shù)值結(jié)果表明了共軛梯度交替算法的可行性與優(yōu)越性.