[R-sig-ME] Fwd: same old question - lme4 and p-values

Andy Fugard a.fugard at ed.ac.uk
Sat Apr 5 17:19:55 CEST 2008 

I made (more) sense of mixed effects models when I went back to standard 
regression and ANOVAs with the philosophy: Everything is a Comparison.

So for instance noting different ways of getting the magical numbers 
that result from doing regression:

---8<--------------------------------------------------------------------

 > m.0 = lm(Fertility ~ 1, data = swiss)
 > m.full = lm(Fertility ~ ., data = swiss)
 > summary(m.full)

...

Residual standard error: 7.17 on 41 degrees of freedom
Multiple R-squared: 0.707,      Adjusted R-squared: 0.671
F-statistic: 19.8 on 5 and 41 DF,  p-value: 5.59e-10
              ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^

 > anova(m.0,m.full)
Analysis of Variance Table

Model 1: Fertility ~ 1
Model 2: Fertility ~ Agriculture + Examination + Education + Catholic +
     Infant.Mortality
   Res.Df  RSS Df Sum of Sq    F  Pr(>F)
1     46 7178
2     41 2105  5      5073 19.8 5.6e-10 ***
       ^^       ^           ^^^^^^^^^^^^
----------------------------------------------------------------->8------

Noting what happens when you compare nested models (incidentally, a 
confusing term when used the context of "multilevel" models):

---8<--------------------------------------------------------------------

 > m.full = lm(Fertility ~ . , data = swiss)
 > summary(m.full)

...

Coefficients:
                  Estimate Std. Error t value Pr(>|t|)
(Intercept)       66.9152    10.7060    6.25  1.9e-07 ***
Agriculture       -0.1721     0.0703   -2.45   0.0187 *
Examination       -0.2580     0.2539   -1.02   0.3155
Education         -0.8709     0.1830   -4.76  2.4e-05 ***
Catholic           0.1041     0.0353    2.95   0.0052 **
Infant.Mortality   1.0770     0.3817    2.82   0.0073 **
                                                ^^^^^^
...

 >
 > m1 = update(full.model, ~. -Infant.Mortality)
 > anova(m.full,m1)
Analysis of Variance Table

Model 1: Fertility ~ Agriculture + Examination + Education + Catholic +
     Infant.Mortality
Model 2: Fertility ~ Agriculture + Examination + Education + Catholic
   Res.Df  RSS Df Sum of Sq    F Pr(>F)
1     41 2105
2     42 2514 -1      -409 7.96 0.0073 **
                                 ^^^^^^

----------------------------------------------------------------->8------

Playing around with this sort of thing made it a lot easier to 
understand why statisticians get annoyed (e.g., 
<http://www.stats.ox.ac.uk/pub/MASS3/Exegeses.pdf> for a good example) 
with us non-statisticians.

On a more practical note, the languageR package has an especially useful 
pvals.fnc function for when you really do need a test of whether a slope 
is significantly different to zero, e.g. when you've got a load of 
categorical predictors and want to see where the difference is.

Hope that wasn't too far off topic.

Cheers,

Andy



Hank Stevens wrote:
> Google:
> p-values lmer wiki
> 
> On Apr 4, 2008, at 9:33 AM, Douglas Bates wrote:
> 
>> ---------- Forwarded message ----------
>> From: Douglas Bates <bates at stat.wisc.edu>
>> Date: Fri, Apr 4, 2008 at 7:54 AM
>> Subject: Re: same old question - lme4 and p-values
>> To: andreas.nord at zooekol.lu.se
>>
>>
>>
>> On Fri, Apr 4, 2008 at 5:24 AM,  <andreas.nord at zooekol.lu.se> wrote:
>>> Dear Prof. Bates,
>>> I've recently switched to using R for my analyses, and I find the
>> lme4 package to be extremely helpful. I have read your explanation
>> (posted on the mailing list) of why you choose not to display
>> p-values. Unfortunately, most of the journals I publish in require
>> that I include p-values, which is why I have to find a way of
>> calculating them from the lmer output. However, not being a trained
>> statistician I have some difficulties following your recommendations
>> given in the explanatory text. In other words, after having fitted my
>> model, I am not at all sure on what to do in order to obtain p-values
>> (or similar).
>>
>>> I am sorry to have to bother you with a question I know you have
>> already answered many times, but perhaps you would be so kind as to
>> give me some hints on how to proceed.
>>
>> I understand your situation.  Statisticians have created the "every
>> question of scientific interest must be answered by a p-value" monster
>> and now it turns on us.  Nevertheless I am reluctant to give advice on
>> p-values in lme4 because apparently I don't know how to do it
>> correctly.
>>
>> May I send a copy of this reply to the
>> R-SIG-Mixed-Models at r-project.org mailing list? ("SIG" == "Special
>> Interest Group")? (I ask your permission to send the copy because I am
>> quoting your original question.)   Some who subscribe to that mailing
>> list may have the courage to wade into this swamp and offer their
>> advice.
>>
>> _______________________________________________
>> R-sig-mixed-models at r-project.org mailing list
>> https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models
> 
> 
> 
> Dr. Hank Stevens, Associate Professor
> 338 Pearson Hall
> Botany Department
> Miami University
> Oxford, OH 45056
> 
> Office: (513) 529-4206
> Lab: (513) 529-4262
> FAX: (513) 529-4243
> http://www.cas.muohio.edu/~stevenmh/
> http://www.cas.muohio.edu/ecology
> http://www.muohio.edu/botany/
> 
> "If the stars should appear one night in a thousand years, how would men
> believe and adore." -Ralph Waldo Emerson, writer and philosopher  
> (1803-1882)
> 
> _______________________________________________
> R-sig-mixed-models at r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models
> 


-- 
Andy Fugard, Postgraduate Research Student
Psychology (Room F3), The University of Edinburgh,
   7 George Square, Edinburgh EH8 9JZ, UK
Mobile: +44 (0)78 123 87190   http://www.possibly.me.uk

The University of Edinburgh is a charitable body, registered in
Scotland, with registration number SC005336.

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