本文主要从下面两个方面展开
我的理解
朴素贝叶斯的本质就是在假设输入特征相互独立的条件下(这个假设是为了方便计算似然函数
P
(
X
∣
Y
)
P(X|Y)
P(X∣Y),因为这里的特征X可能有很多属性,而且取值也很多,先验概率还是后验概率的区分主要是先搞清楚研究的是什么,比如说这里我们想要知道的是类别Y,那么先验概率就是P(Y),后验概率就是P(Y|X)),利用后验概率对样本进行分类的学习方法,后验概率由贝叶斯理论可表示为
P
(
Y
∣
X
)
=
P
(
Y
)
P
(
X
∣
Y
)
P
(
X
)
P(Y|X) = \frac{P(Y)P(X|Y)}{P(X)}
P(Y∣X)=P(X)P(Y)P(X∣Y)
对不同的y求出相应的后验概率,最终输出的y值为后验概率最大的y。在实际问题中,对于同一个输入
X
,
P
(
X
)
X,P(X)
X,P(X)的值是固定的,因此可以将问题转为为
arg max
y
P
(
Y
)
P
(
X
∣
Y
)
\argmax_yP(Y)P(X|Y)
yargmaxP(Y)P(X∣Y),因为X的特征相互独立,所以
P
(
X
∣
Y
)
=
∏
i
=
1
n
P
(
X
(
i
)
∣
Y
)
P(X|Y)=\displaystyle\prod_{i=1}^{n}P(X^{(i)}|Y)
P(X∣Y)=i=1∏nP(X(i)∣Y),n为输入特征维度,先验概率
P
(
Y
)
和
条
件
概
率
P
(
X
(
i
)
∣
Y
)
P(Y)和条件概率P(X^{(i)}|Y)
P(Y)和条件概率P(X(i)∣Y)可以用最大似然估计求得,即可以用样本的频率代替
具体计算方法如下:
首先对公式中的变量进行一些说明,
N
N
N表示样本总数,
n
n
n表示的输入特征X的维度,即特征的个数,
c
k
c_k
ck表示
Y
Y
Y的第k取个值,那么根据最大似然估计可以得到
P
(
Y
=
c
k
)
=
I
(
∑
i
=
1
N
y
i
=
c
k
)
N
P(Y=c_k) = \frac{I(\displaystyle\sum_{i=1}^{N}y_i=c_k)}{N}
P(Y=ck)=NI(i=1∑Nyi=ck)
P
(
X
(
j
)
=
a
j
l
∣
Y
=
c
k
)
=
∑
i
=
1
N
I
(
x
i
(
j
)
=
a
j
l
,
y
i
=
c
k
)
∑
i
=
1
N
I
(
y
i
=
c
k
)
P(X^{(j)}=a_{jl}|Y=c_k)=\frac{\displaystyle\sum_{i=1}^NI(x_i^{(j)}=a_{jl}, y_i=c_k)}{\displaystyle\sum_{i=1}^NI(y_i=c_k)}
P(X(j)=ajl∣Y=ck)=i=1∑NI(yi=ck)i=1∑NI(xi(j)=ajl,yi=ck)
I
I
I为指示函数,
a
j
l
a_{jl}
ajl是特征
X
(
j
)
X^{(j)}
X(j)的第
l
l
l个取值
直接用最大似然估计可能会导致计算条件概率是分母为零的情况,这个时候引入贝叶斯估计,得到先验概率和条件概率的计算公式如下:
P
λ
(
Y
=
c
k
)
=
I
(
∑
i
=
1
N
y
i
=
c
k
)
+
λ
N
+
K
λ
P_\lambda(Y=c_k) = \frac{I(\displaystyle\sum_{i=1}^{N}y_i=c_k)+\lambda}{N+K\lambda}
Pλ(Y=ck)=N+KλI(i=1∑Nyi=ck)+λ
P
λ
(
X
(
j
)
=
a
j
l
∣
Y
=
c
k
)
=
∑
i
=
1
N
I
(
x
i
(
j
)
=
a
j
l
,
y
i
=
c
k
)
+
λ
∑
i
=
1
N
I
(
y
i
=
c
k
)
+
S
j
λ
P_\lambda(X^{(j)}=a_{jl}|Y=c_k)=\frac{\displaystyle\sum_{i=1}^NI(x_i^{(j)}=a_{jl}, y_i=c_k)+\lambda}{\displaystyle\sum_{i=1}^NI(y_i=c_k)+S_j\lambda}
Pλ(X(j)=ajl∣Y=ck)=i=1∑NI(yi=ck)+Sjλi=1∑NI(xi(j)=ajl,yi=ck)+λ
其中K表示Y的种类数,
S
j
S_j
Sj是特征
X
(
j
)
X^{(j)}
X(j)可以取得值的总数
代码实现
下面以手写数字识别为例,将图片的每个像素点作为一个特征,对于28*28的图像一共有784个特征,对像素值进行二值化,大于等于128的为1,小于128的为1,这样手写数字识别问题就转为一个十分类的问题,每个样本一共有784个特征值,每个特征值有两个取值(0和1),具体代码如下,先验概率和条件概率的计算方式都选择贝叶斯估计的方法:
"""
Mnist classification by NaiveBayes
"""
import numpy as np
import time
class NaiveBayes(object):
"""docstring for NaiveBayes"""
def __init__(self, featurenum, classnum):
super(NaiveBayes, self).__init__()
self.featurenum = featurenum
self.classnum = classnum
self.py = np.zeros(classnum)
self.py_x = np.zeros((classnum, featurenum, 2))
def train(self, data, label):
"""
data: np.array (n, m)# (samplenum, featurenum)
label:(n) belong to [0, classnum)
对先验概率和条件概率都取了对数,这样之后的连乘就可以改为加法
这里计算先验概率和条件概率的方法用的都是书中说的贝叶斯估计
"""
n, m = data.shape
#calculate priori probability p(y),用每一类的样本数除以总的样本数
for i in range(self.classnum):
self.py[i] = np.log((sum((label==i).astype('int')) + 1)/ (n + self.classnum))
#计算条件概率p(x|y),循环遍历每一个类别,对每一个类别的data,在样本数量方向求和
#因为每个特征都只有0和1两个取值,求和就得到了特征取1的样本数
for i in range(self.classnum):
self.py_x[i,:,1] = (np.sum(data[label==i], axis=0) + 1) / (sum((label==i).astype('int')) + 2)
self.py_x[:,:,0] = 1 - self.py_x[:,:,1]#用1减去类别为1的概率就得到了零的概率
self.py_x = np.log(self.py_x)
def predict(self, data):
"""
data: np.array (n, m)
"""
n, m = data.shape
prob = np.zeros((n, self.classnum))
featureindex = list(range(self.featurenum))
# py_x = np.tile(self.py_x, (n, 1, 1, 1))
for j in range(n):
#减少了一层循环,速度就提高了大约25倍,原来要50s,现在测试只需要2秒
prob[j] = np.sum(self.py_x[:,featureindex, data[j]], axis=-1) + self.py
# for i in range(self.classnum):
# # prob[j][i] = sum(self.py_x[i,list(range(self.featurenum)),data[j]]) + self.py[i]
# # prob[j][i] = sum(self.py_x[i,featureindex,data[j]]) + self.py[i] #和上面一行相比时间就快了两秒,看来主要是循环造成的
index = prob.argmax(axis=-1) #返回每一行最大值的索引
return index
def getdata(path):
data, label = [], []
with open(path, 'r') as f:
lines = f.readlines()
for line in lines:
items = line.strip().split(',')
label.append(int(items[0]))
data.append([int(int(num) >= 128) for num in items[1:]])
return np.array(data), np.array(label)
if __name__ == '__main__':
traindata, trainlabel = getdata('../mnist_train.csv')
testdata, testlabel = getdata('../mnist_test.csv')
bayes = NaiveBayes(len(traindata[0]), 10)
print('trian start')
start = time.time() #返回1970纪元后经过的浮点秒数
bayes.train(traindata, trainlabel)
print('train model spend {}s'.format(time.time() - start)) #0.93s
print('test start')
start = time.time()
predict = bayes.predict(testdata)
print('test spend {}s'.format(time.time() - start)) #2.59s
accu = sum((predict == testlabel).astype('int')) / len(testlabel)
print('test accuracy is {}%'.format(accu * 100)) #84.3%