机器学习（14）——朴素贝叶斯
前言：在所有的机器学习分类算法中，朴素贝叶斯和其他绝大多数的分类算法都不同。对于大多数的分类算法，比如决策树,KNN,逻辑回归，支持向量机等，他们都是判别方法，也就是直接学习出特征输出Y和特征X之间的关系，要么是决策函数Y=f(X)要么是条件分布P(Y|X)。但是朴素贝叶斯却是生成方法，也就是直接找出特征输出Y和特征X的联合分布然后用P(Y|X)=P(X,Y)/P(X)得出。
朴素贝叶斯很直观，计算量也不大，在很多领域有广泛的应用，
算法思想：基于概率的预测
逻辑回归通过拟合曲线（或者学习超平面）实现分类，决策树通过寻找最佳划分特征进而学习样本路径实现分类，支持向量机通过寻找分类超平面进而最大化类别间隔实现分类。相比之下，朴素贝叶斯独辟蹊径，通过考虑特征概率来预测分类。
贝叶斯思想
那么如何通过概率来进行决策的构建呢？举个可能不太恰当的例子：
我们小时候的看电影，有时候“坏人”的出现非常明显，我们一眼就能看见这个人是“好人”还是“坏人”，我们可以通过这个人的装扮，这个人说话的语气还有这个人的表情和笑容，以笑容来说，假如眼前有100个人，好人和坏人个数差不多，现在要用他们来训练一个“坏蛋识别器”。怎么办呢？咱们不管他们干过啥事，只看他们长啥样（这确实不是个恰当的例子）。也就是说，我们在区分好坏人时，只考虑他们的样貌特征。比如说“笑”这个特征，它的取值可能是“甜美的笑”、“儒雅的笑”、“憨厚的笑”、“没心没肺的笑”、“微微一笑”，等等——这都是“好人的笑”；也可以是“阴险的笑”、“不屑的笑”、“色眯眯的笑”、“任我行似的笑”、“冷笑”、“皮笑肉不笑”，等等——这很可能是“坏人的笑”。单就“笑”这个特征来说，一个好人发出“好人的笑”的概率更大，而且频率更高；而坏人则发出“坏人的笑”的概率更大，频率更高（电视上总能看见作奸犯科的人在暗地里发出挨千刀的笑）。
当然这个例子在如今的影视剧中少多了，好人也有发出坏笑的时候.
贝叶斯公式
再给出贝叶斯公式之前先介绍一下贝叶斯学派。
贝叶斯学派很古老，但是从诞生到一百年前一直不是主流。主流是频率学派。频率学派的权威皮尔逊和费歇尔都对贝叶斯学派不屑一顾，但是贝叶斯学派硬是凭借在现代特定领域的出色应用表现为自己赢得了半壁江山。
贝叶斯学派的思想可以概括为先验概率+数据=后验概率。也就是说我们在实际问题中需要得到的后验概率，可以通过先验概率和数据一起综合得到。数据大家好理解，被频率学派攻击的是先验概率，一般来说先验概率就是我们对于数据所在领域的历史经验，但是这个经验常常难以量化或者模型化，于是贝叶斯学派大胆的假设先验分布的模型，比如正态分布，beta分布等。这个假设一般没有特定的依据，因此一直被频率学派认为很荒谬。虽然难以从严密的数学逻辑里推出贝叶斯学派的逻辑，但是在很多实际应用中，贝叶斯理论很好用，比如垃圾邮件分类，文本分类。
我们先看看条件独立公式，如果X和Y相互独立，则有：
P(X,Y)=P(X)P(Y)
我们接着看看条件概率公式：
P(Y|X)=P(X,Y)/P(X)
P(X|Y)=P(X,Y)/P(Y)
或者说:
P(Y|X)=P(X|Y)P(Y)/P(X)
下面给出几个贝叶斯的基本概念：
先验概率P(A)：在不考虑任何情况下，A事件发生的概率。
注意：1）先验概率是已经给定的；
2）根据历史统计出来的
条件概率P(B|A)：A事件发生的情况下，B事件发生的概率。
后验概率P(A|B)：在B事件发生之后，对A事件发生的概率的重新评估。
全概率：如果A和A'构成样本空间的一个划分，那么事件B的概率为：A和A'的概 率分别乘以B对这两个事件的概率之和。
将后验概率和全概率公式进行组合就可以得到贝叶斯公式，如下：
举个例子：
比如明天晴天的概率是70%（非晴天的概率则为30%），这个概率叫做先验概率P(A)。
然后，晚上观测星星。这个是有频率概率的，历史统计表明，第二天是晴天的情况下，头一天晚上有星星的概率是80%，计做P(B|A)；第二天是非晴天的情况下，头天晚上有星星的概率是30%。晚上观测的结果是：星光灿烂。这个叫做条件P(B)。
最后，推算在该条件下（星星出现）时的明天是晴天的概率，叫做后验概率P(A|B)。
根据贝叶斯公式，P(A|B)=P(A)P(B|A)/P(B)=(0.70.8)/(0.70.8+0.30.3）=86%。
于是，我们可以说，明天下雨的概率，再加入了新的条件后，从70%上升到了86%。
图像表示，其实非常简单。下图左边矩形占全图的比例，即为先验概率。图中横贯东西的一条直线，就是条件。直线下方，左边深灰面积占全部灰色面积（深灰+浅灰）的比例，即为后验概率。假如条件之后使得原本的概率改变了，从面积就可以看得出来。
image.png
朴素贝叶斯算法
朴素贝叶斯(Naive Bayes， NB)是基于“特征之间是独立的”这一朴素假设，应 用贝叶斯定理的监督学习算法
对应给定的样本X的特征向量x1,x2,...,xm；该样本X的类别y的概率可以由贝叶斯公 式得到：
特征属性之间是独立的，所以得到：
在给定样本的情况下，P(x1,x2,...,xm)是常数，所以得到：
我们要求得最终的模型可以写作：
算法流程如下：
通过上述贝叶斯算法的基本公式可以写出贝叶斯算法的流程如下：
1）设x={a1,a2,...,am}为待分类项，其中a为x的一个特征属性
2）类别集合为C={y1,y2,...,yn}
3）分别计算P(y1|x),P(y2|x),....,P(yn|x)的值（贝叶斯公式）
4）如果P(yk|x)=max{P(y1|x),P(y2|x),....,P(yn|x)},那么认为x为yk类型
将以上流程用图表示出来如下图所示：
image.png
朴素贝叶斯按照数据的先验概率的不同可以分为高斯朴素贝叶斯，伯努利朴素贝叶斯，多项式朴素贝叶斯。
1.高斯朴素贝叶斯
因此只需要计算出各个类别中此特征项划分的各个均值和标准差
2.伯努利朴素贝叶斯
Bernoulli Naive Bayes是指当特征属性为连续值时，而且分布服从伯努利分布， 那么在计算P(x|y)的时候可以直接使用伯努利分布的概率公式：
伯努利分布是一种离散分布，只有两种可能的结果。1表示成功，出现的概率为p； 0表示失败，出现的概率为q=1-p；其中均值为E(x)=p，方差为Var(X)=p(1-p)
3.多项式朴素贝叶斯
Multinomial Naive Bayes是指当特征属性服从多项分布，从而，对于每个类别 y，参数为θy=(θy1,θy2,...,θyn)，其中n为特征属性数目，那么P(xi|y)的概率为θyi
示例：文本数据分类
贝叶斯经常用着文本的处理等方面，比如文本的分类和垃圾邮件的过滤等，下面以在新闻中文本的分类为例简单介绍一下贝叶斯的应用。
api介绍：
朴素贝叶斯是一类比较简单的算法，scikit-learn中朴素贝叶斯类库的使用也比较简单。相对于决策树，KNN之类的算法，朴素贝叶斯需要关注的参数是比较少的，这样也比较容易掌握。在scikit-learn中，一共有3个朴素贝叶斯的分类算法类。分别是GaussianNB，MultinomialNB和BernoulliNB。其中GaussianNB就是先验为高斯分布的朴素贝叶斯，MultinomialNB就是先验为多项式分布的朴素贝叶斯，而BernoulliNB就是先验为伯努利分布的朴素贝叶斯。
这三个类适用的分类场景各不相同，一般来说，如果样本特征的分布大部分是连续值，使用GaussianNB会比较好。如果如果样本特征的分大部分是多元离散值，使用MultinomialNB比较合适。而如果样本特征是二元离散值或者很稀疏的多元离散值，应该使用BernoulliNB。
在使用GaussianNB的fit方法拟合数据后，我们可以进行预测。此时预测有三种方法，包括predict，predict_log_proba和predict_proba。
predict方法就是我们最常用的预测方法，直接给出测试集的预测类别输出。
predict_proba则不同，它会给出测试集样本在各个类别上预测的概率。容易理解，predict_proba预测出的各个类别概率里的最大值对应的类别，也就是predict方法得到类别。
predict_log_proba和predict_proba类似，它会给出测试集样本在各个类别上预测的概率的一个对数转化。转化后predict_log_proba预测出的各个类别对数概率里的最大值对应的类别，也就是predict方法得到类别。
下面给一个具体的例子，代码如下：
import numpy as np
X = np.array([[-1, -1], [-2, -1], [-3, -2], [1, 1], [2, 1], [3, 2]])
Y = np.array([1, 1, 1, 2, 2, 2])
from sklearn.naive_bayes import GaussianNB
clf = GaussianNB()
#拟合数据
clf.fit(X, Y)
print("预测的结果是：")
print(clf.predict([[-0.8, -1]]))
print("属于每个类别的概率分别为：")
print(clf.predict_proba([[-0.8, -1]]))
print("属于每个类别的对数转化后的概率分别为：")
print(clf.predict_log_proba([[-0.8, -1]]))
输出的结果如下：
预测的结果是：
[1]
属于每个类别的概率分别为：
[[9.99999949e-01 5.05653254e-08]]
属于每个类别的对数转化后的概率分别为：
[[-5.05653266e-08 -1.67999998e+01]]
分析：从上面的结果可以看出，测试样本[-0.8,-1]的类别预测为类别1。具体的测试样本[-0.8,-1]被预测为1的概率为9.99999949e-01 ，远远大于预测为2的概率5.05653254e-08。这也是为什么最终的预测结果为1的原因了。
新闻分类示例：
我们浏览新闻的时候经常见到各种不同的新闻类别，比如“国际”“国内”“情感”“科技”等，现在需要建立一个模型，根据文章里面的内容，建立分类器，对新闻进行分类。
1.导入模块。
#导入我们要用的包，包括算法数据导入模块，算法评估模块，算法模块，以及画图模块。
这里引入多个分类器，来进行比较，来突出贝叶斯分类器强大的速度。
import numpy as np
from time import time
import matplotlib.pyplot as plt
import matplotlib as mpl
from sklearn.datasets import fetch_20newsgroups#引入新闻数据包
from sklearn.feature_extraction.text import TfidfVectorizer#做tfidf编码
from sklearn.feature_selection import SelectKBest, chi2#卡方检验——特征筛选
from sklearn.linear_model import RidgeClassifier
from sklearn.svm import LinearSVC,SVC
from sklearn.naive_bayes import MultinomialNB, BernoulliNB #引入多项式和伯努利的贝叶斯
from sklearn.neighbors import KNeighborsClassifier
from sklearn.ensemble import RandomForestClassifier
from sklearn.model_selection import GridSearchCV
from sklearn import metricsmpl.rcParams['font.sans-serif'] = [u'SimHei']
mpl.rcParams['axes.unicode_minus'] = False
2.导入数据
#这里的导入数据是直接从scikit-learn中的数据，一共三千多条新闻作为建立贝叶斯分类器的基本。
### 数据加载
print (u'加载数据...')
t_start = time()
## 不要头部信息
remove = ('headers', 'footers', 'quotes')
## 只要这四类数据
categories = 'alt.atheism', 'talk.religion.misc', 'comp.graphics', 'sci.space'
## 分别加载训练数据和测试数据
data_train = fetch_20newsgroups(data_home='./datas/',subset='train', categories=categories, shuffle=True, random_state=0, remove=remove)
data_test = fetch_20newsgroups(data_home='./datas/',subset='test', categories=categories, shuffle=True, random_state=0, remove=remove)
## 完成
print (u"完成数据加载过程.耗时:%.3fs" % (time() - t_start))
查看数据，查看一下总的数据，结果如下：
加载数据...
完成数据加载过程.耗时:2.094s
3.获取加载数据的相关信息
def size_mb(docs):
return sum(len(s.encode('utf-8')) for s in docs) / 1e6
categories = data_train.target_names
data_train_size_mb = size_mb(data_train.data)
data_test_size_mb = size_mb(data_test.data)
print (u'数据类型：', type(data_train.data))
print("%d文本数量 - %0.3fMB (训练数据集)" % (len(data_train.data), data_train_size_mb))
print("%d文本数量 - %0.3fMB (测试数据集)" % (len(data_test.data), data_test_size_mb))
print (u'训练集和测试集使用的%d个类别的名称：' % len(categories))
print(categories)
输出的结果如下：
数据类型： <class 'list'>
2034文本数量 - 2.428MB (训练数据集)
1353文本数量 - 1.800MB (测试数据集)
训练集和测试集使用的4个类别的名称：
['alt.atheism', 'comp.graphics', 'sci.space', 'talk.religion.misc']
4.建划分数据集
x_train = data_train.data
y_train = data_train.target
x_test = data_test.data
y_test = data_test.target
查看一条新闻的数据如下：
image.png
5.数据的处理，先把公章进行分词，然后过滤掉那些文本都会出现的单词，最后再把得到的特征单词进行降维处理
### 文档转换为向量
## 转换
# 要求输入的数据类型必须是list，list中的每条数据，单词是以空格分割开的
vectorizer = TfidfVectorizer(input='content', stop_words='english', max_df=0.5, sublinear_tf=True)
x_train = vectorizer.fit_transform(data_train.data) # x_train是稀疏的，scipy.sparse.csr.csr_matrix
x_test = vectorizer.transform(data_test.data)
print (u'训练集样本个数：%d，特征个数：%d' % x_train.shape)
print (u'停止词:\n')
print(vectorizer.get_stop_words())
## 获取最终的特征属性名称
feature_names = np.asarray(vectorizer.get_feature_names())
通常情况下需要过滤的单词如下：
image.png
特到的特征向量的数目为：
训练集样本个数：2034，特征个数：26576
6.特征选择
## 特征选择
ch2 = SelectKBest(chi2, k=1000)
x_train = ch2.fit_transform(x_train, y_train)
x_test = ch2.transform(x_test)
feature_names = [feature_names[i] for i in ch2.get_support(indices=True)]
把原本两万多个特征向量转变为1000个；通常利用卡方进行特征筛选。
随机查看一下特征属性如下：
7.模型的评估
### 基准模型方法
def benchmark(clf,name):
print (u'分类器：', clf)
## 设置最优参数，并使用5折交叉验证获取最优参数值
alpha_can = np.logspace(-2, 1, 10)
model = GridSearchCV(clf, param_grid={'alpha': alpha_can}, cv=5)
m = alpha_can.size
## 如果模型有一个参数是alpha，进行设置
if hasattr(clf, 'alpha'):
model.set_params(param_grid={'alpha': alpha_can})
m = alpha_can.size
## 如果模型有一个k近邻的参数，进行设置
if hasattr(clf, 'n_neighbors'):
neighbors_can = np.arange(1, 15)
model.set_params(param_grid={'n_neighbors': neighbors_can})
m = neighbors_can.size
## LinearSVC最优参数配置
if hasattr(clf, 'C'):
C_can = np.logspace(1, 3, 3)
model.set_params(param_grid={'C':C_can})
m = C_can.size
## SVM最优参数设置
if hasattr(clf, 'C') & hasattr(clf, 'gamma'):
C_can = np.logspace(1, 3, 3)
gamma_can = np.logspace(-3, 0, 3)
model.set_params(param_grid={'C':C_can, 'gamma':gamma_can})
m = C_can.size * gamma_can.size
## 设置深度相关参数，决策树
if hasattr(clf, 'max_depth'):
max_depth_can = np.arange(4, 10)
model.set_params(param_grid={'max_depth': max_depth_can})
m = max_depth_can.size
## 模型训练
t_start = time()
model.fit(x_train, y_train)
t_end = time()
t_train = (t_end - t_start) / (5*m)
print (u'5折交叉验证的训练时间为：%.3f秒/(5*%d)=%.3f秒' % ((t_end - t_start), m, t_train))
print (u'最优超参数为：', model.best_params_)
## 模型预测
t_start = time()
y_hat = model.predict(x_test)
t_end = time()
t_test = t_end - t_start
print (u'测试时间：%.3f秒' % t_test)
## 模型效果评估
train_acc = metrics.accuracy_score(y_train, model.predict(x_train))
test_acc = metrics.accuracy_score(y_test, y_hat)
print (u'训练集准确率：%.2f%%' % (100 * train_acc))
print (u'测试集准确率：%.2f%%' % (100 * test_acc))
## 返回结果(训练时间耗时，预测数据耗时，训练数据错误率，测试数据错误率, 名称)
return t_train, t_test, 1-train_acc, 1-test_acc, name
8.不同分类器之间进行比较
### 使用不同的分类器对数据进行比较
print (u'分类器的比较：\n')
clfs = [
[RidgeClassifier(), 'Ridge'],
[KNeighborsClassifier(), 'KNN'],
[MultinomialNB(), 'MultinomialNB'],
[BernoulliNB(), 'BernoulliNB'],
[RandomForestClassifier(n_estimators=200), 'RandomForest'],
[SVC(), 'SVM'],
[LinearSVC(loss='squared_hinge', penalty='l1', dual=False, tol=1e-4), 'LinearSVC-l1'],
[LinearSVC(loss='squared_hinge', penalty='l2', dual=False, tol=1e-4), 'LinearSVC-l2']
]
## 将训练数据保存到一个列表中
result = []
for clf,name in clfs:
# 计算算法结果
a = benchmark(clf,name)
# 追加到一个列表中，方便进行展示操作
result.append(a)
print ('\n')
## 将列表转换为数组
result = np.array(result)
结果如下：
image.png
9.将上述结果在直方图上表示，代码如下：
### 获取需要画图的数据
result = [[x[i] for x in result] for i in range(5)]
training_time, test_time, training_err, test_err, clf_names = result
training_time = np.array(training_time).astype(np.float)
test_time = np.array(test_time).astype(np.float)
training_err = np.array(training_err).astype(np.float)
test_err = np.array(test_err).astype(np.float)
### 画图
x = np.arange(len(training_time))
plt.figure(figsize=(10, 7), facecolor='w')
ax = plt.axes()
b0 = ax.bar(x+0.1, training_err, width=0.2, color='#77E0A0')
b1 = ax.bar(x+0.3, test_err, width=0.2, color='#8800FF')
ax2 = ax.twinx()
b2 = ax2.bar(x+0.5, training_time, width=0.2, color='#FFA0A0')
b3 = ax2.bar(x+0.7, test_time, width=0.2, color='#FF8080')
plt.xticks(x+0.5, clf_names)
plt.legend([b0[0], b1[0], b2[0], b3[0]], (u'训练集错误率', u'测试集错误率', u'训练时间', u'测试时间'), loc='upper left', shadow=True)
plt.title(u'新闻组文本数据分类及不同分类器效果比较', fontsize=18)
plt.xlabel(u'分类器名称')
plt.grid(True)
plt.tight_layout(2)
plt.show()
输出的结果如下：
