[R] lme vs. SAS proc mixed. Point estimates and SEs are the same, DFs are different

Peter Dalgaard p.dalgaard at biostat.ku.dk
Tue Jun 5 10:44:24 CEST 2007


John Sorkin wrote:
> R 2.3
> Windows XP
>
> I am trying to understand lme. My aim is to run a random effects regression in which the intercept and jweek are random effects. I am comparing output from SAS PROC MIXED with output from R. The point estimates and the SEs are the same, however the DFs and the p values are different. I am clearly doing something wrong in my R code. I would appreciate any suggestions of how I can change the R code to get the same DFs as are provided by SAS.
>   
This has been hashed over a number of times before. In short:

1) You're not necessarily doing anything wrong
2) SAS PROC MIXED is not necessarily doing it right
3) lme() is _definitely_ not doing it right in some cases
4) both work reasonably in large sample cases (but beware that this is 
not equivalent to having many observation points)

SAS has an implementation of the method by Kenward and Rogers, which 
could be the most reliable general DF-calculation method around (I don't 
trust their Satterthwaite option, though). Getting this or equivalent 
into lme() has been on the wish list for a while, but it is not a 
trivial thing to do.

> SAS code:
> proc mixed data=lipids2;
>   model ldl=jweek/solution;
>   random int jweek/type=un subject=patient;
>   where lastvisit ge 4;
> run;
>
> SAS output:
>                    Solution for Fixed Effects
>
>                          Standard
> Effect       Estimate       Error      DF    t Value    Pr > |t|
>
> Intercept      113.48      7.4539      25      15.22      <.0001
> jweek         -1.7164      0.5153      24      -3.33      0.0028
>
>         Type 3 Tests of Fixed Effects
>
>               Num     Den
> Effect         DF      DF    F Value    Pr > F
> jweek           1      24      11.09    0.0028
>
>
> R code:
> LesNew3 <- groupedData(LDL~jweek | Patient, data=as.data.frame(LesData3), FUN=mean)
> fit3    <- lme(LDL~jweek, data=LesNew3[LesNew3[,"lastvisit"]>=4,], random=~1+jweek)
> summary(fit3) 
>
> R output:
> Random effects:
>  Formula: ~1 + jweek | Patient
>  Structure: General positive-definite, Log-Cholesky parametrization
>  
>
> Fixed effects: LDL ~ jweek 
>                 Value Std.Error DF   t-value p-value
> (Intercept) 113.47957  7.453921 65 15.224144  0.0000
> jweek        -1.71643  0.515361 65 -3.330535  0.0014
>
> John Sorkin M.D., Ph.D.
> Chief, Biostatistics and Informatics
> University of Maryland School of Medicine Division of Gerontology
> Baltimore VA Medical Center
> 10 North Greene Street
> GRECC (BT/18/GR)
> Baltimore, MD 21201-1524
> (Phone) 410-605-7119
> (Fax) 410-605-7913 (Please call phone number above prior to faxing)
>
> Confidentiality Statement:
> This email message, including any attachments, is for the so...{{dropped}}
>
> ______________________________________________
> 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
> and provide commented, minimal, self-contained, reproducible code.
>



More information about the R-help mailing list