6

I'm a new user of mixed models and through the material I've been reading there are always probability values (p>t) or (p>z) that estimate the importance of a level of a factor in the model. However, when using the lmer() function in R, which supposedly gives you those probabilities, I simply don't find them. Here is the output: 

Linear mixed model fit by REML 
Formula: Temp ~ depth + (1 | locality) 
   Data: qminmatrix 
   AIC   BIC logLik deviance REMLdev
 561.3 581.3 -273.7    551.5   547.3
Random effects:
 Groups   Name        Variance Std.Dev.
 locality (Intercept) 4.7998   2.1909  
 Residual             4.0433   2.0108  
Number of obs: 128, groups: locality, 4

Fixed effects:
            Estimate Std. Error t value
(Intercept)  22.0103     1.1500  19.140
depth1        1.9564     0.6832   2.864
depth10       2.6624     0.5756   4.625
depth5        3.0209     0.4932   6.125
depthWS      -2.2585     0.5444  -4.149

Correlation of Fixed Effects:
        (Intr) depth1 dpth10 depth5
depth1  -0.157                     
depth10 -0.175  0.189              
depth5  -0.213  0.313  0.458       
depthWS -0.191  0.334  0.373  0.441
whuber
  • 281,159
  • 54
  • 637
  • 1,101
Agus Camacho
  • 540
  • 2
  • 4
  • 14
  • 6
    Since this question appears to be mostly about interpreting R's output, it may belong on [Stack Overflow](http://stackoverflow.com/) rather than here. OTOH, however, IIRC, the reason the package authors decided not to display the p-values has to do w/ the difficulty of determining the right number of degrees of freedom to use in assessing the t-values. This is a very subtle statistical issue, & requires some nuanced explanations, so it might be worth leaving this Q here to give CVers a chance to attempt to explain it first. – gung - Reinstate Monica Sep 28 '12 at 14:20
  • @gung: Your comment actually answers the question. I've put a link to the author's explanation in my answer, below. – Wayne Sep 28 '12 at 14:47
  • 4
    I think this one should stay here, not go to Stack Overflow. This is about statistics, not about R. – Peter Flom Sep 28 '12 at 21:48
  • The function pvals.fnc does not work with model with random intercept and slope though. Do you have any other idea how I can get the p-values? PS: I do not understand why most journals want p-values despite many articles point out that the confidence intervals are more important than p-values. I cannot find papers with only confidence intervals, but no p-values. – user3288202 Feb 26 '14 at 09:19

3 Answers3

16

The lmer package's author made a conscious choice not to create p-values for the fixed effects. Some packages do, but he feels that they are doing simplistic calculations that are misleading. (Many statisticians feel that there's a p-value obsession that causes confusion in and of itself, but that's a separate matter.)

He addresses the question in: this post. I believe the summary paragraph is:

Most of the research on tests for the fixed-effects specification in a mixed model begin with the assumption that these statistics will have an F distribution with a known numerator degrees of freedom and the only purpose of the research is to decide how to obtain an approximate denominator degrees of freedom. I don't agree.

I don't understand the issue well enough to paraphrase it.

amoeba
  • 93,463
  • 28
  • 275
  • 317
Wayne
  • 19,981
  • 4
  • 50
  • 99
4

use pvals.fnc() function the pMCMC here works like p-value which should be less than 0.05 to reject the null hypothesis.

Hoda
  • 56
  • 1
0

Install the coda and languageR package and run the pvals.fnc as follows for p-value:

Model.pval<-pvals.fnc(Model, nsim = n, withMCMC = TRUE)

Note that this will not work for level 3 or above in nested random effects models.

KarthikS
  • 1,066
  • 1
  • 9
  • 18