[R-sig-ME] extract level 2 residuals of merMod from lme() and test for spatial autocorrelation.

[Dexter Locke]

Thanks very much.

Yes, that was a typo. I'm working with lme4 and the lmer function to make
merMod objects. Sorry for any confusion.

Here is a reproducible example of my attempt to extract residuals at the
population level (ie random effects) by subtracting predictions from
observations.

library(lme4)
fm1 <- lmer(Reaction ~ Days + (1 | Subject), sleepstudy) # load library,
fit merMod object

test <- sleepstudy$Reaction - (predict(fm1, re.form= ~ (1| Subject))) #
subtract predictions from observations

summary(fm1) # examine residuals of random effects, which provides:
# "Residual              960.5   30.99 " for variance and std dev,
respectively
var(test); sd(test) # which provides 869.8171 and 29.49266

# test <- sleepstudy$Reaction - (predict(fm1, re.form= NULL)) #
unsurprisingly has same results as above

960.5 <> 869.8171 and 30.99 <> 29.49266 although they are the same
direction and order of magnitude. I'm finding similar results with my
actual data and model fits.

Any advice is appreciated, thanks again!
Dexter

On Sat, May 7, 2016 at 2:03 PM, Ben Bolker <bbolker at gmail.com> wrote:

> Dexter Locke <dexter.locke at ...> writes:
>
> >
> > Apologies for cross-posting, another version of this was posted on
> > r-sig-geo at ...
> >
> > I want to extract the level 2 residuals from a merMod object created with
> > lme() and am struggling.
>
>   Not sure what you mean here -- this is a little bit inconsistent.
> lmer (from lme4) makes merMod objects, lme (from nlme) makes lme objects.
> I'm going to assume you're talking about lmer (since that's a more
> likely typo).
>
> > How can I pull out a vector that is the residuals
> > at level 2? Everything I see on stack exchange and the mixed models list
> > serve points towards methods like VarCorr(). I'm not interested in a
> > summary of the L2 residuals so much as a vector of the same length as the
> > number of groups I have. Ultimately I want to test them for spatial
> > autocorrelation. Apologies if this is documented somewhere, but I have
> not
> > been able to find help online.
>
>  I always get confused about what levels mean in the context of
> mixed models (i.e. whether one counts down from the population level
> or up from the individual level), but the basic answer here is to
> get the *predictions* at the desired level and subtract them from
> the observed values, see ?predict.merMod -- the re.form argument
> is what you need.  re.form=NA or re.form=~0 gives you predictions
> at the population level (i.e. neglecting all random effects),
> re.form=NULL gives you predictions at the individual level (i.e.
> including all random effects), intermediate cases can be gotten
> by using a formula.
>
> The results will be ordered the same as the original observation
> vector.
>
>
> >
> > While I have found methods of incorporating spatial effects into a mixed
> > model using corStruct, I am interested in first evaluating if that is an
> > appropriate model for a given dataset by examining the level 2 residuals'
> > spatial patterning - or lack thereof.
>
>   Seems reasonable.
>
> >
> > Which slot contains the residuals for level 2 and how are they ordered?
> >
>
> _______________________________________________
> R-sig-mixed-models at r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models
>

	[[alternative HTML version deleted]]