gls vs lme covariance structures
Hi Charles, Well you could post a subset of it, or make up some data that is sharable (whether the data make any sense is not important to us, just nice to have runable code, for example your previous thread about contrasts could have been solved in one email if we could have shown you how to set the contrasts on your data and then it matched your SAS output). In any case, whether you use lme or gls really depends on your question and goals, I think. Generalized least squares is not the same as a random effects model. If you want a random effect, you cannot use gls. If you just want correlated errors, gls is fine. This part of your code strikes me as atypical though I cannot promise it is wrong/not what you want: corr=corAR1(ID) Cheers, Josh
On Thu, May 3, 2012 at 1:44 PM, Charles Determan Jr <deter088 at umn.edu> wrote:
Hi Joshua, Thanks for your response.? It is probably best that I don't post the data as some of it is not yet published.? My main question is whether UN and AR(1) should be done with gls or if I have done the syntax incorrectly with lme. Since AR(1) is replicated perfectly if I put the correct dendf, I can work with it.? And UN is close, so I just want to be sure my use and syntax are correct, not necessarily modifying the data. Regards, Charles On Thu, May 3, 2012 at 3:21 PM, Joshua Wiley <jwiley.psych at gmail.com> wrote:
Hi Charles,
Could you upload the dataset you are using somewhere and post the
link? ?Something like:
##########
dat34 <- read.csv("http:/wherever/you/uploaded/yourdata.csv", header =
TRUE)
## code to convert to factors anything that needs to be etc.
##########
Then it is easier for us to try things that way.
corAR1 and corSymm seem appropriate. ?Have you checked the examples in
their documents? ?I found them helpful.
Cheers,
Josh
On Thu, May 3, 2012 at 12:58 PM, Charles Determan Jr <deter088 at umn.edu>
wrote:
Greetings R users, I have been attempting to replicate various covariance structures from SAS's analysis of mixed models. ?I have successfully been able to replicate compound symmetry, however it becomes confusing with autoregression and unstructured. ?As such, there are two questions regarding this issue. Autoregression SAS output (Type III fixed effects) for covariance structure AR(1) *Type 3 Tests of Fixed Effects* *Effect ? ? ? ? ? ? ? ? ? ?NumDF DenDF F Value Pr > F* *group ? ? ? ? ? ? ? ? ? ? ? ? ?*1 ? ? ? ? 23 ? ? ? ? ?0.99 ? ? ?0.3293 *Event_name ? ? ? ? ? ? ? *5 ? ? ? ? 91 ? ? ? ? 16.23 ? ?<.0001 *Died ? ? ? ? ? ? ? ? ? ? ? ? ? ?*1 ? ? ? ? 23 ? ? ? ? ?1.70 ? ? ?0.2047 *group*Event_name ?*5 ? ? ? ?91 ? ? ? ? ?3.04 ? ? 0.0140 R output (corAR1=AR(1)?) I can replicate these results if I run the following: fit.18=gls(var~group+Event_name+Died+group*Event_name, ? ?data=dat34, ? ?corr=corAR1(, ~1|ID), ? ?weight=varIdent(~1|Event_name)) anova(fit.18, type="marginal", adjustSigma=F) #the DenDF are off with gls, so use the 'correct' ones 1-pf(.9935, 1, 23) 1-pf(16.2323, 5, 91) 1-pf(1.7041, 1, 23) 1-pf(3.0367, 5, 91) #and the output matches exactly However, I can not get the lme function to run the autoregression. ?The output is very different: fit.11=lme(var~group+Event_name+Died+group*Event_name, ? ?data=dat34, ? ?random=~1|ID, ? ?corr=corAR1(ID), ? ?weight=varIdent(~1|Event_name)) anova.lme(fit.11, type="marginal", adjustSigma=F) ? ? ? ? ? ? ? ? ? ? ? numDF denDF ?F-value ? ? ?p-value (Intercept) ? ? ? ? ? ? ? 1 ? ?96 ? ? ? ?9.816419 ?0.0023 group ? ? ? ? ? ? ? ? ? ? ?1 ? ?23 ? ? ? 0.131950 ?0.7197 Event_name ? ? ? ? ? ?5 ? ?96 ? ? ? ?1.081785 ?0.3754 Died ? ? ? ? ? ? ? ? ? ? ? 1 ? ?23 ? ? ? ?0.074428 ?0.7874 Is this type of covariance structure only done with gls and I should continue with the analysis as such or am I doing something silly with lme? My second question is with regards to the unstructured covariance. SAS output (UN) *Type 3 Tests of Fixed Effects* *Effect ? ? ? ?NumDF DenDF F Value Pr > F* *Pig_group ? ?*1 ? ? ? ? ? ?23 ? ? ? ? 2.73 ? ?0.1120 *Event_name *5 ? ? ? ? ? ?23 ? ? ? ? 1.11 ? ?0.3806 *Died ? ? ? ? ? ? ?*1 ? ? ? ? ? ?23 ? ? ? ? 0.51 ? ?0.4833 R output (corSymm = UN?) fit.11=lme(var2~group+Event_name+Died, ? ?data=dat34, ? ?random=~1|ID, ? ?corr=corSymm(, ~1|ID), ? ?weight=varIdent(~1|Event_name)) anova.lme(fit.11, type="marginal", adjustSigma=F) #same as corAR1??? ? ? ? ? ? ? ? ? ? ? ? numDF denDF ?F-value ? ? ?p-value (Intercept) ? ? ? ? ? ? ? 1 ? ?96 ? ? ? ?9.816419 ?0.0023 group ? ? ? ? ? ? ? ? ? ? ?1 ? ?23 ? ? ? 0.131950 ?0.7197 Event_name ? ? ? ? ? ?5 ? ?96 ? ? ? ?1.081785 ?0.3754 Died ? ? ? ? ? ? ? ? ? ? ? 1 ? ?23 ? ? ? ?0.074428 ?0.7874 but with gls fit.18=gls(var~group+Event_name+Died, ? ?data=dat34, ? ?corr=corSymm(~1|ID), ? ?weight=varIdent(~1|Event_name)) anova(fit.18, type="marginal", adjustSigma=F) 1-pf(2.869837, 1, 23) 1-pf(1.126747, 5, 23) 1-pf(.514726, 1, 23) [1] 0.1037549 [1] 0.3742309 [1] 0.4803239 #close but not exact (however I can work with this if it is indeed correct) Overall, I want to clarify the difference between gls and lme and if I am simply making some weird syntax error with lme that I can't seem to get the covariance structures to match. Apologies for lots of information in one go, but hopefully this provides necessary information to point me in the correct direction. Thanks to any and all who give their time to these questions, Regards, Charles ? ? ? ?[[alternative HTML version deleted]]
_______________________________________________ R-sig-mixed-models at r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models
-- Joshua Wiley Ph.D. Student, Health Psychology Programmer Analyst II, Statistical Consulting Group University of California, Los Angeles https://joshuawiley.com/
Joshua Wiley Ph.D. Student, Health Psychology Programmer Analyst II, Statistical Consulting Group University of California, Los Angeles https://joshuawiley.com/