I am trying to define a model that looks like
nvars, ncats = 3, 5
with pm.Model():
α = pm.Exponential('α', lam=λ, shape=(nvars, ncats))
b = pm.Dirichlet('b', a=α)
trace = pm.sample()
This fails, because the Dirichlet distribution expects an array, not a matrix, in PyMC3
(and in numpy
, and in scipy
).
Has anyone had any success in modelling something like this? I guess I can use loops…