基于不完全數(shù)據(jù)的最大似然估計方法
——EM算法

李 順 靜

(重慶工商大學(xué) 數(shù)學(xué)與統(tǒng)計學(xué)院，重慶 400067)

數(shù)理統(tǒng)計的研究總是著眼于總體，著手于樣本。當已知總體分布函數(shù)的形式，但其中有一個或多個參數(shù)未知時，可以根據(jù)樣本來推斷這個或這些未知參數(shù)，這就是參數(shù)估計?？傮w參數(shù)估計法的3個常見方法有順序統(tǒng)計量法、數(shù)字特征法和最大似然法[1]。最大似然估計(MLE)最早是由德國數(shù)學(xué)家C.F.Gauss在1821年針對正態(tài)分布所提出來的，但由于R.A.Fisher在1922年再次提出了這一想法并證明了其一些性質(zhì)，從而使得最大似然估計得到了廣泛的運用，所以一般人們將功勞歸于Fisher。最大似然估計是一種重要的、經(jīng)典的參數(shù)估計方法，人們對其研究一直方興未艾。

EM算法(expectation maximization)，即期望最大化算法。是1977 年由 Dempster、Laird 和 Rudin所提出的用于極大似然估計的迭代算法，同時為缺失數(shù)據(jù)的情況下迭代求似然函數(shù)極值提供了統(tǒng)一的框架。由于 EM 算法具有很多優(yōu)良的性質(zhì)，例如似然函數(shù)值是單調(diào)增加的和對于任意給定的初始條件都穩(wěn)定收斂等，算法很快引來了廣大學(xué)者的關(guān)注和研究。

1 最大似然估計相關(guān)理論

1.1 最大似然估計的基本概念

設(shè)總體的概率函數(shù)為p(x；θ)(當總體是離散型隨機變量時，概率函數(shù)就為分布律；當總體是連續(xù)型隨機變量時，概率函數(shù)就為密度函數(shù))，θ∈Θ，其中θ是一個未知參數(shù)或幾個未知參數(shù)所組成的參數(shù)向量，Θ是參數(shù)空間，x1，…，xn是來自總體的樣本，將樣本的聯(lián)合概率看成θ的函數(shù)，用L(θ；x1，…xn)表示，簡記為L(θ)

L(θ)=L(θ；x1，…xn)=p(x1；θ)p(x2；θ)…p(xn；θ)

(1)

L(θ)稱為樣本的似然函數(shù)。

由于L(x1，…xn；θ)是連乘式，?？紤]將似然函數(shù)兩邊取對數(shù)，從而將乘法轉(zhuǎn)化為加法，這樣會大大簡化計算量。且對數(shù)似然函數(shù)lnL(x)是x的單調(diào)增函數(shù)，lnL(x)與L(x)在同一點達到最大值，所以通常將求似然函數(shù)的極大值問題轉(zhuǎn)化為求對數(shù)似然函數(shù)的極大值。

1.2 最大似然估計失效的情形[3]

例1 設(shè)總體分布有密度函數(shù)

(2)

解分布包含一個參數(shù)θ，θ可以取任何實數(shù)值。通常稱該分布為柯西分布，其密度函數(shù)為x的函數(shù)，關(guān)于θ對稱。故θ是這個分布的中位數(shù)。

現(xiàn)在設(shè)X1，…，Xn為來自該總體的樣本，要估計θ。由于

柯西分布的一階矩不存在，更不用說更高階的矩了。因此矩法估計無法使用，若用最大似然估計，則將得出方程

(3)

在式(3)中有很多根且不容易求。因此對本例而言，最大似然法也不是為理想的方法。

在大多數(shù)情況下，還會遇到樣本數(shù)據(jù)不完整的情況，針對這些情況，本能的尋找更適合更加可行的方法來解決問題。下面所提及的EM算法就是一種專門針對缺失數(shù)據(jù)的參數(shù)估計方法[4]。

2 缺失數(shù)據(jù)模式的相關(guān)知識

缺失數(shù)據(jù)[5]的模式描述了數(shù)據(jù)矩陣中哪些數(shù)值觀測到了，哪些沒有觀測到，以及各個觀測變量間的關(guān)系。目前針對缺失數(shù)據(jù)的分析提出了很多方法，但是不同的方法針對的缺失數(shù)據(jù)模型不一樣。但總的歸納起來有4種：

2.1 基于完全記錄單元的方法

當一些變量對某些單元沒有記錄時，若這些單元數(shù)目不是占很大比例，那么就可以丟棄未完全記錄的單元，只分析有完全記錄的單元。這種方法比較容易實現(xiàn)，且在少量數(shù)據(jù)缺失的時候是可行，但若缺失的數(shù)據(jù)較多，則會導(dǎo)致嚴重的估計偏差。

2.2 加權(quán)方法

通常用預(yù)先設(shè)計好的權(quán)重來加權(quán)抽取的單元，設(shè)計權(quán)是反比于它們的選取概率。設(shè)yi是變量Y在總體單元i的值，總體均值常用如下估計量進行估計

2.3 基于借補的方法

借補值通常是指缺失值預(yù)測分布的平均值，借補方法就是要以觀測數(shù)據(jù)為基礎(chǔ)，為缺失值創(chuàng)建一個預(yù)測分布的方法。存在兩種產(chǎn)生分布的通用途徑：明確建模和模糊建模。

2.4 基于模型的方法

這類方法在大多數(shù)場合均有效。通常對觀測數(shù)據(jù)定義一個模型，然后在模型下根據(jù)適當?shù)姆植甲鐾茢?。這個方法的優(yōu)勢就是靈活，回避特殊情況。

3 EM算法的相關(guān)理論及其應(yīng)用

EM算法將數(shù)據(jù)缺失機制[6]這個概念公式化，然后用估計值代替缺失值，估計參數(shù)；假定新的參數(shù)是正確的，再估計缺失值；再估計參數(shù)，如此迭代直至收斂。

3.1 EM算法的原理

在實際中，假如完全樣本x中有數(shù)據(jù)缺失，只能觀測到缺失數(shù)據(jù)后的不完全樣本y，由于y是不完全的樣本點，其密度函數(shù)比較復(fù)雜，這也就使得它的似然函數(shù)L(θ，y)很復(fù)雜。假如用前述的MLE的方法，一般是得不到解析解，為了使問題簡單化，引入了一種新的算法——EM算法[7]來計算。

EM算法是一種迭代算法，其特點是每一步迭代都由兩步組成：

第一步：E步——求期望；

第二步：M步——求極大值。

這也是EM算法這個名字的由來。EM算法具體的實現(xiàn)如下[8]：

EM算法使用完全樣本點的似然函數(shù)Lc(θ，x)，其對數(shù)似然函數(shù)為logLc(θ，x)，但由于x缺失了部分數(shù)據(jù)而觀測不到，于是就用條件期望來代替。

假定θ的初始值為θ(0)，Θ為參數(shù)空間。

E步：Q(θ，θ(0))=Eθ(0){logLc(θ，x)}

M步：找θ(1)∈Θ，使得Q(θ(1)，θ(0))≥Q(θ，θ(0))，對于?θ∈Θ均成立。

然后將θ(1)代替E步中的θ(0)，不斷的重復(fù)E步和M步，當k+1時，有形式：

E步：Q(θ，θ(k))=Eθ(k){logLc(θ，x)}；

M步：找θ(k+1)∈Θ，使得Q(θ(k+1)，θ(k))≥Q(θ，θ(k))，對?θ∈Θ成立。

對E步與M步的不斷迭代產(chǎn)生的參數(shù)值序列{θ(k)}(不完全數(shù)據(jù)的)似然函數(shù)L(θ，y)是單調(diào)不減的，即

L(θ(k+1))≥L(αθ(k))k=0，1，2，…，

(4)

所以，當L(θ(k+1))-L(θ(k))<ε(?ε>0)時，停止迭代，此時(θ(k+1))就是要尋找的θ的MLE。

3.2 相關(guān)例子

解用y1，y2，y3，y4表示4種結(jié)果發(fā)生的次數(shù)，此時總體分布為多項分布，似然函數(shù)

(5)

∝(2-θ)y1(1-θ)y2(1+θ)y3θy4

∝θz2+y4(1-θ)z1+y2

兩邊取對數(shù)，得其對數(shù)似然函數(shù)為

l(θ；y，z)=(z2+y4)lnθ+(z1+y2)ln (1-θ)

E步：在已有的觀測數(shù)據(jù)y及第i步的估計值θ=θ(i)條件下，求基于完全數(shù)據(jù)的對數(shù)似然函數(shù)的期望(這一步就是把潛在變量z積分掉)：

Q(θ|y，θ(i))=Ezl(θ；y，z)，

M步：求Q(θ|y，θ(i))關(guān)于θ的最大值θ(i+1)，即尋找這樣的θ(i+1)，使得：

這樣就完成了θ(i)到θ(i+1)的一次迭代。如此重復(fù)上述的E步與M步，直至收斂即可得θ的MLE。

在題中，E步為

Q(θ|y，θ(i))=E(z2|y，θ(i)+y4)lnθ+E(z1|y，θ(i)+y2)ln (1+θ)=

M步：在上式兩邊關(guān)于θ求導(dǎo)，并令其等于0，即為

解之，得如下的迭代公式：

開始時，任意取一個初始值進行迭代，如θ(0)=0.5，按照上述的EM算法的E步M步，最后經(jīng)過13次迭代可得θ的MLE為0.606 747。

4 結(jié)束語

根據(jù)應(yīng)用者的實際需求出發(fā)，先介紹了極大似然估計，當所收集的數(shù)據(jù)是缺失數(shù)據(jù)時，又引入了一種新的迭代算法——EM算法。EM算法雖然具有穩(wěn)定收斂的性質(zhì)，但收斂的速度緩慢。當缺失數(shù)據(jù)比例較大時候，它的收斂速率也極其緩慢。將在進一步的工作中，對上述問題作進一步的研究。

參考文獻：

[1] 陳文清.極大似然法的基本思想及其應(yīng)用[J].韶關(guān)大學(xué)學(xué)報：自然科學(xué)版，1996，12(4)：24-27

[2] 宋明珠.關(guān)于條件概率及其應(yīng)用教學(xué)方法的研究[J]. 重慶工商大學(xué)學(xué)報：自然科學(xué)版，2013，30(1)：31-34

[3] 費紹金. 對極大似然估計中幾個問題的探討[J]. 四川文理學(xué)院學(xué)報：自然科學(xué)版，2008，9(5)：1-3

[4] 林鴻. EM算法的改進及其在基因序列分析中的應(yīng)用[D]. 福州大學(xué)碩士學(xué)位論文，2005

[5] 曾傳璜，張鑫，張晶晶，王宏淵. EM算法的研究[J]. 軟件導(dǎo)刊，2008，9(7)：97-98

[6] 孟麗新. 基于EM算法的約束條件下參數(shù)的估計[D]. 東北師范大學(xué)碩士學(xué)位論文，2006

[7] 茆詩松，程依明，濮曉龍. 概率論與數(shù)理統(tǒng)計[M]. 北京：高等教育出版社，2011

[8] 蘇良軍. 高等數(shù)理統(tǒng)計[M]. 北京：北京大學(xué)出版社，2007

[9] 王兆軍，劉民千，鄒長亮，等. 計算統(tǒng)計[M]. 北京：人民郵電出版社，2008