Pystan 后验不确定区间

Question

我在 PyStan 中看到 on another forum that PyStan doesn’t have the same function as RStan where they use posterior_interval(), but we can use numpy.percentile() instead. I’m currently using the pystan.StanModel.optimizing() 函数来获取使后验似然最大化的参数集。我现在还想获得后验结果的外部 95% 置信区间，所以我想知道 numpy.percentile() 函数是否会与 optimizing 函数一起使用？

我尝试找到参数分布的 95% 区间，但这并没有给出围绕结果的良好置信区间。特别是，我不认为它很好，因为当我期望后验呈现多峰分布时，我使用 numpy.percentile() 进行的置信区间在后验二维高斯补丁内。

我觉得95%的区间应该是从后方取的。 我会使用百分位数函数和优化函数来获得 95% 置信度的后验结果吗？

Answer 1

要获得后验估计的界限，需要对后验进行采样，而 pystan.StanModel.optimizing 则不需要。相反，使用 pystan.StanModel.sampling 方法从后验生成 MCMC 绘图。

如果您只需要标准置信区间的读数，那么 pystan.StanFit.stansummary() 方法可能就足够了，因为这将为每个打印出 2.5%、25%、50%、75% 和 97.5% 的分位数范围。例如，

fit = sm.sampling(...) # eight schools model
print(fit.stansummary())

Inference for Stan model: anon_model_19a09b474d1901f191444eaf8a6b8ce2.
4 chains, each with iter=10000; warmup=5000; thin=1;  post-warmup
draws per chain=5000, total post-warmup draws=20000.

           mean se_mean     sd   2.5%    25%    50%   75%  97.5%   n_eff   Rhat 
mu         7.98    0.05   5.04   -2.0   4.76   7.91  11.2   18.2   10614    1.0 
tau        6.54    0.08   5.65   0.24   2.49   5.25   8.98  20.65   4552    1.0 
eta[0]     0.39  6.7e-3   0.94  -1.53  -0.23   0.42   1.02   2.18  20000    1.0 
eta[1]   3.3e-4  6.2e-3   0.88  -1.74  -0.58-2.5e-3   0.57   1.75  20000    1.0 
eta[2]     -0.2  6.6e-3   0.93  -2.01  -0.84  -0.22   0.41   1.68  20000    1.0 
eta[3]    -0.03  6.3e-3   0.89   -1.8  -0.61  -0.03   0.56   1.75  20000    1.0 
eta[4]    -0.35  6.7e-3   0.88  -2.04  -0.94  -0.36   0.22   1.44  17344    1.0 
eta[5]    -0.22  6.6e-3    0.9  -1.96  -0.81  -0.24   0.35   1.59  18298    1.0 
eta[6]     0.34  6.8e-3   0.88  -1.43  -0.23   0.36   0.93   2.04  16644    1.0 
eta[7]     0.05  6.6e-3   0.93  -1.77  -0.58   0.04   0.66   1.88  20000    1.0 
theta[0]   11.4    0.07   8.23  -1.83   6.04  10.22  15.42  31.52  13699    1.0 
theta[1]   7.93    0.04   6.21  -4.58   4.09   7.89  11.79  20.48  20000    1.0 
theta[2]   6.17    0.06   7.72 -11.43   2.06   6.65  10.85  20.53  16367    1.0 
theta[3]   7.72    0.05   6.53  -5.29   3.68    7.7  11.75  21.04  20000    1.0 
theta[4]   5.14    0.04   6.35  -9.35   1.44   5.64   9.38  16.49  20000    1.0 
theta[5]   6.11    0.05   6.66  -8.44   2.22   6.44  10.41  18.52  20000    1.0 
theta[6]  10.63    0.05   6.76  -1.25   6.04  10.08  14.51  25.66  20000    1.0 
theta[7]    8.4    0.06   7.85  -7.56   3.89   8.17  12.76   25.3  16584    1.0 
lp__      -4.89    0.04   2.63 -10.79  -6.47  -4.66  -3.01  -0.43   5355    1.0

或者如果您需要特定的分位数，您可以使用您提到的 numpy.percentile。

但是，正如您正确观察到的那样，这对于多峰分布是不合适的。这种情况在中得到解决，但请注意，如果人们期望多个模式先验，则通常使用混合模型将模式分离为不同的单峰随机变量。

Answer 2

您可以直接从 pystan.stansummary 检索您想要的百分位数：

percentiles = (0.025, 0.25, 0.5, 0.75, 0.975)              # edit these at will
pystan.stansummary(fit=your_fit, probs=percentiles, digits_summary=2)

这应该可以正常工作。

Pystan 后验不确定区间

Pystan Posterior Uncertainty Intervals

python

optimization

confidence-interval

uncertainty

pystan