I am developing a hierarchical model that fits a non-linear regression to individuals nested within species. One of the parameters in the regression controls the convexity of the curve. If I estimate a single convexity value shared across all individuals and all species, the model runs and converges nicely. If I try to estimate convexity values for each individual (Since it’s just a tuning parameter, I don’t bother with nesting them within species), the model fails to initialize.
It’s not a model specification problem, because if I specify inits for this parameter (based loosely on the posterior mean and standard deviation of the parameter estimated from analysis of individual curves), the model runs and converges nicely. To generate the inits, though, I have to use a bit of a hack.
Specifically, I am working with RStan, and I set up an environment to hold the number of individuals:
## define new environment ## init <- new.env() ## set n_indiv to default (illegal) value ## n_indiv <- -999.0 ## return to parent environment ## parent.env(init)
Then just before the call to stan() when I know the number of individuals in the data set I set n_indiv:
## set n_indiv for use in init_theta ## init$n_indiv <- nrow(unique_table)
That all works just fine, but it would be a lot cleaner if I could just specify
init = init_theta(n_indiv)
(or something similar) in the call to stan(). Am I missing something?