How to implement a GLM with a personalized link function different from the Identity function?

IMO, your best bet is to just implement the inverse of your link function using Theano tensor operations. You were on the right track earlier, but you don’t need to actually define a tensor and a Theano function. For your case, just writing

   likelihood = Normal('y', mu=1/(1+tt.exp( - (intercept + x_coeff * x) ) ), 
                        sd=sigma, observed=y)

should work. Not sure, but in your earlier attempt either you were overloading the variable x or PyMC3 wasn’t playing nice with the tensor you defined outside the model.

The above should solve your problem, but I’ll also address your approach using DensityDist. You got that error since givens is not a keyword argument to DensityDist. I think you should be able to just swap givens with observed to get it to work. See this thread for some examples. However, think twice before you do this! If you’re calling scipy.special.expit() when you defined logl(), PyMC3 won’t be able to use the NUTS sampler since it can’t “see” the gradient of an outside black-box function and will fall back on using Metropolis. However, if you use Theano ops when defining logl() then I think it can get the gradients and you wouldn’t have this issue.

In either case, your link function is simple enough that you really just should use the one-liner

likelihood = Normal('y', mu=1/(1+tt.exp( - (intercept + x_coeff * x) ) ), sd=sigma, observed=y)

instead of a DensityDist.

2 Likes