利用朴素贝叶斯算法分析测试集（贝叶斯线性回归）

利用朴素贝叶斯算法分析测试集（贝叶斯线性回归）贝叶斯线性建模的最终结果不是模型参数的单一估计值，而是我们可以用来推断新观测值的分布。这种分布使我们能够证明模型中的不确定性，并且是贝叶斯建模方法的好处之一。随着数据点数量的增加，不确定性应该降低，这表明我们的估计有更高的确定性。在实践中，计算精确的后验分布对于连续值是计算上难以处理的，因此我们转而采用诸如马尔科夫链蒙特卡洛（MCMC）的抽样方法从后面抽取样本以近似后验。蒙特卡罗是指绘制随机样本的一般技术，马尔可夫链意味着下一个绘制的样本仅基于前一个样本值。这个概念是，当我们绘制更多的样本时，后验的近似值将最终收敛于模型参数的真实后验分布。相反，贝叶斯线性回归假定响应是从概率分布（例如正态（高斯）分布）采样的：高斯平均值是参数β和输入X的乘积，标准差是σ。在贝叶斯模型中，不仅假设响应是从分布中采样的，而且参数也是如此。目标是确定给定输入X和输出y的模型参数的后验概率分布：后验等于数据乘以模

在这个贝叶斯机器学习项目的第一部分中，我们概述了我们的问题，进行了全面的探索性数据分析，选择了我们的功能并建立了基准。这里我们将在Python中实现贝叶斯线性回归来构建模型。在我们训练完模型后，我们将解释模型参数并使用模型进行预测。

提醒一下，我们正在研究监督式回归机器学习问题。使用学生成绩数据集，我们希望建立一个模型，可以根据学生的个人和学术特点预测最终学生的成绩。特征选择后的最终数据集是：

我们有6个特征（解释变量）用于预测目标（响应变量），在这种情况下是等级。训练集中有474名学生，测试集中有159名学生。为了了解变量分布（因为我真的很喜欢这个图），这里是一个显示散点图，直方图，密度图和相关系数的变量对图。

其中响应y由模型参数β产生，乘以输入矩阵X，加上由随机采样噪声或潜在变量引起的误差。在普通最小平方（OLS）方法中，模型参数β通过找到使训练数据的误差平方和最小的参数来计算。来自OLS的输出是给定训练数据的“最佳”模型参数的单点估计。这些参数可以用来预测新的数据点。

相反，贝叶斯线性回归假定响应是从概率分布（例如正态（高斯）分布）采样的：

高斯平均值是参数β和输入X的乘积，标准差是σ。在贝叶斯模型中，不仅假设响应是从分布中采样的，而且参数也是如此。目标是确定给定输入X和输出y的模型参数的后验概率分布：

后验等于数据乘以模型参数的先验值除以归一化常数的可能性。如果我们有一些领域知识，我们可以使用它来为模型参数指定先验，或者我们可以使用非信息先验：具有大标准偏差的分布不假定任何有关变量的分布。使用非信息性的先验手段我们“让数据说话”。一个普通的优先选择是使用β的正态分布和σ的半柯西分布。

在实践中，计算精确的后验分布对于连续值是计算上难以处理的，因此我们转而采用诸如马尔科夫链蒙特卡洛（MCMC）的抽样方法从后面抽取样本以近似后验。蒙特卡罗是指绘制随机样本的一般技术，马尔可夫链意味着下一个绘制的样本仅基于前一个样本值。这个概念是，当我们绘制更多的样本时，后验的近似值将最终收敛于模型参数的真实后验分布。

贝叶斯线性建模的最终结果不是模型参数的单一估计值，而是我们可以用来推断新观测值的分布。这种分布使我们能够证明模型中的不确定性，并且是贝叶斯建模方法的好处之一。随着数据点数量的增加，不确定性应该降低，这表明我们的估计有更高的确定性。

在Python中实现贝叶斯线性建模

用于概率编程和Python中贝叶斯推理的最佳库目前是PyMC3。它包括许多用于构建贝叶斯模型和使用MCMC方法来推断模型参数的实用程序。我们将使用PyMC3 的广义线性模型（GLM）模块，特别GLM.from_formula是构建贝叶斯线性模型的功能非常简单。

我们只需要执行两个步骤来执行此模块的贝叶斯线性回归：

• 建立一个与特征相关的公式，并决定数据可能性的先验分布

• 使用MCMC从参数后验分布中抽样

Formula

我们不必为每个模型参数分别定义概率分布，而是传递一个将特征（输入）与目标（输出）相关联的R型公式。以下是关于成绩与学生特点的公式：

Grade ~ failures higher_edu mother_edu studytime father_edu absences

在此语法中，〜被读作“是...的函数”。我们正在告诉模型：Grade是代字号右侧六个特征的线性组合。

该模型使用with语句在上下文中构建。在呼吁GLM.from_formula我们通过公式，数据和数据可能性系列（这实际上是可选的，默认为正态分布）。该函数解析公式，为每个要素添加随机变量（以及标准偏差），为数据添加可能性，并将参数初始化为合理的初始估计。默认情况下，模型参数先验模型为正态分布。

一旦建立了GLM模型，我们使用MCMC算法从后验进行采样。如果我们没有指定哪种方法，PyMC3会自动为我们选择最好的。在下面的代码中，我让PyMC3选择采样器并指定样本数量，2000，链数，2和调整步数500。

import pymc3 as pm

# Context for the model

with pm.Model() as normal_model:

# The prior for the data likelihood is a Normal Distribution

family = pm.glm.families.Normal()

# Creating the model requires a formula and data (and optionally a family)

pm.GLM.from_formula(formula data = X_train family = family)

# Perform Markov Chain Monte Carlo sampling letting PyMC3 choose the algorithm

normal_trace = pm.sample(draws=2000 chains = 2 tune = 500 njobs=-1)

在这种情况下，PyMC3选择了No-U-Turn采样器，并使用jitter adapt_diag初始化采样器。

采样器运行几分钟，我们的结果存储在normal_trace。这包含每个模型参数的所有样本（除了被丢弃的调整样本）。跟踪本质上是我们的模型，因为它包含了我们执行推理所需的所有信息。要了解贝叶斯线性回归的用途，我们可以使用PyMC3中的内置函数来检查跟踪。

跟踪图显示左侧模型参数的后验分布以及右侧变量轨迹中绘制的样本的进展。两种颜色代表两个不同的采样链。

pm.traceplot（normal_trace）

在这里我们可以看到我们的模型参数不是点估计，而是分布。每个分布的均值可以被认为是最可能的估计值，但是我们也使用整个值的范围来表明我们对于真值的不确定性。

查看后验分布的另一种方法是直方图：

pm.plot_posterior（normal_trace）

在这里，我们可以看到我们可以用作最可能估计的均值，以及整个分布。95％HPD代表95％最高后验密度，是我们参数的可靠间隔。甲置信区间是在频率统计置信区间的贝叶斯当量（尽管具有不同的解释）。

我们还可以看到所有模型参数的摘要：

pm.df_summary（normal_trace）

我们可以用与OLS线性回归相同的方式来解释这些权重。例如在模型中：

• 以前的班级失败和缺勤有重要的负面影响

• 高等教育计划和学习时间具有积极的重要性

• 母亲和父亲的教育有一个积极的重量（尽管母亲更积极）

标准差列和hpd限制让我们对模型参数的自信心充满了信心。例如，该father_edu 特征的95％hpd从-0.22到0.27，这意味着我们不能完全确定模型中的效果是负面还是正面的！sd数据可能性也有很大的标准偏差（行），表明目标有很大的不确定性。总的来说，我们发现模型存在相当大的不确定性，因为我们正在处理少量样本。只有几百名学生，我们没有足够的数据来准确地确定模型参数。

解释变量效应

为了看到单个变量对等级的影响，我们可以改变这个变量的值，同时保持其他变量不变，并且看看估计等级如何变化。为了做到这一点，我们使用plot_posterior_predictive函数并假定除了感兴趣的变量（查询变量）之外的所有变量都处于中值。我们为查询变量生成一系列值，并且该函数通过从后验分布中绘制模型参数来估计该范围内的等级。代码如下：

# Examines the effect of changing a single variable

# Takes in the name of the variable the trace and the data

def model_effect(query_var trace X):

# Variables that do not change

steady_vars = list(X.columns)

steady_vars.remove(query_var)

# Linear Model that estimates a grade based on the value of the query variable

# and one sample from the trace

def lm(value sample):

# Prediction is the estimate given a value of the query variable

prediction = sample['Intercept'] sample[query_var] * value

# Each non-query variable is assumed to be at the median value

for var in steady_vars:

# Multiply the weight by the median value of the variable

prediction = sample[var] * X[var].median()

return prediction

figsize(6 6)

# Find the minimum and maximum values for the range of the query var

var_min = X[query_var].min()

var_max = X[query_var].max()

# Plot the estimated grade versus the range of query variable

pm.plot_posterior_predictive_glm(trace eval=np.linspace(var_min var_max 100)

lm=lm samples=100 color='blue'

alpha = 0.4 lw = 2)

# Plot formatting

plt.xlabel('%s' % query_var size = 16)

plt.ylabel('Grade' size = 16)

plt.title("Posterior of Grade vs %s" % query_var size = 18)

plt.show()

结果显示来自后部的100个样本的估计等级对查询变量的范围：

通过从后面的轨迹中选取一组模型参数并评估整个查询变量范围内的预测分数，绘制每条线（每个图中有100条）。线条的分布表明模型参数的不确定性：线条越分散，模型对该变量的影响越不确定。

对于一个变量，父亲的教育，我们的模型甚至不确定增加变量的效果是正面的还是负面的！

如果我们使用这个模型来做出决定，我们可能会想在不首先收集更多数据以形成更多确定的估计值的情况下进行部署。仅有数百名学生，模型参数存在相当大的不确定性。例如，我们不应该提出诸如“父亲的教育水平对成绩有积极影响”这样的说法，因为结果显示这个结论没有什么确定。

如果我们使用Frequentist方法，只看到一个点估计值，那么由于数据量有限，我们可能会做出错误的决定。在数据集有限的情况下，贝叶斯模型是显示模型中不确定性的最佳选择。

做预测

说到预测，贝叶斯模型可以用来估计分布。我们记得贝叶斯线性回归的模型是：

其中β是系数矩阵（模型参数），X是数据矩阵，σ是标准偏差。如果我们想使一个新的数据点的预测，我们可以找到一个正常的分布，通过我们的数据点的模型参数相乘找到均值和使用从模型参数的标准偏差估计输出。

在这种情况下，我们将采用迹线中每个模型参数的均值作为参数的最佳估计值。如果我们采用迹线中的参数的均值，则预测的分布变为：

Grade ~ N(9.20 * Intercept - 1.32 * failures 1.85 * higher_edu 0.26 * mother_edu 0.58 * studytime 0.03 * father_edu - 0.07 * absences 2.28^2)

对于一个新的数据点，我们用变量的值代替，并构建等级的概率密度函数。作为一个例子，下面是来自测试集的一个观察结果以及概率密度函数：

Test Observation:

failures = 0 higher_edu = 1 mother_edu = 2 studytime = 1

father_edu = 2 absences = 8

True Grade = 12

Average Estimate = 11.6763

5% Estimate = 7.7618 95% Estimate = 15.5931

对于这个数据点，平均估计值与实际等级一致，但也有很大的估计间隔。如果我们有更多的学生，估计的不确定性应该更低。

我们还可以对任何不在测试集中的新点进行预测：

New Observation:

absences = 1 failures = 0 father_edu = 1

higher_edu = 1 mother_edu = 4 studytime = 3

Average Estimate = 13.8009

5% Estimate = 10.0696 95% Estimate = 17.4629

与标准机器学习模型的比较

在本系列的第一部分中，我们计算了一些标准机器学习模型的基准以及一个天真的基准。要计算MAE和RMSE指标，我们需要对测试集中的所有数据点进行单点估计。我们可以使用估计分布的平均值做出“最有可能”的预测。结果指标以及基准指标如下所示：

贝叶斯线性回归与最佳标准模型的性能几乎相同！然而，贝叶斯线性建模的主要优点并不在于准确性，而在于我们不确定性的可解释性和量化。任何模型只是对现实世界的估计，在这里我们已经看到了我们在有限数据的训练模型中应该具有的信心。

结论

在这一系列文章中，我们走过了用于解决数据科学问题的完整机器学习过程。我们从探索性数据分析开始，转向建立基线，尝试了几种不同的模型，实现了我们的选择模型，解释了结果，并使用模型做出了新的预测。尽管模型实现细节可能会发生变化，但这种通用结构对于大多数数据科学项目来说都很好。