Dynamic number of Dirichlet distributions?

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…

It should work, have you try specifying the shape of b?

1 Like

I could have sworn I tried that… works fine once the shape is specified!

Once the shape pull request is settled, I might try to make sure that is detected automatically…

Shape issue is one of our daily struggle :joy: