[R-sig-ME] Another case of -1.0 correlation of random effects
Kevin E. Thorpe
kevin.thorpe at utoronto.ca
Fri Apr 9 14:35:29 CEST 2010
Ben Bolker wrote:
> Kevin E. Thorpe wrote:
>> I know this has come up a couple times recently, but I'm still not sure
>> what to do about it in my data. Note that my sessionInfo() will be at
>> the bottom.
>> My data come from a crossover trial and are balanced.
>> > str(gluc)
>> 'data.frame': 96 obs. of 4 variables:
>> $ Subject : int 1 2 3 5 6 7 10 11 12 13 ...
>> $ Treatment: Factor w/ 2 levels "Barley","Oat": 1 1 1 1 1 1 1 1 1 1 ...
>> $ Dose : int 8 8 8 8 8 8 8 8 8 8 ...
>> $ iAUC : num 110 256 129 207 244 ...
>> > xtabs(~Treatment+Dose,data=gluc)
>> Treatment 0 2 4 8
>> Barley 12 12 12 12
>> Oat 12 12 12 12
>> I plot the data (attached as gluc.pdf, if it comes through).
>> From the plot, I think I want to fit the model as:
>> It could possibly be argued that the (Treatment|Subject) part is not
>> needed. When I fit this, I got -1.0 correlation within the Dose random
>> effects. To simplify, I will fit a simpler model, since the issue persists.
>> > lmer(iAUC~Dose+(Dose|Subject),data=gluc,subset=Treatment=="Oat")
>> Linear mixed model fit by REML
>> Formula: iAUC ~ Dose + (Dose | Subject)
>> Data: gluc
>> Subset: Treatment == "Oat"
>> AIC BIC logLik deviance REMLdev
>> 562.6 573.9 -275.3 563.1 550.6
>> Random effects:
>> Groups Name Variance Std.Dev. Corr
>> Subject (Intercept) 8274.324 90.9633
>> Dose 16.214 4.0266 -1.000
>> Residual 4862.319 69.7303
>> Number of obs: 48, groups: Subject, 12
>> Fixed effects:
>> Estimate Std. Error t value
>> (Intercept) 309.352 30.539 10.130
>> Dose -14.424 3.596 -4.012
>> Correlation of Fixed Effects:
>> Dose -0.647
>> Now, a plot created by (and attached as lmlist.pdf):
>> shows (I think) a strong negative correlation between the intercept and
>> slope random effects for Dose.
>> So, I would appreciate some advice on how I might specify these random
>> effects correctly.
> It's not entirely clear whether Subjects are unique within Treatments
> (OK, maybe we should say "nested") or whether the same Subject can have
> both Treatments -- what would xtabs(~Subject+Dose+Treatment) look like?
> *If* each Subject appears in only one Treatment (i.e., 5 doses per
> Subject, each Subject is either Oat or Barley, 12 Subjects per
> Treatment), which is what I initially took for granted, then it makes
> perfect sense to me that you can't specify Treatment|Subject, because
> Treatment does not vary within Subject -- there doesn't seem to be any
> way you can recover information about how the effect of Treatment varies
> among Subject.
> Since your plots do look reasonably linear within Subject,
> seems best to me.
> For what it's worth, this example closely parallels the Orthodont
> example in the nlme package ...
It is a true crossover. There are 12 patients and each patient gets
both treatments and all four doses.
I will look at the Orthodont example too.
Kevin E. Thorpe
Biostatistician/Trialist, Knowledge Translation Program
Assistant Professor, Dalla Lana School of Public Health
University of Toronto
email: kevin.thorpe at utoronto.ca Tel: 416.864.5776 Fax: 416.864.3016
More information about the R-sig-mixed-models