[R] mixed model resuts from SAS and R

aldi at dsgmail.wustl.edu aldi at dsgmail.wustl.edu
Thu May 22 14:08:40 CEST 2008


Hi,

I was wondering if there is a way to figure out why in SAS random beta
coefficients are 0 vs. in R the beta-s are non zero.
The variables of the data are nidl, time, and sub (for subject). Time and
nidl are continuous variables. I am applying random coefficients model.
Any input is greatly appreciated, Thanks, Aldi

1. mixed model in SAS:
======================
ods output SolutionR = out1.randomnidltest2;
proc mixed data = a1 ;
  class sub ;
  model nidl = time /  solution  ;
  random int time /  sub = sub solution;
run;
ods output close;

2. mixed model in R:
====================
a1<-read.table(file="c:\\aldi\\a1.txt",sep=",",header=T)
library(nlme)
fm1nidl.lme<-lme(nidl~time,data=a1,random=~time | sub)
plot(coef(fm1nidl.lme))


3. SAS output:
==============
   Plot of nidl*time.  Symbol used is '*'.

  40 ˆ
     ‚
     ‚*
     ‚
     ‚*
     ‚*           *
     ‚*          *
     ‚*          *      *
  20 ˆ*          *   *
     ‚*          *****
     ‚*          ** **
     ‚*         *** ****          *
     ‚*         *** **  *
     ‚*        * ****** *
     ‚*         ****** *    *
     ‚*       * *******     *          *
   0 ˆ*       * ********   *
     -ˆ-----------ˆ------------ˆ------------ˆ
      0           25           50           75

                        time
NOTE: 684 obs hidden.
            Dimensions
Covariance Parameters             3
Columns in X                      2
Columns in Z Per Subject          2
Subjects                        425
Max Obs Per Subject               2
          Number of Observations
Number of Observations Read             763
Number of Observations Used             763
Number of Observations Not Used           0
 Covariance Parameter Estimates
Cov Parm      Subject    Estimate

Intercept     sub         17.1773
time          sub               0
Residual                  27.0023
The Mixed Procedure

           Fit Statistics

-2 Res Log Likelihood          5005.5
AIC (smaller is better)        5009.5
AICC (smaller is better)       5009.5
BIC (smaller is better)        5017.6


                   Solution for Fixed Effects

                         Standard
Effect       Estimate       Error      DF    t Value    Pr > |t|

Intercept      7.7608      0.3214     424      24.15      <.0001
time         -0.08148     0.01605     337      -5.08      <.0001


                      Solution for Random Effects

                                 Std Err
Effect       sub    Estimate        Pred      DF    t Value    Pr > |t|

Intercept      1      5.6722      3.2426       0       1.75       .
time           1           0           .       .        .         .
Intercept      2      3.6722      3.2426       0       1.13       .
time           2           0           .       .        .         .
Intercept      3     -2.0774      2.7539       0      -0.75       .
time           3           0           .       .        .         .

...      ...          ....                ...           ....    ...

R output:
    (Intercept)          time
1     17.432680 -0.3155496745
2     14.024527 -0.2345787274
3      3.105323  0.0469759240
4     23.047062 -0.5565796200
5     10.708267 -0.1557909941
... ... ...
> summary(fm1nidl.lme)
Linear mixed-effects model fit by REML
 Data: a1
      AIC      BIC    logLik
  4982.15 5009.958 -2485.075

Random effects:
 Formula: ~time | sub
 Structure: General positive-definite, Log-Cholesky parametrization
            StdDev    Corr
(Intercept) 6.0090637 (Intr)
time        0.1717771 -0.831
Residual    4.2885993

Fixed effects: nidl ~ time
                Value Std.Error  DF  t-value p-value
(Intercept)  7.779379 0.3582945 424 21.71225       0
time        -0.086206 0.0158615 337 -5.43494       0
 Correlation:
     (Intr)
time -0.675

Standardized Within-Group Residuals:
       Min         Q1        Med         Q3        Max
-2.0234047 -0.5132952 -0.2246854  0.4249250  3.5611259

Number of Observations: 763
Number of Groups: 425

3. data:
========
see the text file (comma delimited) attached.


-------------- next part --------------
An embedded and charset-unspecified text was scrubbed...
Name: a1.txt
URL: <https://stat.ethz.ch/pipermail/r-help/attachments/20080522/db9c1399/attachment.txt>


More information about the R-help mailing list