Hi,

I have a question about fitting a random coefficients model with different
error structures in lme() function. The data I have in the long format as
the following:

     id gender Occasion Strength
1.0   1   Male        0      161
1.1   1   Male        1      210
1.2   1   Male        2      230
2.0   2   Male        0      215
2.1   2   Male        1      245
2.2   2   Male        2      265
3.0   3   Male        0      134
3.1   3   Male        1      215
................................................................
................................................................

When I fit the following model by imposing a compound symmetry structure on
the error var-cov matrix:

fit <- lme(Strength ~ 1 + Occasion,
           random = ~ 1 + Occasion | id,
           data = grip.L,
           correlation=corCompSymm(),
           control=list(maxIter=1000,msMaxIter=1000, niterEM=1000)
           )

Linear mixed-effects model fit by REML
  Data: grip.L
  Log-restricted-likelihood: -873.9892
  Fixed: Strength ~ 1 + Occasion
(Intercept)    Occasion
144.8304598  -0.6637931

Random effects:
 Formula: ~1 + Occasion | id
 Structure: General positive-definite, Log-Cholesky parametrization
            StdDev   Corr
(Intercept) 72.69465 (Intr)
Occasion    16.14353 -0.366
Residual    16.45378

Correlation Structure: Compound symmetry
 Formula: ~1 | id
 Parameter estimate(s):
Rho
  0
Number of Observations: 174
Number of Groups: 58

The "Rho" is estimated as zero. When I specify corCompSymm(.3), the output
is:

Linear mixed-effects model fit by REML
  Data: grip.L
  Log-restricted-likelihood: -873.9892
  Fixed: Strength ~ 1 + Occasion
(Intercept)    Occasion
144.8304598  -0.6637931

Random effects:
 Formula: ~1 + Occasion | id
 Structure: General positive-definite, Log-Cholesky parametrization
            StdDev   Corr
(Intercept) 72.02091 (Intr)
Occasion    16.14343 -0.369
Residual    19.18853

Correlation Structure: Compound symmetry
 Formula: ~1 | id
 Parameter estimate(s):
      Rho
0.2647256
Number of Observations: 174
Number of Groups: 58

The rho is estimated .264. In this case, I am wondering what really .3 within
corCompSymm(.3) represents. Or, when I say corCompSymm(.3) , what does it me
an?

Also, when I specify corAR1(), I am getting an error message about
convergence.

fit <- lme(Strength ~ 1 + Occasion,
           random = ~ 1 + Occasion | id,
           data = grip.L,
           correlation=corAR1(),
           control=list(maxIter=1000,msMaxIter=1000, niterEM=1000)
           )

Error in lme.formula(Strength ~ 1 + Occasion, random = ~1 + Occasion |  :
  nlminb problem, convergence error code = 1
  message = singular convergence (7)

After a trial-error process, I was able to fit the model by specifying
corAR1(.05):

Linear mixed-effects model fit by REML
 Data: grip.L
       AIC      BIC    logLik
  1761.052 1783.084 -873.5258

Random effects:
 Formula: ~1 + Occasion | id
 Structure: General positive-definite, Log-Cholesky parametrization
            StdDev    Corr
(Intercept) 39.491786 (Intr)
Occasion     5.554041 -0.999
Residual    62.765655

Correlation Structure: AR(1)
 Formula: ~1 | id
 Parameter estimate(s):
      Phi
0.9021773
Fixed effects: Strength ~ 1 + Occasion
                Value Std.Error  DF   t-value p-value
(Intercept) 143.52300  9.736396 115 14.740875  0.0000
Occasion     -0.66379  2.617492 115 -0.253599  0.8003
 Correlation:
         (Intr)
Occasion -0.396

Standardized Within-Group Residuals:
       Min         Q1        Med         Q3        Max
-1.4632691 -0.6942647 -0.1709512  0.4677042  2.0371593

Number of Observations: 174
Number of Groups: 58

Here, the estimated Phi is .9021. I have the same question again. What
(.05)represents within
corAR1(.05)? How is it related to the estimated phi .9021? Should I always do
a trial-error process to figure out what value I should put in paranthesis
to be able fit a model?

I also tried other error structures for the same dataset, but none of them
converged. I used a different larger dataset, but I always get convergence
error when I try to fit different error structures. Any suggestions?

This is for didactic purposes. We would like to teach students how to fit
different error structures for mixed effects models using lme() function, and
I can't get it worked until now for none of the datasets I have.

Thanks in advance.

-- 

Cengiz Zopluoglu

	[[alternative HTML version deleted]]

