[R] lme versus proc mixed in SAS

Dimitris Rizopoulos dimitris.rizopoulos at med.kuleuven.ac.be
Wed May 4 15:42:52 CEST 2005


check this:

library(nlme)
B <- 1000
N <- 100
n <- 5
x <- rep(runif(N, -4, 4), each=n)
gr <- rep(1:N, each=n)
####################
T <- numeric(B)
for(i in 1:B){
    y <- rnorm(N*n, 1 + 1.5*x)
    L0 <- lm(y~x)
    L1 <- lme(y~x, random=~1|gr, method="ML")
    T[i] <- anova(L1, L0)$L.Ratio[2]
}
hist(T, prob=TRUE, breaks=100)


I hope it helps.

Best,
Dimitris

----
Dimitris Rizopoulos
Ph.D. Student
Biostatistical Centre
School of Public Health
Catholic University of Leuven

Address: Kapucijnenvoer 35, Leuven, Belgium
Tel: +32/16/336899
Fax: +32/16/337015
Web: http://www.med.kuleuven.ac.be/biostat/
     http://www.student.kuleuven.ac.be/~m0390867/dimitris.htm



----- Original Message ----- 
From: "Beatrijs Moerkerke" <Beatrijs.Moerkerke at UGent.be>
To: <r-help at stat.math.ethz.ch>
Sent: Wednesday, May 04, 2005 12:03 PM
Subject: [R] lme versus proc mixed in SAS


> Dear all,
>
> I am trying to simulate the null distribution for the likelihood 
> ratio test statistic for testing 1 random effect versus no random 
> effect.  The asymptotic null distribution should be a mixture of a 
> chi-squared distribution with 0 degrees of freedom and a chi-squared 
> distribution with 1 degree of freedom.  This means that I expect a 
> point mass of 50% on 0 for the likelihood ratio test statistic.
> However, when I generate data using no random effects and when I 
> calculate the test statistics for these data, I never obtain exactly 
> zero.  I think this might be due to rounding errors but in fact, 70% 
> of the calculated test statistics are negative.  I have compared a 
> few of these results with the results in proc MIXED and I found that 
> SAS does give test statistics that are exactly zero and gives no 
> negative results.
>
> The code I use for calculating the likelihood ratio test statistics 
> is as follows:
>
> a1<-summary(lme(y~x,random=~1|gr,method="ML"))$logLik
> a2<-logLik(lm(y~x))
> (-2*(a2-a1))
>
> I don't know how I can simulate the null distribution in R using 
> lme.
>
> Thanks for your help,
>
> Kind regards,
> Beatrijs Moerkerke
>
> -- 
> Beatrijs Moerkerke
> Department of Applied Mathematics and Computer Science
> Ghent University
> Krijgslaan 281 - S9
> B-9000 GENT
> Tel: +32-(0)9-264.47.56      Fax: +32-(0)9-264.49.95
> E-mail: Beatrijs.Moerkerke at UGent.be
>
> ______________________________________________
> 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
>




More information about the R-help mailing list