Skip to main content
slightly more detail
Source Link
Ben Bolker
  • 47.4k
  • 3
  • 130
  • 182
  1. the difference between t-tests and Z-tests (as pointed out by @vkehayas); t-tests account for the uncertainty in the estimate of the standard error, so should be preferred to Z-tests where available.
  2. the fact that summary.lme by default adjusts the residual standard error for ML estimates (while glht doesn't); ML estimation in general gives a slightly downward-biased estimate of the standard error (by a factor $\sqrt{(n-p)/n}$), so this adjustment should be preferred where available. This is the adjustSigma parameter of summary.lme:

Both of these adjustments should in general make little difference unless your sample is small, but both adjustSigma=TRUE and t-tests rather than Z-test are technically more correct, so in a pinch you should probably accept the results of summary(.) rather than those of glht().

If you have a factor with more than two levels (so that you need to summarize the joint significance of multiple parameters), you can use anova(), which uses F tests (the analog of t-tests) and includes an adjustSigma option: if you want to do more complicated post hoc testing (e.g. Tukey pairwise comparisons), you will probably need to use glht() and accept that your answers will be slightly anticonservative/optimistic.

  1. the difference between t-tests and Z-tests (as pointed out by @vkehayas)
  2. the fact that summary.lme by default adjusts the residual standard error for ML estimates (while glht doesn't):

Both of these adjustments should in general make little difference unless your sample is small, but both adjustSigma=TRUE and t-tests rather than Z-test are technically more correct, so in a pinch you should probably accept the results of summary(.) rather than those of glht().

  1. the difference between t-tests and Z-tests (as pointed out by @vkehayas); t-tests account for the uncertainty in the estimate of the standard error, so should be preferred to Z-tests where available.
  2. the fact that summary.lme by default adjusts the residual standard error for ML estimates (while glht doesn't); ML estimation in general gives a slightly downward-biased estimate of the standard error (by a factor $\sqrt{(n-p)/n}$), so this adjustment should be preferred where available. This is the adjustSigma parameter of summary.lme:

Both of these adjustments should in general make little difference unless your sample is small, but both adjustSigma=TRUE and t-tests rather than Z-test are technically more correct, so in a pinch you should probably accept the results of summary(.) rather than those of glht().

If you have a factor with more than two levels (so that you need to summarize the joint significance of multiple parameters), you can use anova(), which uses F tests (the analog of t-tests) and includes an adjustSigma option: if you want to do more complicated post hoc testing (e.g. Tukey pairwise comparisons), you will probably need to use glht() and accept that your answers will be slightly anticonservative/optimistic.

added 416 characters in body
Source Link
Ben Bolker
  • 47.4k
  • 3
  • 130
  • 182

Both of these adjustments should in general make little difference unless your sample is small, but both adjustSigmaadjustSigma=TRUE and t-tests rather than Z-test are technically more correct, so you in a pinch you should probably accept the results of summary(.) rather than those of glht().

library(multcomp)
library(nlme)
data("sleepstudy",package="lme4")
                        
m2 <- lme(Reaction~Days, random = ~Days|Subject,
          data=sleepstudy, method="ML")

Results (fancy code with printCoefmat() etc. is just to isolate the information we want from summary(m2)):

printCoefmat(summary(m2)$tTab["Days",,drop=FALSE])
##         Value Std.Error       DF t-value   p-value    
## Days  10.4673    1.5106 161.0000   6.929 9.651e-11 ***
printCoefmat(summary(m2,adjustSigma=FALSE)$tTab["Days"$tTab["Days",,drop=FALSE])
##         Value Std.Error       DF t-value   p-value    
## Days  10.4673    1.50225106 161.0000   6.9678929 79.811e651e-11 ***

With adjustSigma=FALSE, the standard error ischanges from 1.5022 rather than5106 to 1.51065022:

printCoefmat(summary(m2,adjustSigma=FALSE)$tTab["Days",,drop=FALSE])
##         Value Std.Error       DF t-value   p-value    
## Days  10.4673    1.5022 161.0000  6.9678 7.811e-11 ***

In your case most of the difference is from the t- vs Z-test distinction; 2*pt(2.206,df=15,lower.tail=FALSE) (i.e. using the unadjusted standard error with a t test) gives $p=0.043$, most of the way from $p=0.027$ (summary(.) result) to $p=0.049$ (glht(.) result).

Both of these adjustments should in general make little difference unless your sample is small, but both adjustSigma and t-tests are technically more correct, so you in a pinch you should probably accept the results of summary(.) rather than those of glht().

library(multcomp)
library(nlme)
data("sleepstudy",package="lme4")
                        
m2 <- lme(Reaction~Days, random = ~Days|Subject,
          data=sleepstudy, method="ML")
printCoefmat(summary(m2)$tTab["Days",,drop=FALSE])
##         Value Std.Error       DF t-value   p-value    
## Days  10.4673    1.5106 161.0000   6.929 9.651e-11 ***
printCoefmat(summary(m2,adjustSigma=FALSE)$tTab["Days",,drop=FALSE])
##         Value Std.Error       DF t-value   p-value    
## Days  10.4673    1.5022 161.0000  6.9678 7.811e-11 ***

With adjustSigma=FALSE, the standard error is 1.5022 rather than 1.5106

Both of these adjustments should in general make little difference unless your sample is small, but both adjustSigma=TRUE and t-tests rather than Z-test are technically more correct, so in a pinch you should probably accept the results of summary(.) rather than those of glht().

library(multcomp)
library(nlme)
data("sleepstudy",package="lme4")
                        
m2 <- lme(Reaction~Days, random = ~Days|Subject,
          data=sleepstudy, method="ML")

Results (fancy code with printCoefmat() etc. is just to isolate the information we want from summary(m2)):

printCoefmat(summary(m2)$tTab["Days",,drop=FALSE])
##         Value Std.Error       DF t-value   p-value    
## Days  10.4673    1.5106 161.0000   6.929 9.651e-11 ***

With adjustSigma=FALSE, the standard error changes from 1.5106 to 1.5022:

printCoefmat(summary(m2,adjustSigma=FALSE)$tTab["Days",,drop=FALSE])
##         Value Std.Error       DF t-value   p-value    
## Days  10.4673    1.5022 161.0000  6.9678 7.811e-11 ***

In your case most of the difference is from the t- vs Z-test distinction; 2*pt(2.206,df=15,lower.tail=FALSE) (i.e. using the unadjusted standard error with a t test) gives $p=0.043$, most of the way from $p=0.027$ (summary(.) result) to $p=0.049$ (glht(.) result).

Source Link
Ben Bolker
  • 47.4k
  • 3
  • 130
  • 182

There are two things going on here:

  1. the difference between t-tests and Z-tests (as pointed out by @vkehayas)
  2. the fact that summary.lme by default adjusts the residual standard error for ML estimates (while glht doesn't):

adjustSigma: an optional logical value. If ‘TRUE’ and the estimation method used to obtain ‘object’ was maximum likelihood, the residual standard error is multiplied by sqrt(nobs/(nobs - npar)), converting it to a REML-like estimate. ... Default is ‘TRUE’.

Both of these adjustments should in general make little difference unless your sample is small, but both adjustSigma and t-tests are technically more correct, so you in a pinch you should probably accept the results of summary(.) rather than those of glht().

Try to keep in mind that $p=0.0274$ and $p=0.0493$ (from your example) are not very different from each other; in practice people behave as if there's a magic line at $p=0.05$, but there isn't.

Here's an example:

library(multcomp)
library(nlme)
data("sleepstudy",package="lme4")
                        
m2 <- lme(Reaction~Days, random = ~Days|Subject,
          data=sleepstudy, method="ML")
printCoefmat(summary(m2)$tTab["Days",,drop=FALSE])
##         Value Std.Error       DF t-value   p-value    
## Days  10.4673    1.5106 161.0000   6.929 9.651e-11 ***
printCoefmat(summary(m2,adjustSigma=FALSE)$tTab["Days",,drop=FALSE])
##         Value Std.Error       DF t-value   p-value    
## Days  10.4673    1.5022 161.0000  6.9678 7.811e-11 ***

With adjustSigma=FALSE, the standard error is 1.5022 rather than 1.5106

A direct calculation of the p-value using the unadjusted sigma:

2*pt(6.9678,161,lower.tail=FALSE)
## [1]  7.811903e-11

If we instead use a Z-test:

2*pnorm(6.9678,lower.tail=FALSE)
## [1] 3.219354e-12

This agrees with the answer we get from glht:

    summary(glht(m2, linfct=c("Days=0")))
##           Estimate Std. Error z value Pr(>|z|)    
## Days == 0   10.467      1.502   6.968 3.22e-12 ***