Need help modeling a Binomial likelihood with custom p

Questions
1

Hi, I am new to PyMC3 and need some help modeling a custom likelihood function. What I am trying to model is essentially the likelihood of getting zj failures out of a total of nj buildings:
image
But here we replace the p in the binomial with the function of a vulnerability curve. The log likelihood, then, is:
image
The goal is to try to estimate the parameters theta and beta using observations of xj, nj, and zj.

It seems to me that I probably need to define a custom likelihood function using DensityDist (i.e., hard code the entire expression for the likelihood function) for this but I am wondering if there is a more elegant way to do this using the pm.Binomial likelihood? For example, I would like to do something like:

image

but am not sure how I would model the input (custom) p. This is the standard normal CDF.

2

This should work as written, if you implement my_func. I would suggest a logit link instead of a probit link, though, which tends to give answers that are easier to interpret, and I’d reparametrize – the equation should be easier to read/interpret if it’s just written as a typical logistic regression, i.e.
p = pm.invlogit(beta + theta * log(xj))
In that case, β is easy to interpret: It’s the probability that an average building will fail, assuming log(xj) has been standardized (you subtracted out the mean from each observation). That way, it’s easier to set priors for.