[R-sig-ME] Fwd: Re: extract specific values from lmer summary

[Ben Bolker]

 * yes, in general these sorts of question should go directly to
r-sig-mixed-models at r-project.org

 * it's mildly surprising but not shocking that this format doesn't work
any more.

 * without testing, here are a few idioms that should work better:


for(i in 1:nlevels(vital)) {
   lmer.sheep <- glmer(y ~ (1|year),
      family = binomial,
      data = subset(sheep.dat,vital==levels(vital)[i])
}

or

for(i in 1:nlevels(vital)) {
   lmer.sheep <- glmer(y ~ (1|year),
      family = binomial,
      data = sheep.dat,
      subset = vital==levels(vital)[i])
}

 or

model.list <- vector(nlevels(vital),mode="list")
model.list[[1]] <- glmer(y ~ (1|year),
      family = binomial,
      data = sheep.dat,
      subset = vital==levels(vital)[1])
for (i in 2:nlevels(vital)) {
  model.list[[i]] <- update(model.list[[1]],
       subset = vital==levels(vital)[i]
}

  etc.





-------- Original Message --------
Subject: Re: extract specific values from lmer summary
Date: Sun, 04 May 2014 15:25:42 -0800
From: Laura Prugh <lprugh at alaska.edu>
To: Ben Bolker <bbolker at gmail.com>

Hi Ben,
Sorry to bother you again (please let me know if you prefer I post these
questions to a list?)--I updated my version of R to the newest and reran
the code, which had run fine on the older version. It's now giving me an
error referring to how I'm trying to get it to loop through my
dataset--my code is:
for(i in 1:nlevels(vital)) {
lmer.sheep <- glmer(y[vital==levels(vital)[i]] ~
(1|year[vital==levels(vital)[i]]), family = binomial, data = sheep.dat)
}
and it gives the error:
Error in FUN(X[[1L]], ...) :
  Invalid grouping factor specification, year[vital == levels(vital)[i]]

When I put "year[vital == levels(vital)[i]]" directly into R, it appears
to work properly, and this syntax worked yesterday on the older R
version. Is a new syntax required or could this perhaps be a bug?

thanks,
Laura

On 5/2/2014 1:18 PM, Ben Bolker wrote:
> On 14-05-02 03:53 PM, Laura Prugh wrote:
>> Hi Ben,
>> I am trying to save specific outputs from the lmer summary, and I was
>> hoping you might be able to tell me the code to do this. I'm running a
>> random effects model with only the intercept as a fixed effect. I want
>> to extract the fixed effect intercept estimate and the variance (or
>> standard deviation) of the random effect. I seem to have figured out how
>> to get the fixed effect (model at fixedef) but have had no luck getting the
>> random effect. I guessed a variety of things and model at ranef returned
>> something, but it didn't match anything in the summary output. Here is
>> the output of my model, and I would like to extract and save the st dev
>> of the random effect (0.21287). Any help would be greatly appreciated.
>> If you generally know of a resource for figuring out the syntax to
>> extract specific values from summary outputs, that would be great
>> because this is a recurring problem for me! By the way, I plan to have
>> as many of my students at possible attend the workshop you're offering
>> at Chena HS near Fairbanks in August.
>> cheers,
>> Laura
>     This is really more of an r-sig-mixed-models at r-project.org question:
> I'm forwarding the answer there.
>
>    methods(class="merMod") is a good way to check the available accessor
> methods for (g)lmer fits, as is ?getME (an lme4-specific function for
> getting at lower-level components of a (g)lmer fit in a safe/stable way).
>
>
>>> lmer.sheep <- lmer(y ~ (1|year),family=binomial, data = coy,
>> + control=list(msVerbose=F,maxIter=1000,maxFN=3000))
>>> lmer.sheep
>> Generalized linear mixed model fit by the Laplace approximation
>> Formula: y ~ (1 | year)
>> Data: coy
>> Fixed effects:
>> Estimate Std. Error z value Pr(>|z|)
>> (Intercept) -0.8083 0.2168 -3.729 0.000192 ***
>> ---
>> Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
>>> mn <- lmer.sheep at fixef
>>> mn
>> (Intercept)
>> -0.8083453
>>> lmer.sheep at ranef
>> [1] 0.09490052 0.09490052 0.03168227 -0.01783168 -0.14832171 -0.04832885
>     A reproducible example would be slightly better, but here I believe
> you're looking for
>
>    fixef(lmer.sheep)[1]  ## intercept
>    sqrt(unlist(VarCorr(lmer.sheep))  ## std dev of year effect
>
>    It looks like you're still using an older version of lme4: I believe
> the incantations I gave above are backward compatible, but there are a
> fair number of convenience methods in the newer (>1.0) version that
> might be helpful (e.g. as.data.frame(VarCorr(lmer.sheep)))
>
>    Ben Bolker
>
>