[R] lme in R (WinXP) vs. Splus (HP UNIX)

From: Greg Tarpinian <sasprog474474_at_yahoo.com>
Date: Tue 31 Jan 2006 - 12:30:15 EST


R2.2 for WinXP, Splus 6.2.1 for HP 9000 Series, HP-UX 11.0.

I am trying to get a handle on why the same lme( ) code gives such different answers. My output makes me wonder if the fact that the UNIX box is 64 bits is the reason. The estimated random effects are identical, but the fixed effects are very different. Here is my R code and output, with some columns and rows deleted for space considerations:

FOO.lme1 <- lme(fixed = Y ~ X1*X2,
random = ~ X2 | X3, data = FOO,
method = "REML",
control = list(maxIter = 200, msMaxIter = 200, tolerance = 1e-8, niterEM = 200,
msTol = 1e-8, msVerbose = TRUE,
optimMethod = "Nelder-Mead"))

  0      203.991: 0.924323 -0.0884338 0.493897
  1      203.991: 0.924323 -0.0884338 0.493897

> summary(FOO.lme1)

Linear mixed-effects model fit by REML
 Data: FOO

      AIC BIC logLik
  357.484 373.6868 -170.7420

Random effects:
 Formula: ~ X2 | X3
 Structure: General positive-definite,
 Log-Cholesky parametrization

               StdDev Corr

(Intercept)    2.025868 (Intr)
X2             4.908517 -0.475
Residual       4.493171       

Fixed effects: Y ~ X1*X2
                              Value  p-value
(Intercept)                39.08067  0.0000
X1                         12.08700  0.0000
X2                          5.85417  0.1370
X1:X2                       0.80833  0.8299 

Number of Observations: 60
Number of Groups: 3

Here is my Splus code and output, similarly edited for space considerations:

FOO.lme1 <- lme(fixed = Y ~ X1*X2,
random = ~ X2 | X3, data = FOO,
method = "REML",
control = list(maxIter = 200, msMaxIter = 200, tolerance = 1e-8, niterEM = 200,
msTol = 1e-8, msVerbose = TRUE,
optimMethod = "Nelder-Mead"))

Iteration: 0 , 1 function calls, F= 205.3776 Parameters: [1] 0.89525077 0.22549223 -0.05936197

Iteration: 1 , 2 function calls, F= 205.3776 Parameters: [1] 0.89525078 0.22549224 -0.05936187

> summary(FOO.lme1)

Linear mixed-effects model fit by REML
 Data: FOO

       AIC BIC logLik
  360.2566 376.4594 -172.1283

Random effects:
 Formula: ~ X2 | X3
 Structure: General positive-definite

                 StdDev   Corr 
   (Intercept) 2.025861 (Inter
            X2 4.908618 -0.475
      Residual 4.493172       

Fixed effects: Y ~ X1* X2
                        Value  p-value 
         (Intercept) 45.12417  <.0001
               X1     6.04350  <.0001
               X2     6.25833  0.0709
            X1:X2     0.40417  0.8299

Number of Observations: 60
Number of Groups: 3

I have examined the help.start() HTML pages and the lmeScale( ) function looks helpful, but I don't know how to use it to tell R to use a particular set of starting values for the params. It would have been instructive to see what happens if the same starting values are specified for both Splus and R....

Below is a sample dataset that is fairly close to the data that I have actually been using (not specifically the data used to obtain the output above).

Any help would be greatly appreciated.

Kind Regards,

     Greg

Sample Data:

A,X3,X2,X1,Y
1,A,1.6,A,43
1,A,1.6,B,56
2,A,1.6,A,43
2,A,1.6,B,50
3,A,1.6,A,46
3,A,1.6,B,64
4,A,2.0,A,42
4,A,2.0,B,42
5,A,2.0,A,36
5,A,2.0,B,44
6,A,2.0,A,31
6,A,2.0,B,44
7,A,2.0,A,39
7,A,2.0,B,46
8,A,2.4,A,47
8,A,2.4,B,52
9,A,2.4,A,50
9,A,2.4,B,52
10,A,2.4,A,43
10,A,2.4,B,56
11,B,1.6,A,35
11,B,1.6,B,48
12,B,1.6,A,32
12,B,1.6,B,45
13,B,1.6,A,37
13,B,1.6,B,50
14,B,2.0,A,40
14,B,2.0,B,55
15,B,2.0,A,33
15,B,2.0,B,46
16,B,2.0,A,38
16,B,2.0,B,53
17,B,2.0,A,42
17,B,2.0,B,55
18,B,2.4,A,39
18,B,2.4,B,56
19,B,2.4,A,43
19,B,2.4,B,56
20,B,2.4,A,45
20,B,2.4,B,62
21,C,1.6,A,33
21,C,1.6,B,41
22,C,1.6,A,33
22,C,1.6,B,44
23,C,1.6,A,31
23,C,1.6,B,46
24,C,2.0,A,38
24,C,2.0,B,48
25,C,2.0,A,34
25,C,2.0,B,45
26,C,2.0,A,36
26,C,2.0,B,51
27,C,2.0,A,40
27,C,2.0,B,57
28,C,2.4,A,39
28,C,2.4,B,53
29,C,2.4,A,34
29,C,2.4,B,51

30,C,2.4,A,36
30,C,2.4,B,52

R-help@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 Received on Tue Jan 31 12:43:02 2006

This archive was generated by hypermail 2.1.8 : Fri 03 Mar 2006 - 03:42:15 EST