[R] extracting p-values from lmer outputs
Lngmyers
lngmyers at ccu.edu.tw
Sat Mar 4 03:45:55 CET 2006
Thanks to Martin Maechler for the helpful information below.
I guess I didn't make myself clear enough. I had already discovered all
them slots in lmer objects, but these are apparently mostly input things
or things created in the course of the analysis, not the output itself
-- that is, it still takes work to derive from them the text printed on
the screen when lmer is used. I am a mere mortal, not a statistician.
The output of anova(lmerout) (where "lmerout" is a lmer object) is an
ordinary list, so the lines it prints are easy to access, but this is
not true of lmerout itself. I didn't want anova(lmerout) because, as MM
said, certain key values have been suppressed to avoid controversial.
However, the version of lmer currently available on CRAN (lme4 0.995-2
2006-01-17 and Matrix 0.995-5 2006-02-07) still gives p-values for GLMM.
I assumed that they weren't banned there because they're calculated
using z values, making DF irrelevant. The datasets I'll be using always
have over 600 observations, so z values should work fine.
So two questions:
(1) Are the fixed-effects p-values produced by lmer for GLMM still
considered to be valid, or were they left in the last revision by
mistake? Should I give up on these p-values and test for significance of
each fixed effect by using anova(full-model,
model-minus-one-of-the-fixed-effects)?
(2) Assuming the p-values given by lmer for GLMM are valid, is there any
easy way to derive them from a lmer object? Thanks to MM's tip, I looked
more closely and found the estimates for the fixed effects, but I don't
see standard errors or z values. They must be readily calculable from
what's in the lmer object, but I'm too ignorant to figure it out myself.
--
James Myers
Graduate Institute of Linguistics
National Chung Cheng University
168 University Road, Min-Hsiung
Chia-Yi 62102
TAIWAN
Email: Lngmyers at ccu.edu.tw
Web: http://www.ccunix.ccu.edu.tw/~lngmyers/
Phone: 886-5-242-8251
Fax: 886-5-272-1654
Martin Maechler wrote:
>
> >>>>> "lngmyers" == lngmyers <lngmyers at ccu.edu.tw>
> >>>>> on Fri, 3 Mar 2006 21:54:52 +0800 (CST) writes:
>
> lngmyers> I would like to write a function that runs GLMM using lmer
> lngmyers> on a user-input model containing interactions, but if the
> lngmyers> model doesn't produce significant results for the interaction,
> lngmyers> a reduced model will be run without the interaction.
> lngmyers> Doing this seems to require getting the p-values out of an
> lngmyers> lmer object, but I don't know how to do this. (The grand
> lngmyers> DF debate seems to be irrelevant since the number of observations
> lngmyers> for what I want to do with this will always be over 600 or so.)
> lngmyers> The problem is that the output of lmer is a list of length 0.
> lngmyers> My "brilliant" idea, shown below, is to divert the display to
> lngmyers> a file, then read the file as a string. I guess it's useful
> lngmyers> to save the full summary somewhere, but is there really
> lngmyers> no more elegant way to do this?
>
> I'm not answering your question directly,
> just the (implicit)
> ``how can I extract parts from an 'lmer' object ? ''
>
> If you used str(lmerout)
> to inspect the STRucture of the lmerout object,
> you'd seen something like
>
> > str(lmerout)
> Formal class 'lmer' [package "Matrix"] with 35 slots
> ..@ assign : int(0)
> ..@ frame :`data.frame': 180 obs. of 3 variables:
> .. ..$ Reaction: num [1:180] 250 259 251 321 357 ...
> .. ..$ Days : num [1:180] 0 1 2 3 4 5 6 7 8 9 ...
> .. ..$ Subject : Factor w/ 18 levels "308","309","310",..: 1 1 1 1 1 1 1 1 1 1 ...
>
> ............
> ............
>
> which tells you two things (and more) :
>
> 1) you have an object of a *formal* class,
> ``i.e.'' an S4 object
>
> 2) you can (and typically will) use " lmerout @ <slotname> "
> to extract the slots of your lmerout result object.
>
> Note that the above is a low-level / technical approach and you
> should rather work with the official documentation:
>
> ?lmer points you to "lmer-class" and
>
> class ? lmer or
> ?lmer-class {that one is needed for ESS}
>
> contains a description of the "lmer" object and its (most
> important) slots.
>
> Note BTW that in the most uptodate versions of 'lme4'/'Matrix'
> you won't get any P-values anymore directly
> {just because of the never-ending "DF war"}.
>
> Martin Maechler, ETH Zurich
>
> lngmyers> - James Myers
> lngmyers> National Chung Cheng University
> lngmyers> Minhsiung Chiayi 621 Taiwan
>
> lngmyers> lmerout = lmer(Y ~ X1 * X2 + (1|Subj), data = dat, family = "binomial")
> lngmyers> sink("lmerout.txt")
> lngmyers> lmerout
> lngmyers> sink()
> lngmyers> lmerout.string = readChar("lmerout.txt",20000)
> lngmyers> splitstring = strsplit(lmerout.string,"\r\n")
> lngmyers> # If I counted right, the line with the interaction is [19]:
> lngmyers> interstring = strsplit(splitstring[[1]][19]," ")
> lngmyers> # Lines contain TABs, so the p-value is the eighth element along:
> lngmyers> pvalue = as.numeric(interstring[[1]][8])
> lngmyers> if(pvalue < .05) print("The interaction is significant!")
>
> lngmyers> ______________________________________________
> lngmyers> R-help at stat.math.ethz.ch mailing list
> lngmyers> https://stat.ethz.ch/mailman/listinfo/r-help
> lngmyers> PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
More information about the R-help
mailing list