拓端tecdat|R語言有限混合模型(FMM,finite mixture model)及其EM算法聚類分析間歇泉噴
原文鏈接：http://tecdat.cn/?p=22609?
原文出處：拓端數(shù)據(jù)部落公眾號
摘要
本文提供了一套用于分析各種有限混合模型的方法。既包括傳統(tǒng)的方法，如單變量和多變量正態(tài)混合的EM算法，也包括反映有限混合模型的一些最新研究的方法。許多算法都是EM算法或基于類似EM的思想，因此本文包括有限混合模型的EM算法的概述。
1.有限混合模型介紹
人群中的個體往往可以被劃分為群。然而，即使我們觀察到這些個體的特征，我們也可能沒有真正觀察到這些成員的群體。這項任務在文獻中有時被稱為 "無監(jiān)督聚類"，事實上，混合模型一般可以被認為是由被稱為 "基于模型的聚類 "的聚類方法的子集組成。
有限混合模型也可用于那些對個體聚類感興趣的情況之外。首先，有限混合模型給出了整個子群的描述，而不是將個體分配到這些子群中。有時，有限混合模型只是提供了一種充分描述特定分布的手段，例如線性回歸模型中存在異常值的殘差分布。
無論建模者在采用混合模型時的目標是什么，這些模型的大部分理論都涉及到一個假設，即子群是按照一個特定的參數(shù)形式分布的--而這個形式往往是單變量或多變量正態(tài)。
最近的研究目標是放寬或修改多變量正態(tài)假設，有限混合模型分析的計算技術，其中的成分是回歸、多變量數(shù)據(jù)離散化產生的向量，甚至是完全未指定的分布。
2. 有限混合模型的EM算法
EM算法迭代最大化，而不是觀察到的對數(shù)似然Lx(θ)，算式為
1. E步：計算Q(θ|θ(t))
2. M步驟：設定θ(t+1)=argmaxθ∈Φ Q(θ|θ(t))
對于有限混合模型，E步驟不依賴于F的結構，因為缺失數(shù)據(jù)部分只與Z有關。
?
Z是離散的，它們的分布是通過貝葉斯定理給出的。M步驟本身可以分成兩部分，與λ有關的最大化，它不依賴于F，與φ有關的最大化，它必須為每個模型專門處理（例如，參數(shù)化、半?yún)?shù)化或非參數(shù)化）。因此，模型的EM算法有以下共同特點。
11. E步。計算成分包含的 "后驗 "概率（以數(shù)據(jù)和θ（t）為條件）。
對于所有i = 1, . . . ，n和j = 1， . . . 從數(shù)值上看，完全按照公式（2）的寫法來實現(xiàn)是很危險的，因為在xi離任何一個成分都很遠的情況下，所有的φ（t）j 0（xi）值都會導致數(shù)值下溢為零，所以可能會出現(xiàn)不確定的形式0/0。因此，許多例程實際上使用的是等價表達式
或其某種變體。
2. ?λ的M步驟。設
2.3. 一個EM算法的例子
作為一個例子，我們考慮對圖1中描述的間歇泉噴發(fā)間隔時間等待數(shù)據(jù)進行單變量正態(tài)混合分析。這種完全參數(shù)化的情況對應于第1節(jié)中描述的單變量高斯家族的混合分布，其中（1）中的第j個分量密度φj（x）為正態(tài)，均值為μj，方差為σ 2 j。
對于參數(shù)(μj , σ2 j )的M步，j = 1, . . 這個EM算法對這種單變量混合分布的M步驟是很簡單的，例如可以在McLachlan和Peel（2000）中找到。
mixEM(waiting, lambda = .5)
上面的代碼將擬合一個二成分的混合分布（因為mu是一個長度為2的向量），其中標準偏差被假定為相等（因為sigma是一個標量而不是一個向量）。
圖1:對數(shù)似然值的序列，Lx(θ (t))
圖2：用參數(shù)化EM算法擬合間歇泉等待數(shù)據(jù)。擬合的高斯成分。
R> plot(wait1, density = TRUE, cex.axis = 1.4, cex.lab = 1.4, cex.main = 1.8,+ main2 = "Time between Old Faithful eruptions", xlab2 = "Minutes")
兩個圖：觀察到的對數(shù)似然值的序列t 7→Lx(θ (t))和數(shù)據(jù)的直方圖，其中有N(?μj , σ? 2 j)的m(這里m=2)個擬合的高斯分量密度，j=1, . . . ，m，疊加在一起。估計θ?
另外，使用summary也可以得到同樣的輸出。
summary(wait1)
3. Cutpoint methods切割點方法
傳統(tǒng)上，大多數(shù)關于有限混合模型的文獻都假設方程（1）的密度函數(shù)φj（x）來自一個已知的參數(shù)族。然而，一些作者最近考慮了這樣的問題：除了確保模型中參數(shù)的可識別性所需的一些條件外，φj（x）是不指定的。我們使用Elmore等人（2004）的切割點方法。
我們參考Elmore等人從-63開始，一直到63大約以10.5的間隔采用切點。然后從原始數(shù)據(jù)中創(chuàng)建一個多指標數(shù)據(jù)集，如下所示。
R> cutpts <- 10.5*(-6:6)R> mult(data, cuts = cutpts)
一旦創(chuàng)建了多指標數(shù)據(jù)，我們可以應用EM算法估計多指標參數(shù)。最后，計算并繪制出方程的估計分布函數(shù)。圖3給出了3分量和4分量解決方案的圖表；這些圖表與Elmore等人（2004）的圖1和圖2中的相應圖表非常相似。
R> plot(data, posterior, lwd = 2,+ main = "三分量解")
圖3(a)?
圖3(b)
同樣可以用summary來總結EM輸出。
單變量對稱、位置偏移的半?yún)?shù)例子
在φ(-)相對于Lebesgue度量是絕對連續(xù)的額外假設下，Bordes等人(2007)提出了一種估計模型參數(shù)的隨機算法，即(λ, μ, φ)。一個特例
R> plot(wait1, which = 2 )R> wait2 <-EM(waiting)R> plot(wait2, lty = 2)
圖4(a)?
?
圖4(b)
因為半?yún)?shù)版本依賴于核密度估計步驟（8），所以有必要為這個步驟選擇一個帶寬。默認情況下，使用"Silverman的經(jīng)驗法則"（Silverman 1986）應用于整個數(shù)據(jù)集。
R> bw.nrd0(wait)
但帶寬的選擇會產生很大的不同，如圖4(b)所示。
> wait2a <- EM(wait, bw = 1)> plot(wait2a> plot(wait2b
我們發(fā)現(xiàn)，在帶寬接近2的情況下，半?yún)?shù)解看起來非常接近圖2的正態(tài)混合分布解。進一步降低帶寬會導致圖4(b)中的實線所表現(xiàn)出的 "凹凸不平"。另一方面，在帶寬為8的情況下，半?yún)?shù)解效果很差，因為算法試圖使每個成分看起來與整個混合分布相似。
最受歡迎的見解
1.R語言k-Shape算法股票價格時間序列聚類
2.R語言中不同類型的聚類方法比較
3.R語言對用電負荷時間序列數(shù)據(jù)進行K-medoids聚類建模和GAM回歸
4.r語言鳶尾花iris數(shù)據(jù)集的層次聚類
5.Python Monte Carlo K-Means聚類實戰(zhàn)
6.用R進行網(wǎng)站評論文本挖掘聚類
7.用于NLP的Python：使用Keras的多標簽文本LSTM神經(jīng)網(wǎng)絡
8.R語言對MNIST數(shù)據(jù)集分析 探索手寫數(shù)字分類數(shù)據(jù)
9.R語言基于Keras的小數(shù)據(jù)集深度學習圖像分類
