基于自適應(yīng)懲罰的潛變量高斯圖模型結(jié)構(gòu)學(xué)習(xí)

鄭倩貞, 徐平峰

(1. 長春工業(yè)大學(xué) 數(shù)學(xué)與統(tǒng)計學(xué)院, 長春 130012; 2. 東北師范大學(xué) 前沿交叉研究院, 長春 130024)

圖模型的結(jié)構(gòu)學(xué)習(xí)問題, 是指從給定的數(shù)據(jù)集中估計出反映隨機變量間獨立性結(jié)構(gòu)的圖, 也稱為模型選擇問題[1]. 高斯圖模型的參數(shù)估計和模型選擇問題等價于協(xié)方差逆陣的估計問題及協(xié)方差逆陣中零元素的識別問題, 協(xié)方差逆陣的元素為零表示對應(yīng)兩個隨機變量在給定其余隨機變量時具有條件獨立性[2]. 隨著高維復(fù)雜數(shù)據(jù)的不斷增多, 研究稀疏圖模型的估計問題尤為重要. 為估計稀疏協(xié)方差逆陣, 常用的方法為懲罰似然方法, 例如: 文獻[3-5]基于懲罰似然方法研究了變量完全觀測時的高斯圖模型選擇問題. 但在實際生活中, 觀測變量通常與一些潛在的不可觀測變量(即潛變量)相關(guān). 此時, 觀測變量間的圖模型結(jié)構(gòu)不一定具有稀疏性, 需考慮潛變量對觀測變量的影響, 在給定潛變量時探討觀測變量間的條件獨立性. 通常情況下, 潛變量會對圖模型選擇問題造成很大困難, 因為潛變量的個數(shù)以及潛變量與觀測變量之間的關(guān)系可能都是未知的.

Chandrasekaran等[6]首先對含潛變量的圖模型選擇問題進行了研究, 提出了當(dāng)潛變量和觀測變量服從聯(lián)合高斯分布時, 觀測變量邊緣協(xié)方差陣的逆陣可分解為稀疏陣和低秩陣之和, 其中稀疏陣對應(yīng)給定潛變量時觀測變量間的條件獨立性. 基于稀疏低秩分解, 文獻[6]提出了正則化極大似然分解框架進行潛變量圖模型的結(jié)構(gòu)學(xué)習(xí), 其中正則項包含對稀疏陣施加的LASSO(least absolute shrinkage and selection operator)懲罰以及對低秩陣施加的核范數(shù)懲罰. 由于這兩種懲罰均不含自適應(yīng)權(quán)重, 因此稱這兩種懲罰為非自適應(yīng)懲罰, 由此得到的稀疏陣和低秩陣估計稱為非自適應(yīng)懲罰似然估計. 文獻[6]將優(yōu)化問題視為對數(shù)-行列式半正定規(guī)劃問題, 采用Wang等[7]提出的對數(shù)-行列式近端點算法(log-determinant proximal point algorithm, LogdetPPA)進行求解. 針對文獻[6]中的優(yōu)化問題, Ma等[8]提出了兩種交替方向的求解方法: 一種是將目標(biāo)問題視為一致性優(yōu)化問題, 利用傳統(tǒng)的交替方向乘子法(alternating direction method of multipliers, ADMM)[9]進行求解; 另一種是基于近端梯度的交替方向法(proximal gradient-based alternating direction method, PGADM). 這兩種方法都是把包含3個變量的原始問題轉(zhuǎn)換成包含2個變量塊的優(yōu)化問題. PGADM算法速度比LogdetPPA更快, 并具有全局收斂性. 此外, Ma[10]提出了一種交替近端梯度方法求解文獻[6]中的優(yōu)化問題, 并證明了算法的全局收斂性; Meng等[11]對文獻[6]的方法做了進一步的理論研究; 文獻[12-15]等也對含潛變量圖模型結(jié)構(gòu)學(xué)習(xí)進行了相關(guān)研究.

研究表明, 非自適應(yīng)懲罰得到的估計, 包括LASSO懲罰得到的稀疏協(xié)方差逆陣估計[16]以及核范數(shù)懲罰在降秩回歸問題中得到的低秩系數(shù)陣估計[17], 均存在偏差較大的問題, 而自適應(yīng)懲罰通?？山档凸烙嬈? 因此, 本文采用自適應(yīng)懲罰似然方法處理含潛變量的高斯圖模型結(jié)構(gòu)學(xué)習(xí)問題, 對稀疏陣和低秩陣進行估計, 以得到給定潛變量時觀測變量間的條件獨立關(guān)系. 本文對稀疏陣部分施加自適應(yīng)LASSO懲罰[18], 對低秩陣部分施加自適應(yīng)核范數(shù)懲罰[17]. 與非自適應(yīng)懲罰的優(yōu)化問題類似, 采用ADMM算法優(yōu)化求解自適應(yīng)懲罰似然的最小化問題, 并且在求解過程中仍具有顯式表達式, 以確保算法的計算效率. 本文模擬比較了自適應(yīng)懲罰與非自適應(yīng)懲罰在潛變量高斯圖模型結(jié)構(gòu)學(xué)習(xí)和參數(shù)估計上的性能. 結(jié)果表明, 自適應(yīng)懲罰顯著優(yōu)于非自適應(yīng)懲罰, 有效降低了稀疏陣和低秩陣的估計偏差, 能更準(zhǔn)確地學(xué)習(xí)觀測變量間的條件獨立關(guān)系.

1 含潛變量的無向高斯圖模型及其自適應(yīng)懲罰似然

無向高斯圖模型是指與無向圖G=(V,E)相關(guān)的多元正態(tài)分布模型, 其中頂點集V中的每個頂點v∈V表示一個高斯隨機變量Xv, 邊集E中的每條邊(u,v)∈E蘊含了給定其余隨機變量Xv′(v′∈V{u,v})時隨機變量Xu與Xv之間的條件相關(guān)性.當(dāng)采用無向高斯圖模型對實際問題進行分析時, 可能會出現(xiàn)隨機變量為潛變量的情況, 即變量不存在觀測值.因此, 本文考慮含潛變量的無向高斯圖模型結(jié)構(gòu)學(xué)習(xí)問題, 考察給定潛變量時觀測變量間的條件獨立性關(guān)系.

基于懲罰似然方法, 考慮對稀疏陣S施加自適應(yīng)LASSO懲罰[18], 對低秩陣L施加自適應(yīng)核范數(shù)懲罰[17].自適應(yīng)LASSO懲罰和自適應(yīng)核范數(shù)懲罰分別為加權(quán)版的LASSO懲罰和核范數(shù)懲罰. 最初, 自適應(yīng)LASSO懲罰用于線性回歸的變量選擇問題, 自適應(yīng)核范數(shù)懲罰用于高維多元降秩回歸問題. 對于潛變量無向高斯圖模型問題, 自適應(yīng)懲罰似然估計為

其中:S-L?0表示矩陣S-L為正定陣,L0表示矩陣L為半正定陣;λ>0和β>0均為調(diào)整參數(shù), 參數(shù)λ控制懲罰的強度, 參數(shù)β用于權(quán)衡稀疏陣和低秩陣兩項懲罰;σ1(L)≥…≥σp(L)≥0為矩陣L的奇異值;為自適應(yīng)LASSO懲罰,為自適應(yīng)核范數(shù)懲罰,Sij為矩陣S中第(i,j)位置的元素.當(dāng)自適應(yīng)懲罰項中的權(quán)重W=(Wij)p×p及w=(w1,w2,…,wp)T滿足Wij=1,wi=1(i,j=1,2,…,p)時, 式(1)中的估計退化為非自適應(yīng)懲罰似然估計.

自適應(yīng)懲罰要求權(quán)重Wij(或wi)隨著|Sij|(或σi(L))的增大而減小, 以降低估計偏差.為使權(quán)重具有上下界, 本文采用類似于文獻[19]的權(quán)重形式, 即對任意的i,j=1,2,…,p,

(2)

2 優(yōu)化方法

下面采用ADMM算法求解式(1)中的最小化問題, 并給出調(diào)整參數(shù)λ和β的選取方法.

2.1 ADMM算法

ADMM算法是一種迭代算法, 廣泛應(yīng)用于求解線性約束下的優(yōu)化問題中, 如稀疏低秩分解問題[20]. 該算法通過分解-協(xié)調(diào)過程將一個大的全局問題分解成多個小的局部子問題, 從而通過整合各子問題的解得到原問題的最優(yōu)解[9]. 考慮約束條件M=S-L, 將式(1)中的目標(biāo)問題轉(zhuǎn)化成ADMM形式下的優(yōu)化問題:

于是, 增廣Lagrange函數(shù)為

其中Y∈p×p為Lagrange乘子,μ>0為懲罰參數(shù), 〈·,·〉表示矩陣內(nèi)積, ‖·‖F(xiàn)表示矩陣的Frobenius范數(shù).給定第t次ADMM迭代的估計值(Mt,St,Lt,Yt,μt), 則第(t+1)次迭代的更新步驟為

(3)

其中ρ>1為放大因子.下面給出關(guān)于矩陣M,S,L子問題的解.

1) 更新M.

式(3)中矩陣M的優(yōu)化子問題可表示為

A=Udiag(e1(A),…,ep(A))UT,

2) 更新S.

上述優(yōu)化問題完全可分, 即對任意i,j=1,2,…,p, 均有

其中sign(·)為符號函數(shù), (·)+=max{·,0}.

3) 更新L.

令B=St+1-Mt+1-Yt/μt, 則式(3)中稀疏陣L的優(yōu)化子問題可寫成

對矩陣B做特征值分解

B=Qdiag(e1(B),…,ep(B))QT,

其中e1(B)≥…≥ep(B)為矩陣B的特征值,Q為正交陣,QQT=QTQ=I.從而可推得

Lt+1=Qdiag((e1(B)-λw1/μt)+,…,(ep(B)-λwp/μt)+)QT.

2.2 調(diào)整參數(shù)的選取

調(diào)整參數(shù)λ和β控制模型的復(fù)雜度, 不同的參數(shù)設(shè)定可得到不同稀疏度的S及不同秩的L.對于優(yōu)化問題(1), 可采用K折交叉驗證的方式選擇最優(yōu)參數(shù)組合(λ*,β*).將樣本分成K折互不相交的子集, 記為Tk(k=1,2,…,K).定義K折交叉驗證的得分函數(shù)為

3 數(shù)值模擬

在模擬實驗中, 對自適應(yīng)懲罰(γ1≠0,γ2≠0)和非自適應(yīng)懲罰(γ1=0,γ2=0)在含潛變量的無向高斯圖模型上的模型選擇和參數(shù)估計性能進行比較.對于自適應(yīng)懲罰情況, 僅考慮γ1=γ2=1.

為評價模型選擇性能, 比較真陽率(TPR)、 陽性預(yù)測率(PPV)及馬修斯相關(guān)系數(shù)(MCC):

其中TP為真陽類個數(shù), TN為真陰類個數(shù), FP為假陽類個數(shù), FN為假陰類個數(shù). 表1列出了不同情形下TPR,PPV和MCC的均值及標(biāo)準(zhǔn)差. 由表1可見, 在所有情形下由自適應(yīng)懲罰得到的PPV和MCC均顯著優(yōu)于非自適應(yīng)懲罰. 對于TPR, 當(dāng)n=500時非自適應(yīng)懲罰性能更好, 但當(dāng)樣本量增大(即n=1 000)時, 兩種懲罰的TPR幾乎一樣好, 均接近于1. 因此, 基于自適應(yīng)懲罰似然的潛變量圖模型結(jié)構(gòu)學(xué)習(xí)性能更好.

表1 不同情形下TPR,PPV,MCC的均值及標(biāo)準(zhǔn)差Table 1 Means and standard deviations of TPR,PPV,MCC in different situations

為比較參數(shù)估計結(jié)果, 考慮矩陣M,S,L的估計誤差, 分別為

表2列出了不同情形下參數(shù)估計誤差的均值及標(biāo)準(zhǔn)差.由表2可見, 所有情形下自適應(yīng)懲罰的參數(shù)估計誤差均小于非自適應(yīng)懲罰.由自適應(yīng)懲罰得到的參數(shù)估計更接近于真實值.

表2 不同情形下矩陣M,S,L的估計誤差均值及標(biāo)準(zhǔn)差Table 2 Means and standard deviations of estimation errors of matrices M,S,L in different situations

4 實例分析

下面采用自適應(yīng)懲罰似然方法對枯草芽孢桿菌核黃素(維生素B2)生產(chǎn)數(shù)據(jù)集(該數(shù)據(jù)集可在文獻[21]的補充材料中下載)進行圖模型結(jié)構(gòu)學(xué)習(xí). 該數(shù)據(jù)集的樣本量為71, 變量個數(shù)為4 089, 其中有4 088個變量表示不同基因表達水平的對數(shù), 1個變量表示核黃素生產(chǎn)率的對數(shù). 本文僅對文獻[21]中簡化后的數(shù)據(jù)集riboflavinv100.csv做圖模型推斷, 其中包含經(jīng)驗方差最大的100個基因表達水平變量和1個測量核黃素生產(chǎn)率的變量.

為比較含潛變量和不含潛變量(即固定低秩陣L=0)時觀測變量間的圖結(jié)構(gòu), 將樣本分為訓(xùn)練集和測試集, 其中訓(xùn)練集包含57個樣本, 測試集包含14個樣本.首先, 基于訓(xùn)練集樣本得到稀疏陣和低秩陣的估計; 然后, 基于測試集樣本計算負的對數(shù)似然.圖1為由自適應(yīng)懲罰(γ1=γ2=1)得到的圖模型選擇結(jié)果, 其中(A)為給定潛變量時觀測變量間的條件獨立性關(guān)系, (B)為不考慮潛變量時觀測變量間的條件獨立性關(guān)系. 圖1中橫軸從左至右(縱軸從上至下)依次對應(yīng)第1～101個觀測變量, 黑色表示對應(yīng)觀測變量間條件相關(guān), 白色表示對應(yīng)變量間條件獨立. 當(dāng)考慮潛變量時, 估計出的圖模型共包含531條邊(占觀測變量對總數(shù)的10.5%)以及13個潛變量(對應(yīng)低秩陣的秩), 基于測試集的負對數(shù)似然為49.54. 當(dāng)不考慮潛變量時, 估計出的圖模型共包含1104條邊(占觀測變量對總數(shù)的21.9%), 基于測試集的負對數(shù)似然為50.48. 顯然, 考慮潛變量時能得到更稀疏的圖模型, 且負對數(shù)似然更小, 表明考慮潛變量對觀測變量的影響能更好地擬合枯草芽孢桿菌核黃素生產(chǎn)數(shù)據(jù)集.

綜上所述, 本文提出了一種基于自適應(yīng)懲罰似然的潛變量高斯圖模型結(jié)構(gòu)學(xué)習(xí)方法. 首先, 在觀測似然后面加上兩項自適應(yīng)懲罰項, 分別為稀疏陣部分的自適應(yīng)LASSO懲罰以及低秩陣部分的自適應(yīng)核范數(shù)懲罰; 然后, 通過ADMM算法最小化懲罰似然以求解稀疏陣和低秩陣的參數(shù)估計, 從而得到給定潛變量時觀測變量間的條件獨立關(guān)系. 在模擬實驗中, 通過與非自適應(yīng)懲罰的比較, 驗證了自適應(yīng)懲罰似然方法在模型選擇和參數(shù)估計方面性能均更好.