史加榮, 白姍姍
(西安建筑科技大學 理學院, 西安710055)
非負矩陣分解(non-negative matrix factorization, NMF)是尋找非負的基矩陣和系數(shù)矩陣, 即將非負數(shù)據(jù)表示為非負基下的非負線性組合, 以實現(xiàn)維數(shù)約簡和特征提取[1]. 實際應用中很多物理特性都是非負的, 非負約束使模型具有較強的可解釋性[2]. 近年來, NMF受到廣泛關注, 已被成功應用到數(shù)據(jù)挖掘[3]、 計算機視覺[4-5]和模式識別[6]等領域. 乘性更新(multiplicative update, MU)[7]是求解NMF的一種有效方法, 可保證其解收斂于一個平穩(wěn)點. 但MU計算復雜度較高, 故不適用于求解大規(guī)模矩陣分解問題. 隨著數(shù)據(jù)規(guī)模的不斷擴大, 隨機梯度下降(stochastic gradient descent, SGD)算法[8]已成為機器學習, 特別是深度學習研究的熱點[9-10], 其具有參數(shù)更新過程簡單、 收斂速度快且計算復雜度低等特點. 因此, 將SGD算法與MU規(guī)則相結合是求解NMF的一種有效方法.
Kasai[11]將SGD算法應用到NMF的求解中, 提出了隨機乘性更新(stochastic multiplicative update, SMU)和隨機方差縮減乘性更新(stochastic variance reduced multiplicative update, SVRMU)規(guī)則, 分別將SGD算法和隨機方差縮減梯度(stochastic variance reduction gradient, SVRG)[12]算法與乘性更新規(guī)則相結合, 對梯度估計量進行更新; 之后,Kasai[13]又提出了隨機平均梯度乘性更新(stochastic average gradient multiplicative update, SAGMU)規(guī)則. 雖然這些算法取得了較好的效果, 但都未考慮偏差與方差之間的不平衡性差異, 即未對梯度估計量進行校正, 導致迭代效率較低. 文獻[14]提出了一種性能優(yōu)于SVRG的算法: 隨機方差調整梯度(stochastic variance adjusted gradient, SVAG)算法. 該算法在迭代過程中引入?yún)?shù)對梯度估計量進行不斷調整, 使偏差與方差達到平衡. 基于此, 本文提出一種將SVAG算法和MU規(guī)則相結合的隨機乘性更新算法: 隨機方差參數(shù)調整梯度乘性更新(stochastic variance parameter adjusted gradient multiplicative update, SVPAGMU)算法. 該算法結合方差縮減策略和乘性更新規(guī)則, 在梯度向量更新的同時, 引入一個參數(shù)調整梯度估計量, 不斷校正梯度下降方向使其偏差與方差達到平衡, 最終快速逼近局部最優(yōu)解.
(1)
其中‖·‖F(xiàn)為矩陣的Frobenius范數(shù),W為基向量組成的矩陣,H為系數(shù)矩陣. 由于優(yōu)化模型(1)是非凸的, 因此求解該模型的全局最優(yōu)解非常困難, 塊坐標下降算法[15]是解決該問題的一種主要方法.
針對最優(yōu)化模型(1), 文獻[7]提出了一種簡單而有效的乘性更新規(guī)則:
H←H⊙((WTV)./(WTWH)),
(2)
W←W⊙((VHT)./(WHHT)),
(3)
其中⊙表示Hadamard積, 即對應元素之積, ./表示對應元素之商. MU規(guī)則通過交替更新W和H, 使模型(1)中目標函數(shù)在每步迭代中都減小, 最終達到平穩(wěn)點.
許多機器學習任務都可轉化為經(jīng)驗風險最小化模型:
(4)
其中x∈d為模型參數(shù),fi(x):d→是第i個樣本關于參數(shù)x的損失函數(shù). SGD算法是求解大規(guī)模學習問題的一種有效方法, 下面簡單介紹幾種具有代表性的隨機梯度下降算法, 更多算法可參見文獻[16].
SGD算法以其基本的梯度下降形式更新迭代, 在每輪更新參數(shù)時, 僅隨機選取一個或幾個樣本計算其梯度, 并以此梯度作為全局梯度的估計值, 極大降低了計算復雜度. SGD算法的參數(shù)更新公式為
xt+1=xt-ηfi(xt),
(5)
其中參數(shù)η表示迭代步長(或學習率),i∈[N]([N]={1,2,…,N})表示第t輪迭代中隨機選取的樣本序號,fi(xt)表示當前迭代梯度. 由于梯度方差隨迭代次數(shù)的增加而不斷累加, 因此導致SGD算法收斂速度較慢.
為平衡有偏低方差的SAG算法和無偏高方差的SAGA算法, 文獻[14]引入了SVAG算法, 該算法引入一個參數(shù)權值θ對梯度估計量進行校正, 即需判斷隨機選取的單個采樣梯度fit(xt)與該樣例處的最新梯度fit(xtit)之間的相關程度, 不斷調整梯度使其偏差與方差達到平衡, 參數(shù)更新公式為
(6)
其中ti表示第t輪迭代中隨機選取的樣本序號i. 特別地,θ=1和θ=N分別對應于SAG算法和SAGA算法. 數(shù)值算例驗證了SVAG算法優(yōu)于SVRG算法及其他隨機梯度下降算法.
設vn和hn分別表示矩陣V和H的第n列, 且n∈[N]是選取的樣本序號. 最優(yōu)化模型(1)等價于
(7)
圖1 SVPAGMU算法的基本流程Fig.1 Basic flow chart of SVPAGMU algorithm
(8)
(9)
在更新W前, 為使其偏差與方差之間達到平衡, 校正梯度下降方向, 需判斷隨機采樣梯度g1與當前快照點處梯度g2的相關程度, 記為
(10)
(11)
若‖g1-g2‖2<ε, 則選擇較大的參數(shù)θ; 否則, 選擇較小的參數(shù)θ. 通過θ取值的不同調整梯度估計量, 使其沿著最佳的梯度下降方向更新W. 根據(jù)SVAG算法可知: 參數(shù)θ的選取與樣本數(shù)目N有關, 根據(jù)不同的數(shù)據(jù), 判斷梯度g1與g2的相關關系, 選取合適的參數(shù)θ對梯度估計量進行不斷調整.
(13)
(14)
(15)
隨著迭代次數(shù)的增加,αs不斷縮減, 即
αs+1=ραs,
(16)
其中ρ∈(0,1). 綜上可得:
算法1隨機方差調整梯度乘性更新(SVPAGMU)算法.
步驟2) Fors=0,1,…, do
步驟4) Fort=1,2,…,T
步驟6) 根據(jù)式(9)更新hm;
步驟7) 根據(jù)式(10)計算g1;
步驟8) 根據(jù)式(11)計算g2;
步驟9) if ‖g1-g2‖2<ε;θ=θmax; elseθ=θmin; end
步驟14) End For
步驟15) 根據(jù)式(16)更新αs;
步驟17) End For.
下面討論SVPAGMU算法在一次外部迭代過程的計算復雜度. 若內層循環(huán)最大迭代次數(shù)為T, 則更新系數(shù)矩陣H的計算復雜度為O(TMK2), 計算梯度g1和g2的復雜度均為O(TMK), 計算Q和P的復雜度均為O(TMNK), 計算S的復雜度為O(TM2K2), 最終得到更新基矩陣W的計算復雜度為O(T(M2K2+MNK)). 因此, 一次迭代總的計算復雜度為O(T(M2K2+MNK)).
下面在人臉圖像集和氣象數(shù)據(jù)集上進行實驗, 并將SVPAGMU算法與SMU,SVRMU和SAGMU三種算法進行比較. 使用MATLAB R2019a進行編程, 涉及的全部算法程序均在IntelCoreTMi7-8565U CPU @1.80 GHz 1.99 GHz處理器上運行, 操作系統(tǒng)為64位Windows 10系統(tǒng).
(17)
相對誤差越小, 逼近性能越好.
通過對人臉圖像集進行測試, 比較上述4種算法的實際性能. 實驗選擇ORL數(shù)據(jù)集, 該數(shù)據(jù)集來自于UCI機器學習庫, 是由40個不同年齡、 不同性別的人在不同時間獲取的共400幅灰度圖像組成, 即每個人的人臉圖像10幅, 每個像素的灰度大小為0~255, 每幅圖像的分辨率為64×64.
圖2 不同算法對人臉圖像集的迭代效率(A)和時間效率(B)對比Fig.2 Comparison of iterative efficiency (A) and time efficiency (B) of different algorithms on face image set
下面考慮K取不同值時SVPAGMU算法迭代效率和時間效率的比較. 選取K∈{20,30,40,50}, 實驗結果如圖3所示. 由圖3(A)可見, 相對誤差隨迭代次數(shù)的增加而逐漸改善, 并且在同一迭代次數(shù)中K值越大, 相對誤差越小; 由圖3(B)可見, 隨著K值的增加, 其運行時間逐漸延長. 因此, 若將SVPAGMU算法應用于較高維數(shù)據(jù), 需選取合適的K值, 使其在迭代效率和時間效率上均達到相對最優(yōu).
圖3 不同K值下SVPAGMU算法對人臉圖像集的迭代效率(A)和時間效率(B)對比Fig.3 Comparison of iterative efficiency (A) and time efficiency (B) of SVPAGMU algorithm on face image set under different K values
為更直觀地說明SVPAGMU算法引入?yún)?shù)θ的有效性, 圖4展示了目標函數(shù)的震蕩軌跡, 其中橫坐標表示迭代次數(shù), 縱坐標表示當前目標函數(shù)值與最優(yōu)目標函數(shù)值之差, 圖中實心圓表示快照點. 由圖4可見, SVPAGMU算法目標函數(shù)的震蕩范圍隨著迭代次數(shù)的增加而不斷減小, 并且震蕩幅度較小, 同時快照點的變化軌跡一直保持平緩下降, 算法優(yōu)化效率較高. 這是因為SVPAGMU算法在梯度更新時, 引入了調整梯度估計量的參數(shù)θ, 使其方差與偏差達到平衡, 不斷地校正梯度下降方向, 沿著最優(yōu)軌跡逼近局部最優(yōu)解, 從而提升了算法性能.
圖4 SVPAGMU算法目標函數(shù)的震蕩軌跡Fig.4 Oscillation trajectory of objective function for SVPAGMU algorithm
選取中國661個氣象臺站2004—2013年共10年的相對濕度數(shù)據(jù)[19]. 為避免氣象數(shù)據(jù)存在異常波動的影響, 使用每個臺站10年的日平均數(shù)據(jù). 因此, 對氣象變量中的相對濕度, 每個臺站的觀測值由一個d=365維列向量表示, 則所有數(shù)據(jù)可表示為365×661維的非負矩陣V. 比較SMU,SVRMU,SAGMU和SVPAGMU四種算法的性能, 分別選K=15,K=30進行對比, 實驗結果如圖5和圖6所示.
圖5 不同算法在不同K值下對平均相對濕度的迭代效率對比Fig.5 Comparison of iterative efficiency of different algorithms on average relative humidity under different K values
由圖5可見, 各算法的實際性能趨勢與人臉圖像集下的性能基本保持一致. SVPAGMU算法性能明顯優(yōu)于其他3種對比算法, 相對誤差隨迭代次數(shù)的增加而迅速減小. 當K=30時, SVPAGMU算法相對誤差已達0.112 5, 說明該算法矩陣分解結果更準確.
由圖6可見, SAGMU算法所耗費的時間最長, 可能是由于添加了梯度平均, 增加了計算復雜度. 隨著K值的增加, 每個算法的相對誤差都在減少, 但運行時間也逐漸延長. 以SVPAGMU算法為例, 當K=15和K=30時, 相對誤差分別為0.137 9,0.127 5, 運行時間分別為1.63,2.47 s. 因此, 針對不同的數(shù)據(jù)需選取合適的K值, 使其在迭代效率和時間效率上都能夠達到相對最優(yōu).
圖6 不同算法在不同K值下對平均相對濕度的時間效率對比Fig.6 Comparison of time efficiency of different algorithms on average relative humidity under different K values
為檢驗矩陣分解的逼近性能, 先對平均相對濕度數(shù)據(jù)V應用SVPAGMU算法進行非負矩陣分解, 再將其分解得到的矩陣進行重構. 通過實驗對比恢復結果與真實數(shù)據(jù)的差異, 參數(shù)設置如下: 最大迭代次數(shù)為50,K=30, 其他參數(shù)與之前設置保持一致. 部分實驗結果如圖7所示. 圖7中隨機選取的4個氣象臺站(臺站號)分別為: 內蒙古朱日和鎮(zhèn)(53276)、 西藏類烏齊縣(56128)、 浙江麗水市(58646)和新疆巴音布魯克(51542). 由圖7可見, 重構數(shù)據(jù)與原始數(shù)據(jù)非常接近, 且重構數(shù)據(jù)振幅小, 這可能是因為NMF有效地去除了噪聲. 此外, 4個氣象臺站的平均相對誤差分別為0.084 3,0.077 2,0.089 7,0.082 3, 表明該算法逼近性能較好.
圖7 部分氣象臺站相對濕度的原始數(shù)據(jù)與重構數(shù)據(jù)對比Fig.7 Comparison between original data and reconstructed data of relative humidity on some meteorological stations
綜上所述, 針對非負矩陣分解問題, 本文提出了一種隨機方差參數(shù)調整梯度乘性更新算法. 該算法將SVAG算法與MU規(guī)則相結合, 通過引入對梯度估計量進行不斷調整的參數(shù)θ, 使其偏差與方差達到平衡, 從而有效提高了矩陣分解效果. 選取已有的3種算法SMU,SVRMU,SAGMU與本文SVPAGMU算法進行實驗比較, 同時選取不同的K值進行分析, 結果驗證了本文算法的可行性與高效性.