- 线性回归
import numpy as np N=100 x1 = np.random.randn(N) x2 = np.random.randn(N) alpha = 1 sigma = 1 beta = [1,2.5] y = alpha + beta[0]*x1 + beta[1]*x2 + np.random.randn(N)*sigma
np.random.randn()返回指定维度的array，一个位置参数为一个维度，位置参数的值为该维度的数据量，具有标准正态分布N(0,1)
此处为返回一维array，具有100个标准正态数据
其中sigma描述误差，beta为权重，alpha为偏置，此处是已知值
通过随机的输入数据以及随机的误差获得观测的array为y
if __name__ == '__main__': basic_model = pm.Model() with basic_model: alpha = pm.Normal('alpha',mu=0,sd=10) beta = pm.Normal('beta',mu=0,sd=10,shape=2) sigma = pm.HalfNormal('sigma',sd=1) mu = alpha + beta[0]*x1+beta[1]*x2 y_obs = pm.Normal('y_obs',mu=mu,sd=sigma,observed=y)
此处，认为sigma，beta，alpha，均为未知变量，尝试通过观测值为其恢复参数（也即进行后验估计）
在贝叶斯理论中，未知的变量都应指定先验分布
此处认为alpha，beta，sigma均具有如上所示的先验分布（后面会根据观测值对先验进行更新）
mu为随机变量之间的运算结果，注意此时输入仍然是x1和x2，只是sigma，beta，alpha不再是确定值，而是随机变量
获得此时的输出mu可以作为y的均值mu表示（实际上其为确定性随机变量，通过随机变量的运算得到的量）
此时的方差就是之前设置的误差sigma（也就是线性回归误差）
Y_obs为一个观测随机变量，表示模型数据的可能性，通过observed参数来告诉这个观测变量是被观测到的，与其他直接假设先验的随机变量不同，其不会因为拟合算法而发生改变
最后导入观测数据，就可以获得指定先验分布的变量的一些估计值
map_estimate = pm.find_MAP(model=basic_model)
此处是MAP方法，即极大后验估计方法
虽然极大后验估计是一个简单快速的估计未知参数的办法，但是没有对不确定性进行估计是其缺陷
为了使用MCMC采样以获得后验分布的样本，PyMC3需要指定和特定MCMC算法相关的步进方法（采样方法）,如Metropolis, Slice sampling, or the No-U-Turn Sampler (NUTS)。
PyMC3中的 step_methods 子模块提供了如下采样器： NUTS, Metropolis, Slice, HamiltonianMC, and BinaryMetropolis。可以由PyMC3自动指定也可以手动指定。自动指定是根据模型中的变量类型决定的：
* 二值变量：指定 BinaryMetropolis
* 离散变量：指定 Metropolis
* 连续变量：指定 NUTS
trace = pm.sample(2000) pm.traceplot(trace) plt.show() print(pm.summary(trace))
sample 函数用指定的迭代器（采样算法）进行了2000次迭代，收集到的采样值存储在 Trace 对象中，按照采样值获得的先后顺序存储。
Trace对象是一个字典中包含了多个map，比如可以使用trace[‘alpha']来进行使用，其是一个array对象
可以通过 step 参数指定特定的采样器（迭代器）来替换默认的迭代器NUTS
# 用MAP获得初始点 start = pm.find_MAP(fmin=optimize.fmin_powell) # 实例化采样器 step = pm.Slice(vars=[sigma]) # 对后验分布进行5000次采样 trace = pm.sample(5000, step=step,start=start)
可以设置sample的起始点为极大后验点
sigma 变量指定 Slice 采样器（step），那么其他两个变量（ alpha 和 beta）的采样器会自动指定（NUTS）
pm.traceplot(trace) 函数来绘制后验采样的趋势图
左侧是每个随机变量的边际后验的直方图（纵坐标为frequency，横坐标为变量自身的取值），使用核密度估计进行了平滑处理
右侧是马尔可夫链采样值按顺序绘制（纵坐标为采样值，横坐标为采样数量）
对于向量参数 beta 会有两条后验分布直方图和后验采样值
本例中达到的目的是
在拥有观测值集合y（里面包含了误差等外界因素），输入数据（x1和x2），输入和输出的函数关系，
去求解sigma（线性回归的误差，表示为方差），beta（权重），alpha（偏置）的估计值（也就是线性回归参数）
通过nuts采样，获得待求解变量的后验分布数据，从而获得其最高概率时的取值
这样就获得了参数值，完成了一个线性回归的拟合
- 抛硬币问题
硬币扔了n次，正面朝上是h次
想知道 p 的可能性。给定 n 扔的次数和 h 正面朝上次数，p 的值很可能接近 0.5，比如说在 [0.48，0.52]？