ml課程：最大熵與EM演算法及應用（含程式碼實現）

阿新 • • 發佈：2018-11-25

以下是我的學習筆記，以及總結，如有錯誤之處請不吝賜教。

本文主要介紹最大熵模型與EM演算法相關內容及相關程式碼案例。

關於熵之前的文章中已經學習過，具體可以檢視：ml課程：決策樹、隨機森林、GBDT、XGBoost相關（含程式碼實現），補充一些

基本概念：

資訊量：資訊的度量，即一件事情發生的概率 $P(x_{i})$ 。那麼熵既可以表示為資訊量的期望，也就是 $H(x_{i})$ 。

聯合熵（joint entropy）：是聯合概率分佈或者多個變數的熵。

條件熵（conditionam）：公式如下：

變換一下：

互資訊（mutual information）：是衡量兩個隨機變數之間的相關性，常用來特徵選擇或者降維；他與 $\rho =\frac{COV(X,Y)}{\sqrt{Var(x)}\sqrt{Var(y)}}$

的區別是可以衡量非線性的相關性，具體公式表示：

推倒一下得到：

我們可以根據VN圖來表示以上三個概念的

相對熵（relative entropy），即KL散度（kullback-leibler divergence）:是衡量兩個概率分佈的差異。具體表示為：

性質：

是一個非對稱的近似度量；
值越大分佈之間的差異越大；
始終存在 $D(p||q)\geqslant 0$

關於性質3具體證明如下：

證明：假設 $y_{i}=\frac{q_{i}}{p_{i}}$ ， $f(y_{i})=log\frac{q_{i}}{p_{i}}$ ,同時由於log函式為凹函式，根據琴生不等式得到：E(f(x)) $\leq$ f(E(X)):

將上述幾個式子代入得到： $D(p||q)\geqslant 0$

交叉熵（cross entropy）

：是衡量兩個概率分佈間的差異性資訊，公式如下：

$CH(p,q)=-\sum p(x_{i})logq(x_{i})$

變化一下得到：

$CH(p,q)=H(p)+D(p||q)$

隱變數：表示不可觀測變數Y={y1,y2,y3...,yn}；

顯變數：表示可觀測變數Z={z1,z2,z3...,zn}，兩個變數的性質有：

$P(Y)=\sum_{Z}^{N}P(Z)P(Y|Z)=\sum_{Z}^{N}P(Y,Z)$
$P(Y|Z)=\frac{P(Y,Z)}{P(Z)}=\frac{P(Y,Z)}{\sum P(Y)P(Z|Y)}$

最大似然（極大似然）估計：N個誤差的聯合概率密度為：p(y1，y2，...，yN)，由於它們相互獨立，故有：

$likelihold=\prod_{j}^{N}p(y _{i})$

那麼條件概率的最大似然可以表示為：

$p(y|\Theta )=\prod_{j}^{N}p(y _{i}|\Theta )$

這個概率反映了，在概率密度函式的引數是θ時，得到y這組樣本的概率。需要找到一個引數θ，其對應的似然函式L(θ)最大，這個叫做θ的最大似然估計量，記為:

由於似然函式連乘求導無法求解，因此對其取對數：

$H(\Theta )=ln(p(y|\Theta ))=ln\prod_{j}^{N}p(y _{i}|\Theta )=\sum_{j}^{N}lnp(y _{i}|\Theta)$

最後求導令其為0，得到似然方程引數。

最大熵模型（max entropy model）：

在學習概率模型的時候，在所有可能的概率模型中，熵最大的模型是最好的模型；通俗說就是承認已知事物，對未知事物沒有任何偏見。

目標函式：

$maxH(Y|X)=-\sum \sum p(x_{i},y_{i})logp(y_{i},x_{i})$

然後通過最大熵求解得到 $p^{*}$ 以此來分類：

$p^{*}=argmaxH(Y|X)$

我們要求解以上的目標函式，需要構建一個特徵函式：

同時通過概率的條件知：

$\sum p(y_{i},x_{i})=1$

同時我們定義理想狀態的概率為 $\tilde{p}$ ，那麼可以得到：

總結以上得到最大熵模型：

因此我們構造拉格朗日乘子：

其中： $q=p(y|x)$ .

對其求偏導數：

令其為0，得到：

因此最終得到：

$p^{*}(y|x)=\frac{1}{z_{\lambda }(x)}e^{\sum \lambda _{i}f(x,y)}$

其中：

但是上式很難求得解析解 $\lambda$ ，因此我們使用近似逼近解：

IIS（improved iterative scaling）方法：該方法類似座標上升法，是假設最大熵模型下不斷迭代 $\lambda +\delta$ 使得：

$L(\lambda +\delta )-L(\lambda )\geqslant 0$

EM演算法（Expection Maximization）：

用迭代方式交替求解優化值，這個演算法可以從兩個方向求解，如下圖：

下面我們先從右邊爬山，從最大似然的方向推到EM演算法：

根據前面基本概念，將含有隱變數的極大似然估計的方法寫為：

$L(\Theta )=\prod_{j}^{N}p_{\Theta }(y _{i})=\prod_{j}^{N}\sum_{z}^{N}p_{\Theta }(y _{i}|z)p_{\Theta }(z)$

對其兩邊取對數：

$l(\Theta )=lnL(\Theta )=ln\prod_{j}^{N}p_{\Theta }(y _{i})=\sum_{j}^{N}ln\sum_{z}^{N}p_{\Theta }(y _{i}|z)p_{\Theta }(z)$

我們迭代進行求解,使得 $l(\Theta _{n+1})> l(\Theta _{n})$ ,即：

$l(\Theta )-l(\Theta _{n})=\sum_{j}^{N}ln\sum_{z}^{N}p_{\Theta }(y _{i}|z)p_{\Theta }(z)-\sum_{j}^{N}ln\sum_{z}^{N}p_{\Theta n}(y _{i}|z)p_{\Theta n}(z)$

即：

$l(\Theta )-l(\Theta _{n})=\sum_{j}^{N}[ln\sum_{z}^{N}p_{\Theta }(y _{i}|z)p_{\Theta }(z)-ln\sum_{z}^{N}p_{\Theta n}(y _{i}|z)p_{\Theta n}(z)]$

分別對前後兩項進行湊項和琴生不等式得到：

$l(\Theta )-l(\Theta _{n})\geqslant \sum_{z}^{N}p_{\Theta n}(z|y)ln\frac{p_{\Theta }(y|z)p_{\Theta }(z)}{p_{\Theta n}(z|y)p_{\Theta n}(y)}$

$l(\Theta )\geqslant l(\Theta _{n})+ \sum_{z}^{N}p_{\Theta n}(z|y)ln\frac{p_{\Theta }(y|z)p_{\Theta }(z)}{p_{\Theta n}(z|y)p_{\Theta n}(y)}$

令 $Q(\Theta |\Theta _{n})$ 等於後面的式子，則：

$l(\Theta )\geqslant Q(\Theta |\Theta _{n})$

具體推推倒過程如下圖：

因此我們得到EM演算法流程為：

隨機生成 $\Theta$ ，以此來迭代 $\Theta _{n}\rightarrow \Theta _{n+1}$
E-step：令 $\Theta _{_{n}}$ 為已經求得的數，根據前面得到的式子求 $Q(\Theta |\Theta _{n})$ ，也就是其下邊界函式，也是熵值。
M-step：求導=0，最大化 $\Theta ^{*}$
重複2,3兩步，直至 $\Theta _{n+1}-\Theta _{n}\leqslant \varepsilon$ 為止。

下面我們從左邊上山，從K-means（參考：ml課程：聚類概述及K-means講解（含程式碼實現））方向推到EM演算法：

K-means演算法可以表示為：

隨機選擇k個點作為 $\mu_{k}^{0}$ ,即簇類的中心點；
E-step： $\alpha _{ik}^{1}=\sum_{k=1}^{k}\sum_{i=1}^{N}\alpha _{ik}||x_{i}-\mu_{k}^{0}||^{2}$ ,其中 $\alpha _{ik}=1, or ,0$ 表示一個隱變數概率為0或者1的硬分類；
M-step：對上式求導=0，最小化 $\alpha _{ik}^{1}$ 。
重複2,3兩步，直至 $||\mu _{k}^{n+1}-\mu _{k}^{n}||\leq \varepsilon$ 為止。

EM演算法應用：

高斯混合模型（GMM）:

核心程式碼：

######################################################
# E 步：計算每個模型對樣本的響應度
# Y 為樣本矩陣，每個樣本一行，只有一個特徵時為列向量
# mu 為均值多維陣列，每行表示一個樣本各個特徵的均值
# cov 為協方差矩陣的陣列，alpha 為模型響應度陣列
######################################################
def getExpectation(Y, mu, cov, alpha):
    # 樣本數
    N = Y.shape[0]
    # 模型數
    K = alpha.shape[0]

    # 為避免使用單個高斯模型或樣本，導致返回結果的型別不一致
    # 因此要求樣本數和模型個數必須大於1
    assert N > 1, "There must be more than one sample!"
    assert K > 1, "There must be more than one gaussian model!"

    # 響應度矩陣，行對應樣本，列對應響應度
    gamma = np.mat(np.zeros((N, K)))

    # 計算各模型中所有樣本出現的概率，行對應樣本，列對應模型
    prob = np.zeros((N, K))
    for k in range(K):
        prob[:, k] = phi(Y, mu[k], cov[k])
    prob = np.mat(prob)

    # 計算每個模型對每個樣本的響應度
    for k in range(K):
        gamma[:, k] = alpha[k] * prob[:, k]
    for i in range(N):
        gamma[i, :] /= np.sum(gamma[i, :])
    return gamma

######################################################
# M 步：迭代模型引數
# Y 為樣本矩陣，gamma 為響應度矩陣
######################################################
def maximize(Y, gamma):
    # 樣本數和特徵數
    N, D = Y.shape
    # 模型數
    K = gamma.shape[1]

    #初始化引數值
    mu = np.zeros((K, D))
    cov = []
    alpha = np.zeros(K)

    # 更新每個模型的引數
    for k in range(K):
        # 第 k 個模型對所有樣本的響應度之和
        Nk = np.sum(gamma[:, k])
        # 更新 mu
        # 對每個特徵求均值
        for d in range(D):
            mu[k, d] = np.sum(np.multiply(gamma[:, k], Y[:, d])) / Nk
        # 更新 cov
        cov_k = np.mat(np.zeros((D, D)))
        for i in range(N):
            cov_k += gamma[i, k] * (Y[i] - mu[k]).T * (Y[i] - mu[k]) / Nk
        cov.append(cov_k)
        # 更新 alpha
        alpha[k] = Nk / N
    cov = np.array(cov)
    return mu, cov, alpha

######################################################
# 高斯混合模型 EM 演算法
# 給定樣本矩陣 Y，計算模型引數
# K 為模型個數
# times 為迭代次數
######################################################
def GMM_EM(Y, K, times):
    Y = scale_data(Y)  #定義scale_data()資料預處理，將其縮放至0和1之間
    mu, cov, alpha = init_params(Y.shape, K)  #定義init_params()初始化模型引數
    for i in range(times):
        gamma = getExpectation(Y, mu, cov, alpha)  
        mu, cov, alpha = maximize(Y, gamma)
    debug("{sep} Result {sep}".format(sep="-" * 20))  #定義debug()表示除錯輸出變數
    debug("mu:", mu, "cov:", cov, "alpha:", alpha, sep="\n")
    return mu, cov, alpha

具體實現歡迎關注：我的github

後續將繼續介紹EM演算法相關的應用，To be continue...

ml課程：最大熵與EM演算法及應用（含程式碼實現）

基本概念：

最大熵模型（max entropy model）：

EM演算法（Expection Maximization）：

ml課程：最大熵與EM演算法及應用（含程式碼實現）

ml課程：概率圖模型—貝葉斯網路、隱馬爾可夫模型相關（含程式碼實現）

ml課程：SVM相關（含程式碼實現）

ml課程：決策樹、隨機森林、GBDT、XGBoost相關（含程式碼實現）

ml課程：線性迴歸、邏輯迴歸入門（含程式碼實現）

ml課程：聚類概述及K-means講解（含程式碼實現）

機器學習筆記：最大熵（模型，推導，與似然函式關係的推導，求解）

C語言：最大值與最小值對調

最大熵與邏輯迴歸的等價性

OpenCV程式設計：最大熵閾值分割演算法實現(程式碼可執行)

自然語言處理之一：最大熵模型

影象閾值分割：最大熵法

ml課程：XGBoost和lightGBM工具庫學習及相關案例程式碼

最大熵模型及其演算法

Java併發程式設計（3）：執行緒掛起、恢復與終止的正確方法（含程式碼）

10 大深度學習架構：計算機視覺優秀從業者必備（附程式碼實現）

吳恩達機器學習課程：完全用Python完成，可以的！（附程式碼）

乾貨 | 十大經典排序演算法最強總結（內含程式碼實現）

機器學習筆記（十）EM演算法及實踐（以混合高斯模型（GMM）為例來次完整的EM）

【Java併發程式設計】之二十：併發新特性—Lock鎖和條件變數（含程式碼）

ml課程：最大熵與EM演算法及應用（含程式碼實現）

基本概念：

最大熵模型（max entropy model）：

EM演算法（Expection Maximization）：

相關推薦