python - 使用 pymc 模拟从瓮中抽骰子

Question

(来自 https://stats.stackexchange.com/ 的交叉发帖，因为问题被标记为该网站的偏离主题)

尝试做http://greenteapress.com/thinkbayes/html/thinkbayes004.html#toc23在 pymc 中(本质上只是从瓮中拉出一个骰子并滚动几次)

import pymc as mc

die_faces = [4,6,8,12,20]
die_tosses = [x-1 for x in [6,8,7,7,5,4]] # 3 is a 4 on 0th indexed dice

which_die = mc.Categorical('which_die', [1./len(die_faces) for _ in range(len(die_faces))])

#imp = 0 # impossible
imp = 0.00000000000001  # impossible
die_probs = [[(1./face_count if i<face_count else imp) for i in range(max(die_faces))] 
            for face_count in die_faces]

@mc.deterministic
def die_results(which_die=which_die):
    return die_probs[which_die]

observation = mc.Categorical('obs', die_results, value=die_tosses, observed=True, size=len(die_tosses))

model = mc.Model([which_die, die_results, observation])
map_ = mc.MAP( model )
map_.fit() #stores the fitted variables' values in foo.value
mcmc = mc.MCMC( model )
mcmc.sample( 50000, 25000 )


die_data = mcmc.trace('which_die')[:]

for i in range(len(die_faces)):
    print 1.*list(die_data).count(i)/len(die_data)

这段代码似乎基本上得到了预期的答案(与已知的确切答案略有不同)。

问题在于，使用 0 值来表示不可能似乎更符合逻辑/正确。当我使用它时，我收到以下错误(间歇性):

Traceback (most recent call last):
  File "manydice.py", line 19, in <module>
    observation = mc.Categorical('obs', die_results, value=die_tosses, observed=True, size=len(die_tosses))
  File "/usr/local/lib/python2.7/dist-packages/pymc/distributions.py", line 3179, in __init__
    verbose=verbose, **kwds)
  File "/usr/local/lib/python2.7/dist-packages/pymc/PyMCObjects.py", line 772, in __init__
    if not isinstance(self.logp, float):
  File "/usr/local/lib/python2.7/dist-packages/pymc/PyMCObjects.py", line 929, in get_logp
    raise ZeroProbability(self.errmsg)
pymc.Node.ZeroProbability: Stochastic obs's value is outside its support, or it forbids its parents' current values.

是否有更惯用和/或更正确的方法来表示这个问题？

Answer 1

发布了类似的问题here其中解决方案是使用多项式而不是分类变量。

还有一个similar thread on Google Groups解决方案正是您所做的:使概率非常接近但等于零。

最后，在处理概率时，通常无法得到零概率的东西，因此我认为没有比将概率设置为非常接近 0 更好的方法来解决该问题。请记住概率无论如何，这么小实际上是“不可能的”。

但是，这个问题并没有明确的答案。如果您感兴趣，这里有一些关于 zero probability or impossibilities 意味着什么的冗长对话和哲学辩论。在数学网站上。 (另一个类似的问题:0除以0的值是无穷大还是零？)