random slopes model specification
On Wed, Apr 18, 2018 at 3:04 PM, Karista Hudelson <karistaeh at gmail.com> wrote:
Hello Mixed Modelers, I was hoping to get a bit of feedback on the random effects structure of the model below. I wanted to allow two of the fixed effects (my effects of interest) to have random slopes and also random intercepts for each watershed (there are 5 watersheds, 27 observations of Spring_MST & Summer_Rain, and 790 fish (so there are 790 measurements of Hg and FishLength). model6<-lmer(Hg~FishLength+Spring_MST+Summer_Rain+(1+Spring_MST+Summer_Rain|watershed),data=FSV,REML=FALSE) this model converges, and the output looks "like I want it to look"; ie, slopes are different for the two effects which logically should have different slopes:
coef(model6)
$watershed
(Intercept) FishLength Spring_MST Summer_Rain
877911 -1.7734565 0.02147013 0.03245297 -0.04841394
3319324 -1.2429440 0.02147013 0.03037712 -0.04472220
11597701 -1.6688121 0.02147013 0.03204351 -0.04768574
31121323 -0.9214729 0.02147013 0.02911923 -0.04248513
97032352 -1.4819562 0.02147013 0.03131236 -0.04638544
Notice that there is not much variability in your coefficients across watersheds (zero variability in FishLength effects because you left it out of your random effects on purpose, about 10% variability in Spring_MST and Summer_Rain). Normally it's not recommended to fit a model with three correlated RE terms (intercept, spring_MST, summer_rain) to only a grouping factor with only 5 levels -- this is equivalent to trying to estimate a 3x3 variance-covariance matrix with 5 observations (I'm a little surprised this worked at all). Are all the values in getME(fitted_model,"theta") 'reasonable', i.e. larger than (say) 1e-4? There is a huge, ongoing debate about what to do in this case, where the RE model you want is theoretically identifiable (you do have multiple measurements of spring_MST and summer_RAIN in each watershed), but practically very poorly constrained, suggestions ranging from "keep it maximal" (Barr et al): fit the full model, reduce when the model is singular try to constrain complexity _a priori_, reduce by model selection (AIC or p-value based) use Bayesian methods to add priors/regularize the problem substitute fixed effects instead.
Here is a bit of the summary() for model6:
Fixed effects:
Estimate Std. Error df t value Pr(>|t|)
(Intercept) -1.417728 0.144246 6.000000 -9.829 6.18e-05 ***
Length 0.021470 0.001085 783.800000 19.794 < 2e-16 ***
Spring_MST 0.031061 0.007518 312.000000 4.132 4.63e-05 ***
Summer_Rain -0.045938 0.007276 132.100000 -6.314 3.82e-09 ***
It seems to be doing what I intended, but the random effects structure
looks like a bit of a run-on sentence compared to most of the other models
I've seen. Does this appear correct?
Also....if I'm allowing for the slopes of the effects of interest to be
different for each watershed and allowing intercepts to be different, in
your opinion(s) would it be easier to understand if I just run separate
models for each lake? Here how that looks for the Lake with watershed 877911
(corresponding to top line of coef() table above):
clim<-Hg~FishLength+Spring_MST+Summer_Rain climS<-lm(clim,data=Smid) summary(climS)
Estimate Std. Error t value Pr(>|t|)
(Intercept) -1.3531180 0.2152754 -6.286 2.19e-09 ***
Length 0.0276015 0.0027818 9.922 < 2e-16 ***
Spring_MST 0.0247300 0.0081347 3.040 0.002699 **
Summer_Rain -0.0025433 0.0007216 -3.525 0.000531 ***
Any thoughts are appreciated. Thanks in advance for your time and brain
power.
Fitting separate models will be almost equivalent to fitting a single Hg~(FishLength+Spring_MST+Summer_Rain):watershed model -- the only difference will be whether the residual variance is pooled/assumed to be the same across watersheds.
Much Obliged,
Karista Hudelson
[[alternative HTML version deleted]]
_______________________________________________ R-sig-mixed-models at r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models