I’m trying to to estimate a structural model with a simple Bernoulli likelihood with a probability depending on the solution of an implicit function with no closed form solution. Therefore, I’m trying algebra_solver() for the first time.
As I learn about how doable this is I’m starting with some very simple fake data and I simplified the implicit function I need to solve to something I can actually solve without algebra_solver() (I dropped the nonlinear component). In other words, I can fit the model fine without using algebra_solver(), so I tried to fit the model with and without it. In fact, when I tried calling algebra_solver() I ignored the results and used the closed form solution I can calculate. I expected both runs to produce the same outcome (perhaps a lot slower when I use algebra_solver()). What is surprising is that when I comment out the line where algebra_solver() is called the model is estimated nicely (no divergent transitions, R hat < 1.1 and parameters close to the true value), but when I leave the algebra_solver() call in there (even though I am not actually using its output) the model fit has problems (30% divergent transitions, 25% of parameters have R hat > 1.1, but the parameters are still well estimated).
Perhaps I don’t understand how algebra_solver() works, but is it normal for it to have a side effect on the MCMC estimation when I don’t even update the log probability with its output?
I’m attaching my .stan file. The statement at line 133 is the one that changes things.
structural_dist_model.stan (5.9 KB)
I appreciate any help with this!