李航《统计学习方法》——第六章最大熵模型

超、凢脫俗 2022-05-30 10:28 220阅读 0赞

由于网上资料很多，这里就不再对算法原理进行推导，仅给出博主用Python实现的代码，供大家参考

**适用问题**：多类分类  
  
下面用改进的迭代尺度法（IIS）学习最大熵模型，将特征函数定义为：  
  
[![f(x,y)=\\left\\\{\\begin\{matrix\} 1 & (x,y)\\in train set \\\\ 0 & else \\end\{matrix\}\\right.][f_x_y_left_begin_matrix_ 1 _ _x_y_in train set _ 0 _ else _end_matrix_right.]][f_x_y_left_begin_matrix_ 1 _ _x_y_in train set _ 0 _ else _end_matrix_right._f_x_y_left_begin_matrix_ 1 _ _x_y_in train set _ 0 _ else _end_matrix_right.]  
  
与其他分类器不同的是，最大熵模型中的f(x,y)中的x是单独的一个特征，不是一个n维特征向量，因此我们需要对每个维度特征加一个区分标签;如X=(x0,x1,x2,…)变为X=(0\_x0,1\_x1,2\_x2,…)

**测试数据集**：[train.csv][]  
  
**实现代码**：

# encoding=utf-8
    
    import pandas as pd
    import time
    import math
    
    from collections import defaultdict
    
    from sklearn.cross_validation import train_test_split
    from sklearn.metrics import accuracy_score
    
    
    class MaxEnt(object):
    
        def init_params(self, X, Y):
            self.X_ = X
            self.Y_ = set()
    
            self.cal_Vxy(X, Y)
    
            self.N = len(X)                 # 训练集大小,如P59例子中为15
            self.n = len(self.Vxy)          # 数据集中(x,y)对数,如P59例子中为6+3+3+5=17对
            self.M = 10000.0                # 设置P91中的M，可认为是学习速率
    
            self.build_dict()
            self.cal_Pxy()
    
        def cal_Vxy(self, X, Y):
            ''' 计算v(X=x,Y=y),P82 '''
            self.Vxy = defaultdict(int)
    
            for i in range(len(X)):
                x_, y = X[i], Y[i]
                self.Y_.add(y)
    
                for x in x_:
                    self.Vxy[(x, y)] += 1
    
        def build_dict(self):
            self.id2xy = {}
            self.xy2id = {}
    
            for i, (x, y) in enumerate(self.Vxy):
                self.id2xy[i] = (x, y)
                self.xy2id[(x, y)] = i
    
        def cal_Pxy(self):
            ''' 计算P(X=x,Y=y),P82 '''
            self.Pxy = defaultdict(float)
            for id in range(self.n):
                (x, y) = self.id2xy[id]
                self.Pxy[id] = float(self.Vxy[(x, y)]) / float(self.N)
    
    
        def cal_Zx(self, X, y):
            ''' 计算Zw(x/yi)，根据P85公式6.23，Zw(x)未相加前的单项 '''
            result = 0.0
            for x in X:
                if (x,y) in self.xy2id:
                    id = self.xy2id[(x, y)]
                    result += self.w[id]
            return (math.exp(result), y)
    
        def cal_Pyx(self, X):
            ''' 计算P(y|x),根据P85公式6.22 '''
            Pyxs = [(self.cal_Zx(X, y)) for y in self.Y_]
            Zwx = sum([prob for prob, y in Pyxs])
            return [(prob / Zwx, y) for prob, y in Pyxs]
    
        def cal_Epfi(self):
            ''' 计算Ep(fi),根据P83最上面的公式 '''
            self.Epfi = [0.0 for i in range(self.n)]
    
            for i, X in enumerate(self.X_):
                Pyxs = self.cal_Pyx(X)
    
                for x in X:
                    for Pyx, y in Pyxs:
                        if (x,y) in self.xy2id:
                            id = self.xy2id[(x, y)]
    
                            self.Epfi[id] += Pyx * (1.0 / self.N)
    
    
        def train(self, X, Y):
            ''' IIS学习算法 '''
            self.init_params(X, Y)
    
            # 第一步： 初始化参数值wi为0
            self.w = [0.0 for i in range(self.n)]
    
            max_iteration = 500  # 设置最大迭代次数
            for times in range(max_iteration):
                print("the number of iterater : %d " % times)
    
                # 第二步：求δi
                detas = []
                self.cal_Epfi()
                for i in range(self.n):
                    deta = 1 / self.M * math.log(self.Pxy[i] / self.Epfi[i])  # 指定的特征函数为指示函数，因此E~p(fi)等于Pxy
                    detas.append(deta)
    
                # if len(filter(lambda x: abs(x) >= 0.01, detas)) == 0:
                # break
    
                # 第三步：更新Wi
                self.w = [self.w[i] + detas[i] for i in range(self.n)]
    
        def predict(self, testset):
            results = []
            for test in testset:
                result = self.cal_Pyx(test)
                results.append(max(result, key=lambda x: x[0])[1])
            return results
    
    
    def rebuild_features(features):
        ''' 最大熵模型中的f(x,y)中的x是单独的一个特征,不是一个n维特征向量，因此我们需要对每个维度特征加一个区分标签 具体地：将原feature的（a0,a1,a2,a3,a4,...） 变成 (0_a0,1_a1,2_a2,3_a3,4_a4,...)形式 '''
        new_features = []
        for feature in features:
            new_feature = []
            for i, f in enumerate(feature):
                new_feature.append(str(i) + '_' + str(f))
            new_features.append(new_feature)
        return new_features
    
    
    if __name__ == '__main__':
    
        print("Start read data...")
    
        time_1 = time.time()
    
        raw_data = pd.read_csv('../data/train.csv', header=0)  # 读取csv数据
        data = raw_data.values
    
        features = data[:5000:, 1::]
        labels = data[:5000:, 0]
    
        # 避免过拟合，采用交叉验证，随机选取33%数据作为测试集，剩余为训练集
        train_features, test_features, train_labels, test_labels = train_test_split(features, labels, test_size=0.33, random_state=0)
    
        train_features = rebuild_features(train_features)
        test_features = rebuild_features(test_features)
    
        time_2 = time.time()
        print('read data cost %f seconds' % (time_2 - time_1))
    
        print('Start training...')
        met = MaxEnt()
        met.train(train_features, train_labels)
    
        time_3 = time.time()
        print('training cost %f seconds' % (time_3 - time_2))
    
        print('Start predicting...')
        test_predict = met.predict(test_features)
        time_4 = time.time()
        print('predicting cost %f seconds' % (time_4 - time_3))
    
        score = accuracy_score(test_labels, test_predict)
    print("The accruacy score is %f" % score)

代码可从这里[maxEnt/maxEnt.py][maxEnt_maxEnt.py]获取  
  
**运行结果**：  
  
![maxEnt_result.png][]

[f_x_y_left_begin_matrix_ 1 _ _x_y_in train set _ 0 _ else _end_matrix_right.]: https://latex.codecogs.com/gif.latex?f%28x,y%29=\left\%7b\begin%7bmatrix%7d&space;1&space;&&space;%28x,y%29\in&space;train&space;set&space;\\&space;0&space;&&space;else&space;\end%7bmatrix%7d\right.
[f_x_y_left_begin_matrix_ 1 _ _x_y_in train set _ 0 _ else _end_matrix_right._f_x_y_left_begin_matrix_ 1 _ _x_y_in train set _ 0 _ else _end_matrix_right.]: https://private.codecogs.com/eqnedit.php?latex=f%28x,y%29=\left\%7b\begin%7bmatrix%7d&space;1&space;&&space;%28x,y%29\in&space;train&space;set&space;\\&space;0&space;&&space;else&space;\end%7bmatrix%7d\right.
[train.csv]: https://github.com/fuqiuai/lihang_algorithms/blob/master/data/train.csv
[maxEnt_maxEnt.py]: https://github.com/fuqiuai/lihang_algorithms/blob/master/maxEnt/maxEnt.py
[maxEnt_result.png]: /images/20220530/e00a49358ee44d20b9e44374f35ac744.png