[R-sig-ME] Fwd: Linear Mixed-Effects Models - lme() command
ONKELINX, Thierry
Thierry.ONKELINX at inbo.be
Thu Aug 8 15:28:37 CEST 2013
Dear Domenico,
I think you want lme(Values~0 + time:Cat1+time:Cat2, random=~1|ID, data=D)
ir. Thierry Onkelinx
Instituut voor natuur- en bosonderzoek / Research Institute for Nature and Forest
team Biometrie & Kwaliteitszorg / team Biometrics & Quality Assurance
Kliniekstraat 25
1070 Anderlecht
Belgium
+ 32 2 525 02 51
+ 32 54 43 61 85
Thierry.Onkelinx op 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-sig-mixed-models-bounces op r-project.org [mailto:r-sig-mixed-models-bounces op r-project.org] Namens Domenico Di Carlo
Verzonden: dinsdag 6 augustus 2013 18:44
Aan: R-SIG-Mixed-Models op r-project.org
Onderwerp: [R-sig-ME] Fwd: Linear Mixed-Effects Models - lme() command
---------- Forwarded message ----------
From: Domenico Di Carlo <pvt.math82 op gmail.com>
Date: 2013/8/6
Subject: Linear Mixed-Effects Models - lme() command
To: bates op stat.wisc.edu
Dear Mr. Douglas Bates,
I have some questions about Linear Fixed Effects Models and about *lme()*command of
*nlme* library of *R*.
This is dataset (called D) structure
*'data.frame': 187 obs. of 5 variables:*
* $ ID : num 2492 2492 2492 2492 2492 ...*
* $ Cat1 : Factor w/ 3 levels "One","three",..: 2 2 2 2 2 2 2 2 2 2 ...*
* $ Cat2 : Factor w/ 3 levels "first","second",..: 1 1 1 1 1 1 1 1 1 1 ...*
* $ Values: num 253 320 424 432 476 434 527 476 468 544 ...*
* $ time : num 0 1 1 1 1 1 2 2 2 3 ...*
The aim is to undestand *Values* changes in *Cat1* and *Cat2* levels across the *time*. The line command I have used is:
*> library(nlme)*
*> model_1<-lme(Values~time:Cat1+time:Cat2,*
*+ random=~1|ID, data=D)*
*> summary(model_1)*
*Linear mixed-effects model fit by REML*
* Data: D *
* AIC BIC logLik*
* 2295.15 2320.738 -1139.575*
*
*
*Random effects:*
* Formula: ~1 | ID*
* (Intercept) Residual*
*StdDev: 107.022 114.4641*
*
*
*Fixed effects: Values ~ time:Cat1 + time:Cat2 *
* Value Std.Error DF t-value p-value*
*(Intercept) 324.0124 43.49774 175 7.448947 0.0000*
*time:Cat1One 51.4522 17.16804 175 2.996976 0.0031*
*time:Cat1three 47.0102 8.56644 175 5.487712 0.0000*
*time:Cat1two 69.5621 15.68291 175 4.435532 0.0000*
*time:Cat2[T.second] -24.4386 13.75771 175 -1.776360 0.0774* *time:Cat2[T.third] -21.3869 19.35976 175 -1.104710 0.2708*
* Correlation: *
* (Intr) tm:C1O tm:Ct1th tm:Ct1tw tm:Ct2[T.s]*
*time:Cat1One -0.084 *
*time:Cat1three -0.191 0.497 *
*time:Cat1two -0.104 0.711 0.546 *
*time:Cat2[T.second] 0.026 -0.803 -0.605 -0.879 *
*time:Cat2[T.third] 0.018 -0.882 -0.430 -0.625 0.711 *
*
*
*Standardized Within-Group Residuals:*
* Min Q1 Med Q3 Max *
*-3.0730790 -0.4341496 0.0301899 0.4023551 5.4252466 *
*
*
*Number of Observations: 187*
*Number of Groups: 7 *
About *time:Cat1One*, *time:Cat1three* and *time:Cat1two*, I guess (maybe you can correct me) they are slopes of *Cat1* levels, I mean the average increase of *Value* year after year. About *time:Cat2[T.second]* and * time:Cat[T.third]*, I guess they are differences from slope of dummy of * Cat2*, I mean *first* level of the variable. But to understand the real meaning of this difference, I need to see the real value of * time:Cat[T.first]*. Is it possible to see it? I was thinking that maybe this value is the *intercept*, but its value *(324.0124)* is not the estimate I was waiting for.
I have fitted a second model:
*> model_2<-lme(Values~time*Cat1+time*Cat2,*
*+ random=~1|ID, data=D)*
*> summary(model_2)*
*Linear mixed-effects model fit by REML*
* Data: D *
* AIC BIC logLik*
* 2253.048 2291.162 -1114.524*
*
*
*Random effects:*
* Formula: ~1 | ID*
* (Intercept) Residual*
*StdDev: 100.3245 114.4606*
*
*
*Fixed effects: Values ~ time * Cat1 + time * Cat2 *
* Value Std.Error DF t-value p-value*
*(Intercept) 200.91095 171.01061 175 1.1748450 0.2417*
*time 54.11604 17.93540 175 3.0172746 0.0029*
*Cat1[T.three] 150.09485 151.92072 2 0.9879814 0.4273*
*Cat1[T.two] 270.12004 131.59306 2 2.0526921 0.1765*
*Cat2[T.second] -44.75005 133.55838 2 -0.3350599 0.7695*
*Cat2[T.third] 58.30254 202.25915 2 0.2882566 0.8003*
*time:Cat1[T.three] -8.11917 15.53307 175 -0.5227025 0.6018*
*time:Cat1[T.two] 10.91819 13.08954 175 0.8341161 0.4054*
*time:Cat2[T.second] -22.90242 14.49702 175 -1.5798018 0.1160* *time:Cat2[T.third] -22.43406 20.29314 175 -1.1054997 0.2705*
* Correlation: *
* (Intr) time Ct1[T.th] Ct1[T.tw] Ct2[T.s] Ct2[T.t]*
*time -0.314 *
*Cat1[T.three] -0.888 0.259 *
*Cat1[T.two] -0.507 0.136 0.571 *
*Cat2[T.second] -0.781 0.268 0.575 0.000 *
*Cat2[T.third] -0.846 0.265 0.751 0.429 0.660 *
*time:Cat1[T.three] 0.265 -0.866 -0.299 -0.157 -0.185 -0.224 *
*time:Cat1[T.two] 0.143 -0.475 -0.161 -0.298 0.000 -0.121 *
*time:Cat2[T.second] 0.259 -0.808 -0.175 0.000 -0.332 -0.219 *
*time:Cat2[T.third] 0.277 -0.884 -0.229 -0.120 -0.237 -0.314 *
* tm:Ct1[T.th] tm:Ct1[T.tw] tm:Ct2[T.s]*
*time *
*Cat1[T.three] *
*Cat1[T.two] *
*Cat2[T.second] *
*Cat2[T.third] *
*time:Cat1[T.three] *
*time:Cat1[T.two] 0.548 *
*time:Cat2[T.second] 0.576 0.000 *
*time:Cat2[T.third] 0.765 0.420 0.714 *
*
*
*Standardized Within-Group Residuals:*
* Min Q1 Med Q3 Max *
*-3.19403849 -0.41767624 0.03823682 0.44202369 5.35682026 *
*
*
*Number of Observations: 187*
*Number of Groups: 7 *
In this last model, I guess interactions terms differences are related to the term *time*, but I do not understand if this term is the dummy of *Cat1*or *Cat2*.
I am sorry if I have written something wrong, I hope you can help me to understand my doubts. I am sending to you all the materials.
Best regards
Domenico Di Carlo
* * * * * * * * * * * * * D I S C L A I M E R * * * * * * * * * * * * *
Dit bericht en eventuele bijlagen geven enkel de visie van de schrijver weer en binden het INBO onder geen enkel beding, zolang dit bericht niet bevestigd is door een geldig ondertekend document.
The views expressed in this message and any annex are purely those of the writer and may not be regarded as stating an official position of INBO, as long as the message is not confirmed by a duly signed document.
More information about the R-sig-mixed-models
mailing list