基于GEMS算法的潛變量高斯圖模型結(jié)構(gòu)學(xué)習(xí)

鄭倩貞，徐平峰，曹 蕾

(長春工業(yè)大學(xué)數(shù)學(xué)與統(tǒng)計學(xué)院，吉林 長春 130012)

高斯圖模型能夠清晰直觀地反應(yīng)變量間的相互關(guān)系，被廣泛應(yīng)用于高維情形.在對實際問題進(jìn)行圖模型結(jié)構(gòu)學(xué)習(xí)時，僅考慮觀測變量有時并不能正確反應(yīng)變量間的相互關(guān)系，因此需考慮潛變量對可觀測變量的影響，在給定潛變量時探討可觀測變量間的條件關(guān)系.Chandrasekaran等[1]將可觀測變量的邊緣協(xié)方差陣的逆陣分解為一個稀疏陣和一個低秩陣，提出了懲罰似然的方法，對稀疏陣和低秩陣分別施加1范數(shù)及核范數(shù)懲罰，并結(jié)合凸優(yōu)化和代數(shù)幾何對潛變量圖模型選擇問題進(jìn)行了研究.Yuan[2]基于Chandrasekaran等人的研究，將懲罰似然的核范數(shù)懲罰項替換為對低秩陣的秩的約束條件，提出了潛變量GLasso(LVglasso)方法，并結(jié)合EM算法對高維情形下的潛變量圖模型選擇問題進(jìn)行了模擬研究.Lauritzen等[3]對懲罰似然做了與Yuan相似的處理，采用插補(bǔ)的方法，結(jié)合EM算法和GLasso算法進(jìn)行模型選擇.

但上述與EM相結(jié)合的方法需要先給定1組正則化參數(shù)，然后對每個正則化參數(shù)利用EM算法求懲罰似然的最小值點.如果正則化參數(shù)選取不當(dāng)，將會導(dǎo)致每次迭代的模型離真模型越來越遠(yuǎn)，而且增加計算時間.本文基于期望模型選擇(EMS)算法[4]的思想，在每次迭代時從候選模型中選取期望信息準(zhǔn)則最小的模型作為下一步的當(dāng)前模型，下一次迭代時在當(dāng)前模型下求候選模型的期望信息準(zhǔn)則的值.但由于可能的模型太多，在模型選擇時遍歷全部模型不可行，因此只選出部分模型作為候選模型.這里的候選模型也可以通過1組正則化參數(shù)來確定，但每次的正則化參數(shù)不一定相同.稱這種方法為廣義期望模型選擇(GEMS)算法.模擬實驗顯示，基于GEMS的LVglasso方法收斂速度快，計算時間短.

1 高斯圖模型

Yuan等[5]提出通過最小化負(fù)1懲罰對數(shù)似然的方法去估計高斯圖模型的協(xié)方差逆陣Ω，懲罰似然為

2 LVglasso方法

其中：S-L?0表示S-L為正定矩陣；L0表示L為非負(fù)定矩陣；為ΩO的估計，為的估計；為可觀測樣本的對數(shù)似然函數(shù)，即

基于Chandrasekaran等[1]提出的懲罰似然，Yuan[2]提出了計算更加方便的LVglasso方法：

其中：0≤r≤p，S?=S-diag(S).限制條件rank(L)≤r相當(dāng)于假設(shè)存在r個潛變量.

考慮完全數(shù)據(jù)x=(x1，…，xn)T=(xO，xH)，xi=(xO，i，xH，i)T，其中xO，i為第i個樣本的可觀測數(shù)據(jù)，xH，i為第i個樣本的不可觀測數(shù)據(jù)，i=1，…，n.Ω的LVglasso估計為

3 GEMS算法

EMS算法[4]是一種迭代算法，用于處理缺失數(shù)據(jù)情形下的模型選擇問題.該算法的每次迭代都需給定當(dāng)前模型Mc和當(dāng)前模型下的參數(shù)θc∈ΘMc，并依次進(jìn)行期望步(E步)和模型選擇步(MS步)，直至滿足停止準(zhǔn)則得到最優(yōu)的模型估計M*和參數(shù)估計θ*∈ΘM*.本文的GEMS算法與EMS算法類似，不同之處在于：GEMS算法的MS步不遍歷全部模型，而是通過GLasso算法找出候選模型，候選模型可由1組正則化參數(shù)來確定，且每次迭代的正則化參數(shù)可能不同.在這些候選模型中選擇期望BIC最小的模型，然后將該模型及其對應(yīng)的參數(shù)作為下一次迭代的當(dāng)前模型和當(dāng)前參數(shù).從部分而非全部模型中選擇最優(yōu)模型可大大減小計算成本，提高計算效率，尤其是在高維問題中.考慮潛變量高斯圖模型的結(jié)構(gòu)學(xué)習(xí)問題，現(xiàn)有當(dāng)前模型G(t)及當(dāng)前參數(shù)Ω(t)，則第(t+1)次迭代如下：

(1) E步

進(jìn)而，有Q函數(shù)

Q(G，Ω|G(t)，Ω(t))=Et(-2(G，Ω))+log(n)dfG= -nlogdetΩ+ntr[ΩEt(Σn)]+log(n)dfG.

(2) MS步

因此，可得到1組候選模型G={G(t)，Gλm，m=1，…，k}.針對每個候選模型G∈G，計算Q函數(shù)，從而得到Ω的估計Ω(t+1)=argminΩQ(G，Ω|G(t)，Ω(t))及其對應(yīng)的圖模型G(t+1)=argminG∈GQ(G，Ω(t+1)|G(t)，Ω(t)).

將G(t+1)，Ω(t+1)作為下一次迭代的當(dāng)前模型和當(dāng)前參數(shù)，重復(fù)以上步驟直至滿足停止準(zhǔn)則.

4 模擬實驗及實例分析

在模擬實驗中，考慮不同情形下潛變量圖模型的結(jié)構(gòu)學(xué)習(xí)問題，對GEMS算法和EM算法在LVglasso估計求解問題上的模擬結(jié)果進(jìn)行了比較.模擬實驗覆蓋p=48，98，148，198，h=2，r=2，5，n=500，1 000共16種情形，每種情形各模擬50次.真模型產(chǎn)生機(jī)制與Yuan[2]類似，兩者的不同之處體現(xiàn)在對潛變量的設(shè)定上.在本文真模型中，每個潛變量至少和2個可觀測變量、至多和(p-1)個可觀測變量有關(guān).值得注意的是，當(dāng)在進(jìn)行n=1 000，p=198，r=2情形設(shè)定下的第38次EM算法模擬時，由R中l(wèi)vglasso函數(shù)產(chǎn)生的Ω迭代初值為非對稱陣，所以該種情形只模擬了37次.16種情形的CPU平均運行時間如表1所示.

表1 不同情形下CPU平均運行時間 s

從表1中可看出，對于任意一種情況，EM算法的運行時長都要遠(yuǎn)大于GEMS算法，達(dá)到5倍、10倍，甚至是15倍的差距.GEMS算法大大提升了潛變量圖模型選擇的速度.本文用于評價算法性能的指標(biāo)為：

其中tp，tn，fp，fn分別為真陽類、真陰類、假陽類、假陰類的個數(shù).圖1給出了所有情形下tpr，ppv和mcc的箱線圖.總體上看，GEMS較EM有更優(yōu)的表現(xiàn)，但在極個別情況如n=500，p=198，r=2或5時EM的tpr值較大.同時可看出，樣本量越大，潛變量個數(shù)的假設(shè)越接近真實模型，模型推斷就越準(zhǔn)確.

白色箱子代表EM算法，灰色箱子代表GEMS算法，橫軸為真模型可觀測變量的個數(shù).

基于GEMS算法，對Wille等[7]論文中擬南芥植物類異戊二烯生物合成相關(guān)基因的數(shù)據(jù)進(jìn)行了潛變量高斯圖模型結(jié)構(gòu)學(xué)習(xí)，估計了各基因間的條件相關(guān)性.該數(shù)據(jù)的數(shù)據(jù)來源為https：∥static-content.springer.com/esm/art%3A10.1186%2Fgb-2004-5-11-r92/MediaObjects/13059_2004_896_MOESM1_ESM.txt，數(shù)據(jù)中共有118個樣本，每個樣本包含39個基因表達(dá).若假設(shè)的潛變量個數(shù)不同，則推斷出的各基因間的條件相關(guān)性也不同.假設(shè)潛變量個數(shù)為r=1或r=3時的估計結(jié)果如圖2所示.當(dāng)r=1時，共估計出174條邊，算法運行時間約為43 s；當(dāng)r=3時，共估計出38條邊，算法運行時間約為10 s.

圖2 r=1和r=3時的基因圖模型

5 結(jié)語

本文簡要介紹了高斯圖模型及潛變量高斯圖模型下的LVglasso方法，給出了GEMS算法結(jié)合LVglasso下潛變量圖模型結(jié)構(gòu)學(xué)習(xí)的算法步驟，并從模擬實驗的角度比較了GEMS算法和EM算法在潛變量圖模型選擇問題上的優(yōu)劣.通過多種不同情形下的模擬實驗，可以發(fā)現(xiàn)，樣本量越大，潛變量個數(shù)的假設(shè)越接近真實模型，模型推斷就越準(zhǔn)確.結(jié)合tpr，ppv，mcc以及CPU平均運行時間，無論在何種模擬情形，GEMS算法在模型選擇上的表現(xiàn)較EM算法優(yōu)越.