复习笔记1：朴素贝叶斯（面向机器学习入门 公式推导+代码实现）

/**

• Author: Julian Lee
• date: 23:45 01/06/2019
  */

写在前面

我的第一篇博客。临近考试月，静下心来把机器学习这门课的课件，代码和reading list温习一遍，颇有心得。写在这里主要是帮自己加深记忆，如果让翻到这篇的你少踩一些坑，不胜荣幸。

朴素贝叶斯的推导很简单，几乎没有需要计算的内容，我认为重点是理解每个条件概率的“条件“指的是什么。

如果有不对的地方欢迎指正。

祝我拿distinction。

0. 假定你已经具备的知识

1. 矩阵+向量
2. 贝叶斯定理
3. 最大似然估计

1. 两个假设

1. 每条数据独立同分布
2. 各特征条件独立

先记住这两条，至于为什么后面会说。

2. 解构贝叶斯定理

$$\text{posterior（后验）}=\frac{\text{likelihood（似然）}\times \text{prior（先验）}}{\text{marginal likelihood（边际似然）}}$$

目标为求出后验，后验表示为：$P(y^{\ast }|\mathbf{y},\mathbf{X},{\mathbf{x}}^{\ast },\mathbf{\theta })$, 很直观，训练集$\mathbf{X},\mathbf{y}$训练模型，训练完成后，输入测试数据的特征 ${\mathbf{x}}^{\ast }$ ，得到$y^{\ast }$

然而，为了更具解释性，我们先从如下开始，再推导出标准的朴素贝叶斯后验计算公式：
$$P(y^{\ast }|\mathbf{y},\mathbf{X},{\mathbf{x}}^{\ast },\mathbf{\theta })=\frac{p(y\ast ,\mathbf{y},\mathbf{X},{\mathbf{x}}^{\ast }|\mathbf{\theta })}{p(\mathbf{y},\mathbf{X},{\mathbf{x}}^{\ast }|\mathbf{\theta })}$$
${\mathbf{x}}^{\ast },y^{\ast }$指一条测试数据的特征与标签
$\mathbf{X},\mathbf{y}$指训练数据的特征（矩阵）与标签（向量）
$\mathbf{\theta }$指模型参数，我倾向于这样描述，特征与标签的产生机制，若某一特征服从高斯分布，则为均值和方差向量，若为二项分布，则为伯努利分布的参数概率$\pi$

2.1 分子部分

由于独立假设的存在，联合概率可有边际概率简单相乘得到
应用假设1，训练集与测试数据的联合概率可表达为：$$p(y^{\ast },{\mathbf{x}}^{\ast },\mathbf{y},\mathbf{X}|\mathbf{\theta })=p(y^{\ast },{\mathbf{x}}^{\ast }|\mathbf{\theta })p(\mathbf{X},\mathbf{y}|\mathbf{\theta })=p(y^{\ast },{\mathbf{x}}^{\ast }|\mathbf{\theta })\prod _{i=1}^{n}p(y_i,{\mathbf{x}}_i|\mathbf{\theta }).$$
假设2强调条件独立，因此每条数据给定标签的条件概率，可以表示为该数据点的每个特征给定标签的条件概率：$$p({\mathbf{x}}_i|y_i,\mathbf{\theta })=\prod _{j=1}^{p}p(x_{i,j}|y_i,\mathbf{\theta })$$

用上式表达训练集的联合概率，根据贝叶斯定理，还需要乘上$p(y_i)$：
$$p(\mathbf{y},\mathbf{X}|\mathbf{\theta },\pi )=\prod _{i=1}^n\prod _{j=1}^{p}p(x_{i,j}|y_i,\mathbf{\theta })p(y_i|\pi )$$
朴素贝叶斯是一种分类算法，故标签y为定性数据。本文以最简单的二元分类（binary）问题为例，y的取值为{0, 1}或{False, True}。 此时y是服从伯努利分布的，概率为$\pi$,写作$p(y_i|\pi )$.上式具有两个模型参数, $\theta$是特征$\mathbf{X}$所服从的分布,

2.1.1 最大似然

通常我们需要最小化目标函数，所以对上式取负对数（这里的log省略了底数e），得到：
$$E(\mathbf{\theta },\pi )=-\log p(\mathbf{y},\mathbf{X}|\mathbf{\theta },\pi )=-\sum _{i=1}^n\sum _{j=1}^p\log p(x_{i,j}|y_i,\mathbf{\theta })-\sum _{i=1}^n\log p(y_i|\pi ),$$
此时对数似然含有两项，分别只含$\mathbf{\theta }$和$\pi$，可以分开求。
$$E(\mathbf{\theta },\pi )=E(\mathbf{\theta })+E(\pi )=\sum _{i=1}^{2p}E({\theta }_i)+E(\pi )$$
注意，$\mathbf{\theta }$对于每个特征，每种标签均有1个$\theta$，需要逐列去求，根据不同分布去算似然，本文共有p个特征，2类标签，故有2p个$\theta$。而$\pi 只有一个$，先给出结论，$\pi$的取值为所有标签中1的占比。$$\pi =\frac{{\sum }_{i=1}^n{y}_i}{n}$$ 计算伯努利分布solution，推导附在最后。

2.1.2 整合分子部分

为方便读者对照和自己复习，应用假设1，分解测试数据和训练集部分的分子如下：
$$p(y^{\ast },{\mathbf{x}}^{\ast },\mathbf{y},\mathbf{X}|\mathbf{\theta })=p(y^{\ast },{\mathbf{x}}^{\ast }|\mathbf{\theta })p(\mathbf{X},\mathbf{y}|\mathbf{\theta })=p(y^{\ast },{\mathbf{x}}^{\ast }|\mathbf{\theta })\prod _{i=1}^{n}p(y_i,{\mathbf{x}}_i|\mathbf{\theta }).$$
而应用假设2与贝叶斯定理后，分解为
$$p(y^{\ast },\mathbf{y},\mathbf{X},{\mathbf{x}}^{\ast }|\mathbf{\theta })=\prod _{j=1}^{p}p(x_j^{\ast }|y^{\ast },\mathbf{\theta })p(y^{\ast }|\pi )\prod _{i=1}^n\prod _{j=1}^{p}p(x_{i,j}|y_i,\mathbf{\theta })p(y_i|\pi )$$

2.2 分母部分

原式如下：
$$P(y^{\ast }|\mathbf{y},\mathbf{X},{\mathbf{x}}^{\ast },\mathbf{\theta })=\frac{p(y\ast ,\mathbf{y},\mathbf{X},{\mathbf{x}}^{\ast }|\mathbf{\theta })}{p(\mathbf{y},\mathbf{X},{\mathbf{x}}^{\ast }|\mathbf{\theta })}$$
根据全概率公式分母可由分子部分积分或求和得到，先照抄分子，再加上求和符号：
$$P(y^{\ast }|\mathbf{y},\mathbf{X},{\mathbf{x}}^{\ast },\mathbf{\theta })=\frac{{\prod }_{j=1}^{p}p(x_j^{\ast }|y^{\ast },\mathbf{\theta })p(y^{\ast }|\pi ){\prod }_{i=1}^n{\prod }_{j=1}^{p}p(x_{i,j}|y_i,\mathbf{\theta })p(y_i|\pi )}{{\sum }_{y^{\ast }=0}^1{\prod }_{j=1}^{p}p(x_j^{\ast }|y^{\ast },\mathbf{\theta })p(y^{\ast }|\pi ){\prod }_{i=1}^n{\prod }_{j=1}^{p}p(x_{i,j}|y_i,\mathbf{\theta })p(y_i|\pi )}$$

2.3化简

化简后得到：
$$P(y^{\ast }|\mathbf{y},\mathbf{X},{\mathbf{x}}^{\ast },\mathbf{\theta })=\frac{{\prod }_{j=1}^{p}p(x_j^{\ast }|y^{\ast },\mathbf{\theta })p(y^{\ast }|\pi )}{{\sum }_{y^{\ast }=0}^1{\prod }_{j=1}^{p}p(x_j^{\ast }|y^{\ast },\mathbf{\theta })p(y^{\ast }|\pi )}$$
代码实现篇将按照上面式子进行，会在代码实现篇对每个元素进行详细解释，先挖个坑。

这时就能整理成贝叶斯定理的形式：
$$P(y^{\ast }|{\mathbf{x}}^{\ast },\mathbf{\theta })=\frac{p({\mathbf{x}}^{\ast }|y^{\ast },\mathbf{\theta })p(y^{\ast }|\pi )}{p({\mathbf{x}}^{\ast }|\pi ,\mathbf{\theta })}$$
注:并非说明后验与训练集无关，而是训练集的分布模型用参数表示了

3 代码实现

3.1 公式解析

推导出的公式看似略微复杂，其实作为经典算法，具有很好的可解释性。

$P(y^{\ast }|{\mathbf{x}}^{\ast },\mathbf{\theta })$：后验，这里确定了模型的输入和输出，为了与推导过程一致，省略了标签参数$\pi$（其实是存在的），其中$\mathbf{\theta }$指特征参数，输入训练特征，得到预测标签。

$p({\mathbf{x}}^{\ast }|y^{\ast },\mathbf{\theta })$：似然，需要注意的是条件$y^{\ast }$，例如当$y^{\ast }=1$时，选取的会是p个$\theta$，$y^{\ast }=0$时选取的是另外p个$\theta$。

$p(y^{\ast }|\pi )$：先验，所谓先验是没看过特征的时候，你所估猜的标签，前已推导过，先验为该标签的占比。垃圾邮件占全部邮件的1%，在没有邮件过滤系统的情况下，收到一封邮件，你自然会不会认为是垃圾邮件。

$p({\mathbf{x}}^{\ast }|\pi ,\mathbf{\theta })$：边际似然，前已提过这是分子对y求和，所以是对y的不同情况（二元分类问题只有True和False两类情况）下的似然求和。

3.2计算思路

1. 确定数据的特征服从几种分布，若部分特征服从高斯分布，其他服从二项分布，则需要有两个表来记录参数。本文假定所有特征服从高斯分布，参数表如下，行索引下标指代y的真假情况。

	feature1	feature2	feature3
$w_0$
${\sigma }_0^2$
$w_1$
${\sigma }_1^2$

NormalArgs = pd.DataFrame(data=np.zeros((4,len(quant_columns))), columns=quant_columns, index=['w_0', 'sigma2_0', 'w_1', 'sigma2_1'])

2. 计算$\theta$，填充表格

# y的格式为boolean向量，ddof设为0，因为是有偏估计
for column in X_train:
	NormalArgs[column]['mu_0'] = X_train[column][-y].mean()
	NormalArgs[column]['mu_1'] = X_train[column][y].mean()
	NormalArgs[column]['sigma2_0'] = X_train[column][-y].var(ddof=0)
	NormalArgs[column]['sigma2_1'] = X_train[column][y].var(ddof=0)

3. 计算先验
   prior = float(y_train.sum())/len(y_train)#这里计算的是1的先验
4. 逐列求似然并求和

#计算似然
def log_likelihood(x, mu, sigma2):
	return -0.5* np.log(2*np.pi*sigma2)-((x-mu)**2)/(2*sigma2)
#计算后验（可以理解为返回标签是1的概率）
def predict(X_test, NormalArgs, prior):
	#新增两列，储存结果为1的似然和结果为2的似然
    log_positive = pd.Series(data = np.zeros(X_test.shape[0]), index=X_test.index)
    log_negative = pd.Series(data = np.zeros(X_test.shape[0]), index=X_test.index)
    #开始循环计算
    for column in X_test.columns:
        log_positive += log_likelihood(X_test[column], NormalArgs[column]['mu_1'], NormalArgs[column]['sigma2_1'])
        log_negative += log_likelihood(X_test[column], NormalArgs[column]['mu_0'], NormalArgs[column]['sigma2_0'])
    #按照后验公式返回似然（由于求了对数似然 需要指数逆运算回去）
    return np.exp(log_positive + np.log(prior))/(np.exp(log_positive + np.log(prior)) + np.exp(log_negative + np.log(1-prior)))

5.此时返回的是标签为1的概率，阈值通常设为0.5，再转化为标签即可进行模型评估，不再赘述。

伯努利似然推导

对于任意一次伯努利试验，其概率密度为
$$P(y)={\pi }^y(1-\pi {)}^{(1-y)}$$
对于给定的训练集，已经进行过n次试验，根据前面提到的假设2，每个y互相独立，这里给出联合密度：
$$p(\mathbf{y}|\pi )=\prod _{i=1}^n{\pi }^{y_i}(1-\pi {)}^{1-y_i}$$
通常我们需要最小化目标函数，所以对上式取负对数（这里的log省略了底数e），得到：
$$E(\pi )=-\log p(\mathbf{y}|\pi )=-\sum _{i=1}^n{y}_i\log \pi -\sum _{i=1}^n(1-y_i)\log (1-\pi )$$
稍微提下，最大似然估计是通过调整参数（视参数为变量），使模型最符合实际情况，因此需要对上式求导：
$$\frac{\text{d}E(\pi )}{\text{d}\pi }=-\frac{{\sum }_{i=1}^n{y}_i}{\pi }+\frac{{\sum }_{i=1}^n(1-y_i)}{1-\pi }$$
令$\frac{\text{d}E(\pi )}{\text{d}\pi }=0$，即可求出$$\pi =\frac{{\sum }_{i=1}^n{y}_i}{n}$$
上面导数表达式的两项均递增，所以此时确为全局最小值。

Reference:
Rogers, S. and Girolami, M. (n.d.). A first course in machine learning.