[R] lme with/without varPower - can I use AIC?

Dieter Menne dieter.menne at menne-biomed.de
Wed Sep 11 20:27:27 CEST 2002


I want to compare the following two models in AIC
(Treat, Spotter are categorial, p is pressure, Pain is 
continuous)

PainW.lme<-lme(Pain~p+Treat*Spotter,data=saw,random=~p|Pat,
             weights=varPower(form=~Pain))
# AIC= -448

Pain.lme<-lme(Pain~p+Treat*Spotter,data=saw,random=~p|Pat)
#AIC = -19.7

Note the huge differences in AIC, and the estimated power of 6.
A plot of the residual does not show an unusual patterns for
both models.

I do not trust the varPower result, but don't have any rationale
for it.

1) Can I use the AIC at all to compare two weightings?
2) Can I trust such a high power estimate? There is definitively 
some slight dependency of variance on Pain, but it is not extreme.
3) I tried fitted(.) instead of Pain, but it did no converge after
5 Minutes.

---- Detail ---
> PainW.lme<-lme(Pain~p+Treat*Spotter,data=saw,random=~p|Pat,
+             weights=varPower(form=~Pain))
> summary(PainW.lme)
Linear mixed-effects model fit by REML
 Data: saw 
   AIC  BIC logLik
  -448 -413    234

Random effects:
 Formula: ~p | Pat
 Structure: General positive-definite, Log-Cholesky parametrization
            StdDev   Corr  
(Intercept) 6.48e-02 (Intr)
p           4.19e-06 0     
Residual    1.65e-01       

Variance function:
 Structure: Power of variance covariate
 Formula: ~Pain 
 Parameter estimates:

power 
 5.89 

Fixed effects: Pain ~ p + Treat * Spotter 
                     Value Std.Error  DF t-value p-value
(Intercept)          0.747   0.01753 214    42.6  <.0001
p                    0.000   0.00017 214     1.4   0.164
TreatXen            -0.003   0.00555 214    -0.5   0.587
SpotterSp+          -0.015   0.02492  29    -0.6   0.551
TreatXen:SpotterSp+  0.011   0.00847 214     1.4   0.177
 Correlation: 
                    (Intr) p      TretXn SpttS+
p                   -0.151                     
TreatXen            -0.163 -0.088              
SpotterSp+          -0.685 -0.020  0.125       
TreatXen:SpotterSp+  0.106  0.062 -0.655 -0.194

Standardized Within-Group Residuals:
   Min     Q1    Med     Q3    Max 
-2.473 -0.180  0.229  0.839  2.443 

Number of Observations: 248
Number of Groups: 31 

> Pain.lme<-lme(Pain~p+Treat*Spotter,data=saw,random=~p|Pat)
> summary(Pain.lme)
Linear mixed-effects model fit by REML
 Data: saw 
    AIC  BIC logLik
  -19.7 11.8   18.8

Random effects:
 Formula: ~p | Pat
 Structure: General positive-definite, Log-Cholesky parametrization
            StdDev  Corr  
(Intercept) 0.00748 (Intr)
p           0.00746 0.144 
Residual    0.18544       

Fixed effects: Pain ~ p + Treat * Spotter 
                     Value Std.Error  DF t-value p-value
(Intercept)          0.751    0.0360 214   20.87  <.0001
p                    0.008    0.0017 214    4.62  <.0001
TreatXen            -0.047    0.0328 214   -1.45  0.1497
SpotterSp+           0.067    0.0489  29    1.36  0.1839
TreatXen:SpotterSp+ -0.119    0.0472 214   -2.51  0.0127
 Correlation: 
                    (Intr) p      TretXn SpttS+
p                   -0.358                     
TreatXen            -0.444 -0.007              
SpotterSp+          -0.637 -0.013  0.329       
TreatXen:SpotterSp+  0.305  0.017 -0.696 -0.490

Standardized Within-Group Residuals:
   Min     Q1    Med     Q3    Max 
-3.001 -0.586 -0.112  0.380  3.591 

Number of Observations: 248
Number of Groups: 31 
> 
> 
-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-
r-help mailing list -- Read http://www.ci.tuwien.ac.at/~hornik/R/R-FAQ.html
Send "info", "help", or "[un]subscribe"
(in the "body", not the subject !)  To: r-help-request at stat.math.ethz.ch
_._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._



More information about the R-help mailing list