移动间隔时 pymc 中 DiscreteUniform 的意外行为

Unexpected behaviour with DiscreteUniform in pymc when shifting interval

我正在尝试模拟 100 次掷骰子,其中我的 数据 是所有掷骰的总和(类似于 Jaynes 的 Brandeis 骰子的最大熵原理)。这是我后来第一次尝试接近已加载的骰子。

我正在使用 pymc 2.3

如果我用 DiscreteUniform('dice', 1, 6, size=N) 设置从 1 到 6 的骰子值,并设置一个总和值等于平均总和值 100*3.5=350,那么我得到均匀的后验分布,正如预期的那样。

但是如果我设置骰子值从0到5,总和等于100*2.5=250,分布不均匀。 0 值的采样率要低得多!因为我只是将值移动 1 个单位,所以我期望得到相同的结果。知道为什么它们不同吗?我做错了什么?

这是完整代码:

import numpy as np
import pymc as pm
import matplotlib.pyplot as plt

N = 100
shifts = (0, -1)
for shift in shifts:
    obs_mean = 3.5+shift
    obs_total = int(N*obs_mean)
    sigma = 0.01*N
    dice = pm.DiscreteUniform('dice', 1+shift, 6+shift, size=N)

    @pm.deterministic
    def calc_total(d=dice):
        return np.sum(d)

    total = pm.Normal('total', mu=calc_total, tau=1./sigma, observed=True, value=obs_total)

    # package the full model in a dictionary
    model1 = dict(dice=dice, calc_total=calc_total, total=total)

    # run the basic MCMC:
    S = pm.MCMC(model1)
    S.sample(iter=100000, burn=10000)

    dice_trace = S.trace('dice')[:]-shift
    plt.hist(dice_trace.flat, bins=(0.5, 1.5, 2.5, 3.5, 4.5, 5.5, 6.5), normed=True, alpha=0.5)

plt.show()

编辑:根据评论,我做了一个更简单的模型:两个均匀分布,一个从 1 到 6,另一个从 0 到 5,然后是确定性函数 dice2 将其加 1,因此两个模型中的先验 dice2 相同,并且可能性仅取决于 dice2,但是它们的后验分布不同。

另一个有趣的情况是当 shift 设置为 -7 时,结果只是反转骰子的符号,但结果不同。

import numpy as np
import pymc as pm
import matplotlib.pyplot as plt

N = 100
shifts = (0, -1)
for shift in shifts:
    obs_mean = 3.5
    obs_total = int(N*obs_mean)
    sigma = 0.01*N
    dice = pm.DiscreteUniform('dice', 1+shift, 6+shift, size=N)

    @pm.deterministic
    def dice2(d=dice):
        return d-shift

    @pm.deterministic
    def calc_total(d=dice2):
        return np.sum(d)

    total = pm.Normal('total', mu=calc_total, tau=1./sigma, observed=True, value=obs_total)

    # package the full model in a dictionary
    model1 = dict(dice=dice, dice2=dice2, calc_total=calc_total, total=total)

    # run the basic MCMC:
    S = pm.MCMC(model1)
    S.sample(iter=100000, burn=10000)

    dice_trace = S.trace('dice2')[:]
    plt.hist(dice_trace.flat, bins=(0.5, 1.5, 2.5, 3.5, 4.5, 5.5, 6.5), normed=True, alpha=0.5)

plt.show()

不清楚为什么您一定会期望均匀分布。离散制服只是你的先验。该模型的所有信息是 shift=0 的 350 和 shift=-1 的 250 的总和,并将生成具有此期望的参数估计。当我 运行 每个偏移值下的模型并查看轨迹时,我得到 shift=0 的以下分布(仅按唯一值查看摘要):

>>> pd.Series(shift_0.flat).value_counts()

2    1526136
4    1526011
3    1511494
5    1503698
6    1471922
1    1460739

对应于以下期望:

>>> pd.Series(shift_0.flat).mean() * 100

350.02311111111112

shift=-1

>>> pd.Series(shift_1.flat).value_counts()

1    1894489
2    1724072
3    1577420
4    1457896
5    1320425
0    1025698
dtype: int64


>>> pd.Series(shift_1.flat).mean() * 100

250.08703333333332

因此,该模型的行为似乎符合我的预期。

这是 PyMC 2.3.6 版采样步骤的问题。它在版本 2.3.2 中按预期工作。我在 github 中与 Chris Fonnesbeck 讨论了这个问题,他在 PyMC 的开发版本中更正了它。将来,请检查您拥有的版本及其行为方式。