4.3 高斯混合與期望最大化

動機

位於兩個群組之間的數據，未必明顯屬於其中一組。硬分配只保留選中的群組；混合模型則結合多個概率密度，為每筆觀測保留屬於各分量的概率，再用這些概率加權重新估計各分量。

困難在於似然的每一項都包含「總和的對數」。即使單個高斯分量很簡單，不同分量的參數仍會透過總和互相牽連。期望最大化法，簡稱 EM，引入責任度，把參數更新拆成加權估計問題。它的兩個步驟共同改善一個在目前參數處與觀測對數似然相切的下界。

本節會先推導權重與更新式，再完整計算一次迭代。軟概率與硬標籤之間的分別是關鍵：最大化選中分量的密度，與最大化混合密度有不同目標。加權協方差的分母亦必須由所最大化的似然推導。

幾何軟權重

軟權重亦可由幾何定義。設互異站點 $S=\{s_1,\ldots,s_K\}$ 的舊 Voronoi 單元為 $V_j$，距離相同時固定分給其中一個單元。插入查詢點 $x\notin S$，以 $V_x$ 表示它的新單元，並假設該單元有界，體積有限且為正。定義

$$0<\operatorname{vol}(V_x)<\infty,\qquad w_j=\frac{\operatorname{vol}(V_x\cap V_j)} {\operatorname{vol}(V_x)},\qquad \sum_{j=1}^K w_j=1.$$

舊單元把新單元分割，所以重疊體積之和等於新單元體積。各權重非負，沒有正體積重疊的站點得到零權重；這些比例量度新區域原本分屬各舊站點的份額。

一維體積就是長度。取兩個站點及一個新查詢點：

$$S=\{0,4\},\quad x=1,\quad V_x=\left[\frac12,\frac52\right],\quad (w_1,w_2)=\left(\frac{3/2}{2},\frac{1/2}{2}\right) =\left(\frac34,\frac14\right).$$

舊邊界為 $2$，新區間在其左右的重疊長度分別為 $3/2$ 與 $1/2$，總長為 $2$。這是幾何分配權重；下文的高斯責任度則是混合模型下的條件概率。兩者都跨分量規範，但定義所用的資訊不同。

定義

設 $x_1,\ldots,x_n\in\mathbb R^d$ 為獨立觀測，固定高斯分量數目 $K$。分量 $k$ 有權重 $\pi_k>0$、平均 $\mu_k\in\mathbb R^d$ 及正定協方差 $\Sigma_k$，並滿足 $\sum_{k=1}^K\pi_k=1$。以 $\theta$ 表示全部參數，以下對數均為自然對數。

定義

高斯混合與觀測對數似然

高斯密度及混合密度分別是

$$\phi_d(x;\mu,\Sigma) =\frac{\exp\!\left(-\frac12(x-\mu)^\top\Sigma^{-1}(x-\mu)\right)} {(2\pi)^{d/2}\det(\Sigma)^{1/2}}, \qquad p_\theta(x)=\sum_{k=1}^K\pi_k\phi_d(x;\mu_k,\Sigma_k).$$

觀測對數似然為

$$\ell(\theta)=\sum_{i=1}^n\log p_\theta(x_i).$$

潛在分量標籤 $z_i$ 以概率 $\pi_k$ 取值 $k$；分量密度描述給定該標籤時 $x_i$ 的分佈。

混合權重本身是模型的一部分，不是選定分量後才加上的規範。即使兩個分量在某筆觀測的密度相同，權重仍會影響結果。對全部可能標籤求和，才得到觀測密度。

若連標籤也已觀測，對數似然便是 $\log\pi_{z_i}+\log\phi_d(x_i;\mu_{z_i},\Sigma_{z_i})$ 的總和，亦可用每個分量的零一指示量選出相應項。未知標籤令我們無法直接按這種方式分組；EM 改為取這些指示量的條件期望，得到 $Q$ 中的責任度。因此分數權重有明確的概率意義，是目前模型下潛在指示量的期望。

這也區分了推導中的三個函數：觀測對數似然在每個對數內對標籤求和；完整數據對數似然使用已知標籤；$Q$ 則在目前後驗分佈下平均完整數據表達式。分清觀測量、潛在量及暫時固定的量，便能理解加權高斯擬合為何可以改善較複雜的混合似然。

定義

責任度與分量的有效樣本數

在目前參數 $\theta^{\mathrm{old}}$ 下，Bayes 公式給出

$$r_{ik} =\Pr_{\theta^{\mathrm{old}}}(z_i=k\mid x_i) =\frac{\pi_k^{\mathrm{old}} \phi_d(x_i;\mu_k^{\mathrm{old}},\Sigma_k^{\mathrm{old}})} {\sum_{j=1}^K\pi_j^{\mathrm{old}} \phi_d(x_i;\mu_j^{\mathrm{old}},\Sigma_j^{\mathrm{old}})}, \qquad N_k=\sum_{i=1}^n r_{ik}.$$

每一行滿足 $\sum_k r_{ik}=1$。混合權重為正、協方差正定時，每個 $r_{ik}$ 都為正，因此 $N_k>0$，並且 $\sum_kN_k=n$。

$N_k$ 是在目前後驗標籤下，預期歸於分量 $k$ 的觀測數，可以不是整數。責任度是給定目前模型與觀測後的概率，不是隨數據提供的已知類別標籤。

這裏有兩種不同的歸一化。E 步對每筆觀測跨分量規範，故每行總和為一；平均更新則對某個分量跨觀測規範，分母是列總和 $N_k$。若平均更新除以 $n$，便會把分給其他分量的責任度也算入本分量的總量。兩個分母回答的是不同條件問題。

某點的高斯密度本身不是成員概率：密度可以大於一，也會隨量度單位改變尺度。責任度則是比較同一觀測的加權密度後得到的無量綱概率。一般協方差下，比較同時涉及指數中的 Mahalanobis 距離及規範中的行列式。較分散的分量有較大體積及較低峰值，$\pi_k$ 則控制其總體混合權重。等權重與共同各向同性協方差下，分量密度排序可化為歐幾里得距離比較。

定理/命題

定理

高斯 M 步以有效樣本數規範

固定責任度，定義完整數據對數似然的條件期望

$$Q(\theta\mid\theta^{\mathrm{old}}) =\sum_{i=1}^n\sum_{k=1}^K r_{ik}\left[\log\pi_k+\log\phi_d(x_i;\mu_k,\Sigma_k)\right].$$

若以下加權協方差均為正定，其最大化更新為

$$\pi_k^{\mathrm{new}}=\frac{N_k}{n},\qquad \mu_k^{\mathrm{new}}=\frac1{N_k}\sum_{i=1}^nr_{ik}x_i,$$$$\Sigma_k^{\mathrm{new}} =\frac1{N_k}\sum_{i=1}^nr_{ik} (x_i-\mu_k^{\mathrm{new}}) (x_i-\mu_k^{\mathrm{new}})^\top.$$

若模型預先指定共同協方差為 $\sigma^2I_d$、$\sigma^2>0$，便保持它固定，只使用相同的權重及平均更新。

當且僅當加權中心化觀測張成全部特徵方向，這個協方差才是正定的；否則所得矩陣奇異，不屬於正定高斯協方差。固定方差模型讓我們在毋須估計協方差矩陣的情況下，逐步計算有效的模型更新。

定理

有效的 EM 迭代不會降低觀測對數似然

先在 $\theta^{\mathrm{old}}$ 計算精確後驗責任度，再選取滿足以下條件的有效新參數：

$$Q(\theta^{\mathrm{new}}\mid\theta^{\mathrm{old}}) \ge Q(\theta^{\mathrm{old}}\mid\theta^{\mathrm{old}}).$$

則

$$\ell(\theta^{\mathrm{new}})\ge\ell(\theta^{\mathrm{old}}).$$

精確 M 步最大化 $Q$；只令它增加的更新稱為廣義 EM 步。結論針對觀測對數似然的值；全域最優及參數向量收斂是另外的數學主張。

證明思路

推導加權更新

權重部分是在 $\sum_k\pi_k=1$ 下最大化 $\sum_kN_k\log\pi_k$。引入乘數 $\lambda$，駐點條件為 $N_k/\pi_k=\lambda$。對 $k$ 求和得到 $\lambda=n$，所以 $\pi_k=N_k/n$。

各平均的相關項是加權 Mahalanobis 距離平方和的負數，微分給出

$$\Sigma_k^{-1}\sum_i r_{ik}(x_i-\mu_k)=0.$$

由 $\Sigma_k$ 可逆及 $N_k>0$，平均必須是命題中的加權平均。預先指定共同各向同性協方差時，這就是加權最小二乘的最小點。每筆觀測按責任度貢獻給各分量，而不是只貢獻給責任度最大的分量。

推導協方差時，令

$$S_k=\sum_i r_{ik}(x_i-\mu_k^{\mathrm{new}}) (x_i-\mu_k^{\mathrm{new}})^\top,\qquad \Lambda_k=\Sigma_k^{-1}.$$

$Q$ 的相關部分為

$$\frac{N_k}{2}\log\det\Lambda_k -\frac12\operatorname{tr}(\Lambda_kS_k).$$

駐點方程是 $N_k\Lambda_k^{-1}=S_k$。$S_k$ 正定時，由對數行列式的凹性，最大點為 $\Sigma_k=S_k/N_k$。這是似然估計，所以分母為 $N_k$，不減一；省略該分母會得到尺度不同的加權散布矩陣。

與目前似然相切的下界

對任意正概率行 $q_i=(q_{i1},\ldots,q_{iK})$，Jensen 不等式給出

$$\log\sum_k\pi_k\phi_d(x_i;\mu_k,\Sigma_k) \ge\sum_kq_{ik}\log \frac{\pi_k\phi_d(x_i;\mu_k,\Sigma_k)}{q_{ik}}.$$

把右方對觀測求和，定義

$$\mathcal F(q,\theta) =\sum_{i,k}q_{ik} \left[\log\pi_k+\log\phi_d(x_i;\mu_k,\Sigma_k)-\log q_{ik}\right].$$

精確差距為

$$\ell(\theta)-\mathcal F(q,\theta) =\sum_i\operatorname{KL}\!\left(q_i\Vert r_i(\theta)\right)\ge0,$$

其中 $\operatorname{KL}(q\Vert r)=\sum_kq_k\log(q_k/r_k)$，而 $r_i(\theta)$ 是參數 $\theta$ 下的後驗行。當 $q_i=r_i(\theta)$，差距為零。

E 步選取 $q=r(\theta^{\mathrm{old}})$，令下界在目前觀測對數似然處取等號。固定 $q$ 後，其熵 $-\sum_{i,k}q_{ik}\log q_{ik}$ 是常數，因此增加 $Q$ 就會增加 $\mathcal F$，得到

$$\ell(\theta^{\mathrm{new}}) \ge\mathcal F(r^{\mathrm{old}},\theta^{\mathrm{new}}) \ge\mathcal F(r^{\mathrm{old}},\theta^{\mathrm{old}}) =\ell(\theta^{\mathrm{old}}).$$

兩步各有作用：E 步令下界在目前模型處精確，M 步則在固定下界下改善模型。若用硬標籤代替後驗概率，便失去了證明所依賴的相切等號。

一次 M 步全程固定責任度。新平均、新協方差及新權重，都使用完整舊參數所算出的同一後驗矩陣。協方差式以新加權平均中心化，因為它最大化的是同一固定 $Q$ 問題。形成新模型後，下一次 E 步才重新計算後驗；不等式鏈中間的比較，就是保留一個概率矩陣而改變參數。

以 $0\log0=0$ 及連續延伸解讀零權重項，下界與 KL 恆等式亦適用於含零元素的概率行。因此可以比較精確後驗的相切下界，以及硬分配得到的較鬆下界。

例題

例題

1. 分量密度相同，混合權重仍可不同

考慮一維混合模型，平均為 $\mu_1=0$、$\mu_2=2$，共同方差 $\sigma^2=1$，權重為 $(\pi_1,\pi_2)=(3/4,1/4)$。在 $x=1$，到兩個平均的距離相同，故高斯密度相同。代入 Bayes 公式後，共同密度抵消：

$$r_1=\frac{(3/4)\phi_1(1;0,1)} {(3/4)\phi_1(1;0,1)+(1/4)\phi_1(1;2,1)} =\frac34,\qquad r_2=\frac14.$$

中點與兩個分量形狀同樣相容，但模型賦予第一分量三倍的先驗權重。若省略混合權重，只規範兩個密度，就會得到 $1/2$ 與 $1/2$；那屬於另一個混合模型的後驗。因此軟成員概率同時取決於密度與權重。

例題

2. 固定方差下的完整一次 E 步與 M 步

取觀測 $x_1=0$、$x_2=2$，初始平均 $(0,2)$，混合權重相同，並指定共同方差為 $1$。令

$$a=\frac1{1+e^{-2}}\approx0.880797078.$$

每筆觀測到另一平均的距離平方，比到相應平均大 $4$，所以密度比為 $e^{-2}$。因此

$$R=\begin{bmatrix}a&1-a\\1-a&a\end{bmatrix}, \qquad N_1=N_2=1.$$

M 步得到

$$\pi^{\mathrm{new}}=(1/2,1/2),\qquad \mu^{\mathrm{new}}=(2(1-a),2a) \approx(0.238405844,1.761594156).$$

每筆觀測最初雖然恰好位於一個分量平均，它對另一分量仍有正責任度，因而把兩個新平均往內拉。共同方差保持為 $1$。使用完整高斯密度計算，觀測對數似然由 $-2.970315405$ 升至 $-2.886562696$，約增加 $0.083752710$。

以下 NumPy 程式重現同一次迭代。對數規範同時使用全部分量計算後驗，並不先挑選某個分量再形成分母。

import numpy as np

X = np.array([0.0, 2.0])
means = np.array([0.0, 2.0])
weights = np.array([0.5, 0.5])
variance = 1.0

def log_components(means, weights):
    return (np.log(weights)[None, :]
            - 0.5 * np.log(2.0 * np.pi * variance)
            - 0.5 * (X[:, None] - means[None, :]) ** 2 / variance)

joint = log_components(means, weights)
normalizer = np.logaddexp.reduce(joint, axis=1)
responsibilities = np.exp(joint - normalizer[:, None])
mass = responsibilities.sum(axis=0)
new_weights = mass / len(X)
new_means = (responsibilities.T @ X) / mass
before = normalizer.sum()
after = np.logaddexp.reduce(
    log_components(new_means, new_weights), axis=1
).sum()
assert np.allclose(responsibilities.sum(axis=1), 1.0)
assert after >= before
print("responsibilities =", np.round(responsibilities, 6).tolist())
print("mass =", np.round(mass, 6).tolist())
print("weights =", np.round(new_weights, 6).tolist())
print("means =", np.round(new_means, 6).tolist())
print(f"log_likelihood: {before:.6f} -> {after:.6f}")

輸出為

responsibilities = [[0.880797, 0.119203], [0.119203, 0.880797]]
mass = [1.0, 1.0]
weights = [0.5, 0.5]
means = [0.238406, 1.761594]
log_likelihood: -2.970315 -> -2.886563

例題

3. 協方差塌縮可令似然無上界

取觀測 $0$ 與 $3$。兩個分量權重各為 $1/2$，平均為 $0$ 與 $3$，方差分別為 $\varepsilon^2$ 與 $1$。每個 $\varepsilon>0$ 都給出有效高斯分量，其混合密度滿足

$$p_\theta(0)\ge\frac1{2\sqrt{2\pi}\varepsilon}, \qquad p_\theta(3)\ge\frac1{2\sqrt{2\pi}}.$$

所以

$$\ell(\theta)\ge-\log\varepsilon-\log(8\pi) \longrightarrow+\infty \qquad(\varepsilon\downarrow0).$$

把一個分量收縮到某筆觀測附近，可以令該處密度任意高；另一分量繼續解釋其餘觀測。這個不受限制的似然沒有有限最大值。問題在統計目標本身，而不是要求計算奇異協方差的高斯密度。像例題 2 般預先指定正的共同方差，是另一個不會出現這種塌縮的模型。

如何理解 EM 的進度

指定共同各向同性方差時，各分量密度至多為 $(2\pi\sigma^2)^{-d/2}$；其加權平均亦有相同上界，因此

$$\ell(\theta)\le-\frac{nd}{2}\log(2\pi\sigma^2).$$

所以這個模型的不下降似然值有有限極限。這是對實數目標值序列的結論，不會指定最佳混合模型，也不會建立每個參數座標的極限。交換分量標籤，亦會描述相同混合密度。

初始模型決定第一次後驗權重，這些權重再決定第一次加權擬合問題。不同初始模型可以沿着不同的相切下界序列前進。似然可檢查所選模型內的最佳化進度，但不會自動選擇分量數目或比較不同協方差假設。

等權重、共同各向同性方差下，責任度是負距離平方的指數再作規範。方差越小，責任度越集中於唯一最近的平均。對每個正方差，後驗仍然是軟分配；硬最近中心更新則會改變最佳化問題。這解釋了加權平均與中心式聚類的相似之處，同時保留兩種目標的區別。

常見錯誤

常見錯誤

用硬 argmax 代替 E 步

從例題 2 的新平均出發，硬分配會把 $0$ 送到第一分量、$2$ 送到第二分量，再以各組觀測平均更新，就會回到 $(0,2)$。共同方差與等權重不變時，觀測對數似然會由約 $-2.886563$ 降至 $-2.970315$。硬分配不是令 EM 下界與目前似然相切的後驗分佈。

常見錯誤

把散布當成協方差，或省略混合權重

高斯協方差更新把加權散布除以 $N_k$；E 步則跨分量規範 $\pi_k\phi_d(x_i;\mu_k,\Sigma_k)$。兩個因子都來自混合似然。改變任何一個規範，都會改變擬合方程，使上述 EM 推導不再成立。

總結

高斯混合把分量密度加權求和。E 步使用規範後的後驗責任度，M 步以有效樣本數更新權重、平均及需要估計的協方差。Jensen 下界在目前後驗處精確，M 步增加 $Q$ 時下界亦增加，因而證明觀測對數似然不下降。正權重、正定協方差及正有效樣本數，令各公式保持在指定模型內。估計協方差或預先指定方差，是目標的一部分；不受限制的協方差可因分量塌縮而令似然無上界。

練習

思考檢查

1. 在 $x=1$，兩個分量的平均為 $0$ 與 $2$、共同方差為 $1$、權重為 $(3/4,1/4)$。求兩個責任度。若只保留最大責任度所對應的分量，會失去甚麼資訊？

分清相同分量密度與不同混合權重。

思考檢查

2. 對觀測 $(0,2)$，固定責任度為 $R=\begin{bmatrix}0.9&0.1\\0.3&0.7\end{bmatrix}$。指定共同方差時，求 M 步的權重及平均。

先計算每個 $N_k$，再用它規範加權總和。

思考檢查

3. 對觀測 $(0,2,4)$，某分量的責任度為 $(4/5,1/2,1/5)$。求有效樣本數、新平均，以及估計協方差時的一維協方差值。

把圍繞新平均的加權距離平方和除以有效樣本數。

思考檢查

4. 對例題 2 更新後的模型作硬分配及平均更新。比較得到的觀測對數似然與目前值；EM 證明中的哪個等號不再成立？

兩個目前平均的中點仍是 $1$。

解答

解答 · 1. 後驗保留兩個比例

責任度為 $(3/4,1/4)$，因為共同密度會抵消。只保留第一分量，便只記錄一個選定標籤，丟掉了餘下的後驗質量 $1/4$。

解答 · 2. 分數樣本量亦會更新混合權重

有效樣本數為 $(N_1,N_2)=(6/5,4/5)$。所以 $\pi^{\mathrm{new}}=(3/5,2/5)$，且 $\mu_1^{\mathrm{new}}=(0.9\cdot0+0.3\cdot2)/(6/5)=1/2$， 而 $\mu_2^{\mathrm{new}}=(0.1\cdot0+0.7\cdot2)/(4/5)=7/4$。

解答 · 3. 以有效樣本數規範加權散布

有效樣本數是 $N=3/2$，平均為 $((4/5)0+(1/2)2+(1/5)4)/(3/2)=6/5$。 加權距離平方和為 $(4/5)(6/5)^2+(1/2)(4/5)^2+(1/5)(14/5)^2=76/25$， 所以協方差估計為 $(76/25)/(3/2)=152/75$。

解答 · 4. 下界不再與目前似然相切

硬分配令每個分量各有一筆觀測，平均回到 $(0,2)$，權重回到 $(1/2,1/2)$。觀測對數似然由 $-2.886562696$ 降至 $-2.970315405$。硬概率行 $q$ 不同於嚴格為正的後驗，所以 $\mathcal F(q,\theta^{\mathrm{old}})=\ell(\theta^{\mathrm{old}})$ 不再成立。改善該硬分配下界，無法用 EM 的不等式鏈保證觀測似然改善。

可接着閱讀最小二乘與殘差， 或查看課程閱讀書目。