[R-sig-ME] How to get estimates of time (Age) as well as Lots in a mixed model
Jesus Frias
Jesus.Frias at dit.ie
Tue Jun 29 20:26:33 CEST 2010
Dear Robert,
I think that your situation is similar to many other mixed effect
problems (e.g. a clinical trial in phase II where there are a lot of
subjects but there is very few information on each individual patient).
In my personal opinion (and I have used R for something very similar to
what you are trying to do) the nlme library is well able to solve your
problem. However you probably will need a lot more data than the
example you have there to have an estimate of the lot random with some
confidence.
In my opinion
fit1<- lme(Conc ~ Age, data=blots, random=~1|Lot)
will do what you want. You have at least two replicates from each lot
with the same Age.
Be mindful however that the lot-to-lot random effect you are
estimating is only at the intercept. All your lots are behaving with
the same time slope. If you want to study if the lot is affecting how
fast the food product degrades you might need a model like:
fit1<- lme(Conc ~ Age, data=blots, random=~Age|Lot)
And that is where your confounding might give you trouble.
Is nice to see that other people in the area of food tech are using
this fantastic tool!
regards,
Jesus
Robert A. LaBudde wrote:
>I have a small dataset that represents measurements across Lots of a
>processed food, where each lot has a particular Age since manufacturer.
>
>The object of the study is twofold: 1) Determine the Lot-Lot random
>effect size, and 2) determine the degradation slope with Age.
>
>Unfortunately, there is only one unique Lot for each particular
>unique Age. (A Lot is all manufacture from a particular date.) So the
>Lot effect is confounded to some extent with the Age effect.
>
>My question is: What is the best way to get estimate of the two
>effects from this dataset? Should I use the Lot effect from a random
>effects model with no Age variable (so the Lot effect will include
>the Age degradation effect), or should I fit a mixed model anyway,
>and let the software sort them out (which it apparently does by
>dropping a level of Lots)?
>
>I happen to be using the "nlme" package, so here is the example (I
>get a little larger Lot standard deviation when I include Age in the
model):
>
> > blots<- read.table('pseudolots.txt', header=TRUE)
> > blots$ID<- factor(blots$ID)
> > blots$Unit<- factor(blots$Unit)
> > blots
> ID Lot Unit Age Conc
>1 1 A 1 3 1.44
>2 2 A 2 3 1.56
>3 3 B 3 41 1.03
>4 4 B 4 41 1.57
>5 5 C 5 229 1.49
>6 6 C 6 229 1.66
>7 7 D 7 238 0.88
>8 8 D 8 238 0.93
>9 9 E 9 349 1.43
>10 10 E 10 349 1.22
>11 11 F 11 391 1.42
>12 12 F 12 391 1.46
> > #Modeling
> > require('nlme')
>Loading required package: nlme
> > fit1<- lme(Conc ~ Age, data=blots, random=~1|Lot)
> > summary(fit1)
>Linear mixed-effects model fit by REML
> Data: blots
> AIC BIC logLik
> 22.74858 23.95893 -7.374292
>
>Random effects:
> Formula: ~1 | Lot
> (Intercept) Residual
>StdDev: 0.2320592 0.1786757
>
>Fixed effects: Conc ~ Age
> Value Std.Error DF t-value p-value
>(Intercept) 1.3710359 0.18970706 6 7.227121 0.0004
>Age -0.0001449 0.00074846 4 -0.193538 0.8560
> Correlation:
> (Intr)
>Age -0.823
>
>Standardized Within-Group Residuals:
> Min Q1 Med Q3 Max
>-1.59441923 -0.45471612 -0.06071212 0.52440223 1.42781630
>
>Number of Observations: 12
>Number of Groups: 6
> > anova(fit1)
> numDF denDF F-value p-value
>(Intercept) 1 6 154.51072 <.0001
>Age 1 4 0.03746 0.856
> >
> > fit2<- lme(Conc ~ 1, data=blots, random=~1|Lot)
> > summary(fit2)
>Linear mixed-effects model fit by REML
> Data: blots
> AIC BIC logLik
> 8.122374 9.31606 -1.061187
>
>Random effects:
> Formula: ~1 | Lot
> (Intercept) Residual
>StdDev: 0.2010265 0.1786757
>
>Fixed effects: Conc ~ 1
> Value Std.Error DF t-value p-value
>(Intercept) 1.340833 0.09693139 6 13.83281 0
>
>Standardized Within-Group Residuals:
> Min Q1 Med Q3 Max
>-1.57582754 -0.56625710 -0.01917496 0.56893404 1.44640790
>
>Number of Observations: 12
>Number of Groups: 6
> > anova(fit2)
> numDF denDF F-value p-value
>(Intercept) 1 6 191.3466 <.0001
>================================================================
>Robert A. LaBudde, PhD, PAS, Dpl. ACAFS e-mail: ral at lcfltd.com
>Least Cost Formulations, Ltd. URL: http://lcfltd.com/
>824 Timberlake Drive Tel: 757-467-0954
>Virginia Beach, VA 23464-3239 Fax: 757-467-2947
>
>"Vere scire est per causas scire"
>
>_______________________________________________
>R-sig-mixed-models at r-project.org mailing list
>https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models
>
--------------------------------------------------
JesÃºs MarÃa FrÃas Celayeta
Head of the Dept of Food Science
School of Food Sci. and Env. Health
DIT
t + 353 1 402 4459
f + 353 1 402 4495
http://fseh.dit.ie/o4/StaffListing/JesusFrias.html
This message has been scanned for content and viruses by the DIT Information Services E-Mail Scanning Service, and is believed to be clean. http://www.dit.ie
More information about the R-sig-mixed-models
mailing list