[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