4

在另一个论坛上看到 PyStan 没有与他们使用的 RStan 相同的功能posterior_interval(),但我们可以使用numpy.percentile()。我目前正在使用pystan.StanModel.optimizing()PyStan 中的函数来获取最大化后验似然的参数集。我现在还想获得后验结果的外部 95% 置信区间,所以我想知道该numpy.percentile()函数是否会与该optimizing函数一起使用?

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

我认为 95% 的间隔必须取自后部。我会使用百分位数函数和优化函数来获得 95% 置信度的后验结果吗?

4

2 回答 2

3

要获得后验估计的界限,需要对后验进行采样,但这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您提到的。

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

于 2018-12-11T23:46:18.173 回答
3

您可以直接从以下位置检索所需的百分位数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)

这应该可以正常工作。

于 2019-01-21T13:20:12.267 回答