Binomial probit regression, sample code

Sharing
2

Thanks for sharing!

You can actually skip the custom probit_phi and use pm.math.invprobit, as in:

N_ar=[1000,10]
k_ar=[500,0]

with pm.Model() as model:
    theta_p = pm.Normal('theta_p', mu=0, sigma =100)
    theta = pm.math.invprobit(theta_p)

    # likelihood
    y = pm.Binomial('y', p=theta, n=N_ar, observed= k_ar)
    idata = pm.sample()
1 Like