最大熵模型詳解與Python實現

一 最大熵原理

首先介紹一下什麼是熵 ，設X是一個取有限個值的離散隨機變數 其概率分佈如下

則隨機變數X的熵可以表示為
並且我們知道熵越大，不確定性就越大

熵的具體講解參考文章資訊熵公式的詳細講解

最大熵原理說的是，在滿足已知的約束條件的情況下，所有可能的概率模型中，熵最大的模型就是最好的模型。

在什麼情況下熵最大，事實上當一個隨機變數符合均勻分佈時，熵最大

通俗的解釋一下

想象一下如下場景

你的朋友問 明天你要做什麼

你回答說 我80%的可能性是在家睡覺 20%的可能性是寫作業，那你的朋友心裡就有底了，這貨明天要在家睡覺（朋友差不多知道你明天要做什麼了,確定性高)

但是如果你回答說 我50%的可能性是睡覺，50%的可能性是寫作業，這兩種情況的概率一樣，你的朋友無法做出權衡，仍是一頭霧水不知道你明天要做什麼 (朋友不清楚你明天要做什麼，確定性低)

在上面的場景中 概率分佈是這樣的 隨機變數X是 明天你要做什麼

睡覺	寫作業
p	1-p

我們可以看出當 p=0.5時，X服從均勻分佈 所有事件發生的概率一樣，不確定性最高,也就是熵最大

熵滿足如下公式
IXI 是 X 的取值個數，當且僅當 X 的分佈是均勻分佈時右邊的等號成立

為什麼要選擇熵最大的模型呢？

舉一個扔骰子的例子

我現在隨機扔一個骰子，問你扔到6的概率是多少，你可能會說是1/6

當我說了 扔到1的概率是2/6時，問你扔到6的概率是多少，你可能會說2/15

實際上我們在潛意識中都是會這樣回答的，我們在什麼也不知道的情況下，會假設，扔到所有數字的可能性是一樣的， 我們在知道某些條件之後，會在這些條件的基礎上，認為扔到其他數字依然是等可能的，這樣回答是最穩妥的，風險最低.

直觀地，最大熵原理認為要選擇地概率模型首先必須滿足已有的事實，即約束條件，在沒有更多資訊的情況下，那些不確定的部分都是"等可能的" 。 最大熵原理通過熵的最大化表示等可能性。

二 最大熵模型

最大熵模型是一個條件概率分佈P(Y|X)

給定一個訓練資料集 T={(X₁,Y₁) , (X₂,Y₂) , … ,(X_N，Y_N)}

學習的目的是用最大熵原理選擇最好的分類模型

下面介紹一些符號

根據給定的訓練資料集，聯合分佈P(X,Y)的經驗分佈和邊緣分佈P(X)的經驗分佈如下表示

特徵函式 f(x,y) 表示輸入x和輸入y之間滿足某一個事實

定義如下

舉個例子

特徵函式 f(x， y) 關於經驗分佈$\tilde{P}$(X， Y)的期望值如下

特徵函式 f(x,y) 關於模型P(Y|X) 與 經驗分佈 $\tilde{P}$(X) 的期望值如下

如果我們的模型可以從訓練資料中獲取資訊，可以假設

這便是我們模型的約束條件 。假如有 n 個特徵函式 f_i(x, y) , i = 1， 2…n，那麼就有 n 個約束條件。

下面給出最大熵模型的定義

三 最大熵模型學習過程

對於給定的訓練資料集 T = {(X₁， Y₁)， (X₂， Y₂)， … ， (X_N， Y_N)} 以及特徵函式 f_i(x,y) ， i=1,2. . .,n，最大熵模型的學習等價於約束最優化問題:

將求最大值問題改寫為等價的求最小值問題:

簡單介紹一下 拉格朗日函式

通過引進拉格朗日乘子w₀,w₁,w₂,…,w_n 定義一個拉格朗日函式L(P|w）

最優化的原始問題是

為什麼這樣表達？ 建議先閱讀文章 對偶問題 作者fkyyly

原始問題的對偶問題如下

由於拉格朗日函式 L(P， ω) 是 P 的凸函式 (因為 L(P,w) 對P(y|x) 的二階偏導數大於0），原始問題的解與對偶問題的解是等價的。這樣，可以通過求解對偶問題來求解原始問題。

具體地 ，求 L(P,w) 對P(y|x) 的偏導數

令偏導數=0 可以求出最大熵模型 P_w(y|x) 下面給出具體求解過程

四 最大熵模型的極大似然估計

所謂最大似然，即最大可能，在“模型已定，引數θ未知”的情況下，通過觀測資料估計引數θ的一種思想或方法，換言之，解決的是取怎樣的引數θ使得產生已得觀測資料的概率最大的問題。

極大似然函式的一般形式如下

關於這個公式是如何來的，可參考文章最大熵模型中的對數似然函式的解釋

對上式兩邊取對數可得

把p(x)換成p(x,y)也是同樣的道理

由於第二項 是一個常數

所以似然函式最終形式為

將我們前面求的最大熵模型代入 極大似然函式中

而把最大熵模型代入對偶函式中得:

可得

於是證明了最大熵模型學習中的對偶函式極大化等價於最大熵模型的極大似然估計這一事實。

五 最大熵模型學習的最優化演算法

改進的迭代尺度法

改進的迭代尺度法 (improved iterative scaling, IIS) 是一種最大熵模型學習的最優化演算法。

L(w)公式如下

則

因為 -log a>= 1-a a>0
則

則

介紹一個不等式 jensen不等式

令不等式右端為

改進的迭代尺度演算法

六 程式碼實現

最大熵模型

class Maximum_Entropy:
    def __init__(self,X_Train,Y_Train,X_Test,Y_Test,iteration):
        self.X_Train=X_Train #訓練集特徵
        self.Y_Train=Y_Train  #訓練集標籤
        self.X_Test=X_Test #測試集特徵
        self.Y_Test=Y_Test #測試集標籤
        self.feature_num = X_Train.shape[1]
        self.n = 0  # 一共有多少個特徵函式
        self.fi = self.cal_fi()
        self.w=[0]*self.n
        self.sigma=[0]*self.n
        self.fi_index,self.index_fi,self.index_feature= self.create_fi_map()
        self.iteration=iteration

    def cal_fi(self):
        #在本方法中為每一個特徵都建一個字典 並把這些字典存放到fi列表中
        #每個特徵的字典 的鍵是這個特徵的某個取值x 以及這個取值所在的那個樣本的標籤y 所組成的元組(x,y) 它的value是在訓練集中(x,y) 出現了多少次
        #一個元組(x,y)代表的就是一個特徵函式 fi
        fi=[]
        for i in range(self.feature_num): # fi 中是 feature_num個 字典
            fi.append({})
        for i in range(self.X_Train.shape[0]):#遍歷所有資料
            for j in range(self.feature_num):
                if (self.X_Train[i][j],self.Y_Train[i][0]) not in fi[j]:
                     fi[j][(self.X_Train[i][j],self.Y_Train[i][0])]=1
                else:
                     fi[j][(self.X_Train[i][j], self.Y_Train[i][0])] += 1
        for dict in fi:
            self.n+=len(dict)

        return fi

    def create_fi_map(self):
        #本方法是給self.fi中的所有"特徵函式" 賦一個索引值的 為了之後方便查詢
        index=0
        fi_index={}
        index_fi=[0]*self.n
        index_feature=[0]*self.n
        for i in range(self.feature_num):
            for (x,y)  in self.fi[i]:
                fi_index[(x,y)]=index
                index_fi[index]=(x,y)
                index_feature[index]=i
                index+=1
        return fi_index,index_fi,index_feature

    def cal_class(self):
        #本方法是計算資料一共被分為幾類 每一類是如何表示的
        list=[]
        for i in range(self.Y_Train.shape[0]):
            if self.Y_Train[i][0] not in list:
                list.append(self.Y_Train[i][0])
        return list

    def cal_p(self,X,Y):
    #本方法是求解最大熵模型的，對應於李航《統計學習方法》 公式6.22
        sum=0
        for i in range(self.feature_num):
            if (X[i],Y) in self.fi[i]:
               index=self.fi_index[(X[i],Y)]
               sum+=self.w[index]
       # print(sum)
        e1=math.exp(sum)
        label=self.cal_class()
        label.remove(Y)
        e2=0
        for y in label:
           sum=0
           for i in range(self.feature_num):
             if (X[i],y) in self.fi[i]:
               index=self.fi_index[(X[i],y)]
               sum+=self.w[index]
              # print(self.w)
              # print(sum)
               e2+=math.exp(sum)
        return e1/(e1+e2)

    def cal_E_P(self):
        #計算特徵函式關於經驗分佈P(X,Y)的期望值
        ep=[0]*self.n
        for i in range(self.feature_num):
            for (x,y) in self.fi[i]:
                index=self.fi_index[(x,y)]
                ep[index] = self.fi[i][(x,y)]/self.X_Train.shape[0]
        return ep
    def fjing(self,X,Y):# f#(x,y)
        sum=0
        for i in range(self.feature_num):
            if (X[i],Y) in self.fi[i]:
                sum+=1
        return sum
    def cal_px(self,feature,x_value):
        X = self.X_Train[:, feature:feature + 1]
       # print(X)
        dict={}
        for i in range(X.shape[0]):
            if X[i][0] not in dict:
                dict[X[i][0]]=1
            else:
                dict[X[i][0]]+=1
        return dict[x_value] / X.shape[0]

    def gsigma(self,index):
        ep=self.cal_E_P()
        p1=ep[index]
        feature=self.index_feature[index]
        sigma=self.sigma[index]
        label=self.cal_class()
        p2=0
        for i in range(self.X_Train.shape[0]):
            inner_sum=0
            for y in label:
                if (self.X_Train[i][feature],y) == self.index_fi[index]: #滿足特徵函式i
                   inner_sum+= self.cal_p(self.X_Train[i],y)*math.exp(sigma*self.fjing(self.X_Train[i],y))
            p2+=self.cal_px(feature,self.X_Train[i][feature])*inner_sum
        return p1-p2

    def cal_gsigma_derivatives(self,index):
        feature = self.index_feature[index]
        sigma = self.sigma[index]
        label = self.cal_class()
        p = 0
        for i in range(self.X_Train.shape[0]):
            inner_sum = 0
            for y in label:
                if (self.X_Train[i][feature], y) == self.index_fi[index]:  # 滿足特徵函式i
                    inner_sum += self.cal_p(self.X_Train[i], y)*self.fjing(self.X_Train[i], y)* math.exp(sigma * self.fjing(self.X_Train[i], y))
            p += self.cal_px(feature, self.X_Train[i][feature]) * inner_sum
        return -p
    def fit(self):#訓練模型 引數的設定
        for i in range(self.iteration):
            print(i)
            for j in range(self.n):
                self.sigma[j] = self.sigma[j]-(self.gsigma(j)/self.cal_gsigma_derivatives(j))
                self.w[j] = self.w[j] + self.sigma[j]
            print(self.predict())
    def predict(self):
        right=0
        label=self.cal_class()
        for i in range(self.X_Test.shape[0]):
            max = 0
            maxy = 0
            for y in label:
                p=self.cal_p(self.X_Test[i],y)
                if max<p:
                   max=p
                   maxy=y
            if maxy==self.Y_Test[i][0]:
                right+=1
        return right/self.X_Test.shape[0]