鄭偉東，李聲豪，涂志輝，胡文玉，喻高航

(贛南師范大學(xué) 數(shù)學(xué)與計(jì)算機(jī)科學(xué)學(xué)院，江西 贛州 341000)

1 引言

矩陣填充(Matrix Completion, MC)問題主要解決得到數(shù)據(jù)矩陣的某一小部分的數(shù)據(jù)，填充那些未知或缺失的數(shù)據(jù).為了解決這個(gè)問題，通常假設(shè)這個(gè)矩陣存在信息冗余，例如數(shù)據(jù)矩陣是低秩的，即其數(shù)據(jù)分布在一個(gè)低維的線性子空間上.目前，矩陣填充問題在推薦系統(tǒng)、圖像恢復(fù)、數(shù)據(jù)分析以及機(jī)器學(xué)習(xí)等領(lǐng)域中有著廣泛的應(yīng)用.許多學(xué)者[1-8]在矩陣填充中做了大量的工作，這些工作主要圍繞如下優(yōu)化問題實(shí)現(xiàn)矩陣填充：

(1)

其中，X∈Rm×n是存在元素缺失的觀測矩陣，A∶Rm×n→Rn是一個(gè)線性算子，b∈Rn.由于矩陣的秩函數(shù)是一個(gè)非凸且不連續(xù)的函數(shù)，因此，問題(1)是一個(gè)NP難問題.Candès和Recht[1]證明了問題(1)可以通過核范數(shù)近似矩陣的秩進(jìn)行求解：

(2)

(3)

其中，δ≥0是噪聲水平.另一方面，問題(3)還可以轉(zhuǎn)化為如下的無約束的最小二乘問題：

(4)

Micchelli等人[13]提出了基于臨近算子的不動(dòng)點(diǎn)算法，該算法在理論上可以解決兩個(gè)凸函數(shù)相加的優(yōu)化問題. 基于這種動(dòng)機(jī), 本文將問題(2)、(3)和(4)統(tǒng)一的歸納為兩個(gè)凸函數(shù)相加的模型，使得求解該類矩陣填充問題時(shí)更加通用簡潔.考慮如下模型：

(5)

(6)

2 算法

為了進(jìn)一步提高算法的效率，本節(jié)使用了臨近算子并結(jié)合了臨近算子和次微分之間的關(guān)系，提出了一類不動(dòng)點(diǎn)算法來解決問題(5).

2.1 相關(guān)理論

首先給出臨近算子和次微分的定義，并給出它們之間的聯(lián)系和后續(xù)需要使用的一些性質(zhì).

定義1[13](臨近算子)若f是Rm上的實(shí)凸函數(shù)，則函數(shù)f在x∈Rn的臨近算子如下：

(7)

(8)

定義2[13](次微分)若f是Rm上的實(shí)凸函數(shù)，則函數(shù)f在x∈Rn的次微分如下：

?f(x)={y∈Rn∶f(z)≥f(x)+〈y,z-x〉,z∈Rn}.

(9)

次微分是函數(shù)在不可微情形下的推廣.若函數(shù)可微，則次微分就是函數(shù)的導(dǎo)數(shù).由定義容易證明：λ?f(x)=?λf(x)，其中λ≥0.

引理1[15]設(shè)f*是f的共軛函數(shù)，即f*(y)∶= sup{〈x,y〉-f(x)}.若x∈domf，y∈domf*，則

y∈?f(x)??f*(y)∈x

(10)

成立.

引理2[13-14]若f是Rn上的一個(gè)凸函數(shù)，且x∈Rn，則對任意的H∈S+，有

Hy∈?f(x)?x=proxf,H(x+y).

(11)

2.2 不動(dòng)點(diǎn)算法

本小結(jié)設(shè)計(jì)求解優(yōu)化問題(5)的不動(dòng)點(diǎn)算法.假設(shè)X∈Rm×n是問題的最優(yōu)解，則根據(jù)費(fèi)馬定理和鏈?zhǔn)椒▌t，可得：

0∈?f(X)+A*?g(AX).

(12)

其中A*表示線性算子A的共軛算子.因此，存在y∈?g(AX)，使得0∈?f(X)+A*y.由引理1和(12)，可得：

AX∈?g(y), -A*y∈?f(X).

(13)

則對任意的矩陣P∈X+和Q∈S+，有

PP-1AX∈?g*(y), -QQ-1A*y∈?f(X).

(14)

根據(jù)引理2，可得如下與問題(5)等價(jià)的不動(dòng)點(diǎn)方程組：

(15)

根據(jù)不動(dòng)點(diǎn)方程組(15)，可以設(shè)計(jì)不同迭代格式的迭代方法，如

(16)

雖然迭代格式(16)簡單，且是顯式的，但是算法可能不收斂[14].因此本文將構(gòu)造一個(gè)半隱式迭代格式.對m個(gè)向量空間E1,E2,…,Em，定義向量空間Ei的內(nèi)積為〈·,·〉Ei，且在笛卡爾積的內(nèi)積空間定義如下[16]：

(17)

proxΦ,R(Z)=(proxg*,P(y),proxf,Q(X))?Zk+1=proxΦ,R(EZk),

(18)

Zk+1=proxΦ,R(E0Zk+1+E1Zk),

(19)

即

(20)

(21)

2.3 子問題求解

在問題的求解過程中，需要求解函數(shù)的臨近算子.接下來將詳細(xì)說明子問題的求解.

子問題yk+1對于有約束和無約束優(yōu)化問題，計(jì)算g(x)的臨近算子，分兩種情況討論.

因此，根據(jù)g(x)的臨近算子容易求得yk+1.

子問題Xk+1給定值為r的矩陣X∈Rm×n，其奇異值分解為：X=Udiag(max{σi}1≤i≤r)VT，對任意的τ≥0，收縮奇異值閾值算子Dτ定義為：Dτ(X)=Udiag(max{σi-τ,0})VT.

(22)

3 收斂性分析

本節(jié)將證明本文提出的不動(dòng)點(diǎn)算法(20)是收斂的.

引理3[18]若f是從Rn到(-∞,+∞]的真下半連續(xù)(Proper lower semi-continuous)凸函數(shù)的集合，則proxf和I-proxf是擴(kuò)張的.

〈Zk+1-Z*,R(E0Zk+1+E1Z*-Zk+1)-R(EZ*-Z*)〉≥0

?-〈Zk+1-Z*,RE1(Zk+1-Z*)〉+〈Zk+1-Z*,R(E-I)(Zk+1-Z*)〉≥0.

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

在本節(jié)中，我們創(chuàng)建隨機(jī)矩陣和樣本集，然后使用我們的算法來解決無噪聲矩陣填充問題、噪聲矩陣填充問題和低秩圖像恢復(fù)問題.同時(shí)將我們的算法與前面提到的ADMM[7]、IADM-CG[8]和IADM-BB[11]算法進(jìn)行比較.所有的實(shí)驗(yàn)代碼在MATLAB R2016上編寫，以及在2.4GHz的主頻CPU和4 GB的內(nèi)存的計(jì)算機(jī)環(huán)境下運(yùn)行.

部分奇異值計(jì)算計(jì)算Xk+1時(shí)，注意到在每次迭代過程中，需要計(jì)算奇異值分解(SVD). 但是，實(shí)際上只需要比閾值1/μ更大的奇異值.因此可以通過使用軟件包PROPACK[19]來計(jì)算矩陣的部分奇異值分解，該軟件用于計(jì)算大于閾值的奇異值和相應(yīng)的奇異值向量.雖然PROPACK可以計(jì)算固定數(shù)量的奇異值，但它無法自動(dòng)確定大于1/μ的奇異值的數(shù)量. 因此，為了執(zhí)行部分SVD，我們需要預(yù)測每次迭代中大于閾值的奇異值的數(shù)量，具體更新規(guī)則如下：

其中svk表示預(yù)測的奇異值和奇異值向量的數(shù)量，sv0=min(m,n)×0.01，svpk是Xk的正奇異值的數(shù)量.

4.1 矩陣填充問題

在矩陣填充問題中，線性映射A是一個(gè)線性投影算子：AX=XΩ,XΩ是矩陣X的成分組成的向量形式.實(shí)驗(yàn)過程中，我們將根據(jù)三元組{m,r,p/dr}隨機(jī)生成方陣M，其中m是矩陣的維數(shù)，r表示矩陣的秩，p為采樣元素的數(shù)目，dr=r×(m+n-r)為m×n矩陣的自由度，同時(shí)定義樣本采樣率sr=p/(mn).首先生成每個(gè)元素為獨(dú)立高斯隨機(jī)變量的矩陣L=randn(m,r)和R=randn(n,r)，然后計(jì)算M=LRT，接著均勻采樣矩陣的p個(gè)元素.

表1給出的是本文提出的模型和算法與ADM和IADM-CG求解無噪聲矩陣填充問題的數(shù)值實(shí)驗(yàn)結(jié)果，和噪聲矩陣填充問題的數(shù)值結(jié)果，表2是使用本文的模型和ADM[7]、IADM-CG[12]算法求解有噪聲矩陣填充問題的實(shí)驗(yàn)結(jié)果.從表中可以看出，本文的算法在精確度和計(jì)算效率方面優(yōu)于IADM-CG算法，與ADMM算法相比計(jì)算效率相當(dāng)，但本文的算法在精確度上具有優(yōu)勢.

表1 本文算法與ADM、IADM-CG算法求解無噪聲的矩陣填充問題(Time：s)

表2 本文算法與ADM、IADM-CG算法求解有噪聲的矩陣填充問題(Time:s)

4.2 低秩圖像恢復(fù)中的應(yīng)用

在本節(jié)中，我們使用最小二乘問題模型來對低秩圖像進(jìn)行恢復(fù).通過測試512×512灰度圖像(lena, pirate, cameraman)來驗(yàn)證本文的方法在實(shí)際應(yīng)用中的有效性.首先，通過奇異值分解將原始圖像壓縮為秩為40的低秩圖像.然后從低秩矩陣中隨機(jī)獲取一部分元素得到數(shù)據(jù)丟失的圖像.最后，使用本文的模型和算法進(jìn)行恢復(fù).當(dāng)RelChg低于tol=10-6或迭代步數(shù)大于1 000時(shí)，迭代過程停止.原始圖像、低秩圖像、數(shù)據(jù)丟失的圖像和恢復(fù)的圖像在圖1中列出.從恢復(fù)的圖像結(jié)果中可以看出，我們的模型和算法對于低秩圖像的恢復(fù)是有效的.

圖1 第一列：原始圖像；第二列：秩為40的圖像；第三列：采樣率為sr=0.4的圖像；最后一列：恢復(fù)圖像

我們通過計(jì)算各方法的峰值信噪比(PSNR)來更好的比較恢復(fù)圖像的效果，定義如下：

表3給出了通過不同的方法恢復(fù)秩為40的圖像lena的時(shí)間、相對誤差和PSNR.從表3可以看出，對于不同的采樣率，本文的方法與ADM方法的PSNR比IADM-CG方法要好，即恢復(fù)的圖像質(zhì)量更好，并且本文方法的優(yōu)勢在于所需時(shí)間比其它兩種方法較快許多.

表3 本文算法與ADM、IADM-CG算法對于低秩圖像恢復(fù)(lena)

通過圖2描述了恢復(fù)圖像秩的估計(jì)與迭代步數(shù)的關(guān)系，以及相對誤差與算法所用時(shí)間的關(guān)系.由左圖可知，本文的方法在35次迭代之后非常接近真實(shí)的秩，由右圖可知，本文的方法收斂速率更快.

圖2 左圖：秩的估計(jì)；右圖：收斂速率

5 結(jié)語

本文歸納了基于核范數(shù)的矩陣填充模型，提出了一個(gè)統(tǒng)一的數(shù)學(xué)模型，并采用不動(dòng)點(diǎn)算法求解該模型，使得迭代格式變得通用和簡單, 且具有簡潔的收斂性證明.在合成數(shù)據(jù)的矩陣填充數(shù)值實(shí)驗(yàn)中，本文的方法提高了算法的性能；在低秩圖像恢復(fù)中，節(jié)省了大量的運(yùn)行時(shí)間；并且對于相同模型不同的線性算子，算法框架也將適用，可應(yīng)用于求解其它問題.