机器学习(2)-朴素贝叶斯的理解和代码实现

本文主要从下面两个方面展开

我的理解

朴素贝叶斯的本质就是在假设输入特征相互独立的条件下（这个假设是为了方便计算似然函数$P(X|Y)$，因为这里的特征X可能有很多属性，而且取值也很多，先验概率还是后验概率的区分主要是先搞清楚研究的是什么，比如说这里我们想要知道的是类别Y，那么先验概率就是P(Y),后验概率就是P(Y|X)），利用后验概率对样本进行分类的学习方法，后验概率由贝叶斯理论可表示为
$$P(Y|X)=\frac{P(Y)P(X|Y)}{P(X)}$$
对不同的y求出相应的后验概率，最终输出的y值为后验概率最大的y。在实际问题中，对于同一个输入$X,P(X)$的值是固定的，因此可以将问题转为为$\underset{y}{\mathrm{argmax}}P(Y)P(X|Y)$,因为X的特征相互独立，所以$P(X|Y)=\prod _{i=1}^{n}P(X^{(i)}|Y)$,n为输入特征维度，先验概率$P(Y)和条件概率P(X^{(i)}|Y)$可以用最大似然估计求得，即可以用样本的频率代替
具体计算方法如下：
首先对公式中的变量进行一些说明，$N$表示样本总数，$n$表示的输入特征X的维度，即特征的个数，$c_k$表示$Y$的第k取个值，那么根据最大似然估计可以得到
$$P(Y=c_k)=\frac{I(\sum _{i=1}^N{y}_i=c_k)}{N}$$
$$P(X^{(j)}=a_{jl}|Y=c_k)=\frac{\sum _{i=1}^{N}I(x_i^{(j)}=a_{jl},y_i=c_k)}{\sum _{i=1}^{N}I(y_i=c_k)}$$
$I$为指示函数，$a_{jl}$是特征$X^{(j)}$的第$l$个取值
直接用最大似然估计可能会导致计算条件概率是分母为零的情况，这个时候引入贝叶斯估计，得到先验概率和条件概率的计算公式如下：
$$P_{\lambda }(Y=c_k)=\frac{I(\sum _{i=1}^N{y}_i=c_k)+\lambda }{N+K\lambda }$$
$$P_{\lambda }(X^{(j)}=a_{jl}|Y=c_k)=\frac{\sum _{i=1}^{N}I(x_i^{(j)}=a_{jl},y_i=c_k)+\lambda }{\sum _{i=1}^{N}I(y_i=c_k)+S_j\lambda }$$
其中K表示Y的种类数，$S_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%