基于概率编程的贝叶斯推断

(华南农业大学 广东 广州 510000)

基于概率编程的贝叶斯推断

林镇溪

(华南农业大学广东广州510000)

本文介绍了基于概率编程的贝叶斯推断的建模过程和利用贝叶斯推断对挑战者号事故进行分析和推断，建立事故发生概率关于温度的逻辑分布模型p(t)=1/(1+eβt+α)，其中α，β采用正态分布，最后采用伯努利分布将事故发生的概率和事故发生的结果联系起来，使用PyMC计算后验分布，分析得出当温度为31度时事故发生的概率几乎是必然的。

概率编程；贝叶斯推断；事故发生概率；PyMC

一、实验背景

贝叶斯推断是一种统计学方法，用来估计统计量的性质。贝叶斯推断与其他统计学推断方法截然不同，他建立在主观判断的基础上，也就是说，你不需要客观证据，先估计一个值，然后根据实际结果不断的修正。概率编程是对问题建模，然后利用计算机采样的方法进行自动的贝叶斯推理，得出未知参数的概率分布。PyMC是一个马尔科夫链蒙特卡洛采样(Markov chain Monte Carlo Sampling)工具包，是一个做贝叶斯分析使用的Python库,包含了画图，拟合优化度和聚类诊断的方法。本文将利用PyMC来进行贝叶斯的推断[1]。

二、实验方法

贝叶斯是通过引入先验的不确定性，我们事实上允许了我们的主观判断可能是错误的。在观察数据，证据或其他信息后，我们不断的更新我们的判断使得它错的不那么离谱。所以在贝叶斯推断的过程中，最重要的就是先验分布的构建和确立。第一步，先确定需要推断的对象以及和对象相关的参数。例如我们要推断一个西瓜是否是好瓜，那么我们需要推断的对象就是西瓜的好坏，而和西瓜好坏相关的参数可以是西瓜的色泽，根蒂，敲声等等[2]。第二步，确定参数和对象的模型。根据对象以及参数的类型以及取值范围构建合理的模型。例如对象的取值为0或者1，那么可以猜测参数和对象之间的模型为逻辑模型。第三步，确定先验分布的模型。根据对象的数据类型我们需要确定一个先验分布，例如如果数据是离散型的，那么对象的分布可能符合Poisson分布；如果数据是连续的，那么对象的分布可能符合指数分布。

三、实验过程

(一)实验背景

1986年1月28号，挑战者号起飞不久后一个火箭推动器发生了爆炸。事故的起因是因为连接在火箭推进器上的O型圈有缺陷，这种缺陷来自于设计的不合理，这种设计使得O型圈对很多因素包括外界温度都是非常敏感的。我们将利用挑战者号的数据来分析最后一次飞行发生事故的概率,即当温度为31度时，事故发生的概率。

(二)数据整理

数据中包括了日期，温度以及事故是否发生。由于数据中包含了缺失值，可能会影响推断的结果，所以我们需要剔除这些包含缺失值的样本，而且根据背景我们知道事故的发生和温度有关，而和日期的关系不太大，所以我们也需要剔除那些对结果没什么影响的因素。

(三)构建模型

1.数据观察

由于直观数据很难展现温度和事故发生的关系，所以我们首先作图来观察它们之间的关系。因为事故发生的情况只有两种，0(不发生)和1(发生)。我们可以清晰地看出，随着外界温度的下降，发生事故的概率变得更高。因为温度和事故发生之间没有一个严格的转折点，所以我们需要对事故发生的概率p建模。

2.构建先验分布

在这个模型模型中，t为温度，p(t)为事故发生的概率，β是个不确定变量。

(2)确定β和α的模型。虽然我们确定了温度和事故发生概率的模型，但是我们引进了两个辅助参数，然而我们对这两个参数没有任何的信息，而且可以看出β和α取值可为正可为负，一般这种情况我们可以想到的正态分布来模拟，因为世间万物很多都满足正态分布。

(3)确定先验分布模型。我们前面分析了两个模型，事故发生的概率和温度的模型，事故发生的概率中参数β和α的模型，相当于如果我们确定了β和α的值，我们只能算出一个概率，但是最后的推断结果为0或者1。所以我们需要构建概率和推断结果的模型。这里可以很容易想到伯努利分布。

所以我们最终确定的先验分布模型如下：

Di～Ber(p(ti)),i=1,2,…,N

其中，

α～N(μ,τ)

β～N(μ,τ)

3.利用PyMC3计算后验分布

PyMC3可以很方便地帮我们算出模型的后验分布，只需要我们确定参数以及模型即可。

(四)结果分析

根据MCMC采样结果，下面我们作出了事故发生的概率关于温度的曲线图。

上面的图描绘了事故发生概率的后验分布图，由于受到温度，α，β的影响，我们不能从图中判断在哪个温度时概率发生明显的变化。但是从图中的曲线斜率的变化趋势也可以看出在温度在[60,70]区间变化最快，所以我们猜测温度的影响值在[60,70]之间。我们作出当温度t=31时事故发生的后验分布如下。

从上面的直方图可以看出当温度t=31°时事故发生的概率基本为99.9%以上，说明1986年1月28号O型圈发生事故几乎是必然的了。除此之外，我们还可以预测其他温度发生事故的概率。

[1]贝叶斯方法.概率编程与贝叶斯推断，Cameron Davidson-Pilon

[2]周志华.机器学习

林镇溪，男，汉族，广东廉江，本科，华南农业大学，研究方向机器学习。