How to determine random effects in mixed model

I want to model frequency of a social effect by state over several years as potentially predicted by several demographic factors. I have been reading that one is supposed to determine the random effect structure before trying to determine the fixed effects structure. One thing I want to be sure to test for fixed effects is an interaction with time for each fixed effect.

Do I want to just have intercept|State or do I want to look at slope of each possible fixed factor per state, as well. How would I do the latter if I have not yet built the fixed effects portion of the model? Would I also look at fixed:time for each potential fixed effect as random effects?

I want to avoid dredging. I have only 459 points for each variable (9 years, 50 states + DC). However, I also want to be able to contrast model predictions for the full model vs. a fixed-effect only model.

• REs depend on FEs and FEs on REs. They both need to be recalculated to conduct nested tests. The best way to fit a model is simply to state what the important variables are and defend them for scientific reasons, not statistical ones. This is confirmatory data analysis. If you fit multiple models and compare them using ICs or other criteria, then it is data dredging if you find a "good" model and claim you thought of it without multiple testing. – AdamO Apr 28 '17 at 21:13
• I am doing exploratory modeling and intend to be honest about that. – Bryan Apr 29 '17 at 18:46

You can test if the variance in slopes (and covariance between slope and intercept) is significant by modeling one model with just the random intercept and another model with the random slope and intercept. Then you can do a nested model comparison between the two:

mod1 <- lmer(... + (1|state), ...)
mod2 <- lmer(... + (1+predictor|state), ...)
anova(mod1,mod2)


This will give you a p-value that you need to correct, though. You can do so with this code:

1-(.5*pchisq(anova(mod1,mod2, refit=FALSE)$Chisq[[2]],df=2)+ .5*pchisq(anova(mod1,mod2, refit=FALSE)$Chisq[[2]],df=1))


P-value correction is found here or here.

• Isn't that just "hunting p"? – Bryan May 1 '17 at 11:58
• Not at all. It's not selectively dropping cases or adding covariates or performing multiple comparisons to attain p < .05; it's simply performing a significance test for the variance of slopes. – Mark White May 1 '17 at 23:01
• Can you elaborate on what's the problem with anova call, why the p-value needs to be corrected, and what is the logic behind this correction? – amoeba says Reinstate Monica May 31 '17 at 13:44
• Also, your correction code uses mod3 and mod4 but before you have only defined mod1 and mod2. – amoeba says Reinstate Monica May 31 '17 at 13:44
• I fixed the error to say mod1 and mod2. Sorry, copy-and-paste error. The model comparison is testing two parameters at once: First, the variance of the slope of predictor. Since this is bounded at > 0, a two-tailed test is unnecessary since it literally cannot be below zero. Second, the covariance between the slopes and intercepts is estimated. This is can be a negative or positive number. Because one estimate is bounded at zero and the other isn't, people recommend using a mixture chi-square distribution to test them both simultaneously. – Mark White May 31 '17 at 13:49

Exploratory analyses for describing correlation structures in dependent data include variograms for continuous spatio-temporal data, intraclass correlation coefficients for clustered data, lorelograms for binary outcomes.

Other descriptive statistics include bootstrapped or profile likelihood confidence intervals for variance components in random effects models, goodness-of-fit tests with unstructured covariance structures or saturated specification of random and fixed effects.

Formal inference for the hypothesis of random effects having 0 dispersion can be done testing nested models with likelihood ratio tests. Be sure to change software settings to fit models with maximum likelihood and not REML when conducting these tests.

• What about "informal comparison" with information criteria? Do those only work for nested models? – Bryan May 3 '17 at 12:43
• @Bryan I don't know why you would use ICs if you have nested fully parametric models. The problem with informal comparisons is that "a little bit better" performance with overadjustment may overcome the degree-of-freedom penalty but not achieve statistical significance. Nonetheless, it is an approach you can take. – AdamO May 3 '17 at 15:58
• They are not nested models. Where did I say they were nested? For example, I have two alternate measures of "income inequality", but it would be absurd to have a model with both measures. – Bryan May 3 '17 at 19:51
• @Bryan we can't divinate the actual regression coefficients. I hope they have clear cut purpose for your analyses. If you start with a base model of fixed effects, any additional random effects like a random intercept will yield a nested model. – AdamO May 3 '17 at 20:00
• I never saw lorelogram before. Is it this? biostat.jhsph.edu/~fdominic/teaching/bio655/references/extra/… – kjetil b halvorsen Jun 9 '19 at 9:18