第六章-邏輯回歸與最大熵模型

本文轉載自查看原文 2019-09-06 13:34 281 統計學習方法-李航/ 機器學習/ 學習筆記/ 最大熵模型/ 邏輯回歸

邏輯斯諦回歸是統計學習中的經典分類方法，和最大熵模型相比，具有以下的共同點和區別：

共同點
- 都屬於概率模型，該模型要尋找的是給定一個x，得到輸出變量Y的概率分布P(Y|x)，如果是二分類，Y取值為0或1，如果是多分類，Y有K個不同的類別。
- 都屬於對數線性模型，對概率分布P(Y|x)取對數，可得lnP(Y|x)=w * x關於x的線性函數。
兩個模型之前的區別和聯系是Logistic回歸和最大熵模型都屬於判別模型。在最大熵模型中，不僅x有隨機性，Y也具有隨機性，是一個隨機變量。

Logistic回歸模型（Logistic Regression Model）

給定輸入變量，輸出變量為 $y \in {0,1}$ ，將 y=1 的概率記作 $\pi(x)=P(Y=1|x) \in [0,1]$ ，我們已經知道該模型是一個線性模型，可以用 $w \cdot x$ 線性函數來表示，由於 $w \cdot x \in (-\infty, +\infty)$ ，那么要如何將 $\pi(x)$ 與 $w \cdot x$ 對應起來呢？這就需要用到一個變換，該變換稱為Logit變換。

Logit變換： $\displaystyle \text{logit}(\pi(x))=\log \frac{\pi(x)}{1-\pi(x)} \in (-\infty, +\infty)$
$\displaystyle \therefore \ln \frac{\pi(x)}{1-\pi(x)}=w \cdot x$ ，這個就是Logistic回歸模型的一個形式。
$\displaystyle \therefore \pi(x)=\frac{\exp (w \cdot x)}{1+\exp (w \cdot x)}$ ，其中 $\pi(x)$ 就是給定的條件下， Y=1 的概率。
所以可得下面兩個公式： $\displaystyle \log \frac{P(Y=1 | x)}{1-P(Y=1 | x)}=w \cdot x$
$\displaystyle P(Y=1 | x)=\frac{\exp (w \cdot x)}{1+\exp (w \cdot x)}$
有了這個模型之后，需要求解參數，一旦求出，那么任意給定一個輸入變量，就可以得到 Y=1 的概率，如果該概率值大於0.5，就將該類判定為1，如果小於0.5，將該類判定為0。
求解使用的方法是極大似然估計，給定參數，求樣本的聯合概率密度，最大化該聯合概率，從而求出參數。

已知 $y_i \in {0,1}$ ， $\pi(x)$
可得 $P_w(y|x)=\pi(x)^y[1-\pi(x)]^{1-y}$
似然函數為 $\displaystyle L(w)=\prod_{i=1}^{N} \pi(x)^{y_i}[1-\pi(x)]^{1-y_i}$

多項Logistic回歸模型

在多分類Logistic回歸模型中，取值為 $1,\cdots,K$ ，那么Logit變換是 $\displaystyle \ln \frac{p_(Y=k|x)}{P(Y=K|x)} = w_k \cdot x \in (-\infty, +\infty)$ ，因為的取值為 $1,\cdots, K-1$ ，所以變為 w_k ，一共需要求取 K-1 個參數向量。
所以可得下面兩個公式：
$\displaystyle P(Y=k|x)=\frac{\exp (w_k \cdot x)}{1+\sum_{k=1}^{K-1} \exp (w_k \cdot x)}, k=1,\cdots,K-1$
$\displaystyle P(Y=K|x)=\frac{1}{1+\sum_{k=1}^{K-1} \exp (w_k \cdot x)}$

最大熵模型（Maximum Entropy Model）

最大熵模型主要是在自然語言處理中應用的比較廣泛。最大熵模型是由最大熵原理推導實現的。而最大熵原理指的是在滿足約束條件的模型集合中選擇熵最大的模型。那么我們為何需要選擇熵最大的模型呢？在決策樹中，我們已經提到熵的概念，熵其實表示的是數據的混亂程度。在我們不知道任何信息的情況下，我們就只能假設數據在各個取值上比較平均，比較散亂，即比較數據的混亂程度。這樣描述，可能不太清楚，直接看下面的一個例題：

假設隨機變量X由5個取值{A,B,C,D,E}，要估計取各個值的概率P(A)，P(B)，P(C)，P(D)，P(E)。

對於一個概率問題，我們可以知道：

P(A) + P（B）+ P(C) + P(D) + P(E) =1

根據最大熵原理，取值隨機且平均，那么P(A)=P(B)=P(C)=P(D)=P(E)=1/5。

當然有時候我們也能得到一些先驗條件，比如:P(A)+P(B)=3/10。對於這種問題，我們在問題中就出現了一些約束條件，此時的熵就是條件熵，即：

\[H(y|x)=-\sum_y P(y|x)lnP(y|x) \]

對於一個訓練數據集樣本，我們可以確定聯合分布P(X,Y)的經驗分布和邊緣分布P(X)的經驗分布，但是求解條件熵仍然存在困難，因為在所有數據集中x的取值是不一樣的，而且不固定。此時，我們可以采用經驗分布的期望來進行比較。為何可以這樣說呢？這就需要我們了解一個二值函數——特征函數(Feature Function)。該函數描述的是輸入x和輸出y之間的某一個事實，可以定義為：

\[f(x,y) = \left\{ \begin{aligned} 1 & , & x與y滿足某一個事實 \\ 0 & , & 否則 \end{aligned} \right. \]

直觀的看，比較難以理解。可以舉個例子，比如我們在做一個機器英譯漢的翻譯的時候，如果你輸入一句話中帶有take單詞的英文句子I take a lot of money to buy a house.,此時take的含義表示花費，但是如果你輸入She decided to take a bus to work since yesterday.,此時take的含義表示乘坐。take單詞后面跟着不同單詞，就需要翻譯成不同的意思，此時next word是一個條件，在這種條件下，take翻譯成不同的含義，滿足這一個事實。由此，可以得到特征函數關於經驗分布P(X,Y)和P(Y|X)與P(X)的期望值如下表述：

\[經驗分布\widetilde{P}(X,Y)的期望：\\ E_{\widetilde{P}}(f)= \sum_{x,y}\widetilde{P}(x,y)f(x,y) \\ 經驗分布\widetilde{P}(Y|X)與\widetilde{P}(X)的期望：\\ E_{P}(f)= \sum_{x,y}\widetilde{P}(x)P(y|x)f(x,y) \\ \]

如果模型能夠獲取訓練數據集中的信息，那么可以假設兩個期望相等。得到以下最大熵模型。

\[最大熵模型定義：假設滿足所有約束條件的模型集合為:\\ \zeta = \{P \in Q|E_p(f_i)=E_{\widetilde{P}}(f_i),i=1,2,3...n\}\\ 定義在條件概率分布P(Y|X)上的條件熵為：\\ H(P)=-\sum_{x,y}\widetilde{P}(x)P(y|x)logP(y|x) \]

最大熵模型的數學推導

最大熵模型學習過程中，最重要的是最大熵模型的學習，也就是求解約束條件下的最優化問題。

按照最優化的習慣問題，如果需要求H(P)的最大值，可以轉化成求-H(P)的最小值。

\[max_{P \in C} H(P)=-\sum_{x,y}\widetilde{P}(x)P(y|x)logP(y|x) \\ s.t. E_p(f_i)=E_{\widetilde{P}}(f_i),i=1,2,3...n\\ \sum_yP(y|x)=1 \\ =>轉化成： min_{P \in C} -H(P)=\sum_{x,y}\widetilde{P}(x)P(y|x)logP(y|x) \\ s.t. E_p(f_i)=E_{\widetilde{P}}(f_i),i=1,2,3...n\\ \sum_yP(y|x)=1 \]

此處可以將約束最優化問題轉換成無約束最優化的對偶問題，通過求解對偶問題來求解原始問題。

第一步，首先針對原始問題定義拉格朗日乘子w0,w1...wn,拉格朗日函數為：

\[L(P,w)=-H(p)+w_0(1-\sum_yP(y|x))+\sum_{i=1}^nw_i(E_p(f_i)-E_{\widetilde{P}}(f_i))\\ =\sum_{x,y}\widetilde{P}(x)P(y|x)logP(y|x)+w_0(1-\sum_yP(y|x))+\sum_{i=1}^nw_i(E_p(f_i)-E_{\widetilde{P}}(f_i))\\ 原始問題：min_{p \in C}max_wL(P,w)\\ =>對偶問題：max_wmin_{p \in C}L(P,w) \]
由於拉格朗日函數是凸函數，因此原始問題的解和對偶問題的解是等價的，即求解原始問題的極小極大值可以轉換成極大極小值。對L(P,w)求偏導。

\[令對偶函數為\Psi(w)=min_{P \in C} L(P,w)=L(P_w,w)\\ 因此對偶函數求解得:P_w=argmin_{P \in C}L(P,w)=P_w(y|x)\\ \frac{\partial L(P,w)}{\partial w} = \sum_{x,y}(\widetilde{P}(x)logP(y|x)+1)-\sum_yw_0-\sum_{x,y}(\widetilde{P}(x)\sum_{i=1}^nw_if_i(x,y))\\ = \sum_{x,y}\widetilde{P}(x)\{logP(y|x)+1-\sum_yw_0-\sum_{x,y}(\sum_{i=1}^nw_if_i(x,y)\}\\ 令偏導數為0，\widetilde{P}(x) >0的情況下，可得到： \\ P(y|x)=exp\{\sum_{i=1}^nw_if_i(x,y)+w_0-1\}\\ =\frac{exp(\sum_{i=1}^nw_if_i(x,y))}{exp(1-w_0)}\\ 由於\sum_y(P(y|x))=1,得：\\ P_w(y|x)=\frac{1}{Z_w(x)}exp(\sum_{i=1}^nw_if_i(x,y)),其中Z_w(x)=exp(\sum_{i=1}^nw_if_i(x,y)),\\ 在這里，Z_w(x)被稱為歸一化因子，f_i是特征函數，w_i是特征的權值。 \]
求解對偶問題的外部的極大值問題。

\[w*=argmax_w \Psi(w)\\ 即：P*=P_{w*}=P_{w*}(y|x)是學習到的最優化模型，也就是最大熵模型。 \]

Python實現

自編程實現

## 使用梯度下降法實現邏輯斯蒂回歸算法
import numpy as np
import time
import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import Axes3D
plt.rcParams['font.sans-serif']=['SimHei'] #用來正常顯示中文標簽
plt.rcParams['axes.unicode_minus']=False #用來正常顯示負號

class LogisticRegression:
    def __init__(self,learn_rate=0.1,max_iter=10000,tol=1e-2):
        self.learn_rate=learn_rate  #學習率
        self.max_iter=max_iter      #迭代次數
        self.tol=tol                #迭代停止閾值
        self.w=None                 #權重

    def preprocessing(self,X):
        """將原始X末尾加上一列，該列數值全部為1"""
        row=X.shape[0]
        y=np.ones(row).reshape(row, 1)
        X_prepro =np.hstack((X,y))
        return X_prepro

    def sigmod(self,x):
        return 1/(1+np.exp(-x))

    def fit(self,X_train,y_train):
        X=self.preprocessing(X_train)
        y=y_train.T
        #初始化權重w
        self.w=np.array([[0]*X.shape[1]],dtype=np.float)
        k=0
        for loop in range(self.max_iter):
            # 計算梯度
            z=np.dot(X,self.w.T)
            grad=X*(y-self.sigmod(z))
            grad=grad.sum(axis=0)
            # 利用梯度的絕對值作為迭代中止的條件
            if (np.abs(grad)<=self.tol).all():
                break
            else:
                # 更新權重w 梯度上升——求極大值
                self.w+=self.learn_rate*grad
                k+=1
        print("迭代次數：{}次".format(k))
        print("最終梯度：{}".format(grad))
        print("最終權重：{}".format(self.w[0]))

    def predict(self,x):
        p=self.sigmod(np.dot(self.preprocessing(x),self.w.T))
        print("Y=1的概率被估計為：{:.2%}".format(p[0][0]))  #調用score時，注釋掉
        p[np.where(p>0.5)]=1
        p[np.where(p<0.5)]=0
        return p

    def score(self,X,y):
        y_c=self.predict(X)
        error_rate=np.sum(np.abs(y_c-y.T))/y_c.shape[0]
        return 1-error_rate

    def draw(self,X,y):
        # 分離正負實例點
        y=y[0]
        X_po=X[np.where(y==1)]
        X_ne=X[np.where(y==0)]
        # 繪制數據集散點圖
        ax=plt.axes(projection='3d')
        x_1=X_po[0,:]
        y_1=X_po[1,:]
        z_1=X_po[2,:]
        x_2=X_ne[0,:]
        y_2=X_ne[1,:]
        z_2=X_ne[2,:]
        ax.scatter(x_1,y_1,z_1,c="r",label="正實例")
        ax.scatter(x_2,y_2,z_2,c="b",label="負實例")
        ax.legend(loc='best')
        # 繪制p=0.5的區分平面
        x = np.linspace(-3,3,3)
        y = np.linspace(-3,3,3)
        x_3, y_3 = np.meshgrid(x,y)
        a,b,c,d=self.w[0]
        z_3 =-(a*x_3+b*y_3+d)/c
        ax.plot_surface(x_3,y_3,z_3,alpha=0.5)  #調節透明度
        plt.show()

Sklearn庫實現

from sklearn.linear_model import LogisticRegression
import numpy as np

def main():
    # 訓練數據集
    X_train=np.array([[3,3,3],[4,3,2],[2,1,2],[1,1,1],[-1,0,1],[2,-2,1]])
    y_train=np.array([1,1,1,0,0,0])
    # 選擇不同solver,構建實例，進行訓練、測試
    methodes=["liblinear","newton-cg","lbfgs","sag","saga"]
    res=[]
    X_new = np.array([[1, 2, -2]])
    for method in methodes:
        clf=LogisticRegression(solver=method,intercept_scaling=2,max_iter=1000)
        clf.fit(X_train,y_train)
        # 預測新數據
        y_predict=clf.predict(X_new)
        #利用已有數據對訓練模型進行評價
        X_test=X_train
        y_test=y_train
        correct_rate=clf.score(X_test,y_test)
        res.append((y_predict,correct_rate))

實現結果圖：

免責聲明！

本站轉載的文章為個人學習借鑒使用，本站對版權不負任何法律責任。如果侵犯了您的隱私權益，請聯系本站郵箱yoyou2525@163.com刪除。

猜您在找 邏輯回歸與最大熵模型邏輯回歸和最大熵模型李航-統計學習方法-筆記-6：邏輯斯諦回歸與最大熵模型最大熵模型邏輯斯蒂回歸3 -- 最大熵模型之改進的迭代尺度法(IIS) 第六章最大熵模型 VLSI基礎-- 第六章時序邏輯電路第六章字典最大熵模型推導