Dear Rusers,
I have used R,S-PLUS and SAS to analyze the sample data "bacteria" in
MASS package. Their results are listed below.
I have three questions, anybody can give me possible answers?
Q1:From the results, we see that R get 'NAs'for AIC,BIC and logLik, while
S-PLUS8.0 gave the exact values for them. Why? I had thought that R should
give the same results as SPLUS here.
Q2: The model to analyse the data is logity=b0+u+b1*trt+b2*I(week>2), but
the results for Random effects in R/SPLUS confused me. SAS may be clearer.
Random effects:
Formula: ~1 | ID
(Intercept) Residual
StdDev: 1.410637 0.7800511
Which is the random effect 'sigma'? I think it is "1.410637", but what
does "0.7800511" mean? That is, i want ot know how to explain/use the
above two data for Random effects.
Q3:In SAS and other softwares, we can get *p*-values for the random effect
'sigma', but i donot see the *p*-values in the results of R/SPLUS. I have
used attributes() to look for them, but no *p* values. Anybody knows how to
get *p*-values for the random effect 'sigma',.
Any suggestions or help are greatly appreciated.
#R Results:MASS' version 7.2-44; R version 2.7.2
library(MASS)
summary(glmmPQL(y ~ trt + I(week > 2), random = ~ 1 | ID,family = binomial,
data = bacteria))
Linear mixed-effects model fit by maximum likelihood
Data: bacteria
AIC BIC logLik
NA NA NA
Random effects:
Formula: ~1 | ID
(Intercept) Residual
StdDev: 1.410637 0.7800511
Variance function:
Structure: fixed weights
Formula: ~invwt
Fixed effects: y ~ trt + I(week > 2)
Value Std.Error DF t-value
p-value
(Intercept) 3.412014 0.5185033 169 6.580506 0.0000
trtdrug -1.247355 0.6440635 47 -1.936696 0.0588
trtdrug+ -0.754327 0.6453978 47 -1.168779 0.2484
I(week > 2)TRUE -1.607257 0.3583379 169 -4.485311 0.0000
Correlation:
(Intr) trtdrg trtdr+
trtdrug -0.598
trtdrug+ -0.571 0.460
I(week > 2)TRUE -0.537 0.047 -0.001
#S-PLUS8.0: The results are the same as R except the followings:
AIC BIC logLik
1113.622 1133.984 -550.8111
#SAS9.1.3
proc nlmixed data=b;
parms b0=-1 b1=1 b2=1 sigma=0.4;
yy=b0+u+b1*trt+b2*week;
p=1/(1+exp(-yy));
Model response~binary(p);
Random u~normal(0,sigma) subject=id;
Run;
-2 Log Likelihood = 192.2
AIC (smaller is better)=200.2
AICC (smaller is better) =200.3
BIC (smaller is better)= 207.8
Parameter Estimates
Standard
Parameter Estimate Error DF t Value Pr > |t| Alpha
Lower Upper Gradient
b0 3.4966 0.6512 49 5.37 <.0001 0.05
2.1880 4.8052 -4.69E-6
trt -0.6763 0.3352 49 -2.02 0.0491 0.05
-1.3500 -0.00266 -0.00001
I(week>2) -1.6132 0.4785 49 -3.37 0.0015 0.05 -2.5747
-0.6516 -9.35E-7
sigma 1.5301 0.9632 49 1.59 0.1186 0.05
-0.4054 3.4656 -2.42E-6
--
With Kind Regards,
oooO:::::::::
(..):::::::::
:\.(:::Oooo::
::\_)::(..)::
:::::::)./:::
::::::(_/::::
:::::::::::::
[***********************************************************************]
ZhiJie Zhang ,PhD
Dept.of Epidemiology, School of Public Health,Fudan University
Office:Room 443, Building 8
Office Tel./Fax.:+86-21-54237410
Address:No. 138 Yi Xue Yuan Road,Shanghai,China
Postcode:200032
Email:[EMAIL PROTECTED] <[EMAIL PROTECTED]>
Website: www.statABC.com <http://www.statabc.com/>
[***********************************************************************]
oooO:::::::::
(..):::::::::
:\.(:::Oooo::
::\_)::(..)::
:::::::)./:::
::::::(_/::::
:::::::::::::
[[alternative HTML version deleted]]
______________________________________________
[email protected] 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.
