Dear Thierry,
Thanks a lot for pointing me to the right direction! I still have some
questions, really appreciate if you could provide any help:
Is there a relationship between the nested mixed model that I used vs. the model
that you gave (using biop location as random slope of pid), i.e.
lme(y~age + time * trt, random=~1|pid/biopsy.site, data = test,
correlation=corAR1())
vs.
lme(y~age + time * trt, random=~0? + biopsy.site|pid, data = test,?
correlation=corAR1())
In my nested mixed model, a variance of pid and a variance of biopsy.site within
pid will be estimated. In your mixed model, there is no variance estimated at
the pid level , instead variance for each biopsy.site is given. I thought the
statistical model was always the same for both mixed models:
y_ijk=fixed_effect + b_i + b_ij + e_ijk
where b_i is the random effect for pid, and v_ij is the random effect for
biopsy.site within pid. I thought that the difference is in my nested mixed
model, b_ij is independent of each other and has the same variance, whereas in
your mixed model, b_ij is modeled by the pdClasses() chosen. If my thought was
correct, then my model should be the same as
lme(y~age + time * trt, random=list(pid=pdIdent(~0? + biopsy.site)), data =
test,? correlation=corAR1())
But they are not the same. What did I misunderstand here?
Many thanks
John
----- Original Message -----
From: "ONKELINX, Thierry" <Thierry.ONKELINX at inbo.be>
To: array chip <arrayprofile at yahoo.com>; "r-sig-mixed-models at
r-project.org" <r-sig-mixed-models at r-project.org>
Cc:
Sent: Monday, August 8, 2011 1:19 AM
Subject: RE: [R] mixed model fitting between R and SAS
Please don't cross-post.
- corAR1() models the correlation between the residuals of the two time points.
- if you want a specific correlation structure for biopsy locations the you must
use on of the pdClasses() and use biopsy location as random slope of pid rather
than random effect nested in pid.
#basic structure = positive definitive symmetrical variance/covariance matrix
lme(y~age + time * trt, random=~0? + biopsy.site|pid, data = test,?
correlation=corAR1(~time))
#no correlation between biopsy location and different variance
lme(y~age + time * trt, random=list(pid? = pdDiag(~0? + biopsy.site), data =
test,? correlation=corAR1(~time))
#no correlation between biopsy location and equal variance
lme(y~age + time * trt, random=list(pid? = pdIdent(~0? + biopsy.site), data =
test,? correlation=corAR1(~time))
Note that since you have only two biopsy locations there will be no difference
between pdSymm (the default) and pdCompSymm
Best regards,
Thierry
----------------------------------------------------------------------------
ir. Thierry Onkelinx
Instituut voor natuur- en bosonderzoek
team Biometrie & Kwaliteitszorg
Gaverstraat 4
9500 Geraardsbergen
Belgium
Research Institute for Nature and Forest
team Biometrics & Quality Assurance
Gaverstraat 4
9500 Geraardsbergen
Belgium
tel. + 32 54/436 185
Thierry.Onkelinx at inbo.be
www.inbo.be
To call in the statistician after the experiment is done may be no more than
asking him to perform a post-mortem examination: he may be able to say what the
experiment died of.
~ Sir Ronald Aylmer Fisher
The plural of anecdote is not data.
~ Roger Brinner
The combination of some data and an aching desire for an answer does not ensure
that a reasonable answer can be extracted from a given body of data.
~ John Tukey
> -----Oorspronkelijk bericht-----
> Van: r-help-bounces at r-project.org [mailto:r-help-bounces at
r-project.org]
> Namens array chip
> Verzonden: maandag 8 augustus 2011 8:48
> Aan: r-sig-mixed-models at r-project.org
> CC: r-help
> Onderwerp: [R] mixed model fitting between R and SAS
>
> Hi al,
>
> I have a dataset (see attached), which basically involves 4 treatments for
a
> chemotherapy drug. Samples were taken from 2 biopsy locations, and biopsy
> were taken at 2 time points. So each subject has 4 data points (from 2
biopsy
> locations and 2 time points). The objective is to study treatment
difference.
>
> I used lme to fit a mixed model that uses "biopsy.site nested within
pid" as a
> random term, and used corAR1() as the correlation structure for between the
2
> time points:
>
>
> library(nlme)
>
>
test<-read.table("test.txt",sep='\t',header=T,row.names=1)
> fit<-lme(y~age + time * trt, random=~1|pid/biopsy.site, data = test,
> correlation=corAR1())
>
> First, by above model specification, corAR1() is used for the correlation
between
> the 2 time points; what is the correlation structure implicitly used for
between
> biopsy locations? How do I specify a particular correlation structure for
between
> biopsy locations in this situation?
>
> Second, does anyone know how to write the above mixed model in SAS? One of
> my colleagues wrote the following, but it gave me different results:
>
> proc mixed data=test;
>
> class time trt pid biopsysite;
> model y=age time trt time*trt;
> random biopsysite
> repeated pid / type=ar(1)
> run;
>
> Is there anyone familiar with SAS and know if the above SAS code does what
the
> R code does?
>
> Many thanks
>
> John