朴素贝叶斯和半朴素贝叶斯（AODE）分类器Python实现

一、概述

机器学习最后一次实验，要求实现朴素贝叶斯和AODE的半朴素贝叶斯分类器。由于老师说可以调用现成的相关机器学习的库，所以我一开始在做朴素贝叶斯分类器的时候，直接调用了sklearn库，很方便，可是问题来了，在做AODE半朴素贝叶斯分类器的时候，并没有找到集成好的方法。所以就想着自己把半朴素贝叶斯分类器实现了，朴素贝叶斯分类就直接调用库算了。可是让人头大的是，上来就直接实现AODE分类器还是不太科学，还得从基本的贝叶斯原理到朴素贝叶斯开始，所以我又从头看，主要看的这篇博客 西瓜书读书笔记——第七章：贝叶斯分类器，写的非常好，看完之后就基本弄懂了。我感觉实现起来，朴素贝叶斯和半朴素贝叶斯有很多相似之处，索性就先把朴素贝叶斯实现了，正好也能加深一下理解，然后再实现AODE半朴素贝叶斯分类器就会容易些了。

二、贝叶斯分类器

2.1 朴素贝叶斯分类器

（1）贝叶斯分类基本思想简单解释

首先，贝叶斯分类的思想很简单，假设数据一共有$n$种类别，即 L a b e l = { L 1 , L 2 , ⋯   , L n } Label=\{L_1,L_2,\cdots,L_n\} Label={L1​,L2​,⋯,Ln​}，给定一个样本数据 x = { x 1 , x 2 , ⋯   , x m } x=\{x_1,x_2,\cdots,x_m\} x={x1​,x2​,⋯,xm​}，注意，这里的$x$是一个样本数据，$x_1,x_2,\cdots ,x_m$表示这个数据有$m$维特征，当知道了样本数据$x$，根据$x$计算出来这个数据属于每一种类别的概率 P = { P L 1 , P L 2 , ⋯   , P L n } P=\{P_{L_1},P_{L_2},\cdots,P_{L_n}\} P={PL1​​,PL2​​,⋯,PLn​​}，最后将这个数据$x$划分为最大概率对应的类别。比如，如果 a r g m a x { P L 1 , P L 2 , ⋯   , P L n } = L 2 argmax\{P_{L_1},P_{L_2},\cdots,P_{L_n}\}=L_2 argmax{PL1​​,PL2​​,⋯,PLn​​}=L2​，那么$x$就被划分为$L_2$。

（2）贝叶斯分类基本原理

其次，贝叶斯分类的实现也很简单。现在已经知道了基本原理，设$x$所属的类别为$c$，$x_i$代表样本$x$的第$i$个属性值（第$i$个维度的值），那么上面所说的要求样本$x$属于每种类别的概率就记为$P(c|x)$，根据贝叶斯模型和极大似然估计原理，那么就有：
$$P(c|x)=\frac{P(x|c)P(c)}{P(x)}=\frac{P(c)}{P(x)}\prod _{i=1}^{m}P(x_i|c)$$
其中$P(x)={\prod }_{i=1}^{m}P(x_i)$，对于数据$x$来说，计算对应的每一种类别的概率时，$P(x)$是相同的，所以为了计算方便，可以省略掉，记为：
$$P(c|x)\propto P(c)\prod _{i=1}^{m}P(x_i|c)$$（注：在实现时，为了防止概率连乘导致趋近于0，将${\prod }_{i=1}^{m}P(x_i|c)$取对数变成连加）所以最终计算$x$所对应的类别就有：
$$L(x)=argmaxP(c)\prod _{i=1}^{m}P(x_i|c)$$这里c是所要求的值。根据这个公式就知道，要求$x$所属的类别，只要求出$P(c)$和$P(x_i|c)$就行了。
设整个样本数据集为$D$，当$|D|$足够大时（即样本数量足够多），就可以利用频率估计概率（大数定律）来计算出先验概率$P(c)$，
$$P(c)=\frac{|D_c|}{|D|}$$$|D|$代表所有数据的数量，$|D_c|$表示类别为$c$的样本的数量。现在$P(c)$很容易的求出来了，然后就是$P(x_i|c)$了。

然而，对于样本属性$x_i$来说，其可分为连续性样本属性和离散型样本属性。先理解一下什么是“离散型样本属性”，那么“连续型”就比较容易理解了。

“离散型样本属性”：比如，对于西瓜$A$来说，$A$就是整个西瓜家族里面的一个样本，那么$A$的属性就会有 { 外 皮 颜 色 ， 敲 击 声 音 ， 触 摸 手 感 ， ⋯   } \{外皮颜色，敲击声音，触摸手感，\cdots\} {外皮颜色，敲击声音，触摸手感，⋯}，而对于“外皮颜色”这个属性来说，它的取值可能有 { 黄 色 ， 绿 色 ， 青 绿 色 ， ⋯   } \{黄色，绿色，青绿色，\cdots\} {黄色，绿色，青绿色，⋯}，这个属性的取值是离散的，有限的，那么这个属性就是“离散型”属性了，另外两个属性也是一样。对于离散型样本属性的条件概率计算公式也是根据频率估计概率得到：
$$P(x_i|c)=\frac{|D_{c,x_i}|}{|D_c|}$$$|D_{c,x_i}|$为类别$c$中的所有数据的第$i$个属性值为$x_i$的样本的数量。

“连续型样本属性”：还是对西瓜$A$来说，$A$除了上面提到的那些属性之外，还有 { 含 糖 量 ， 水 分 } \{含糖量，水分\} {含糖量，水分}这些属性，这两个属性如果以定量方法（精确地测量数值）来表示，比如含糖量为$45.678\%$，这样就叫做“连续型”属性了，但是如果以定性方法来表示，比如将含糖量划分为 { 低 ， 中 ， 高 } \{低，中，高\} {低，中，高}三个等级，那么就是“离散型”属性了。对于连续型样本属性的条件概率计算公式使用高斯核密度估计(应该是这么叫吧，希望数理统计没白学)得到：
$$P(x_i|c)=\frac{1}{\sqrt{2\pi }{\sigma }_{c,i}}exp\left(-\frac{(x_i-{\mu }_{c,i}{)}^2}{2{\sigma }_{c,i}^2}\right)$$其中${\mu }_{c,i}$和${\sigma }_{c,i}$分别为第$c$类样本在第$i$个属性取值的均值和标准差。

为了避免数据集不充分而导致估计概率为0的情况，需要在实现过程中引入拉普拉斯修正（具体拉普拉斯修正方法也很简单），但我为了图方便，直接在先验概率和条件概率（针对离散型属性）的分子和分母都 +1 了（实际测试中，并没有发现有太大的差别）。

2.2 AODE半朴素贝叶斯分类器

首先，要知道什么是AODE(Averaged One-Dependent Estimator)，就要先了解什么是SPODE(Super-Parent One-Dependent Estimator)。在上面所讲的朴素贝叶斯分类都是基于样本的所有属性是相互独立的，那么如果某些属性存在依赖关系怎么办呢。SPODE就是假设样本的每个属性都依赖于某一个属性，这个属性就叫做“超父”。比如，将$x_1$属性设置为“超父”，那么计算$x$属于第$c$类数据的概率公式为：
$$P(c|x)\propto P(c)\prod _{i=1}^{m}P(x_i|c,x_1)$$而AODE就是在此基础上，将样本的属性依次作为超父来计算概率,最后求和，同时，假设类别$c$也依赖于样本的属性，那么计算$x$属于类别$c$的概率公式为：
$$P(c|x)\propto \sum _{i=1}^{m}P(c,x_i)\prod _{j=1}^{m}P(x_j|c,x_i)$$同样的，利用频率估计概率的方法可以得到：
$$P(x_j|c,x_i)=\frac{|D_{c,x_i,x_j}|}{|D_{c,x_i}|}$$与朴素贝叶斯分类相同，实现时也要引入拉普拉斯修正，我还是分子分母都 +1 了。

通过观察$x$属于类别$c$的概率公式就能看出，这个方法的时间复杂度很高，实际上，对于每一个样本，计算出分别属于哪一个类的概率一共有三层循环（不知是否有优化方法，但是时间复杂度几乎不可能降低了）。实现之后，这个方法并没有进行实验结果的验证，因为时间代价太大，所以我猜测，这也是为什么sklearn中有朴素贝叶斯却没有AODE半朴素贝叶斯方法的原因了。

三、实验结果与代码

3.1 实验结果

对于朴素贝叶斯分类器的实现，我只实现了离散型样本属性的分类（连续性属性也比较容易，只要把条件概率函数换成高斯核即可），使用了MNIST、Yale和COIL20这三个数据对其进行了实验，使用ACC评价标准对其进行评价。（注：由于AODE时间成本过高，不太适合属性维度较多数据，我没有进行数据实验，所以我也不知道对不对，不过按照朴素贝叶斯的思路来应该也不会错吧。。。）

方法\数据	MNIST	Yale	COIL20
McQueen_NBC	0.95	1	1
GaussianNB	0.55	0.8	1

这里我就没有将数据集划分成测试集和训练集了，这三组数据测试集是在训练集中每个类别数据分别抽取前0.5%，10%，1%得到的。McQueen_NBC方法是我自己实现的朴素贝叶斯分类器（实际上是多项式贝叶斯），适用于离散属性样本，而GaussianNB是调用的sklearn库的方法，这个方法是基于连续型属性的，由结果可以看出，在MNIST和Yale这两个数据上，连续型准确率不如离散型（因为这两个数据是离散型样本数据）。

3.2 完整代码

datadvi.py

from scipy.io import loadmat
import numpy as np

def divdata():
    filename = 'C:/Users/ALIENWARE/Documents/作业/机器学习/datasets/' + input("input name of data file: ")
    data = loadmat(filename)
#    print(data['X'])


    if filename == 'C:/Users/ALIENWARE/Documents/作业/机器学习/datasets/COIL20.mat':
        dataX = data['fea']
        dataY = data['gnd'][0]

    else:
        dataX = data['X']
        dataY = data['Y'].T[0]

    print(len(dataX[0]))

    divideornot = input("divide data or not?(Yes/No): ")
    if divideornot == 'Yes':
        dataX_train = []
        dataX_predict = []
        dataY_train = []
        dataY_predict = []
        num_Y = np.unique(dataY).astype(int)
        for i in range(len(num_Y)):
            temp = dataY == num_Y[i]
            temp.astype(float)
            num_Y[i] = np.sum(temp)
            flag = 0
            for j in range(len(dataY)):
                if temp[j] == 1:
                    if flag < int(round(0.9 * num_Y[i])):
                        dataX_train.append(dataX[j])
                        dataY_train.append(dataY[j])
                        flag += 1
                    else:
                        dataX_predict.append(dataX[j])
                        dataY_predict.append(dataY[j])

        dataX_train = np.array(dataX_train)
        dataX_predict = np.array(dataX_predict)
        dataY_train = np.array(dataY_train)
        dataY_predict = np.array(dataY_predict)
        return dataX_train,dataX_predict,dataY_train,dataY_predict
    else:
        return dataX,dataX,dataY,dataY

def decreaseData(dataX,dataY):
    dataX_train = []
    dataY_train = []
    num_Y = np.unique(dataY).astype(int)
    print("this data has {} samples".format(len(dataX)))
    ratio = float(input("input the ratio you want to decrease: "))
    for i in range(len(num_Y)):
        temp = dataY == num_Y[i]
        temp.astype(float)
        num_Y[i] = np.sum(temp)
        flag = 0
        for j in range(len(dataY)):
            if temp[j] == 1:
                if flag < round(ratio * num_Y[i]):
                    dataX_train.append(dataX[j])
                    dataY_train.append(dataY[j])
                    flag += 1

    dataX_train = np.array(dataX_train)
    dataY_train = np.array(dataY_train)
    print(dataX_train)

    return dataX_train,dataY_train

Acc.py

import numpy as np

def acc(L1, L2):
    sum = np.sum(L1[:]==L2[:])
    return sum/len(L2)

NBC.py

import math
import numpy as np
import datadvi
import Acc

#加载数据
def loadData(filename):
    return datadvi.divdata()

#按标签类别生成不同标签样本组成的集合，返回值为每种类别样本的索引
def divSamples(dataY):
    label = np.unique(dataY)
    D = []
    for i in label:
        D.append(np.argwhere(dataY==i).T[0])
    return np.array(D)

# 计算第c类样本在第i个属性上取值的均值和标准差,smaple_cIndx是第c类样本的索引(用于连续型属性，此次未用到)
def calcMuSig(sample_cIndx,i,D):
    mu =  np.average(D[sample_cIndx][:,i])
    sigma = np.std(D[sample_cIndx][:,i])
    return mu,sigma

#计算类先验概率P(c),
def beforeProb(sample_cIndx,D):
    return float(len(sample_cIndx)+1)/(D.shape[0]+1)

#计算离散型条件概率P(xi|c)
def condProb_disp(i,xi,sample_cIndx,D):
    numerator = np.sum(D[sample_cIndx][:,i]==xi)+1
    denominator = len(sample_cIndx)+1
    return float(numerator)/denominator

#计算连续型条件概率P(xi|c)
def condProb_cont(i,xi,sample_cIndx,D):
    mu,sigma = calcMuSig(sample_cIndx,i,D)
    prob = 1/(math.sqrt(2*3.14)*sigma)*math.exp(-float((xi-mu)*(xi-mu))/(2*sigma*sigma))
    return prob

#计算类后验概率P(c|x)
def afterProb(sample_x,c,dataX,dataY):
    sample_c = divSamples(dataY)
    p = beforeProb(sample_c[c],dataX)
    #p = beforeProb(sample_c[c],dataX)
    p1 = 0
    for i in range(len(sample_x)):
        p1 += math.log10(condProb_disp(i,sample_x[i],sample_c[c],dataX))
        #p1 *= condProb_cont(i,sample_x[i],sample_c[c],dataX) #会下溢
    return p*p1

#计算最大概率对应的类
def argMaxProb_c(sample_x,dataX,dataY):
    label = np.unique(dataY)
    argProb1 = []
    for c in label:
        temp_prob = afterProb(sample_x,c-1,dataX,dataY)
        argProb1.append(temp_prob)
    argProb = np.array(argProb1)
    return label[np.argmax(argProb)]


#将所有数据分类
def bayesClassifier(dataPredict,dataX,dataY):
    pred = []
    for sample_x in dataPredict:
        pred.append(argMaxProb_c(sample_x,dataX,dataY))
        print(len(pred))
    return pred

dataX_train, dataX_predict, dataY_train, dataY_predict = datadvi.divdata()
dataX_predict,dataY_predict = datadvi.decreaseData(dataX_predict,dataY_predict)
print(len(dataX_predict))
sample_c = divSamples(dataY_train)
pred = bayesClassifier(dataX_predict,dataX_train,dataY_train)
print(pred)
print(Acc.acc(pred,dataY_predict))

AODE.py

import math
import numpy as np
import datadvi
import Acc
import RBFNN

#加载数据
def loadData(filename):
    return datadvi.divdata()

#按标签类别生成不同标签样本组成的集合，返回值为每种类别样本的索引
def divSamples(dataY):
    label = np.unique(dataY)
    D = []
    for i in label:
        D.append(np.argwhere(dataY==i).T[0])
    return np.array(D)


#计算先验概率P(c,xi)
def beforeProb(sample_cIndx,i,xi,D):
    numerator = len(np.argwhere(D[sample_cIndx][:,i]==xi).T[0])+1 #索引改变了，但是数量没变
    denominator = D.shape[0]+1 #此处1需被替换为N*Ni
    return float(numerator)/denominator


#计算条件概率P(xj|c,xi)
def condProb(i,xi,j,xj,sample_cIndx,D):
    D_c = D[sample_cIndx]
    D_c_xi = D_c[np.argwhere(D_c[:,i]==xi).T[0]]
    D_c_xi_xj = D_c_xi[np.argwhere(D_c_xi[:,j]==xj).T[0]]
    numerator = len(D_c_xi_xj)+1
    denominator = len(D_c_xi)+1
    return float(numerator)/denominator

#计算后验概率P(c|x)
def afterProb(sample_x,c,dataX,dataY):
    sample_c = divSamples(dataY)
    prob = 0
    for i in range(len(sample_x)):
        p1 = 0
        p = beforeProb(sample_c[c],i,sample_x[i],dataX)
        for j in range(len(sample_x)):
            p1 += math.log10(condProb(i,sample_x[i],j, sample_x[j],sample_c[c], dataX))  #防止下溢
        prob += p*p1
    return prob


#计算最大概率对应的类
def argMaxProb_c(sample_x,dataX,dataY):
    label = np.unique(dataY)
    argProb1 = []

    for c in label:
        temp_prob = afterProb(sample_x, c - 1, dataX, dataY)
        argProb1.append(temp_prob)

    argProb = np.array(argProb1)
    return label[np.argmax(argProb)]


#将所有数据分类
def bayesClassifier(dataPredict,dataX,dataY):
    pred = []
    for sample_x in dataPredict:
        label_pred = argMaxProb_c(sample_x,dataX,dataY)
        pred.append(label_pred)
        print(len(pred))
    return pred


dataX_train, dataX_predict, dataY_train, dataY_predict = datadvi.divdata()
dataX_predict,dataY_predict = datadvi.decreaseData(dataX_predict,dataY_predict)


print(len(dataX_predict))

pred = bayesClassifier(dataX_predict,dataX_train,dataY_train)
print(Acc.acc(pred,dataY_predict))

四、总结

机器学习四次实验全部结束了，学到了很多知识，感觉很充实，体会到了Python的强大，还好之前自学了Python，不然真可能现在被逼的去用MATLAB（不是针对MATLAB，我只想说C++天下第一，手动狗头）了。虽然机器学习课时短，加上理论课没有好好听，浪费了很多时间，不过实验好好做了感觉把知识又全都补回来了。接下来又得去做计算机视觉的实验了，我哭了。。。