How would you go about marginalizing over discrete parameters in PyMC?

In STAN you cannot sample from discrete variables. You have to marginalize over them with a custom for loop, and then you can draw samples from the posterior by calling the rng function. From what I understand this is a mix of engine limitation as well as intentional design, as samplers for discrete variables are not very efficient / robust.

Is there any easy way to do the same with PyMC? Adding the probability term seems simple with a pm.Potential and a for loop (although I have not tried any example). But I don’t see how to generate back samples from the discrete parameters? Can we somehow use the pm.Deterministic for this?

The recommended approach is still rewrite your model into a marginalized version, usually in a form of a mixture model (Frequently Asked Questions) There are many automatic marginalization effort going on that detects the discrete variable and marginalized them, but we dont have any solid implementation yet in PyMC3