PyMC3版本差异致az.summary输出HPD/HDI命名不一致咨询
问题描述
我尝试复现多个使用PyMC3实现的示例并对比运行结果,以下为用于估计HPD的示例代码:
import pymc3 import arviz as az import numpy as np import warnings warnings.filterwarnings("ignore") import matplotlib.pyplot as plt import graphviz import pymc3 as pm from pymc3 import Model, Normal, HalfNormal from pymc3 import find_MAP basic_model = Model() with basic_model: # Priors for unknown model parameters alpha = Normal('alpha', mu=0, sd=5) beta = Normal('beta', mu=0, sd=5, shape=2) sigma = HalfNormal('sigma', sd=4) # Expected value of outcome mu = alpha + beta[0]*X1 + beta[1]*X2 # Deterministic variable, to have PyMC3 store mu as a value in the trace use # mu = pm.Deterministic('mu', alpha + beta[0]*X1 + beta[1]*X2) # Likelihood (sampling distribution) of observations Y_obs = Normal('Y_obs', mu=mu, sd=sigma, observed=Y) pm.model_to_graphviz(basic_model) with basic_model: # obtain starting values via MAP start = find_MAP(fmin=optimize.fmin_powell) # instantiate sampler - not really a good practice step = NUTS(scaling=start) # draw 2000 posterior samples trace = sample(2000, step, start=start)
运行摘要统计的命令如下:
az.summary(trace)
示例网站的运行结果截图中,区间列名为hpd(即后验HDI):
本地使用PyMC3 3.11.5运行相同命令时,输出区间列名为HDI:
目前不确定两类输出的参数含义是否一致,想确认是命名规则发生了变更,还是新版本需要调整其他配置才能得到真实的HPD值。
解答
这是依赖库ArviZ的命名规则调整导致的,两者计算逻辑完全一致,不需要额外配置就能得到正确的区间结果:
- 版本对应关系:PyMC3 3.10.0依赖的ArviZ版本低于0.12,当时最高密度区间的缩写用的是HPD(最高后验密度,Highest Posterior Density);PyMC3 3.11.5依赖的ArviZ版本在0.12及以上,官方统一将该指标命名改为HDI(最高密度区间,Highest Density Interval),适配更通用的贝叶斯统计术语表述。
- 数值完全等价:两个命名对应的区间计算方法没有任何改动,都是返回覆盖指定比例后验样本的最窄连续区间,默认覆盖比例都是94%,本地输出的
hdi_3%、hdi_97%和旧示例里的hpd_3%、hpd_97%含义完全相同,不存在计算偏差。 - 如果需要和旧示例输出格式完全对齐,可以在调用
az.summary()时手动重命名列,没有其他配置需要调整。 - 额外注意:贴出的示例代码缺了几个导入和变量定义,包括
import scipy.optimize as optimize、from pymc3 import NUTS, sample以及预测变量X1、X2和观测值Y的定义,直接运行会触发NameError。
内容的提问来源于stack exchange,提问作者Niko Gamulin
相关产品推荐
相关产品推荐

