stan-dev
stan_nlmer formula and fixed effects
Summary:
Problem in specifying fixed effects in stan_nlmer formula.
Description:
Following the documentation of stan_nlmer the formula needs to be in the same form as the nlmer one: a three-part “nonlinear mixed model” formula, of the form resp ~ Nonlin(...) ~ fixed + random (from the nlmer help page).
I can't find a way to specify the fixed effect variable. When I do it, after a long list of rejected values I get:
Chain 1:
Chain 1: Initialization between (-2, 2) failed after 100 attempts.
Chain 1: Try specifying initial values, reducing ranges of constrained values, or reparameterizing the model.
[1] "Error in sampler$call_sampler(args_list[[i]]) : Initialization failed."
error occurred during calling the sampler; sampling not done
Error in check_stanfit(stanfit) :
Invalid stanfit object produced please report bug
Also if I try to specify my own function with the suggested deriv option (check example 3 of the nlmer help page), I'm not able to fit the model with stan_nlmer.
I suspect is a long standing problem of nlmer..
Reproducible Steps:
Here an example of the problem with the code you have in your help page
data("Orange", package = "datasets")
Orange$circumference <- Orange$circumference / 100
Orange$age <- Orange$age / 100
#define a new to use as fixed effect:
Orange$fixed.eff <- ifelse(Orange$Tree=="3"|Orange$Tree=="5",1,2)
#try to include it following the formula definition:
fit <- stan_nlmer(
circumference ~ SSlogis(age, Asym, xmid, scal) ~ fixed.eff + (Asym|Tree),
data = Orange,
#for speed only
chains = 1,
iter = 1000
)
RStanARM Version:
packageVersion("rstanarm") [1] ‘2.21.2’ packageVersion("lme4") [1] ‘1.1.23’ packageVersion("rstan") [1] ‘2.21.2’
R Version:
[1] ‘4.0.2’
Operating System:
OS X 10.15.6
I think the problem here is that for some reason it wants (Asym|Tree) to be Asym|Tree (no parentheses). If I change that then this runs for me. @gioiadc Does that work for you too?
@bgoodri Is this supposed to work with the parentheses too? If so then we have a bug. If not then I can add something to catch this and throw a better error message. Right now it goes all the way through to Stan and the user just sees initialization errors.
Thanks @jgabry, it works without parentheses but then fixed.effect is not considered as a fixed effect, am I wrong?
@gioiadc You are right, sorry. I saw that it ran ok, but I didn't check that the results made sense!
@bgoodri Now I definitely think there's a bug. It doesn't let you use parentheses around the random effects but because of that (presumably) if you use a fixed effect it thinks it's part of the random part.
Printing @gioiadc's model (after removing parentheses so it runs) gives
Median MAD_SD
Asym 1.9 0.1
xmid 7.2 0.4
scal 3.5 0.4
Auxiliary parameter(s):
Median MAD_SD
sigma 0.1 0.0
Error terms:
Groups Name Std.Dev. Corr
Tree fixed.eff 0.116
Asym 0.343 -0.29
Residual 0.087
Num. levels: Tree 5
which shows that the predictor fixed.eff is being treated as part of the 'random' part not 'fixed' part. (It's not just the print output either, it's really estimating fixed.eff as if it varies by level of Tree.)
@bgoodri Any idea what's going on here? I thought we were using lme4's own functions to parse the formula, in which case I'm not sure how this would happen. I'm not very familiar with stan_nlmer() but if you don't have time to look into this let me know and I can probably spend a little while on it.
If the formula parses in lme4::nlmer(), then it is a bug in rstanarm::stan_nlmer().
On Mon, Aug 10, 2020 at 6:03 PM Jonah Gabry [email protected] wrote:
@gioiadc https://github.com/gioiadc You are right, sorry. I saw that it ran ok, but I didn't check that the results made sense!
@bgoodri https://github.com/bgoodri Now I definitely think there's a bug. It doesn't let you use parentheses around the random effects but because of that (presumably) if you use a fixed effect it thinks it's part of the random part.
Printing @gioiadc https://github.com/gioiadc's model (after removing parentheses so it runs) gives
Median MAD_SDAsym 1.9 0.1 xmid 7.2 0.4 scal 3.5 0.4
Auxiliary parameter(s): Median MAD_SD sigma 0.1 0.0
Error terms: Groups Name Std.Dev. Corr Tree fixed.eff 0.116 Asym 0.343 -0.29 Residual 0.087 Num. levels: Tree 5
which shows that the predictor fixed.eff is being treated as part of the 'random' part not 'fixed' part. (It's not just the print output either, it's really estimating fixed.eff as if it varies by level of Tree.)
@bgoodri https://github.com/bgoodri Any idea what's going on here? I thought we were using lme4's own functions to parse the formula, in which case I'm not sure how this would happen. I'm not very familiar with stan_nlmer() but if you don't have time to look into this let me know and I can probably spend a little while on it.
— You are receiving this because you were mentioned. Reply to this email directly, view it on GitHub https://github.com/stan-dev/rstanarm/issues/463#issuecomment-671612604, or unsubscribe https://github.com/notifications/unsubscribe-auth/AAZ2XKRONKMCTWQANWPY653SABVETANCNFSM4P2KXEFA .
Definitely a bug, but I think maybe in the Stan program, not the R code. Using formula
circumference ~ SSlogis(age, Asym, xmid, scal) ~ X + (Asym|Tree)
we get an initialization failure plus this exception from Stan:
Chain 1: Exception: Exception: x is the wrong length (in '/functions/SSfunctions.stan' at line 142; included from 'model_continuous' at line 8)
(in 'model_continuous' at line 178)
So I think it's the Stan program that can't handle so-called fixed-effects terms.
Line 178 in the Stan program is
matrix[len_y, K] P = reshape_vec(eta, len_y, K);
so there's some dimension mismatch or something coming from reshape_vec from functions/SSfunctions.stan. According to that function rows(eta) is supposed to equal len_y * K but apparently doesn't.
@bgoodri Can you take a look if you have a chance? I'm not at all familiar with the Stan code for the nlmer stuff, so you can probably figure this out much quicker.