13

I was under the impression that the function lmer() in the lme4 package didn't produce p-values (see lmer, p-values and all that).

I've been using MCMC generated p values instead as per this question: Significant effect in lme4 mixed model and this question: Can't find p-values in the output from lmer() in the lm4 package in R.

Recently I tried a package called memisc and its getSummary.mer() to get the fixed effects of my model into a csv file. As if by magic, a column called p appears which matches my MCMC p-values extremely closely (and doesn't suffer the processing time that comes with using pvals.fnc()).

I've tentatively had a look at the code in getSummary.mer and have spotted the line that generates the p-value:

p <- (1 - pnorm(abs(smry@coefs[, 3]))) * 2

Does this mean p values can be generated directly from lmer's output rather than running pvals.fnc? I realise this will no doubt start the 'p-value fetishism' debate but I'm interested to know. I've not heard memisc mentioned before when it comes to lmer.

To be more succinct: What is the benefit (if any) of using MCMC p-values over those generated by getSummary.mer()?

Community
  • 1
Tumbledown
  • 899
  • 9
    I am the original author of the getSummary.mer function. The reported $p$-values should only be used as a quick check. If I recall, I actually only included the $p$-values to make it work within the framework provided by memisc. But this should really be provided with an appropriate warning to the user, and I will contact the package maintainer to see about getting this added. My advice is to follow that provided by Doug Bates: MCMC is the safe bet (assuming others don't have better options). – Jason Morgan Sep 03 '13 at 00:55
  • 4
    @JasonMorgan This looks to me to be a pretty reasonable answer to the question. – Glen_b Sep 03 '13 at 01:22
  • 2
    @JasonMorgan I agree with what you say but currently mcmcsamp() is not available because of a number of issues (One can check the Status of mcmcsamp section in glmm.wikidot.com/faq for more details). I feel that at the moment time, probably (parametric?) bootstrapping is a viable -and not too hard to implement- alternative; the bootMer() functiom can be of service. – usεr11852 Sep 03 '13 at 03:10
  • @Glen_b Should I just copy-paste my comment as an answer? – Jason Morgan Sep 03 '13 at 03:21
  • 1
    @JasonMorgan I think a copy-paste would be fine, since to my eyes your valuable comment does actually answer the question. (However, if you feel you can expand even a little on your explanation of why the p-values should really only be used as a quick check rather than as p-values, or when they may be more likely to be poor approximations, or why MCMC is a safe bet, so much the better.) – Glen_b Sep 03 '13 at 06:15
  • @JasonMorgan Thanks for your comment, happy to accept it as an answer if you paste it as such. One thing I did note was that pvals.fnc does appear to generate identical figures to getSummary.mer() in its Pr(>|t|) output. I've seen in another tutorial (the location escapes me) that these values should be compared to the MCMC p values as a sense check, just as you mentioned. – Tumbledown Sep 03 '13 at 06:39
  • 3
    One other technical/terminological note to the OP. The p-values obtained using the method in memisc are the p-values from treating the observed test statistics as Wald statistics (treating the t as a Wald z in this case). Such a test relies on the "large sample" assumption and so is more and more trustworthy as your sample sizes grow larger. The MCMC-based value, to my knowledge, does not rely on such an assumption. So anyway, reading a little about Wald tests and alternatives to them could help to shed further light on your question. – Jake Westfall Sep 03 '13 at 06:47
  • 1
    @JakeWestfall Is there any document demonstrating that the MCMC-based method has satisfactory frequentist properties ? – Stéphane Laurent Sep 03 '13 at 08:02
  • 1
    @StéphaneLaurent It is a good question and unfortunately I don't know the answer. I don't know too much about MCMC to be honest. Maybe someone else who knows more about the usage here could chime in, or maybe you could even expand your query into a whole separate question. – Jake Westfall Sep 03 '13 at 08:08

1 Answers1

11

I am the original author of the getSummary.mer function. The reported $p$-values should only be used as a quick check. If I recall, I actually only included the $p$-values to make it work within the framework provided by memisc. But this should really be provided with an appropriate warning to the user, and I will contact the package maintainer to see about getting this added. My advice is to follow that provided by Doug Bates: MCMC is the safe bet (assuming others don't have better options).

Jason Morgan
  • 1,214
  • 11
  • 15