Thanks Thierry for the reply. I think I now have a better understanding for the specification of the random effects when using lme function. Are my interpretations below correct?
random=~ 1 | individual (same random intercept no random slope) random=~ 1 +method| individual (same random intercept and same random slope) random=~ 1 +method:time| individual (same random intercept and different random slope for different method) random=~ 1 +method + method:time| individual (different random intercept and different random slope for different method The summary results from the lme function shows whether the slopes for the three methods are equal (parallelism). I also wanted to test the hypotheses that each of the fixed slopes (corresponding to the three methods) equals 0, can I use multicomp package for that purpose? I am confused on how to make correct specifications in glht function to test these hypotheses. Hanna > summary(mod1) Linear mixed-effects model fit by REML Data: one AIC BIC logLik 304.4703 330.1879 -140.2352 Random effects: Formula: ~1 + time | individual Structure: General positive-definite, Log-Cholesky parametrization StdDev Corr (Intercept) 0.2487869075 (Intr) time 0.0001841179 -0.056 Residual 0.3718305953 Variance function: Structure: Different standard deviations per stratum Formula: ~1 | method Parameter estimates: 3 1 2 1.00000 26.59750 24.74476 Fixed effects: reponse ~ method * time Value Std.Error DF t-value p-value(Intercept) 96.65395 3.528586 57 27.391694 0.0000 method2 1.17851 4.856026 57 0.242689 0.8091 method3 5.87505 3.528617 57 1.664973 0.1014 *time 0.07010 0.250983 57 0.279301 0.7810 method2:time -0.12616 0.360585 57 -0.349877 0.7277 method3:time -0.08010 0.251105 57 -0.318999 0.7509* Correlation: (Intr) methd2 methd3 time mthd2: method2 -0.726 method3 -0.999 0.726 time -0.779 0.566 0.779 method2:time 0.542 -0.712 -0.542 -0.696 method3:time 0.778 -0.566 -0.779 -0.999 0.696 Standardized Within-Group Residuals: Min Q1 Med Q3 Max -2.67575293 -0.51633192 0.06742723 0.59706762 2.81061874 Number of Observations: 69 Number of Groups: 7 > ---------- Forwarded message ---------- From: Thierry Onkelinx <thierry.onkel...@inbo.be> Date: 2016-05-30 4:40 GMT-04:00 Subject: Re: [R] model specification using lme To: li li <hannah....@gmail.com> Cc: r-help <r-help@r-project.org> Dear Hanna, None of the models are correct is you want the same random intercept for the different methods but different random slope per method. You can random = ~ 1 + time:method | individual The easiest way to get alpha_0 and tau_i is to apply post-hoc contrasts. That is fairly easy to do with the multcomp package. alpha_0 = (m1 + m2 + m3) / 3 m1 = intercept m2 = intercept + method2 m3 = intercept + method3 hence alpha_0 = intercept + method2/3 + method3/3 m1 = alpha_0 + tau_1 tau_1 = intercept - method2/3 - method3/3 Best regards, 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 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 2016-05-29 21:23 GMT+02:00 li li <hannah....@gmail.com>: > Hi all, > For the following data, I consider the following random intercept and > random slope model. Denote as y_ijk the response value from *j*th > individual within *i*th method at time point *k*. Assume the following > model for y_ijk: > > y_ijk= (alpha_0+ tau_i +a_j(i))+(beta_i+b_j(i)) T_k + e_ijk > > > Here alpha_0 is the grand mean; > tau_i is the fixed effect for ith method; > a_j(i) is random intercept corresponding to the *j*th individual > within *i*th method, assumed to be common for all three methods; > beta_i is the fixed slope corresponding to the ith method; > b_j(i) is the random slope corresponding to jth individual for > the ith method, assumed to be different for different methods; > T_k is the time corresponding to y_ijk; > e_ijk is the residual. > > For this model, I consider the three specification using the lme function > as follows: > > > mod1 <- lme(fixed= reponse ~ method*time, random=~ 1 +time | individual, > data=one, weights= varIdent(form=~1|method), > control = lmeControl(opt = "optim")) > > mod2 <- lme(fixed= reponse ~ method*time, random=~ 0 +time | individual, > data=one, weights= varIdent(form=~1|method), > control = lmeControl(opt = "optim")) > > mod3 <- lme(fixed= reponse ~ method*time, random=~ method +time | > individual, data=one, weights= varIdent(form=~1|method), > control = lmeControl(opt = "optim")) > > I think mod1 is the correct one. However, I am kind of confused with the > right usage of lme function. Can someone familiar with this give some help > here? > > Another question is regarding the fixed effect tau_1, tau_2 and tau_3 > (corresponding to the three methods). One main question I am interested in > is whether each of them are statistically different from zero. In the > summary results below (shaded part), it looks that the result for method 2 > and 3 are given with reference to method 1). Is there a way to obtain > specific result separately for alpha_0 (the overall mean) and also tau_1, > tau_2 and tau3? > > Thanks very much for the help! > Hanna > > > > summary(mod1)Linear mixed-effects model fit by REML > Data: one > AIC BIC logLik > 304.4703 330.1879 -140.2352 > > Random effects: > Formula: ~1 + time | individual > Structure: General positive-definite, Log-Cholesky parametrization > StdDev Corr > (Intercept) 0.2487869075 (Intr) > time 0.0001841179 -0.056 > Residual 0.3718305953 > > Variance function: > Structure: Different standard deviations per stratum > Formula: ~1 | method > Parameter estimates: > 3 1 2 > 1.00000 26.59750 24.74476 > Fixed effects: reponse ~ method * time > Value Std.Error DF t-value p-value(Intercept) > 96.65395 3.528586 57 27.391694 0.0000 > method2 1.17851 4.856026 57 0.242689 0.8091 > method3 5.87505 3.528617 57 1.664973 0.1014time > 0.07010 0.250983 57 0.279301 0.7810 > method2:time -0.12616 0.360585 57 -0.349877 0.7277 > method3:time -0.08010 0.251105 57 -0.318999 0.7509 > Correlation: > (Intr) methd2 methd3 time mthd2: > method2 -0.726 > method3 -0.999 0.726 > time -0.779 0.566 0.779 > method2:time 0.542 -0.712 -0.542 -0.696 > method3:time 0.778 -0.566 -0.779 -0.999 0.696 > > Standardized Within-Group Residuals: > Min Q1 Med Q3 Max > -2.67575293 -0.51633192 0.06742723 0.59706762 2.81061874 > > Number of Observations: 69 > Number of Groups: 7 > > > > > > > > one response individual time method > 1 102.9 3 0 3 > 2 103.0 3 3 3 > 3 103.0 3 6 3 > 4 102.8 3 9 3 > 5 102.2 3 12 3 > 6 102.5 3 15 3 > 7 103.0 3 18 3 > 8 102.0 3 24 3 > 9 102.8 1 0 3 > 10 102.7 1 3 3 > 11 103.0 1 6 3 > 12 102.2 1 9 3 > 13 103.0 1 12 3 > 14 102.8 1 15 3 > 15 102.8 1 18 3 > 16 102.9 1 24 3 > 17 102.2 2 0 3 > 18 102.6 2 3 3 > 19 103.4 2 6 3 > 20 102.3 2 9 3 > 21 101.3 2 12 3 > 22 102.1 2 15 3 > 23 102.1 2 18 3 > 24 102.2 2 24 3 > 25 102.7 4 0 3 > 26 102.3 4 3 3 > 27 102.6 4 6 3 > 28 102.7 4 9 3 > 29 102.8 4 12 3 > 30 102.5 5 0 3 > 31 102.4 5 3 3 > 32 102.1 5 6 3 > 33 102.3 6 0 3 > 34 102.3 6 3 3 > 35 101.9 7 0 3 > 36 102.0 7 3 3 > 37 107.4 3 0 1 > 38 101.3 3 12 1 > 39 92.8 3 15 1 > 40 73.7 3 18 1 > 41 104.7 3 24 1 > 42 92.6 1 0 1 > 43 101.9 1 12 1 > 44 106.3 1 15 1 > 45 104.1 1 18 1 > 46 95.6 1 24 1 > 47 79.8 2 0 1 > 48 89.7 2 12 1 > 49 97.0 2 15 1 > 50 108.4 2 18 1 > 51 103.5 2 24 1 > 52 96.4 4 0 1 > 53 89.3 4 12 1 > 54 112.6 5 0 1 > 55 93.3 6 0 1 > 56 99.6 7 0 1 > 57 109.5 3 0 2 > 58 98.5 3 12 2 > 59 103.5 3 24 2 > 60 113.5 1 0 2 > 61 94.5 1 12 2 > 62 88.5 1 24 2 > 63 99.5 2 0 2 > 64 97.5 2 12 2 > 65 98.5 2 24 2 > 66 103.5 4 0 2 > 67 89.5 5 0 2 > 68 87.5 6 0 2 > 69 82.5 7 0 2 > > > [[alternative HTML version deleted]] > > ______________________________________________ > R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see > https://stat.ethz.ch/mailman/listinfo/r-help > PLEASE do read the posting guide > http://www.R-project.org/posting-guide.html > and provide commented, minimal, self-contained, reproducible code. > > [[alternative HTML version deleted]] ______________________________________________ R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.