I just wanted to point out that for simple repeated measures, longitudinal data, it --is-- possible to mimic some simple reasonable lme models with gls. If lme uses the random effects specification random = ~1 | Subject Then gls can mimic that by using correlation = corCompSymm(form = ~1 | Subject) The resulting parameterizations are different, but they give the same predictions and SEs for the "fixed"/model effects. ----- I often take this gls paradigm a step further and make it a habit to use correlation = corrExp(form = LongitudinalTimeVariable | Subject, nugget = TRUE) In the limit of the range parameter going to +infinity, the model reduces to the above corCompSymm model equivalent to an lme model. In the limit as the nugget parameter going to zero, the model reduces to a CAR1 model. (see ?corrExp and other corXXX) I see cases of both of these. But what astonishes me is that how often this two parameter correlation model fits my data so well: it mirrors nearly exactly what a full unrestricted correlation structure, using corSymm, gives. (of course, that is my data, YMMV.) Furthermore, the product of the nugget parameter and the overall variance gives an estimate of cross sectional replication variance, a nicety since I don't believe one can use replicated data with the corXXX functions. John John Szumiloski, Ph.D. Senior Biometrician Biometrics Research WP53B-120 Merck Research Laboratories P.O. Box 0004 West Point, PA 19486-0004
(215) 652-7346 (PH) (215) 993-1835 (FAX)
-----Original Message----- From: r-sig-mixed-models-bounces at r-project.org [mailto:r-sig-mixed-models-bounces at r-project.org] On Behalf Of r-sig-mixed-models-request at r-project.org Sent: Thursday, May 03, 2012 5:17 PM To: r-sig-mixed-models at r-project.org Subject: R-sig-mixed-models Digest, Vol 65, Issue 12 Send R-sig-mixed-models mailing list submissions to r-sig-mixed-models at r-project.org To subscribe or unsubscribe via the World Wide Web, visit https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models or, via email, send a message with subject or body 'help' to r-sig-mixed-models-request at r-project.org You can reach the person managing the list at r-sig-mixed-models-owner at r-project.org When replying, please edit your Subject line so it is more specific than "Re: Contents of R-sig-mixed-models digest..." Today's Topics: 1. Re: gls vs lme covariance structures (Joshua Wiley) 2. Re: Nested subject-longitudinal logit design (arun) 3. Re: Extracting variances of the estimated variance components in lme4 (Joshua Wiley) ---------------------------------------------------------------------- Message: 1 Date: Thu, 3 May 2012 13:59:15 -0700 From: Joshua Wiley <jwiley.psych at gmail.com> To: Charles Determan Jr <deter088 at umn.edu> Cc: r-sig-mixed-models at r-project.org Subject: Re: [R-sig-ME] gls vs lme covariance structures Message-ID: <CANz9Z_JYcSput9=kB8ibkvooKcABZWfTRDm4QZU9tmYYcxnOHQ at mail.gmail.com> Content-Type: text/plain; charset=UTF-8 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/ [...snip...] ------------------------------ _______________________________________________ R-sig-mixed-models mailing list R-sig-mixed-models at r-project.org https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models End of R-sig-mixed-models Digest, Vol 65, Issue 12 ************************************************** Notice: This e-mail message, together with any attachme...{{dropped:11}}