[R-sig-ME] Covariance models with lme and gls

Ben Bolker bbolker at gmail.com
Tue Mar 6 17:44:03 CET 2012


Student Stats <statsstud8 at ...> writes:


> I'm trying to fit several covariance models using gls and lme. My aim is to
> identify which covariance model fits my data better. I'm afraid, however,
> that I'm not specifying the code properly. Could someone take a look on my
> code and help me figuring out whether I'm pursuing everything correctly?
> 

 You just posted this to StackOverflow.  It isn't explicitly proscribed
(as cross-posting among R help lists is), but I don't like it because it
diffuses information and potentially wastes effort.  I told you there that
this seems like a bit of a "please debug my code for me question" -- what
have you tried already in order to figure out whether you're doing it
right or not?  Do you have reasons to suspect you're doing something
wrong?  It's easier to answer a specific question than to look over
a big hunk of code.

   I suggested to you on Stack Overflow that using update() would
make it easier to look at your code.  In fact you can structure it
a bit further:

vs <- varIdent(form=~1|time)
cs <- corSymm(form=~1|id)
ccs <- corCompSymm(form=~1|id)
car1 <- corAR1(form=~1|id)

> # Independence covariance matrix
> IN <- gls(y ~ ses + time, data, corr=NULL, weights=NULL, method="REML",
> control=lmeControl(msMaxIter = 500, msVerbose = TRUE))
> 

UN <- update(IN, corr=cs, weights=vs)

> # Fit Random Intercept Model (RI)
> RI <- lme(y ~ ses + time, data, na.action=na.omit, method="REML",
> random=~1|id, control=lmeControl(msMaxIter = 200, msVerbose = TRUE))
> 

RIAS <- update(RI, random=~time|id)
CS <- update(IN,correlation=ccs)
CSH <- update(CS,weights=vs)
AR1 <- update(IN,correlation=car1)

  et cetera.




More information about the R-sig-mixed-models mailing list