[R-sig-ME] random slopes model specification

Karista Hudelson karistaeh at gmail.com
Thu Apr 19 17:13:25 CEST 2018


Greetings again Modelers,

Thanks for the reply Dr. Bolker.  You brought up some excellent points: the
theta values were ~ 0; I was trying to "keep it maximal" (Schielzeth &
Forstmeier 2008); I was also surprised this converged/ "worked".  Your
point about the residual pooled variance vs allowing the variance to be
different speaks directly about my hypothesis that "one model should rule
them all", if climate is truly driving variance, it should influence all
the lakes, but in fact the models I've made thus far do not support this
very well.  Also noted from your reply: I have too many fixed effects and
not enough variance to explain with them...yes.

Nonetheless, I have one follow up question, I would very much appreciate an
experienced point of view:
If I'm mostly interested in the effect of the two climate effects (and only
including watershed and FishLength because they're influential), then
should I try to allow the slopes to vary, but the intercept, which is set
mainly due to watershed and FishLength, not to vary (because the intercept
is not informative for my hypothesis, it's more about the slopes, I
think)?

For example:

model7<-lmer(Hg~FishLength+Spring_MST+Summer_Rain+(Spring_MST+Summer_Rain-1|watershed),data=FSV,REML=FALSE)

> coef( model7 )
$watershed
         (Intercept)  FishLength    Spring_MST Summer_Rain
877911     -2.111044 0.03653274 -0.0001905891 -0.04581202
3319324    -2.111044 0.03653274  0.0446875053 -0.03808708
11597701   -2.111044 0.03653274  0.0337457453 -0.03997050
31121323   -2.111044 0.03653274 -0.2581203000 -0.09020993
97032352   -2.111044 0.03653274  0.0492101881 -0.03730858


> getME(model7,"theta")
            WA.Spring_MST WA.Summer_Rain.Spring_MST
*WA.Summer_Rain *
               0.37626815                0.06476772               *
0.00000000  #:[*

from summary(model7):

Random effects:
 Groups   Name        Variance  Std.Dev. Corr
 wat...   Spring_MST  0.0144248 0.12010
          Summer_Rain 0.0004274 0.02067  1.00
 Residual             0.1018861 0.31920
Number of obs: 790, groups:  WA, 5

Fixed effects:
              Estimate Std. Error         df t value Pr(>|t|)
(Intercept)  -2.111044   0.070860 786.000000 -29.792   <2e-16 ***
FishLength    0.036533   0.001764 785.500000  20.711   <2e-16 ***
Spring_MST   -0.026133   0.055546   5.100000  -0.470   0.6575
Summer_Rain  -0.050278   0.015880   8.600000  -3.166   0.0121 *
---

Thanks modelers for any advice...except it the advice is along the lines of
square peg, round hole, don't force it.  I see that too.

Karista


On Wed, Apr 18, 2018 at 3:54 PM, Ben Bolker <bbolker at gmail.com> wrote:

> 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
>



-- 
Karista

	[[alternative HTML version deleted]]



More information about the R-sig-mixed-models mailing list