Help with Poisson model

In another thread I wrote Recommendations for what to do when k exceeds 0.5 in the loo package? - #11 by avehtari
that I think this is a problem with initialization. I copy that comment here and add a bit more in the end.

see ?stan and option seed. Try runinng stan(...., seed=1337) with different seed values.

see ?stan and option init. The default is to randomly generate initial values between -2 and 2 on the unconstrained support. If you transform -2 and 2 from unconstrained to constrained space, are these initial values sensible?

If you are unlucky initial values are all 2, and if n=1 and x=1, then
log(n) + x*2 + 2 + exp(2)*2 + epx(2)*2*4
is about 79 which when exponentiated to Poisson parameter is about 2e34, which is a lot. You probably need to initialize log(a_std) and log(b_std) with much smaller range than [-2, 2] or scale them in Stan code by dividing them with a big number. For example, try dividing ‘a’ and ‘b’ by 100. Then maximum Poisson parameter value (with n=1 and x=1) is about 311, which should work fine.

Since initial values are randomly chosen, one chain can work, but running several can fail. Adding _rng to generated quantities changes the behavior as _rng is called between sampling iterations and it changes the state of the random number generator.

I recommend to test with different seeds, and with more specific inits or scaled a and b. For checking the effect of reasnoable inits/scaling ,you could add a print statement in Stan code which would print the value of
log(n) + x_beta + gamma[index_time] + a[practice] + b[practice] .* time

I’ve seen this same problem before (I just didn’t realize that this one is likely to be the same), and the common thing is a distribution with positive constrained parameter which is a function of positive constrained parameters so that effectively there is something like exp(a + b*exp(c)) and if b and c are both initialized to be close to 2, then that expression has very high value.

2 Likes