孫聚波，徐平峰

(長春工業(yè)大學(xué) 基礎(chǔ)科學(xué)學(xué)院，吉林 長春130012)

0 引 言

在統(tǒng)計(jì)分析中，列聯(lián)表是分類數(shù)據(jù)的一種常用、直觀的表示方法［1］.很多高維分類數(shù)據(jù)都具有某種特殊結(jié)構(gòu)，通?？梢杂脠D模型表示［2］.其中，模型選擇是統(tǒng)計(jì)推斷的重要方面，許多學(xué)者提出了各種各樣的方法，例如，Akaike 提出了基于似然的AIC 準(zhǔn)則［3］，Schwarz 提出了BIC 準(zhǔn)則［4］，Edwards和Havranek［5］提出了從大規(guī)模模型族中選擇出可接受的一組模型等等.當(dāng)數(shù)據(jù)有些值缺失時(shí)，模型選擇的最新成果是，Jiang，Nauyen 和Rao 于2014年提出的E－MS 算法［6］.本文對(duì)于含有缺失值的列聯(lián)表數(shù)據(jù)，在E－MS 算法的框架下，進(jìn)行圖模型的模型選擇，其中廣義信息準(zhǔn)則采用BIC 準(zhǔn)則.對(duì)于含3、4 和5 個(gè)分類變量的情形進(jìn)行了模擬實(shí)驗(yàn).

1 列聯(lián)表及圖模型

1.1 列聯(lián)表的概率結(jié)構(gòu)

一般來說，觀測數(shù)據(jù)按兩個(gè)或多個(gè)屬性分類時(shí)所列出的頻數(shù)表稱為列聯(lián)表.本文中，令V 表示分類變量的集合.對(duì)分類變量γ ∈V，Yγ表示γ 的有限的水平集，YV=×γ∈VYγ表示分類集V 的水平集，則列聯(lián)表中的一個(gè)格子對(duì)應(yīng)集合YV中的一個(gè)元素y=(yγ)γ∈V.把n 個(gè)觀測個(gè)體按V 進(jìn)行分類，令計(jì)數(shù)xy為落入格子y 的觀測頻數(shù)，則這些計(jì)數(shù)xy的集合便構(gòu)成了一個(gè)列聯(lián)表，記作X={xy|y ∈YV}，且設(shè)一個(gè)個(gè)體落入格子y 的概率為p(y)，滿足p(y)3≥0，且于是整個(gè)表的聯(lián)合概率分布是一個(gè)多項(xiàng)分布:

在高維列聯(lián)表中，飽和模型的參數(shù)個(gè)數(shù)一般大于樣本量，不僅統(tǒng)計(jì)上無法處理，計(jì)算上也不可行.但事實(shí)上，很多高維數(shù)據(jù)都具有某種特殊結(jié)構(gòu)，可用圖模型表示，以減少參數(shù)個(gè)數(shù).

1.2 圖模型及其概率分解

在圖模型中，隨機(jī)變量由圖的頂點(diǎn)表示，有直接關(guān)聯(lián)隨機(jī)變量，對(duì)應(yīng)的頂點(diǎn)用邊相連，這樣構(gòu)成一個(gè)圖G(V，E)，這里V 為頂點(diǎn)集，E 為邊集.相對(duì)于圖G 滿足馬爾科夫性的概率分布族，即為圖模型［3］.在圖模型中，可利用圖的語言表示概率統(tǒng)計(jì)相關(guān)問題，并依據(jù)圖論的理論和算法幫助進(jìn)行概率統(tǒng)計(jì)推斷［2］.

在圖G 中，如果子集c ?V 中任意兩個(gè)頂點(diǎn)都是相連的，則稱c 是完全的.如果完全子集c 相對(duì)于包含運(yùn)算而言是最大的，稱它為團(tuán)，K(G)為圖G的所有團(tuán)的集合.對(duì)于團(tuán)c，稱yc=(yγ)γ∈C為邊緣格子，相應(yīng)的邊緣計(jì)數(shù)為，邊緣概率為許多學(xué)者指出，馬爾科夫性與概率的基于團(tuán)的對(duì)數(shù)展開形式之間存在等價(jià)關(guān)系.其中，Lauritzen［2］給出了嚴(yán)格的數(shù)學(xué)敘述:概率p 相對(duì)于圖G 滿足馬爾科夫性當(dāng)且僅當(dāng)存在函數(shù)φc(yc)滿足將該式變形得p(y)=于是p 的參數(shù)為θG={φc(yc)，c ∈K(G)}.分解后參數(shù)個(gè)數(shù)大大減少.例如，每個(gè)變量取2 個(gè)值，則p(y)分解后需要個(gè)參數(shù)來描述，遠(yuǎn)遠(yuǎn)小于飽和模型下的參數(shù)個(gè)數(shù)2|v|.

1.3 圖模型的極大似然估計(jì)

圖模型的統(tǒng)計(jì)推斷中，求參數(shù)的極大似然估計(jì)是非常重要的一個(gè)方面.對(duì)于列聯(lián)表X，對(duì)數(shù)似然為:

Lauritzen［3］證明了上式中極大似然估計(jì)^p 的唯一性，給出了相應(yīng)的似然方程組:

上述似然方程組往往可采用迭代比例擬合(IPS)算法［3，4］來求解.在投往CSDA 的論文中給出了基于團(tuán)分劃改進(jìn)的IPS 算法［7］.

當(dāng)數(shù)據(jù)含缺失值時(shí)，記觀測數(shù)據(jù)為Xobs，此時(shí)可用EM 算法求極大似然估計(jì).設(shè)qtG 為當(dāng)前參數(shù)值，則E 步為求，M步對(duì)Q 函數(shù)關(guān)于團(tuán)求導(dǎo)得，對(duì)所有的yc∈Yc，c ∈K(G)，

2 基于E－MS 算法的缺失數(shù)據(jù)的圖

運(yùn)用E－MS 算法進(jìn)行缺失數(shù)據(jù)的圖模型選擇.首先規(guī)定圖模型選擇空間為所有可能的圖模型的集合G.在E－MS 算法的E 步中，當(dāng)前圖模型Gc以及Gc的當(dāng)前參數(shù)估計(jì)，令θG為候選圖模型G 的參數(shù).定義Q(G，θG)為給定觀測數(shù)據(jù)Xobs下，完全數(shù)據(jù)X 的對(duì)數(shù)似然的條件期望，即:

圖模型G 的改變對(duì)BICG值的影響遠(yuǎn)大于θG的變化對(duì)BICG值的影響，所以在模擬中，若連續(xù)5 次迭代的最優(yōu)圖模型^Gopt都為同一個(gè)模型，就停止迭代.

3 模擬研究

表1 E－MS 算法選擇的模型與真模型相同的比率

在表1 中，隨著變量個(gè)數(shù)增加，經(jīng)E－MS 算法選擇的模型與真模型相同的比率有所下降，這是因?yàn)楹蜻x模型個(gè)數(shù)成指數(shù)型增長.實(shí)際工作中，可根據(jù)其他先驗(yàn)知識(shí)或方法縮小候選模型的空間，提高準(zhǔn)確率.另外，表2 給出了多邊比率和少邊比率.其中多邊比率為選擇的多余的邊數(shù)除以真模型的邊數(shù)，再求平均值;少邊比率為沒有找的真模型中的邊數(shù)除以真模型的邊數(shù)，再求平均值.由表2 可分析出，數(shù)據(jù)缺失概率越大，經(jīng)E－MS 算法選擇的模型與真模型不同的邊數(shù)越大;其次隨著樣本量的增加，相同缺失概率下，E－MS 算法選擇的模型與真模型邊數(shù)不同的比率在減小.總的來說，無論在哪種情形，多邊比率與少邊比率之和不大于0.4，也就是說我們可以至少0.6 的概率正確的判斷邊是否存在.

表2 E－MS 算法選擇的模型與真模型的比較

4 結(jié) 語

通過對(duì)模擬研究的分析可知，E－MS 算法在候選模型個(gè)數(shù)較少時(shí)具有良好的性能，當(dāng)候選模型個(gè)數(shù)非常多時(shí)，E－MS 算法的性能下降較多，這也是現(xiàn)在大多模型選擇方法的共同不足之處.對(duì)此，可結(jié)合其它方法或先驗(yàn)知識(shí)縮小候選模型集合，然后再利用E－MS 算法進(jìn)行模型選擇.

［1］ Agresti，A 著，張淑梅，王睿，曾莉譯.屬性數(shù)據(jù)分析引論［M］，第二版.北京:高等教育出版社，2008.

［2］ Lauritzen，S L.Graphical Models［M］.Oxford:Clarendon Press，1996.

［3］ Akaike，H.An Information Criterion［J］.Math Sci，1976.

［4］ Schwarz，G.Estimating the Dimension of a Model［J］.Annals of Statistics，1978，6(2):461－464.

［5］ Edwards，D，Havranek，T.A Fast Model Selection Procedure for Large Families of Models［J］.Journal of the American Statistical Association，1987，82(397):205－213.

［6］ Jiang，J，Nauyen，T，Rao，J S.The E－MS Algorithm:Model Selection with Incomplete Data［J］.To appear in Journal of the American Statistical Association，2014.

［7］ Xu，P F，Sun，J，Shan，N.Local Computations of the Iterative Proportional Scaling procedure for Hierarchical Models［J］.Revised in Computational Statistics and Data Analysis，2015.