朴素贝叶斯
朴素贝叶斯法是基于贝叶斯定理和特征条件独立假设的分类方法,由于现实中很难满足『特征条件独立假设』,这也是为什么他被叫做『朴素』贝叶斯的原因。
下面将用朴素贝叶斯发来实现手写数字的识别。
问题描述
我用的数据集是sklearn
中的一个数据集,先导入它
from sklearn.datasets import load_digits
digits = load_digits()
简单介绍一下这个数据集,它包含
1797
1797
1797张
8
×
8
8\times 8
8×8大小的图片,每张图片也就是
64
64
64个特征,我们可以把他标记为:
x
=
(
x
1
,
⋯
,
x
64
)
T
x = (x_1, \cdots, x_{64})^T
x=(x1,⋯,x64)T
每张图片都有对应一个输出
y
y
y,范围为
0
0
0到
9
9
9的整数,也就是图片对应的数字。
为了更加直观地感受,可以拿一部分出来可视化一下
fig = plt.figure(figsize=(6, 6))
fig.subplots_adjust(left=0, right=1, bottom=0, top=1, hspace=0.05, wspace=0.05)
for i in range(64):
ax = fig.add_subplot(8, 8, i + 1, xticks=[], yticks=[])
ax.imshow(digits.images[i], cmap=plt.cm.binary, interpolation='nearest')
ax.text(0, 7, str(digits.target[i]))
plt.show()
效果如下所示:
图片的左下角就是真实的数字,也就是每个数据
x
x
x对应的
y
y
y值。
然后我们要做的就是根据输入 x x x,也就是每个像素点的值,来判断它对应的数字 y p r e d i c t y_{predict} ypredict。
理论部分
对于一个给定的输入
x
=
(
x
1
,
⋯
,
x
n
)
x=(x_1, \cdots, x_n)
x=(x1,⋯,xn),假设他的输出一共有
M
M
M种可能,我们把他记为
y
=
{
y
1
,
⋯
,
y
M
}
y=\{y_1, \cdots, y_M\}
y={y1,⋯,yM},然后我们要做的就是,在输入为
x
x
x情况下,计算输出为
y
y
y中每一个值的概率,我们把他记为
p
m
p_m
pm,
p
m
=
P
r
(
Y
=
y
m
∣
X
=
x
)
=
P
r
(
Y
=
y
m
∣
X
1
=
x
1
,
⋯
,
X
n
=
x
n
)
p_m = {\rm Pr}(Y = y_m | X = x) = {\rm Pr}(Y = y_m | X_1 = x_1, \cdots, X_n = x_n)
pm=Pr(Y=ym∣X=x)=Pr(Y=ym∣X1=x1,⋯,Xn=xn)
最后我们返回的就是其中概率最大的一种可能,也就是
y
p
r
e
d
i
c
t
=
a
r
g
m
a
x
y
m
P
r
(
Y
=
y
m
∣
X
1
=
x
1
,
⋯
,
X
n
=
x
n
)
y_{predict} = {\rm arg max}_{y_m} {\rm Pr}(Y = y_m | X_1 = x_1, \cdots, X_n = x_n)
ypredict=argmaxymPr(Y=ym∣X1=x1,⋯,Xn=xn)
下面要处理的就是一个寻找最大值的问题了,根据贝叶斯定理则有
P
r
(
Y
=
y
m
∣
X
=
x
)
=
P
r
(
Y
=
y
m
)
P
r
(
X
=
x
∣
Y
=
y
m
)
P
r
(
X
=
x
)
{\rm Pr}(Y = y_m | X = x) = \frac{{\rm Pr}(Y = y_m) {\rm Pr}(X=x|Y=y_m)}{{\rm Pr}(X = x)}
Pr(Y=ym∣X=x)=Pr(X=x)Pr(Y=ym)Pr(X=x∣Y=ym)
由于朴素贝叶斯法是基于每个变量相对独立假设的,也就是说我们假设
(
x
1
,
⋯
,
x
n
)
(x_1, \cdots, x_n)
(x1,⋯,xn)的分布是相互独立的,上面的式子可以进一步写为
P
r
(
Y
=
y
m
∣
X
=
x
)
=
P
r
(
Y
=
y
m
)
P
r
(
X
1
=
x
1
∣
Y
=
y
m
)
×
⋯
×
P
r
(
X
n
=
x
n
∣
Y
=
y
m
)
P
r
(
X
=
x
)
=
P
r
(
Y
=
y
m
)
∏
i
P
r
(
X
i
=
x
i
∣
Y
=
y
m
)
P
r
(
X
=
x
)
\begin{aligned} {\rm Pr}(Y = y_m | X = x) &= \frac{{\rm Pr}(Y = y_m) {\rm Pr}(X_1=x_1|Y=y_m)\times \cdots\times{\rm Pr}(X_n=x_n|Y=y_m)}{{\rm Pr}(X = x)}\\ &=\frac{{\rm Pr}(Y = y_m) \prod_i{\rm Pr}(X_i=x_i|Y=y_m)}{{\rm Pr}(X = x)} \end{aligned}
Pr(Y=ym∣X=x)=Pr(X=x)Pr(Y=ym)Pr(X1=x1∣Y=ym)×⋯×Pr(Xn=xn∣Y=ym)=Pr(X=x)Pr(Y=ym)∏iPr(Xi=xi∣Y=ym)
因为每个估计值
p
m
p_m
pm求解式子分母都是
P
r
(
X
=
x
)
{\rm Pr}(X=x)
Pr(X=x),所以我们要处理的问题就变成了
y
p
r
e
d
i
c
t
=
a
r
g
m
a
x
y
m
P
r
(
Y
=
y
m
)
∏
i
P
r
(
X
i
=
x
i
∣
Y
=
y
m
)
y_{predict} = {\rm arg max}_{y_m} {\rm Pr}(Y = y_m) \prod_i{\rm Pr}(X_i=x_i|Y=y_m)
ypredict=argmaxymPr(Y=ym)i∏Pr(Xi=xi∣Y=ym)
为了便于计算,我们对上面式子做
l
o
g
{\rm log}
log处理,使得乘法变为加法,便于计算
y
p
r
e
d
i
c
t
=
a
r
g
m
a
x
y
m
l
o
g
P
r
(
Y
=
y
m
)
+
∑
i
l
o
g
P
r
(
X
i
=
x
i
∣
Y
=
y
m
)
y_{predict} = {\rm arg max}_{y_m} {\rm logPr}(Y = y_m) + \sum_i{\rm log Pr}(X_i=x_i|Y=y_m)
ypredict=argmaxymlogPr(Y=ym)+i∑logPr(Xi=xi∣Y=ym)
最后的一个问题就是如何计算
P
r
(
Y
=
y
m
)
{\rm Pr}(Y = y_m)
Pr(Y=ym)和
P
r
(
X
i
=
x
i
∣
Y
=
y
m
)
{\rm Pr}(X_i=x_i|Y=y_m)
Pr(Xi=xi∣Y=ym)。
其中
P
r
(
Y
=
y
m
)
{\rm Pr}(Y = y_m)
Pr(Y=ym)很好处理,我们只需要从测试集中统计每种
y
m
y_m
ym出现的概率,也就是先验概率。
P
r
(
X
i
=
x
i
∣
Y
=
y
m
)
{\rm Pr}(X_i=x_i|Y=y_m)
Pr(Xi=xi∣Y=ym)相对复杂一些,我们需要获取测试集中输出为
y
m
y_m
ym特征值
x
i
x_i
xi的均值
μ
m
i
\mu_{mi}
μmi和方差
σ
m
i
\sigma_{mi}
σmi,然后通过以下公式计算得到
P
r
(
X
i
=
x
i
∣
Y
=
y
m
)
=
1
2
π
σ
m
i
2
e
−
(
x
i
−
μ
x
i
)
2
2
σ
m
i
2
{\rm Pr}(X_i=x_i|Y=y_m) = \frac{1}{\sqrt{2\pi\sigma^2_{mi}}}e^{-\frac{(x_i - \mu_{xi})^2}{2\sigma_{mi}^2}}
Pr(Xi=xi∣Y=ym)=2πσmi21e−2σmi2(xi−μxi)2
小结
我们先整理一下思路,来看看我们下面具体要做些什么。
- 获取一个测试集,统计测试集中每种输出 y m y_m ym对应特征值 x x x的先验概率 P r ( Y = y m ) {\rm Pr}(Y=y_m) Pr(Y=ym)、均值 μ m \mu_m μm、方差 σ m \sigma_m σm。
- 输入任意一个张图片的像素信息 x x x利用第一步中得到的先验概率 P r ( Y = y m ) {\rm Pr}(Y=y_m) Pr(Y=ym)、均值 μ m \mu_m μm、方差 σ m \sigma_m σm,估计出 y p r e d i c t y_{predict} ypredict并与真实值 y t r u e y_{true} ytrue进行对比,查看效果。
Python实现
划分数据集
这里直接调用sklearn
中的api来实现数据集的划分,
from sklearn.datasets import load_digits
from sklearn.model_selection import train_test_split
digits = load_digits()
x_train, x_test, y_train, y_test = train_test_split(digits.data,
digits.target,
test_size=0.3)
其中x_train
和y_train
是我们用来训练的数据,分别对应于特征值和正确输出,完成模型的学习后,我们利用x_test
和y_test
来评估模型的好坏。
模型训练
根据之前的理论,我们知道,要先用训练集来获取一些样本的均值、方差、先验概率,注意这里不需要什么梯度下降之类的优化。我们可以定义一个class
来便于我们的调用。
class NaiveBayes:
def train(self, X: np.ndarray, y: np.ndarray):
"""
训练函数,获取训练集的均值,方差,先验概率
"""
# 获取数据集大小
n_samples, n_features = X.shape
# 获取数据集类别和类别数
self._classes = np.unique(y)
n_classes = len(self._classes)
# print(self._classes)
# 初始化均值,方差,先验概率
self._mean = np.zeros((n_classes, n_features), dtype=float)
self._var = np.zeros((n_classes, n_features), dtype=float)
self._priors = np.zeros(n_classes, dtype=float)
# 计算均值,方差,先验概率
for c in self._classes:
# 找到输出为c的样本
X_c = X[c == y]
# 获取输出为c的均值,方差,先验概率
self._mean[c, :] = X_c.mean(axis=0)
self._var[c, :] = X_c.var(axis=0)
self._priors[c] = X_c.shape[0] / float(n_samples)
模型预测
获取训练集的均值、方差、先验概率后,我们就可以对测试集做出预测了
def predict(self, X) -> List[int]:
"""
获取测试集的预测
:return: 预测结果
"""
pred = []
for x in X:
# 统计每中输出的可能性
probs = []
for y in range(len(self._classes)):
probs.append(self._cal_prob(x, y))
# 找到概率最大的那个类别
pred.append(np.argmax(probs))
return pred
最重要的部分,便是self._cal_prob(x, y)
的编写
def _cal_prob(self, x, y) -> float:
"""
计算 log p(x_1|y) + ... + log p(x_n|y) + log p(y)
注意:这里的 y 为在 self._classes 中的下标
"""
# log p(y)
p_y = np.log(self._priors[y])
# print(p_y)
# log p(x_1|y) + ... + log p(x_n|y)
mean = self._mean[y]
# print("mean are: \n", mean)
var = self._var[y]
# print("var: \n", var)
p_xy = np.log(np.exp(-(x-mean)**2 / (2*var)) / np.sqrt(2*np.pi*var))
# p_xy 的计算结果可能会出现 nan,我们将它取值改为0,让他不对输出造成影响
p_xy[np.isnan(p_xy)] = 0
return p_y + np.sum(p_xy)
模型评估
我们通过计算我们的模型预测的准确性来评估我们的模型,代码如下:
def accuracy(y_pred, y_true):
return np.sum(y_pred == y_true) / len(y_true)
效果展示
走完整个过程的代码如下,
digits= load_digits()
x_train, x_test, y_train, y_test = train_test_split(digits.data,
digits.target,
test_size=0.3)
bayes = NaiveBayes()
bayes.train(x_train, y_train)
y_pred = bayes.predict(x_test)
print(accuracy(y_pred, y_test))
我运行的输出结果为
0.8962962962962963
我们也可以随机选择部分测试集查看效果,测试代码如下
def show_test_result(nb: NaiveBayes):
fig = plt.figure(figsize=(6, 6)) # figure size in inches
fig.subplots_adjust(left=0, right=1, bottom=0, top=1, hspace=0.05, wspace=0.05)
n = digits.data.shape[0]
for i in range(64):
ax = fig.add_subplot(8, 8, i + 1, xticks=[], yticks=[])
# 随机选择一张图片
idx = random.randint(0, n-1)
ax.imshow(digits.images[idx], cmap=plt.cm.binary, interpolation='nearest')
data = np.asarray(digits.data[idx]).reshape(1, -1)
# 左下角为预测值
ax.text(0, 7, str(nb.predict(data)[0]), color="green")
# 左上角为真实值
ax.text(0, 1, str(digits.target[idx]), color="red")
plt.show()
由于我们假设了每个像素点之间为相对独立的,这与实际不符合,故最终也能发现一些错误。
完整测试代码
#!/usr/bin/python3
import numpy as np
import matplotlib.pyplot as plt
import random
from typing import List, Tuple
from sklearn.decomposition import PCA
from sklearn.datasets import load_digits
from sklearn.model_selection import train_test_split
from sklearn.naive_bayes import GaussianNB
import warnings
warnings.filterwarnings('ignore')
def show_dataset():
fig = plt.figure(figsize=(6, 6)) # figure size in inches
fig.subplots_adjust(left=0, right=1, bottom=0, top=1, hspace=0.05, wspace=0.05)
for i in range(64):
ax = fig.add_subplot(8, 8, i + 1, xticks=[], yticks=[])
ax.imshow(digits.images[i], cmap=plt.cm.binary, interpolation='nearest')
# label the image with the target value
ax.text(0, 7, str(digits.target[i]))
plt.show()
class NaiveBayes:
def train(self, X: np.ndarray, y: np.ndarray):
"""
训练函数
:return: 训练集的均值,方差,先验概率
"""
# 获取数据集大小
n_samples, n_features = X.shape
# print(n_features)
# 获取数据集类别和类别数
self._classes = np.unique(y)
n_classes = len(self._classes)
# print(self._classes)
# 初始化均值,方差,先验概率
self._mean = np.zeros((n_classes, n_features), dtype=float)
self._var = np.zeros((n_classes, n_features), dtype=float)
self._priors = np.zeros(n_classes, dtype=float)
# 计算均值,方差,先验概率
for c in self._classes:
X_c = X[c == y]
self._mean[c, :] = X_c.mean(axis=0)
self._var[c, :] = X_c.var(axis=0)
self._priors[c] = X_c.shape[0] / float(n_samples)
def _cal_prob(self, x, y) -> float:
"""
计算 log p(x_1|y) + ... + log p(x_n|y) + log p(y)
注意:这里的 y 为在 self._classes 中的下标
"""
# log p(y)
p_y = np.log(self._priors[y])
# print(p_y)
# log p(x_1|y) + ... + log p(x_n|y)
mean = self._mean[y]
# print("mean are: \n", mean)
var = self._var[y]
# print("var: \n", var)
p_xy = np.log(np.exp(-(x-mean)**2 / (2*var)) / np.sqrt(2*np.pi*var))
p_xy[np.isnan(p_xy)] = 0
return p_y + np.sum(p_xy)
def predict(self, X) -> List[int]:
"""
测试数据
:return: 预测结果
"""
pred = []
for x in X:
probs = []
for y in range(len(self._classes)):
probs.append(self._cal_prob(x, y))
# 找到概率最大的那个类别
pred.append(np.argmax(probs))
return pred
def accuracy(y_pred, y_true):
return np.sum(y_pred == y_true) / len(y_true)
def show_test_result(nb: NaiveBayes):
fig = plt.figure(figsize=(6, 6)) # figure size in inches
fig.subplots_adjust(left=0, right=1, bottom=0, top=1, hspace=0.05, wspace=0.05)
n = digits.data.shape[0]
for i in range(64):
ax = fig.add_subplot(8, 8, i + 1, xticks=[], yticks=[])
# 随机选择一张图片
idx = random.randint(0, n-1)
ax.imshow(digits.images[idx], cmap=plt.cm.binary, interpolation='nearest')
data = np.asarray(digits.data[idx]).reshape(1, -1)
# 左下角为预测值
ax.text(0, 7, str(nb.predict(data)[0]), color="green")
# 左上角为真实值
ax.text(0, 1, str(digits.target[idx]), color="red")
plt.show()
digits = load_digits()
if __name__ == '__main__':
x_train, x_test, y_train, y_test = train_test_split(digits.data,
digits.target,
test_size=0.3)
bayes = NaiveBayes()
bayes.train(x_train, y_train)
y_pred = bayes.predict(x_test)
print(accuracy(y_pred, y_test))
show_test_result(bayes)