低秩Toeplitz 張量的高精度隨機填充算法

溫瑞萍 李文韋

(太原師范學(xué)院數(shù)學(xué)與統(tǒng)計學(xué)院,晉中,030619)

1 引言

由于張量在多模態(tài)結(jié)構(gòu)信息表中具有優(yōu)勢,所以對張量的研究成為近些年來的熱點問題.一個已經(jīng)引起學(xué)界高度關(guān)注的現(xiàn)象是: 在張量數(shù)據(jù)的采集、存儲、傳輸與處理等環(huán)節(jié),常常會發(fā)生數(shù)據(jù)的丟失或腐蝕,因此張量數(shù)據(jù)的填充與恢復(fù)問題的研究就尤為重要.低秩張量填充問題作為低秩矩陣填充問題的推廣,主要是根據(jù)已知的部分張量數(shù)據(jù)信息來填充其未知或缺失的數(shù)據(jù).張量填充方法已經(jīng)廣泛應(yīng)用于圖像處理[1]、數(shù)據(jù)挖掘[2]、信號處理[3]及機器學(xué)習(xí)[4]等領(lǐng)域.低秩張量填充的數(shù)學(xué)模型一般為:

類似于矩陣填充的秩極小化模型,問題(1.1)屬于非凸優(yōu)化模型,是NP 難問題.如果待填充張量C的秩已知或容易估計,一些學(xué)者提出了類似于求解矩陣填充秩極小化模型的解法,但填充效果不甚理想.因而,在實際處理中,通常用逼近張量T的一個最緊凸包‖T‖*來代替rank(T),將問題(1.1)轉(zhuǎn)化為下列凸優(yōu)化模型,即張量填充的核范數(shù)極小化數(shù)學(xué)模型:

其中‖·‖*表示相應(yīng)張量的核范數(shù).

求解(1.2)通常的作法是將張量按照不同的模展開成矩陣,對矩陣進(jìn)行填充操作后再折疊復(fù)原回張量形式.針對矩陣填充的核范數(shù)極小化問題求解,常見的經(jīng)典算法有: Toh 等學(xué)者提出的加速近端梯度(APG)算法[5]?Cai 等提出的奇異值閾值(SVT)算法[6]?Lin 等提出的增廣拉格朗日乘子(ALM)算法等[7].因為張量的分解方式有CP 分解、Tucker 分解、張量鏈分解、張量環(huán)分解和t-SVD 分解等[8],所以張量填充問題的求解要比矩陣填充問題的求解復(fù)雜得多.根據(jù)各種分解,許多學(xué)者提出了求解問題(1.2)的有效方法.比如,Kolda 等提出了高階SVD 算法[9].Gandy 等提出了采用ADM 的框架求解低秩稀疏張量恢復(fù)問題[10].Liu 等基于Tucker 分解提出了三種張量填充方法[11]: 利用松弛技術(shù),采用塊坐標(biāo)下降(BCD)方法求全局最優(yōu)解的簡單低秩張量填充算法(SiLRTC)?將原非光滑問題轉(zhuǎn)化為光滑問題,進(jìn)而求解張量跡范數(shù)極小化問題的快速低秩張量填充算法(FaLRTC)?在SiLRTC 算法的基礎(chǔ)上應(yīng)用交替方向乘子方法的高精度低秩張量填充算法(HaLRTC).

由于待填充矩陣往往具有特殊結(jié)構(gòu),Toeplitz 矩陣填充問題近年來受到許多學(xué)者的關(guān)注.比如,王川龍等提出了Toeplitz 矩陣填充的均值算法[12]、修正的增廣拉格朗日乘子算法[13]、保結(jié)構(gòu)算法[14]、子空間算法[15]等.溫瑞萍等提出了求解Toeplitz 矩陣填充問題的逐步結(jié)構(gòu)化增廣拉格朗日乘子算法[16]、l步結(jié)構(gòu)化增廣拉格朗日乘子算法[17]、均值修正增廣拉格朗日乘子算法[18]及尾端修正增廣拉格朗日乘子算法[19].同樣地,在張量填充的某些應(yīng)用問題中,雖然待填充數(shù)據(jù)表的規(guī)模很大,但是填充問題本身往往具有對稱,Hankel 及Toeplitz 等特殊結(jié)構(gòu),因而,結(jié)構(gòu)化張量的填充和恢復(fù)在實際問題中的應(yīng)用是很有意義且非常必要的.諸如,Badeac 等針對Hankel 張量、對稱張量、Toeplitz 張量填充問題分別提出了奇異值分解快速算法[20]?聶家旺等學(xué)者討論了Hankel 張量分解和秩,證明了Hankel 張量的CP 秩、對稱秩及Vandermonde 秩的一致性[21]?王川龍等提出了一種快速且具有較高精度的Hankel 張量填充算法[22].

本文采用HaLRTC 算法[11]求解Toeplitz 張量填充問題,即當(dāng)(1.2)中的待填充張量C恰好具有Toeplitz 結(jié)構(gòu)時的情形.由于HaLRTC 算法在N次迭代的每一步均包括了需要較大計算花費的奇異值分解(SVD),本研究主要從以下兩個方面改進(jìn)它以降低其每次迭代中計算奇異值分解的代價: 一方面是利用具有特殊結(jié)構(gòu)的Toeplitz 張量的快速奇異值分解方法,縮減奇異值分解的時間?另一方面是應(yīng)用隨機化思想,將在每一步迭代需按全部模展開Toeplitz 張量的操作改進(jìn)為按第n模隨機展開一次,減少奇異值分解的次數(shù),從而達(dá)到提高算法效率的目的.最后通過數(shù)值實驗驗證改進(jìn)算法的有效性.

2 預(yù)備知識

定義2.1(奇異值分解(SVD),[23])如果秩為r的矩陣A∈Rm×n的全部非零奇異值依次為σ1≥σ2≥···≥σr ＞0,則A的奇異值分解定義為:

其中U=[u1,u2,...,ur]∈Rm×r和V=[v1,v2,...,vr]∈Rn×r分別為矩陣A的前r個左、右奇異向量構(gòu)成的列正交矩陣.

的矩陣稱為一個n×n的Toeplitz 矩陣.

顯然,T由其第一行和第一列決定,共有t-n+1,t-n+2,...,tn-1至多2n-1 個互異元素,它們分布于T的2n-1 條對角線上.若用l表示Toeplitz 矩陣對角線的指標(biāo),則有l(wèi) ∈{-n+1,-n+2,...,n-1}.相應(yīng)地,Toeplitz 矩陣填充問題的元素采樣指標(biāo)集??{-n+1,-n+2,...,n-1}.用diag(T,l)表示Toeplitz 矩陣T的第l條對角線元素構(gòu)成的向量,則

定義2.4(Toeplitz 均值結(jié)構(gòu)化算子,[12])對于任意矩陣X=(xij)∈Rn×n,Toeplitz 均值結(jié)構(gòu)化算子T定義如下:

其中,

可見,T(X)是由矩陣X派生的Toeplitz 矩陣,即任意矩陣都可以通過結(jié)構(gòu)化算子T(·)轉(zhuǎn)化為一個Toeplitz 矩陣.

則稱X為Toeplitz 張量.

類似于矩陣的情形,將任意張量X進(jìn)行Toeplitz 結(jié)構(gòu)化的運算表示為Γ(X)=T.

3 Toeplitz 張量的填充算法

本節(jié)將文獻(xiàn)[11]給出的高精度低秩張量填充算法(算法3.1)進(jìn)行修改,應(yīng)用于求解模型(1.2)中的待填充張量C恰好具有Toeplitz 結(jié)構(gòu)的低秩Toeplitz 張量填充問題(算法3.2).

算法3.1 高精度低秩張量填充算法(High-accuracy Low Rank Tensor Completion Algorithm,HaLRTC)[11]

初始化:給定下標(biāo)集合?,參數(shù)?,ρ0,t,X0=PΩ(C),Yi,0=0,k:=0;

第1 步:計算

end for?

第2 步:計算

第3 步:若‖Xk+1-Xk‖F(xiàn)/‖Xk‖F(xiàn) ＜?,則停止迭代,否則,轉(zhuǎn)第4 步?

第4 步:計算Yi,k+1=Yi,k-ρk(Mi,k+1-Xk+1),ρk+1=tρk,k:=k+1,轉(zhuǎn)第1 步.

算法3.2 高精度低秩張量隨機填充算法(High-accuracy Low Rank Tensor Random Completion Algorithm,HaLRTRC)

初始化:給定下標(biāo)集合?,參數(shù)?,ρ0,t,X0=PΩ(C),Yi,0=0,k:=0;

第1 步:隨機生成一整數(shù)i=randperm(N,1),并計算

第2 步:令Ti,k+1=Γ(Mi,k+1);

第3 步:計算

第4 步:若‖Xi,k+1-Xi,k‖F(xiàn)/‖Xi,k‖F(xiàn) ＜?,則停止迭代,否則,轉(zhuǎn) 第4 步?

第5 步:計算Yi,k+1=Yi,k-ρk(Ti,k+1-Xi,k+1),ρk+1=tρk,k=k+1,轉(zhuǎn) 第1 步.

算法3.1 與算法3.2 的主要區(qū)別在于算法3.2 利用張量結(jié)構(gòu)的同時應(yīng)用了隨機思想.將算法3.1 的每步迭代中張量按每個n模展開修改為算法3.2 中的張量隨機地僅按某一k模展開一次.雖然算法的迭代步數(shù)和計算復(fù)雜度有所增加,但是避免了算法3.1 每步迭代中N次展開后帶來的奇異值分解計算量,進(jìn)而體現(xiàn)出算法3.2 在計算時間上具有優(yōu)勢.結(jié)構(gòu)化張量不僅在數(shù)據(jù)傳輸和算法中能夠減少存儲空間,而且在算法中每次迭代均可保持Toeplitz 結(jié)構(gòu).這一方面提高了迭代矩陣的精度,另一方面可以用快速算法進(jìn)行結(jié)構(gòu)化矩陣的奇異值分解.

4 收斂性分析

引理4.1由算法3.2 產(chǎn)生的序列{Yi,k+1}是有界的.

由算法3.2 的第5 步可得

此外,

再由算法3.2 的第2 步可知

再由文獻(xiàn)[7] 中引理1 和引理4 可得Yi,k+1=Yi,k -ρk(Ti,k+1-Xi,k+1)∈?‖Ti,k+1‖*.令Ti,k+1=UΣVT.再由文獻(xiàn)[6]可知

因此,

所以序列{Yi,k+1}是有界的.

令(T *,X*)為問題(1.2)的最優(yōu)解,Y*是文獻(xiàn)[7]中的對偶問題的最優(yōu)解.我們首先證明

所以(4.1)式成立.

類似于文獻(xiàn)[7]中定理2 的證明即可得T *是(1.2)的最優(yōu)解.

根據(jù)文獻(xiàn)[20],可得下列平凡的結(jié)論.

定理4.2設(shè)M∈Γ(M)是Toeplitz 張量.則對于任意的Toeplitz 張量Y,有

定理4.3設(shè)Ti,k+1是在算法3.2 中由Mi,k+1所生成的Toeplitz 張量.則有

其中T *是問題(1.2)的最優(yōu)解.

證明由于

5 數(shù)值實驗

本節(jié)通過數(shù)值實驗驗證算法3.1 與算法3.2 的有效性.實驗包括兩個部分: 分別利用T-HaLRTC與T-HaLRTRC 算法求解Toeplitz 張量填充問題及分別利用MT-HaLTRC 與MT-HaLRTRC 算法求解Toeplitz 均值張量填充問題.所有數(shù)值實驗均在戴爾(DELL)PowerEdge R740xd 高性能2U 機架式并行計算服務(wù)器上,MATLAB(R2019b)環(huán)境下運行實現(xiàn).從算法的迭代步數(shù)、CUP 總的計算時間(單位: 秒(s))及精度三個方面進(jìn)行比較.

在實驗中,p表示采樣率,參數(shù)αi=,ρ0=10-7,其余參數(shù)均采用文獻(xiàn)[11] 的建議.另外,利用隨機方法生成秩不同且維數(shù)分別為50×50×50,100×100×100,150×150×150,200×200×200,250×250×250 和300×300×300 且秩為(10,10,10)的Toeplitz 張量和Toeplitz均值張量.

5.1 Toeplitz 張量

兩種算法在采樣率分別為0.6,0.5,0.4,0.3,0.2 時的實驗結(jié)果分別見表1 至表5.下文中X*表示輸出的Toeplitz 張量,C為原始的Toeplitz 張量?加速比SP 和相對誤差RSE 的計算式如下:

表1 當(dāng)p=0.6 時,算法T-HaLRTC 和T-HaLRTRC 的比較

表2 當(dāng)p=0.5 時,算法T-HaLRTC 和T-HaLRTRC 的比較

表3 當(dāng)p=0.4 時,算法T-HaLRTC 和T-HaLRTRC 的比較

表4 當(dāng)p=0.3 時,算法T-HaLRTC 和T-HaLRTRC 的比較

表5 當(dāng)p=0.2 時,算法T-HaLRTC 和T-HaLRTRC 的比較

包括[11] 在內(nèi)的文獻(xiàn)中的數(shù)值實驗,通常只針對max{I1,I2,I3} 小于或等于100 的張量數(shù)據(jù)進(jìn)行數(shù)值測試,而本文的數(shù)值實驗規(guī)模達(dá)到了I1=I2=I3=300.由表1 至表5 可知,在采樣率p ∈{0.2,0.3,0.4,0.5,0.6}時,T-HaLRTC 算法與算法3.2 均收斂,兩個算法都滿足了規(guī)定的停止準(zhǔn)則.雖然新算法3.2 的迭代步數(shù)總是大于T-HaLRTC 算法的迭代步數(shù),但是由于嵌入Toeplitz 結(jié)構(gòu)化操作,使得精度提高,而且加入隨機思想減少了奇異值分解的計算量,從而總的運行時間顯著地縮短,新算法3.2 比T-HaLRTC 算法花費更少.特別是當(dāng)張量的采樣率較小而秩比較大的時候,算法3.2 更彰顯優(yōu)勢,其加速比高達(dá)約50%.

此外,兩個算法的計算時間比較圖(圖1)更加直觀地體現(xiàn)了新算法的有效性,其中圖1 (a),(b),(c),(d)分別顯示在不同采樣率下,T-HaLRTC 算法與T-HaLRTRC 算法在不同維數(shù)且秩為(10,10,10)時對應(yīng)的CPU 耗時.從圖1 可知,T-HaLRTRC 算法的CPU 時間顯著少于T-HaLRTC 算法的CPU時間,表明T-HaLRTRC 算法比T-HaLRTC 算法在總計算代價方面更具有優(yōu)勢.

圖1 不同采樣率下T-HaLRTC 算法與T-HaLRTRC 算法時間的比較

5.2 Toeplitz 均值張量

本小節(jié)針對Toeplitz 均值張量,在采樣率分別為0.5,0.4,0.3,0.2 時對兩種算法進(jìn)行比較,實驗結(jié)果見表6 至表9.本小節(jié)針對Toeplitz 均值張量也進(jìn)行了大量的數(shù)值實驗.同5.1 節(jié)的情形,實驗規(guī)模達(dá)到了I1=I2=I3=300.

表6 當(dāng)p=0.5 時,算法MT-HaLRTC 和MT-HaLRTRC 的比較

表7 當(dāng)p=0.4 時,算法MT-HaLRTC 和MT-HaLRTRC 的比較

表8 當(dāng)p=0.3 時,算法MT-HaLRTC 和MT-HaLRTRC 的比較

表9 當(dāng)p=0.2 時,算法MT-HaLRTC 和MT-HaLRTRC 的比較

從表6 至表9 可以看出,在采樣率p ∈{0.2,0.3,0.4,0.5}時,MT-HaLRTC 算法與本文的新算法3.2 均收斂,均達(dá)到了規(guī)定的停止準(zhǔn)則.這表明兩種數(shù)值算法都具有較好的穩(wěn)定性.雖然增加了隨機思想的算法3.2 的迭代步數(shù)總是大于MT-HaLRTC 算法的迭代步數(shù),但是在引入隨機思想的同時嵌入了Toeplitz 結(jié)構(gòu)化操作,這不僅使得精度提高,而且降低了計算量,從而總的運行時間明顯縮短,特別是當(dāng)張量的采樣率較小(這意味著原始的已知數(shù)據(jù)很少) 而秩比較大的時候,優(yōu)勢更加顯著,其加速比高達(dá)50%.

圖2 直觀地體現(xiàn)了新算法的有效性,其中,(a),(b),(c),(d) 分別顯示了在不同采樣率下,MT-HaLRTC 算法與MT-HaLRTRC 算法在不同規(guī)模且秩為(10,10,10) 時對應(yīng)的CPU 耗時.圖2 表明MT-HaLRTRC 算法比MT-HaLRTC 算法在總計算花費方面更具有優(yōu)勢.

圖2 不同采樣率下MT-HaLRTC 算法與MT-HaLRTRC 算法時間的比較

6 總結(jié)

本文通過在高精度填充算法中引入隨機思想給出了一種新算法并應(yīng)用于求解低秩Toeplitz 張量與低秩Toeplitz 均值張量填充問題.新算法按預(yù)期達(dá)到了停止準(zhǔn)則.文中給出了算法的收斂性分析,并且通過大量的數(shù)值實驗驗證了新算法的有效性.