Dirichlet-Multinomial Model does not appear to be sampling at all

Questions
4

Actually, upon a closer look it is actually a bug, if you print the model test_point:

with pm.Model() as mod:
    theta=pm.Dirichlet('theta',beta/sparsity)
    transition=pm.Multinomial('transition',1,theta,shape=4)
mod.test_point

which returns:

Out[14]: 
{'theta_stickbreaking__': array([ 0.,  0.,  0.]),
 'transition': array([0, 0, 0, 0])}

The transition is all zeros, which is wrong. I will report the bug to Github, meanwhile you can bypass this error by providing a valid start value. For example you can try the code below:

import numpy as np
import pymc3 as pm
import scipy.stats as st

sparsity = 1 #not zero
beta = np.ones([3,4])/4 #input for dirichlet
n = 10
testval = np.asarray([st.multinomial.rvs(p=a, n=n) for a in beta])
with pm.Model() as mod:
    theta=pm.Dirichlet('theta', beta/sparsity, shape=(3, 4))   
    transition=pm.Multinomial('transition', n, theta, shape=(3, 4), testval=testval)
    trace = pm.sample(5000, step=pm.Metropolis())
    
pm.traceplot(trace);
trace['transition']

Keep in mind that Metropolis acceptance is still quite low in this situation. And I observed that if the n in Multinominal is 1 the trace hardly move when the sparsity is high.

1 Like