基于EM算法的混合t-分布模型參數(shù)估計(jì)

王小英，李迎華，楊雪梅

（華北電力大學(xué) 數(shù)理學(xué)院，北京 102206）

0 引言

混合分布模型是分析復(fù)雜現(xiàn)象的一個(gè)靈活而強(qiáng)有力的建模工具，它提供了用簡(jiǎn)單結(jié)構(gòu)模擬復(fù)雜密度的一個(gè)有效方法，從而將聚類問題轉(zhuǎn)化為統(tǒng)計(jì)分析問題?；旌细咚鼓Ｐ陀捎谟?jì)算上的便利，在聚類方法研究中應(yīng)用較多，但通常我們收集到的很多數(shù)據(jù)并不是嚴(yán)格的服從正態(tài)分布，而是較明顯的服從重尾分布?；旌蟭-分布模型由于其具有較長的尾巴，可對(duì)重尾點(diǎn)和異常點(diǎn)有效地降低權(quán)值，因此，相對(duì)于高斯分布混合模型，可以獲得較強(qiáng)的精度和穩(wěn)健性。

文獻(xiàn)[1]中Dempster等提出的EM算法成為了混合模型參數(shù)估計(jì)極有效的工具。Peel和McLachlan在文獻(xiàn)[2]中指出EM算法可以獲得有限混合模型任意分布的極大似然估計(jì)。對(duì)于單一的t-分布，為了使M步更好求解，Meng和Rubin在文獻(xiàn)[3]中用一種受限制的最大化CM步來替代M步，得到期望條件最大化算法(ECM)；Peel和McLachlan在文獻(xiàn)[2]和文獻(xiàn)[4]中提出混合t-分布模型，用標(biāo)準(zhǔn)EM算法求解混合t-分布模型參數(shù)的極大似然估計(jì)，并給出了ECM算法的一個(gè)應(yīng)用；在此基礎(chǔ)上，Liu和Rubin在文獻(xiàn)[5]中對(duì)ECM算法進(jìn)行兩處修改，得到收斂速度更快的雙期望條件最大化算法(ECME)。隨著計(jì)算機(jī)性能的快速發(fā)展，基于EM算法的混合t-分布模型已越來越廣泛地應(yīng)用到諸多領(lǐng)域，如楊云飛在文獻(xiàn)[6]中提出了自適應(yīng)均值濾波的多元t-分布混合模型，對(duì)醫(yī)學(xué)圖像分割進(jìn)行了研究；熊太松在文獻(xiàn)[7]中對(duì)伯克利圖像數(shù)據(jù)用視覺和量化對(duì)比兩種評(píng)估方式，證明了基于空間平滑的t-分布混合模型在真實(shí)圖像分割中的有效性；朱志娥在文獻(xiàn)[8]中針對(duì)偏t正態(tài)數(shù)據(jù)、異方差和線性回歸提出了偏t正態(tài)數(shù)據(jù)下混合線性聯(lián)合位置與尺度模型，詳細(xì)介紹了該模型下的EM算法并進(jìn)行了有效的模擬驗(yàn)證。在算法初始化方面，冉延平在文獻(xiàn)[9]中用k-means方法確定混合高斯分布的最大混合子分布數(shù)目以及混合比例；史鵬飛在文獻(xiàn)[10]中通過k-means方法先給出混合數(shù)據(jù)的一個(gè)粗糙分組，然后根據(jù)分組數(shù)據(jù)給出參數(shù)的一個(gè)粗略估計(jì)值，作為混合高斯分布EM算法的迭代初始值。

在前人研究的基礎(chǔ)上，本文研究了基于EM算法的一元混合t-分布模型參數(shù)的極大似然估計(jì)，克服了多元混合t-分布模型中協(xié)方差矩陣向一元混合t-分布模型中尺度參數(shù)的轉(zhuǎn)變過程中參數(shù)推導(dǎo)的困難，并首次將k-means方法用于該模型下算法初值的選取。引進(jìn)了混合高斯模型，然后分別在三種不同類型數(shù)據(jù)下進(jìn)行對(duì)比模擬實(shí)驗(yàn)，驗(yàn)證了本文研究的模型和方法的有效性以及其在處理重尾數(shù)據(jù)上的優(yōu)勢(shì)。

1 有限t-分布混合模型

1.1 一元學(xué)生t-分布

設(shè)隨機(jī)變量y服從一元學(xué)生t-分布，記做 y～t(y|μ，σ，ν)，概率密度函數(shù)定義為[11]：

其中參數(shù)μ和σ分別表示t-分布的位置參數(shù)和尺度參數(shù)，Γ(·)表示伽馬函數(shù)。參數(shù)ν稱為t-分布的自由度，當(dāng)ν=1時(shí)，t-分布就成為了柯西分布；如果ν＞1，μ就為t-分布的均值；當(dāng)ν→∞時(shí)，t-分布就以相同的均值μ和方差υ2趨近于高斯分布。

1.2 t-分布有限混合模型

為方便起見，本文只研究?jī)蓚€(gè)子分布的情況，即取m=2，則式（2）化為：

其中，t(y |θk)為第k個(gè)子分布的概率密度函數(shù)，具體形式見式（1）。

2 模型參數(shù)極大似然估計(jì)的EM算法

本文要研究的模型為上文所提到的式（3）。假設(shè)兩個(gè)子分布的自由度相同，即ν1=ν2=ν?；旌夏Ｐ蛥?shù)的求解大多采用EM算法，它提供了一種近似計(jì)算含有隱變量概率模型的極大似然估計(jì)的方法，具有簡(jiǎn)單性和普適性。在EM算法的基本框架下，引入隱變量以得到完整數(shù)據(jù)集。完整數(shù)據(jù)集定義為Yc={Y ，Z，U }，其中，Z 為標(biāo)簽變量 Z={z1，z2，…，zN}，且：

則完整數(shù)據(jù)的對(duì)數(shù)似然函數(shù)：

EM算法是一種迭代求解算法，它主要分兩步進(jìn)行：E步是對(duì)對(duì)數(shù)似然函數(shù)求期望，M步是最大化對(duì)數(shù)似然函數(shù)以獲得新的參數(shù)值。

應(yīng)用EM算法于上式，求解第 j次各參數(shù)的極大似然更新表達(dá)式。

E步：對(duì)對(duì)數(shù)似然函數(shù)求期望確定Q函數(shù)：

首先計(jì)算關(guān)于隱變量Z，U的條件分布的期望：

利用Q函數(shù)對(duì)各參數(shù)求偏導(dǎo)數(shù)并令其等于零，求解得到各參數(shù)的第 j+1次迭代更新表達(dá)式：

自由度ν(j+1)是非線性方程（13）的解：

式（13）是關(guān)于ν的非線性方程，文獻(xiàn)[5]中采用搜索ν的空間求出ν的估計(jì)值，但計(jì)算量大。文獻(xiàn)[12]中給出了一個(gè)計(jì)算量相對(duì)較小的可直接計(jì)算ν近似解的方法。在這里，給出該方法：

3 數(shù)值模擬

為了驗(yàn)證上述參數(shù)估計(jì)方法的有效性，本文共采用三大類數(shù)據(jù)進(jìn)行模擬研究；為了體現(xiàn)混合t-分布處理重尾數(shù)據(jù)的優(yōu)勢(shì)所在，用混合高斯分布模型[13]與之作對(duì)比。算法的初始化均采用k-means方法。此外，為了便于比較，由t-分布的方差與尺度參數(shù)的關(guān)系，將混合t-分布EM算法參數(shù)估計(jì)結(jié)果中的尺度參數(shù)σ轉(zhuǎn)化為標(biāo)準(zhǔn)差υ，再與混合高斯分布EM算法估計(jì)的參數(shù)υ作比較。參數(shù)估計(jì)的精確度采用均方誤差來衡量，如混合比例π1的均方誤差定義為：

其中，π1(0)是π1的真值，n為模擬次數(shù)。

3.1 混合高斯分布數(shù)據(jù)

給定真值 π1(0)=0.3和 0.5，μ1(0)=2、μ2(0)=15、υ1(0)=1、υ2(0)=1，分別取樣本量 N=500，1000，共產(chǎn)生4組混合高斯分布數(shù)據(jù)。對(duì)混合t-分布模型，分別取自由度ν=3[14]，15，30。重復(fù)模擬100次，模擬結(jié)果如表1至表3所示：

表1 ν=3的模擬結(jié)果

表2 ν=15的模擬結(jié)果

表3 ν=30的模擬結(jié)果

由表1至表3可知：ν=3時(shí)，混合高斯模型參數(shù)估計(jì)的均方誤差均比混合t-分布模型參數(shù)估計(jì)的均方誤差小，

給定真值π1(0)=0.3和0.5，μ1(0)=2、μ2(0)=15、σ1(0)=1、σ2(0)=1、ν=3[14]，15，30。分別取樣本量 N=500，1000，共產(chǎn)生12組混合t-分布數(shù)據(jù)。重復(fù)模擬100次，模擬結(jié)果如表4至表6所示。這一點(diǎn)在υ1、υ2上更為明顯；在ν=15，30時(shí)，兩種方法對(duì)各個(gè)參數(shù)估計(jì)的均方誤差，幾乎無差。此外，隨著自由度的增大，混合t-分布模型參數(shù)估計(jì)的均方誤差變小；整體來看，樣本量越大，MSE越??；同一樣本量下，除 μ2和υ2外，兩種方法在混合比例π1=0.5時(shí)的估計(jì)結(jié)果均好于混合比例π1=0.3時(shí)的估計(jì)結(jié)果。

3.2 混合t-分布數(shù)據(jù)

表4 ν=3的模擬結(jié)果

表5 ν=15的模擬結(jié)果

表6 ν=30的模擬結(jié)果

由表4至表6可知：混合t-分布模型可以較好地?cái)M合該數(shù)據(jù)，參數(shù)估計(jì)值與真值十分接近。當(dāng)ν=3時(shí)，對(duì)所有參數(shù)的估計(jì)，混合t-分布模型參數(shù)估計(jì)的均方誤差均比混合高斯分布模型參數(shù)估計(jì)的均方誤差小，這一點(diǎn)在υ1、υ2上更為明顯；ν=15時(shí)，混合t-分布模型參數(shù)估計(jì)的均方誤差比混合高斯分布模型參數(shù)估計(jì)的均方誤差略小，但相差不大；在ν=30時(shí)，兩種方法下參數(shù)估計(jì)的均方誤差相比，幾乎無差。此外，隨著自由度的增大，混合t-分布模型參數(shù)估計(jì)的均方誤差變??；整體來看，樣本量越大，MSE越小，估計(jì)結(jié)果越好；在同一樣本量下，除μ2和υ2外，兩種方法在混合比例π1=0.5時(shí)的估計(jì)結(jié)果均好于混合比例π1=0.3時(shí)的估計(jì)結(jié)果，但相差不大。

3.3 含噪聲的混合高斯數(shù)據(jù)

因?yàn)閠-分布混合模型相對(duì)于高斯混合模型有著較好的穩(wěn)健性，這種穩(wěn)健性尤其體現(xiàn)在對(duì)重尾數(shù)據(jù)(含噪聲點(diǎn)、異常點(diǎn)數(shù)據(jù))的處理。而處理重尾數(shù)據(jù)的另一種方法是在高斯分布的基礎(chǔ)上添加一個(gè)均勻分布的成分[6]。因此，本文在高斯數(shù)據(jù)的基礎(chǔ)上添加一個(gè)均勻分布的部分作為重尾數(shù)據(jù)，然后再分別用混合t-分布模型和混合高斯模型進(jìn)行擬合并作比較。因?yàn)樯衔囊呀?jīng)對(duì)自由度、樣本量和混合比例進(jìn)行了研究比較，并且發(fā)現(xiàn)在自由度取ν=15時(shí)兩種方法的估計(jì)效果已相差不大，因此這里不再考慮此三者的影響。取噪聲所占比例分別為5%和10%，混合比例π1=0.3，自由度ν=15，樣本量 N=1000。重復(fù)模擬100次，模擬結(jié)果如表7和表8所示。

表7 含5%噪聲的混合高斯數(shù)據(jù)下參數(shù)估計(jì)結(jié)果

表8 含10%噪聲的混合高斯數(shù)據(jù)下參數(shù)估計(jì)結(jié)果

由表7和表8知：通過比較兩種模型下參數(shù)的估計(jì)結(jié)果和均方誤差可以得到，混合t-分布模型對(duì)該類型數(shù)據(jù)擬合的較好，均方誤差較小。因此相對(duì)于混合高斯分布，混合t-分布模型可以更好地?cái)M合含噪聲的混合高斯數(shù)據(jù)，這也正說明了混合t-分布模型較于混合高斯模型能夠更好地處理重尾數(shù)據(jù)。

4 結(jié)論

本文主要研究了一元混合t-分布模型，給出了EM算法下該模型參數(shù)的極大似然估計(jì)，并采用k-means方法進(jìn)行算法初始化，最后在模擬的三種類型的數(shù)據(jù)下與混合高斯模型進(jìn)行了對(duì)比分析。在前兩類數(shù)據(jù)的模型參數(shù)估計(jì)結(jié)果中可以看出，每個(gè)子分布的自由度固定且取相同的值的情況下，對(duì)于混合高斯數(shù)據(jù)，當(dāng)自由度的取值足夠大時(shí)，基于混合t-分布模型的EM算法的參數(shù)估計(jì)結(jié)果并不比基于混合高斯模型的EM算法差；對(duì)于混合t-分布數(shù)據(jù)，基于混合t-分布模型的EM算法能夠得到較好的估計(jì)結(jié)果并優(yōu)于基于混合高斯模型的EM算法的估計(jì)結(jié)果，且隨著自由度的增大，效果會(huì)更好；而在第三類含噪聲的混合高斯分布數(shù)據(jù)下，混合t-分布模型比混合高斯分布模型擬合效果更好，說明了混合t-分布模型在處理重尾數(shù)據(jù)上更具優(yōu)勢(shì)。以上結(jié)果驗(yàn)證了本文研究的模型和方法的有效性。