王小英,李迎華,楊雪梅
(華北電力大學 數理學院,北京 102206)
混合分布模型是分析復雜現象的一個靈活而強有力的建模工具,它提供了用簡單結構模擬復雜密度的一個有效方法,從而將聚類問題轉化為統計分析問題。混合高斯模型由于計算上的便利,在聚類方法研究中應用較多,但通常我們收集到的很多數據并不是嚴格的服從正態分布,而是較明顯的服從重尾分布?;旌蟭-分布模型由于其具有較長的尾巴,可對重尾點和異常點有效地降低權值,因此,相對于高斯分布混合模型,可以獲得較強的精度和穩健性。
文獻[1]中Dempster等提出的EM算法成為了混合模型參數估計極有效的工具。Peel和McLachlan在文獻[2]中指出EM算法可以獲得有限混合模型任意分布的極大似然估計。對于單一的t-分布,為了使M步更好求解,Meng和Rubin在文獻[3]中用一種受限制的最大化CM步來替代M步,得到期望條件最大化算法(ECM);Peel和McLachlan在文獻[2]和文獻[4]中提出混合t-分布模型,用標準EM算法求解混合t-分布模型參數的極大似然估計,并給出了ECM算法的一個應用;在此基礎上,Liu和Rubin在文獻[5]中對ECM算法進行兩處修改,得到收斂速度更快的雙期望條件最大化算法(ECME)。隨著計算機性能的快速發展,基于EM算法的混合t-分布模型已越來越廣泛地應用到諸多領域,如楊云飛在文獻[6]中提出了自適應均值濾波的多元t-分布混合模型,對醫學圖像分割進行了研究;熊太松在文獻[7]中對伯克利圖像數據用視覺和量化對比兩種評估方式,證明了基于空間平滑的t-分布混合模型在真實圖像分割中的有效性;朱志娥在文獻[8]中針對偏t正態數據、異方差和線性回歸提出了偏t正態數據下混合線性聯合位置與尺度模型,詳細介紹了該模型下的EM算法并進行了有效的模擬驗證。在算法初始化方面,冉延平在文獻[9]中用k-means方法確定混合高斯分布的最大混合子分布數目以及混合比例;史鵬飛在文獻[10]中通過k-means方法先給出混合數據的一個粗糙分組,然后根據分組數據給出參數的一個粗略估計值,作為混合高斯分布EM算法的迭代初始值。
在前人研究的基礎上,本文研究了基于EM算法的一元混合t-分布模型參數的極大似然估計,克服了多元混合t-分布模型中協方差矩陣向一元混合t-分布模型中尺度參數的轉變過程中參數推導的困難,并首次將k-means方法用于該模型下算法初值的選取。引進了混合高斯模型,然后分別在三種不同類型數據下進行對比模擬實驗,驗證了本文研究的模型和方法的有效性以及其在處理重尾數據上的優勢。
設隨機變量y服從一元學生t-分布,記做 y~t(y|μ,σ,ν),概率密度函數定義為[11]:

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

為方便起見,本文只研究兩個子分布的情況,即取m=2,則式(2)化為:

其中,t(y |θk)為第k個子分布的概率密度函數,具體形式見式(1)。
本文要研究的模型為上文所提到的式(3)。假設兩個子分布的自由度相同,即ν1=ν2=ν?;旌夏P蛥档那蠼獯蠖嗖捎肊M算法,它提供了一種近似計算含有隱變量概率模型的極大似然估計的方法,具有簡單性和普適性。在EM算法的基本框架下,引入隱變量以得到完整數據集。完整數據集定義為Yc={Y ,Z,U },其中,Z 為標簽變量 Z={z1,z2,…,zN},且:

則完整數據的對數似然函數:

EM算法是一種迭代求解算法,它主要分兩步進行:E步是對對數似然函數求期望,M步是最大化對數似然函數以獲得新的參數值。
應用EM算法于上式,求解第 j次各參數的極大似然更新表達式。
E步:對對數似然函數求期望確定Q函數:

首先計算關于隱變量Z,U的條件分布的期望:


利用Q函數對各參數求偏導數并令其等于零,求解得到各參數的第 j+1次迭代更新表達式:

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

式(13)是關于ν的非線性方程,文獻[5]中采用搜索ν的空間求出ν的估計值,但計算量大。文獻[12]中給出了一個計算量相對較小的可直接計算ν近似解的方法。在這里,給出該方法:

為了驗證上述參數估計方法的有效性,本文共采用三大類數據進行模擬研究;為了體現混合t-分布處理重尾數據的優勢所在,用混合高斯分布模型[13]與之作對比。算法的初始化均采用k-means方法。此外,為了便于比較,由t-分布的方差與尺度參數的關系,將混合t-分布EM算法參數估計結果中的尺度參數σ轉化為標準差υ,再與混合高斯分布EM算法估計的參數υ作比較。參數估計的精確度采用均方誤差來衡量,如混合比例π1的均方誤差定義為:

其中,π1(0)是π1的真值,n為模擬次數。
給定真值 π1(0)=0.3和 0.5,μ1(0)=2、μ2(0)=15、υ1(0)=1、υ2(0)=1,分別取樣本量 N=500,1000,共產生4組混合高斯分布數據。對混合t-分布模型,分別取自由度ν=3[14],15,30。重復模擬100次,模擬結果如表1至表3所示:

表1 ν=3的模擬結果

表2 ν=15的模擬結果

表3 ν=30的模擬結果
由表1至表3可知:ν=3時,混合高斯模型參數估計的均方誤差均比混合t-分布模型參數估計的均方誤差小,
給定真值π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,共產生12組混合t-分布數據。重復模擬100次,模擬結果如表4至表6所示。這一點在υ1、υ2上更為明顯;在ν=15,30時,兩種方法對各個參數估計的均方誤差,幾乎無差。此外,隨著自由度的增大,混合t-分布模型參數估計的均方誤差變??;整體來看,樣本量越大,MSE越??;同一樣本量下,除 μ2和υ2外,兩種方法在混合比例π1=0.5時的估計結果均好于混合比例π1=0.3時的估計結果。

表4 ν=3的模擬結果

表5 ν=15的模擬結果

表6 ν=30的模擬結果
由表4至表6可知:混合t-分布模型可以較好地擬合該數據,參數估計值與真值十分接近。當ν=3時,對所有參數的估計,混合t-分布模型參數估計的均方誤差均比混合高斯分布模型參數估計的均方誤差小,這一點在υ1、υ2上更為明顯;ν=15時,混合t-分布模型參數估計的均方誤差比混合高斯分布模型參數估計的均方誤差略小,但相差不大;在ν=30時,兩種方法下參數估計的均方誤差相比,幾乎無差。此外,隨著自由度的增大,混合t-分布模型參數估計的均方誤差變小;整體來看,樣本量越大,MSE越小,估計結果越好;在同一樣本量下,除μ2和υ2外,兩種方法在混合比例π1=0.5時的估計結果均好于混合比例π1=0.3時的估計結果,但相差不大。
因為t-分布混合模型相對于高斯混合模型有著較好的穩健性,這種穩健性尤其體現在對重尾數據(含噪聲點、異常點數據)的處理。而處理重尾數據的另一種方法是在高斯分布的基礎上添加一個均勻分布的成分[6]。因此,本文在高斯數據的基礎上添加一個均勻分布的部分作為重尾數據,然后再分別用混合t-分布模型和混合高斯模型進行擬合并作比較。因為上文已經對自由度、樣本量和混合比例進行了研究比較,并且發現在自由度取ν=15時兩種方法的估計效果已相差不大,因此這里不再考慮此三者的影響。取噪聲所占比例分別為5%和10%,混合比例π1=0.3,自由度ν=15,樣本量 N=1000。重復模擬100次,模擬結果如表7和表8所示。

表7 含5%噪聲的混合高斯數據下參數估計結果

表8 含10%噪聲的混合高斯數據下參數估計結果
由表7和表8知:通過比較兩種模型下參數的估計結果和均方誤差可以得到,混合t-分布模型對該類型數據擬合的較好,均方誤差較小。因此相對于混合高斯分布,混合t-分布模型可以更好地擬合含噪聲的混合高斯數據,這也正說明了混合t-分布模型較于混合高斯模型能夠更好地處理重尾數據。
本文主要研究了一元混合t-分布模型,給出了EM算法下該模型參數的極大似然估計,并采用k-means方法進行算法初始化,最后在模擬的三種類型的數據下與混合高斯模型進行了對比分析。在前兩類數據的模型參數估計結果中可以看出,每個子分布的自由度固定且取相同的值的情況下,對于混合高斯數據,當自由度的取值足夠大時,基于混合t-分布模型的EM算法的參數估計結果并不比基于混合高斯模型的EM算法差;對于混合t-分布數據,基于混合t-分布模型的EM算法能夠得到較好的估計結果并優于基于混合高斯模型的EM算法的估計結果,且隨著自由度的增大,效果會更好;而在第三類含噪聲的混合高斯分布數據下,混合t-分布模型比混合高斯分布模型擬合效果更好,說明了混合t-分布模型在處理重尾數據上更具優勢。以上結果驗證了本文研究的模型和方法的有效性。