《统计学习方法》 朴素贝叶斯 极大使然估计 Python实现

代码可在Github上下载:(https://github.com/FlameCharmander/MachineLearning)

前言

使用贝叶斯法来进行分类，我们的模型就是\(f(x)=\mathop {\arg \max }\limits_{{c_k}} P\left( {Y = {C_k}|X = x} \right)\)（后验概率P(隐变量|观测变量)。

算法理论

先介绍贝叶斯法用于模型的分类。

\[\begin{array}{c}
P\left( {Y = {c_k}|X = x} \right) = \frac{{P\left( {X = x|Y = {c_k}} \right)P\left( {Y = {c_k}} \right)}}{{P\left( {X = x} \right)}}\\
 = \frac{{P\left( {X = x|Y = {c_k}} \right)P\left( {Y = {c_k}} \right)}}{{\sum\nolimits_k {P\left( {X = x|y = {c_k}} \right)P\left( {Y = {c_k}} \right)} }}
\end{array}\]

要求后验概率，\({P\left( {X = x} \right)}\)可以用全概率公式即可，所以只需要求出先验\(P\left( {Y = {c_k}} \right)\)以及条件概率\(P\left( {X = x|Y = {c_k}} \right) = P\left( {{X^{(1)}} = {x^{\left( 1 \right)}},...,{X^{\left( n \right)}} = {x^{\left( n \right)}}|Y = {c_k}} \right)\)即可。其中先验概率计算简单，通过极大似然估计，把频数除以总的次数就可以得到，难的是这个条件概率，一般我们可以通过链式求导得到。

\[\begin{array}{l}
P\left( {X = x|Y = {c_k}} \right) = P\left( {{X^{\left( 1 \right)}} = {x^{\left( 1 \right)}},...,{X^{\left( n \right)}} = {x^{\left( n \right)}}|Y = {c_k}} \right)\\
 = P\left( {{X^{\left( 1 \right)}} = {x^{\left( 1 \right)}}|Y = {c_k}} \right) \cdot P\left( {{X^{\left( 2 \right)}} = {x^{\left( 2 \right)}}|{X^{\left( 1 \right)}} = {x^{\left( 1 \right)}},Y = {c_k}} \right)... \cdot P\left( {{X^{\left( n \right)}} = {x^{\left( n \right)}}|{X^{\left( 1 \right)}} = {x^{\left( 1 \right)}},...,{X^{\left( n \right)}} = {x^{\left( n \right)}},Y = {c_k}} \right)
\end{array}\]

其中，\(P\left( {{X^{\left( 1 \right)}} = {x^{\left( 1 \right)}}|Y = {c_k}} \right)\)，\(P\left( {{X^{\left( 2 \right)}} = {x^{\left( 2 \right)}}|{X^{\left( 1 \right)}} = {x^{\left( 1 \right)}},Y = {c_k}} \right)\)这种叫做参数。共有\(\prod\limits_{j = 1}^n {{S_j}} \)种参数的取法，乘以各个类别的数量，那么其拥有\(K\prod\limits_{j = 1}^n {{S_j}} \)参数，这数量级之大以致难以计算。

这里举个例子，比如X有两个特征，每个特征2种取值，并有两个类别。

$$P\left( {{X^{\left( 1 \right)}} = 0,{X^{\left( 2 \right)}} = 0|Y = 0} \right)...P\left( {{X^{\left( 1 \right)}} = 0,{X^{\left( 2 \right)}} = 1|Y = 0} \right)...P\left( {{X^{\left( 1 \right)}} = 1,{X^{\left( 2 \right)}} = 1|Y = 1} \right)$$

共有8（\(K\prod\limits_{j = 1}^n {{S_j}} \)）种取值，所以这里一个一个计算是不可能的，因此也有了朴素贝叶斯的出现。

朴素贝叶斯主要是对条件概率做了条件独立性的假设，可曾还记得两个随机变量是独立的，\(P\left( {XY} \right) = P\left( X \right)P\left( Y \right)\)。这里类似，利用条件独立性假设可将公式改写为：

\[\begin{array}{c}
P\left( {X = x|Y = {c_k}} \right) = P\left( {{X^{(1)}} = {x^{\left( 1 \right)}},...,{X^{\left( n \right)}} = {x^{\left( n \right)}}|Y = {c_k}} \right)\\
 = \prod\limits_{{\rm{j}} = 1}^n {P\left( {{X^{\left( j \right)}} = {x^{\left( j \right)}}|Y = {c_k}} \right)} 
\end{array}\]

\({P\left( {{X^{\left( j \right)}} = {x^{\left( j \right)}}|Y = {c_k}} \right)}\)这个意思就是在指定的类别\({c_k}\)(比如类别-1或者类别1)下，\({X^{(i)}}\)（比如某个样本的第一维特征）这个维度的特征等于取值\(x^{(i)}\)（比如这个维度特征取值可以有"M"，"S"，"L"）的概率。这样就可以通过极大似然估计非常简单地求出，其次通过这个公式，可求出联合概率分布，可以得到后验概率。

分类器一个直观理解就是在通过以上的后验概率得到每个类别的概率，并输出最高类别的概率为分类结果。

$$y = f\left( x \right) = \mathop {\arg \max }\limits_{{c_k}} \left( {{{P\left( {X = x|Y = {c_k}} \right)P\left( {Y = {c_k}} \right)} \over {\sum\nolimits_k {P\left( {X = x|y = {c_k}} \right)P\left( {Y = {c_k}} \right)} }}} \right)$$

可以看到以上的公式分母一样，所以只要求分子就行了。

\[y = f\left( x \right) = \mathop {\arg \max }\limits_{{c_k}} \left( {P\left( {X = x|Y = {c_k}} \right)P\left( {Y = {c_k}} \right)} \right)\]

这样我们也就只需要求得\(P\left( {X = x|Y = {c_k}} \right)\)，\(P\left( {Y = {c_k}} \right)\)。

\[P\left( {Y = {c_k}} \right){\rm{ = }}\frac{{\sum\limits_{i = 1}^N {I\left( {{y_i} = {c_k}} \right)} }}{N}\]

其中\(I\)是指示函数，当括号的值相同时为1，不同时为0，意思就是把每个类别出现的概率。

\[P\left( {{X^{\left( j \right)}} = {x_{jl}}|Y = {c_k}} \right) = \frac{{\sum\limits_{i = 1}^N {I\left( {{x_i}^j = {a_{jl}},{y_i} = {c_k}} \right)} }}{{\sum\limits_{i = 1}^N {I\left( {{y_i} = {c_k}} \right)} }}\]

 

算法实现

实现的是P50页的例4.1，分类器公式如上。

# coding:utf-8
import numpy as np

class Bayes:
    def create_vocab_list(self, dataSet): #创建词汇表 create vocab list, like[1, 2, 3, 'S', 'M', 'L']
        vocab_set = set()
        for document in dataSet:
            vocab_set = vocab_set | set(document)
        return list(vocab_set)

    def set_of_word2vec(self, vocab_list, input_set):   #词条向量 return feature vector, like a feature vector is [1, 'S']，if vocab_list is [1, 2, 3, 'S', 'M', 'L']，return vocab vector [1, 0, 0, 1, 0, 0]
        vocab_vec = [0] * len(vocab_list)    #vocablist大小的零向量 zero vector
        for word in input_set:   #遍历输入样本的每个特征   iterating every feature
            if word in vocab_list:
                vocab_vec[vocab_list.index(word)] = 1    #如果发现有匹配的值就设置为1
        return vocab_vec

    def train(self, dataSet, labels):   #训练样本 train
        self._vocab_list = self.create_vocab_list(dataSet)  #创建特征词汇表 create vocab list
        train_matrix = []    #多条词条向量的矩阵（一个词条向量代表着一个样本在词条中出现的次数） matrix consists of vocab vector
        for line in dataSet:    #将每个训练样本转换为词条向量 feature vector to vocab vector
            train_matrix.append(self.set_of_word2vec(self.vocab_list, line))
        n = len(self.vocab_list) #词条的特征数   feature num
        negative_feature_num = np.zeros(n)    #在类别为-1时，出现特征的次数向量（n1 means negative 1）,the vector of counting num of every feature when label equal -1
        positve_feature_num = np.zeros(n)    #在类别为1时，出现特征的次数向量（）
        negative_num = 0    #标签中出现-1的次数 counting the number of negative label
        positive_num = 0
        for i in range(len(train_matrix)):
            if labels[i] == 1:
                positive_num += 1
                positve_feature_num += train_matrix[i]
            else:
                negative_feature_num += train_matrix[i] #与词条向量相加
                negative_num += 1
        self._positive_vec = positve_feature_num / positive_num   #类别为1的各个随机向量（特征）的概率分布    the probability of feture num
        self._negative_vec = negative_feature_num / negative_num
        self._p_positive = positive_num / float(len(labels)) #p(y=1)的概率 the probability of positive label
        # return self._positive_vec, self._negative_vec, self._p_positive

    def predict(self, input_data):   #预测函数
        input_vec = self.set_of_word2vec(self.vocab_list, input_data)#测试样本的词条向量
        # np.multiply(self.p1Vect ,inputVec)
        p_positive = self.p_positive   #按照公式需要乘以p(y=1)的值，我们就以此为初始值
        p_negative = (1 - self.p_positive)
        for num in np.multiply(self.positive_vec ,input_vec):  #概率分布和词条向量进行相乘，得出p(x=xi|y=1)的概率，然后相乘
            if (num > 0):
                p_positive *= num
        for num in np.multiply(self.negative_vec ,input_vec):
            if (num > 0):
                p_negative *= num
        print(p_positive, p_negative)
        if (p_positive > p_negative):  #相比，谁大就倾向谁 up to max probability
            return 1
        else:
            return -1

    @property
    def vocab_list(self):
        return self._vocab_list

    @property
    def positive_vec(self):
        return self._positive_vec

    @property
    def negative_vec(self):
        return self._negative_vec

    @property
    def p_positive(self):
        return self._p_positive

Line1~2 导入numpy包和编码集

Line5~9 首先需要创建一个词条集，就是一个包含着样本所有特征取值的词汇表。这里的词条是[1, 2, 3, 'S', 'M', 'L']。为了得到词条，需要对每一个样本进行取特征，但是词条里是不重复的，这里用了set函数（元素不重复的集合）来保存词条。

Line11~16 将特征向量转为词条向量，比如特征向量为[1, 'S']，那么词条是[1, 2, 3, 'S', 'M', 'L']，那么词条向量为[1, 0, 0, 1, 0, 0]。

Line19~22 将所有的样本都转为词条矩阵，词条矩阵由多个向量组成，比如有2个样本[1, 'S']，[2, 'S']，那么矩阵为[[1, 0, 0, 1, 0, 0],[0, 1, 0, 1, 0, 0]]

Line24~25 分别对应正样本（类别为1）和负样本（类别为-1）时，这个是一个用zeros初始化的向量，我们这个例子中的特征数是6个，那么negative_feature_num=[0,0,0,0,0,0]，对应着每个特征出现的次数（做什么用请看下文）。

Line25~26 正样本出现的次数 ，负样本出现的次数。

Line28~34 接下来遍历trainMatrix，判断类别为1，还是-1。因为我们需要求出以及p(x=xi|y=1)和p(x=xi|y=-1)的概率分布，当类别为-1时，negative_feature_num跟词条向量相加，就得出-1类别时，这个样本的特征出现的个数。比如-1类别的样本[1, 'M']对应[1, 0, 0, 1, 0, 0]与negative_feature_num[0, 0, 0, 0, 0, 0]相加=[1, 0, 0, 1, 0, 0],当第二个类别为1的样本时[1, 'S']对应[1, 1, 0, 0, 0, 0]，相加=[2, 1, 0, 1, 0, 0]，那么用这个向量除以类别出现的次数，就可以得到概率分布了。

Line35~37 \(P\left( {{X^{\left( 1 \right)}}{\rm{ = }}1{\rm{|Y = 1}}} \right)P\left( {{X^{\left( 2 \right)}}{\rm{ = }}1{\rm{|Y = 1}}} \right)\)

Line40~55 最后我们求出p(y=1)的概率之后，就可以，这里为什么不求出p(y=-1)的概率，是因为就2个类别，p(y=-1) = 1-p(y=1)。好了，我们有了概率分布和类别的概率之后就可以用这个来进行预测了。我们要判别的样本是[2, "S"]。我们有了概率分布，再跟词条向量点乘，就得出了后面p(x=xi|y=1)的概率了。按照公式，我们需要计算2个类别的概率，然后判断谁大就属于谁。

以上就是朴素贝叶斯的极大似然估计的python实现。

我们也提供了一份训练数据在完整代码里，有需要的去github下载，如果方便的话，麻烦点个赞，谢谢。

下载地址：贝叶斯极大似然估计代码