extract level 2 residuals of merMod from lme() and test for spatial autocorrelation.

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?