将不确定性纳入 pymc3 模型

Question

我有一组数据，每个点都有平均值、标准差和观察次数（即，我知道测量的准确性）。在我只看手段的传统 pymc3 模型中，我可能会按照以下方式做一些事情：

x = data['mean']

with pm.Model() as m:
    a = pm.Normal('a', mu=0, sd=1)
    b = pm.Normal('b', mu=1, sd=1)
    y = a + b*x

    eps= pm.HalfNormal('eps', sd=1)
    likelihood = pm.Normal('likelihood', mu=y, sd=eps, observed=x)

将有关观测值方差的信息合并到模型中的最佳方法是什么？显然，结果应该比高方差（不太确定）的观察更重视低方差观察。

统计学家建议的一种方法是执行以下操作：

x = data['mean']   # mean of observation
x_sd = data['sd']  # sd of observation
x_n = data['n']    # of measures for observation
x_sem = x_sd/np.sqrt(x_n)

with pm.Model() as m:
    a = pm.Normal('a', mu=0, sd=1)
    b = pm.Normal('b', mu=1, sd=1)
    y = a + b*x

    eps = pm.HalfNormal('eps', sd=1)
    obs = mc.Normal('obs', mu=x, sd=x_sem, shape=len(x))
    likelihood = pm.Normal('likelihood', mu=y, eps=eps, observed=obs)

然而，当我运行这个时，我得到：

TypeError: observed needs to be data but got: <class 'pymc3.model.FreeRV'>

我正在运行ning pymc3 的 master 分支（3.0 有一些性能问题导致采样时间非常慢）。

Answer 1

你已经很接近了，你只需要做一些小的改变。主要原因是对于 PyMC3 数据总是不变的。检查以下代码：

with pm.Model() as m:
    a = pm.Normal('a', mu=0, sd=1)
    b = pm.Normal('b', mu=1, sd=1)
    mu = a + b*x

    mu_est = pm.Normal('mu_est', mu, x_sem, shape=len(x))

    likelihood = pm.Normal('likelihood', mu=mu_est, sd=x_sd, observed=x)

请注意，我保持数据不变，并在两点引入观察到的不确定性：mu_est 的估计和可能性。当然，您可以不使用 x_sem 或 x_sd 而是估计它们，就像您在代码中使用变量 eps 所做的那样。

根据历史记录，带有 "random data" 的代码曾经在 PyMC3 上工作（至少对于某些模型），但考虑到它并不是真正设计成那样工作的，开发人员决定阻止用户使用随机数据，这就解释了您收到的消息。

将不确定性纳入 pymc3 模型

Incorporating uncertainty into a pymc3 model

statistics

pymc3