Thanks - @ricardoV94’s solution worked a treat.
A quick follow-up question is: what is the best way to demand that the values in a CustomDist are bounded? I don’t want the sampler to waste time searching in zero-likelihood regions, so I’d like to impose a lower and upper bound on:
n = pm.CustomDist(‘n’, N, p, size = len(obs_n), logp = binom_logp, random = binom_draw, observed = obs_n)
so that the imputed missing values are sampled efficiently. Many of the inbuilt distributions have a ‘lower’ and ‘upper’ keyword argument; does CustomDist have the same?
Thanks!