(本文假设读者懂贝叶斯统计分析的基本框架,主要是为了记录自己的学习,因此可能很多地方比较简略,并不适合当作入门读物。代码部分可以跳过。)
1986年,美国挑战者号航天飞机升空失败,导致数名宇航员死亡。事后调查发现,事故原因是 O 型环的设计问题,导致其对一些外在条件非常敏感,其中一个影响条件是外部温度。在此前的 24 次类似发射中,有 23 次的发射数据保留了下来,其中 7 次发生了O型环出故障。在挑战者号发射前夜,虽然工作人员讨论了 O 型环问题,但是认为这过往的 7 次故障没有体现出明显的规律,因此没有留意。
下图是过往 23 次发射记录中的O型环故障情况图,横轴是外部温度(华氏度),纵轴表示是否发生故障。
从图中来看,一个大致的趋势是温度越低,越可能发生故障(65度左右是一个可能的分水岭)。我们是否可以从中得到一个关于温度和故障几率的函数 p(t),给定一个温度 t,返回在该温度下O型环会故障的概率(或概率分布)?
让我们试着用贝叶斯统计来算出这个函数。
确定模型
很多函数可以做到给定一个值,返回一个 [0, 1] 之间的值,其中最常见的就是 logistic 函数:
其中 β 和 α 是需要通过贝叶斯统计计算出的参数。如果没有 α,则 t=0 会是 p(t) 的转折点,为了让分水岭发生在 t=65 度左右,所以需要一个 α 参数来平移函数图像
计算模型参数
α 和 β 可以取任何值,但是肯定有一个相对比较合理的区间,只是我们不知道这个区间,因此可以使用正态分布来模拟这两个参数。
import pymc3 as pm
import theano as tt
with pm.Model() as model:
beta = pm.Normal(