# [R] help with gls

Zhenqiang Lu lu25 at stat.purdue.edu
Sat Jan 6 07:21:29 CET 2007

```Hello R-users,

I am using gls function in R to fit a model with certain correlation
structure.

The medol as:
fit.a<-gls(y~1,data=test.data,correlation=corAR1(form=~1|aa),method="ML")
mu<-summary(fit.a)\$coefficient

With the toy data I made to test, the estimate of  mu is exactly equal to
the overall mean of y which can not be true.

But, if I make a toy data with y more than two replicates (for each level
of aa we have more than 2 abservations, for example N=4), the estimates of
mu will not be as the same as the overall mean of y.

Would you please tell me why this happens?
The following is my testing code.
Thanks.
James

require(mvtnorm)
rho=-0.5
N=2
sigma.a=2

Rho.a<-matrix(c(1,rho^(1:(N-1)))[outer(X=1:N,Y=1:N,function(x,y)
1+abs(x-y))],ncol=N)
Sigma.a<-sigma.a^2 * Rho/(1-rho.a^2)
y<-rep(0,N*20)
for (ii in 1:20)
{y[((ii-1)*N+1):(ii*N)]<-rmvnorm(1, mean=rep(0,N), sigma=Sigma.a)
}
test.data<-data.frame(y=y,aa=rep(1:20,each=N,len=N*20))
fit.a<-gls(y~1,data=test.data,correlation=corAR1(form=~1|aa),method="ML")
mu.a<-summary(fit.a)\$coefficient
rho.a<-getVarCov(fit.a)[1,2]/getVarCov(fit.a)[1,1]
print(c(mean(y),mu.a))

______________________________________________________
Zhenqiang (James) Lu

Department of Statistics
Purdue University
West Lafayette, IN 47906
TEL: (765) 494-0027
FAX: (765) 494-0558

```