[R] summary.lme: argument "adjustSigma"

[Spencer Graves]

	  It appears that the "adjustSigma" argument of 'summary.lme' does 
nothing with the default method, "REML".  To check, I tried the 
following modification of the 'summary.lme' example:

fm1 <- lme(distance ~ age, Orthodont, random = ~ age | Subject)
fm1sa <- summary(fm1)
fm1s <- summary(fm1, adjustSigma=FALSE)
fm1saa <- summary(fm1, adjustSigma=TRUE)
all.equal(fm1s, fm1sa)
TRUE
all.equal(fm1s, fm1saa)
TRUE
#################

	  When I changed 'method' to "ML" in this example, the result suggested 
that the adjustment to sigma also affected the standard errors and t 
values for the fixed effects.  If the p-values had not been 0, they also 
would have been affected:

fm2 <- lme(distance ~ age, Orthodont, random = ~ age | Subject,
            method="ML")
(fm2sa <- summary(fm2))
Linear mixed-effects model fit by maximum likelihood
  Data: Orthodont
        AIC      BIC    logLik
   451.2116 467.3044 -219.6058

Random effects:
  Formula: ~age | Subject
  Structure: General positive-definite, Log-Cholesky parametrization
             StdDev    Corr
(Intercept) 2.1940998 (Intr)
age         0.2149245 -0.581
Residual    1.3100399

Fixed effects: distance ~ age
                 Value Std.Error DF   t-value p-value
(Intercept) 16.761111 0.7678975 80 21.827277       0
age          0.660185 0.0705779 80  9.353997       0
  Correlation:
     (Intr)
age -0.848

Standardized Within-Group Residuals:
          Min           Q1          Med           Q3          Max
-3.305969237 -0.487429631  0.007597973  0.482237063  3.922789795

Number of Observations: 108
Number of Groups: 27
 > (fm2s <- summary(fm2, adjustSigma=FALSE))
Linear mixed-effects model fit by maximum likelihood
  Data: Orthodont
        AIC      BIC    logLik
   451.2116 467.3044 -219.6058

Random effects:
  Formula: ~age | Subject
  Structure: General positive-definite, Log-Cholesky parametrization
             StdDev    Corr
(Intercept) 2.1940998 (Intr)
age         0.2149245 -0.581
Residual    1.3100399

Fixed effects: distance ~ age
                 Value Std.Error DF  t-value p-value
(Intercept) 16.761111 0.7607541 80 22.03223       0
age          0.660185 0.0699213 80  9.44183       0
  Correlation:
     (Intr)
age -0.848

Standardized Within-Group Residuals:
          Min           Q1          Med           Q3          Max
-3.305969237 -0.487429631  0.007597973  0.482237063  3.922789795

Number of Observations: 108
Number of Groups: 27
 >
	  Does this answer the question?
	  spencer graves
	
Christoph Buser wrote:

> Dear R-list
> 
> I have a question concerning the argument "adjustSigma" in the
> function "lme" of the package "nlme". 
> 
> The help page says:
> 
> "the residual standard error is multiplied by sqrt(nobs/(nobs - 
> npar)), converting it to a REML-like estimate."
> 
> Having a look into the code I found:
> 
> stdFixed <- sqrt(diag(as.matrix(object$varFix)))
> 
> if (object$method == "ML" && adjustSigma == TRUE) {
>     stdFixed <- sqrt(object$dims$N/(object$dims$N - length(stdFixed))) * 
>         stdFixed
> }
> 
> tTable <- data.frame(fixed, stdFixed, object$fixDF[["X"]], 
>     fixed/stdFixed, fixed)
> 
> 
> To my understanding, only the standard error for the fixed
> coefficients is adapted and not the residual standard error. 
> 
> Therefore only the tTable of the output is affected by the
> argument "adjustSigma", but not the estimate for residual
> standard error (see the artificial example below). 
> 
> May someone explain to me if there is an error in my
> understanding of the help page and the R code? 
> Thank you very much.  
> 
> Best regards,
> 
> Christoph Buser
> 
> --------------------------------------------------------------
> Christoph Buser <buser at stat.math.ethz.ch>
> Seminar fuer Statistik, LEO C13
> ETH Zurich	8092 Zurich	 SWITZERLAND
> phone: x-41-44-632-4673		fax: 632-1228
> http://stat.ethz.ch/~buser/
> --------------------------------------------------------------
> 
> 
> Example
> -------
> 
> set.seed(1)
> dat <- data.frame(y = rnorm(16), fac1 = rep(1:4, each = 4),
>                   fac2 = rep(1:2,each = 8))
> 
> telme <- lme(y ~ fac1, data = dat, random = ~ 1 | fac2, method = "ML")
> summary(telme)
> summary(telme, adjustSigma = FALSE)
> 
> ______________________________________________
> R-help at stat.math.ethz.ch mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html