[R] gls() error

toby909 at gmail.com toby909 at gmail.com
Fri May 18 04:02:32 CEST 2007


Hi All

How can I fit a repeated measures analysis using gls? I want to start with a 
unstructured correlation structure, as if the the measures at the occations are 
not longitudinal (no AR) but plainly multivariate (corSymm).

My data (ignore the prox_pup and gender, occ means occasion):
 > head(dta,12)
    teacher occ prox_self prox_pup gender
1        1   0      0.76     0.41      1
2        1   1      1.00     1.05      1
3        1   3      0.81     0.91      1
4        2   0      0.96     0.64      0
5        3   0      1.12     1.13      1
6        3   2      1.34     1.35      1
7        4   1      0.35    -0.40      0
8        4   2      0.25     0.27      0
9        4   3      0.54     0.26      0
10       5   0      0.65     1.02      1
11       5   1      0.68     0.87      1
12       5   2      1.04     0.98      1


x=factor(dta$occ)

Following gives me an error message:

gls(prox_pup~x-1, dta, corSymm(, ~x-1|teacher))
 > gls(prox_pup~x-1, dta, corSymm(, ~x-1|teacher))
Error in Initialize.corSymm(X[[1]], ...) :
         Covariate must have unique values within groups for corSymm objects
In addition: There were 50 or more warnings (use warnings() to see the first 50)

I checked that the covariate, occ, has unique values within each of the teachers.

Aside, lme actually gives me what I want, except that the residual variance is 
not where I want to have it. I want the residuals being part of the covariance 
matrix to be estimated rather than in the 'level one' residual, ie the residuals 
included on the diagonal of "G" displayed below.

 > a4 = lme(prox_pup~x-1, dta, ~x-1|teacher)
Linear mixed-effects model fit by REML
   Data: dta
   Log-restricted-likelihood: -53.91059
   Fixed: prox_pup ~ x - 1
        x0        x1        x2        x3
0.5739072 0.7169963 0.6503699 0.6567064

Random effects:
  Formula: ~x - 1 | teacher
  Structure: General positive-definite, Log-Cholesky parametrization
          StdDev    Corr
x0       0.5424187 x0    x1    x2
x1       0.4326164 0.739
x2       0.3343281 0.611 0.681
x3       0.3638630 0.569 0.752 0.900
Residual 0.0929820

Number of Observations: 153
Number of Groups: 51

 > G = lapply(pdMatrix(a4$modelStruct$reStruct), "*", a4$sigma^2)
$teacher
           x0         x1         x2        x3
x0 0.2942180 0.17330375 0.11089028 0.1123597
x1 0.1733037 0.18715693 0.09847681 0.1183089
x2 0.1108903 0.09847681 0.11177526 0.1094374
x3 0.1123597 0.11830892 0.10943738 0.1323963



Thanks for your help on this.

Toby



More information about the R-help mailing list