鄭倩貞，徐平峰，曹 蕾

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

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

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

1 高斯圖模型

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

2 LVglasso方法

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

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

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

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

3 GEMS算法

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

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

4 模擬實(shí)驗(yàn)及實(shí)例分析

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

表1 不同情形下CPU平均運(yùn)行時(shí)間 s

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

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

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

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

圖2 r=1和r=3時(shí)的基因圖模型

5 結(jié)語(yǔ)

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