Dear all,
I'm trying to make multiple comparisons for an lme-object. The data is for
an experiment on parental work load in birds, in which adults at different
sites were induced to work at one of three levels ('treat'; H, M, L).
The
response is 'feedings', which is a quantitative measure of nest
provisioning
per parent per chick per hour. Site is included as a random effect (one
male/female pair per site).
My final model takes the following form:
feedings ~ treat + year + data^2, random = ~1|site,data=feed.df
For this model, I would like to do multiple comparisons on 'treat',
using
the multcomp package:
summary(glht(m4.feed,linfct=mcp(treat="Tukey")))
However, this does not work, and I get the below error message.
Error in if (is.null(pkg) | pkg == "nlme") terms(formula(x)) else
slot(x, :
argument is of length zero
Error in factor_contrasts(model) :
no ?model.matrix? method for ?model? found!
I suspect this might have quite a straightforward solution, but I'm stuck at
this point.
Any help would be most appreciated. Sample data below.
Kind regards,
Andreas Nord
Sweden
=============feedings sex site treat year date^2
1.8877888 M 838 H 2009 81
1.9102787 M 247 H 2009 81
1.4647229 M 674 H 2010 121
1.4160590 M 7009 M 2009 144
1.3106749 M 863 M 2010 196
1.2718121 M 61 M 2009 225
1.2799263 M 729 L 2009 256
1.5829564 M 629 L 2009 256
1.4847251 M 299 L 2010 324
1.2463151 M 569 L 2010 324
2.1694169 F 838 H 2009 81
1.5966899 F 247 H 2009 81
2.4136983 F 674 H 2010 121
1.7784873 F 7009 M 2009 144
1.6681317 F 863 M 2010 196
2.3691275 F 61 M 2009 225
2.0672192 F 729 L 2009 256
1.6389902 F 629 L 2009 256
0.9307536 F 299 L 2010 324
1.6786767 F 569 L 2010 324
=============
--
View this message in context:
http://r.789695.n4.nabble.com/Problems-with-glht-function-for-lme-object-tp3179128p3179128.html
Sent from the R help mailing list archive at Nabble.com.
anord <andreas.nord <at> zooekol.lu.se> writes:> > > Dear all, > > I'm trying to make multiple comparisons for an lme-object. The data is for > an experiment on parental work load in birds, in which adults at different > sites were induced to work at one of three levels ('treat'; H, M, L). The > response is 'feedings', which is a quantitative measure of nest provisioning > per parent per chick per hour. Site is included as a random effect (one > male/female pair per site). >(For more complicated mixed model questions you might want to try the r-sig-mixed-models list instead.) It works for me: feed.df <- read.table(textConnection(" feedings sex site treat year date 1.8877888 M 838 H 2009 81 1.9102787 M 247 H 2009 81 1.4647229 M 674 H 2010 121 1.4160590 M 7009 M 2009 144 1.3106749 M 863 M 2010 196 1.2718121 M 61 M 2009 225 1.2799263 M 729 L 2009 256 1.5829564 M 629 L 2009 256 1.4847251 M 299 L 2010 324 1.2463151 M 569 L 2010 324 2.1694169 F 838 H 2009 81 1.5966899 F 247 H 2009 81 2.4136983 F 674 H 2010 121 1.7784873 F 7009 M 2009 144 1.6681317 F 863 M 2010 196 2.3691275 F 61 M 2009 225 2.0672192 F 729 L 2009 256 1.6389902 F 629 L 2009 256 0.9307536 F 299 L 2010 324 1.6786767 F 569 L 2010 324"), header=TRUE) library(nlme) m4.feed <- lme(feedings ~ treat + year + date, random = ~1|site,data=feed.df) library(multcomp) gg <- glht(m4.feed,linfct=mcp(treat="Tukey")) plot(gg) summary(gg) It works for me (although it doesn't look like there's anything going on there ...) ======== Simultaneous Tests for General Linear Hypotheses Multiple Comparisons of Means: Tukey Contrasts Fit: lme.formula(fixed = feedings ~ treat + year + date, data = feed.df, random = ~1 | site) Linear Hypotheses: Estimate Std. Error z value Pr(>|z|) L - H == 0 -0.6452 0.7667 -0.842 0.592 M - H == 0 -0.3996 0.4276 -0.934 0.529 M - L == 0 0.2456 0.4229 0.581 0.771 (Adjusted p values reported -- single-step method)>sessionInfo() R version 2.12.1 (2010-12-16) Platform: i486-pc-linux-gnu (32-bit) locale: [1] LC_CTYPE=en_CA.UTF-8 LC_NUMERIC=C [3] LC_TIME=en_CA.UTF-8 LC_COLLATE=en_CA.UTF-8 [5] LC_MONETARY=C LC_MESSAGES=en_CA.UTF-8 [7] LC_PAPER=en_CA.UTF-8 LC_NAME=C [9] LC_ADDRESS=C LC_TELEPHONE=C [11] LC_MEASUREMENT=en_CA.UTF-8 LC_IDENTIFICATION=C attached base packages: [1] splines stats graphics grDevices utils datasets methods [8] base other attached packages: [1] multcomp_1.2-4 survival_2.36-3 mvtnorm_0.9-95 nlme_3.1-97 loaded via a namespace (and not attached): [1] grid_2.12.1 lattice_0.19-18>
On Jan 7, 2011, at 8:28 AM, anord wrote:> > Dear all, > > I'm trying to make multiple comparisons for an lme-object. The data > is for > an experiment on parental work load in birds, in which adults at > different > sites were induced to work at one of three levels ('treat'; H, M, > L). The > response is 'feedings', which is a quantitative measure of nest > provisioning > per parent per chick per hour. Site is included as a random effect > (one > male/female pair per site). > > My final model takes the following form: > feedings ~ treat + year + data^2, random = ~1|site,data=feed.dfI am guessing that problems could easily arise as a result of your variable name containing "^". That is an invalid name in an un-back- quoted state. You have clearly not given us a copy of a console session since your variable is date^2 below and data^2 above. -- David.> > For this model, I would like to do multiple comparisons on 'treat', > using > the multcomp package: > summary(glht(m4.feed,linfct=mcp(treat="Tukey"))) > > However, this does not work, and I get the below error message. > > Error in if (is.null(pkg) | pkg == "nlme") terms(formula(x)) else > slot(x, : > argument is of length zero > Error in factor_contrasts(model) : > no ?model.matrix? method for ?model? found! > > I suspect this might have quite a straightforward solution, but I'm > stuck at > this point. > Any help would be most appreciated. Sample data below. > > Kind regards, > Andreas Nord > Sweden > > =============> feedings sex site treat year date^2 > 1.8877888 M 838 H 2009 81 > 1.9102787 M 247 H 2009 81 > 1.4647229 M 674 H 2010 121 > 1.4160590 M 7009 M 2009 144 > 1.3106749 M 863 M 2010 196 > 1.2718121 M 61 M 2009 225 > 1.2799263 M 729 L 2009 256 > 1.5829564 M 629 L 2009 256 > 1.4847251 M 299 L 2010 324 > 1.2463151 M 569 L 2010 324 > 2.1694169 F 838 H 2009 81 > 1.5966899 F 247 H 2009 81 > 2.4136983 F 674 H 2010 121 > 1.7784873 F 7009 M 2009 144 > 1.6681317 F 863 M 2010 196 > 2.3691275 F 61 M 2009 225 > 2.0672192 F 729 L 2009 256 > 1.6389902 F 629 L 2009 256 > 0.9307536 F 299 L 2010 324 > 1.6786767 F 569 L 2010 324 > =============> > > -- > View this message in context: http://r.789695.n4.nabble.com/Problems-with-glht-function-for-lme-object-tp3179128p3179128.html > Sent from the R help mailing list archive at Nabble.com. > > ______________________________________________ > R-help at r-project.org mailing list > 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.David Winsemius, MD West Hartford, CT
anord wrote:> > I'm trying to make multiple comparisons for an lme-object. > ... > > My final model takes the following form: > feedings ~ treat + year + data^2, random = ~1|site,data=feed.df > > ... > =============> feedings sex site treat year date^2 > 1.8877888 M 838 H 2009 81 > 1.9102787 M 247 H 2009 81 >1) Try I(date^2) 2) Check the version you are using. For some older versions, multcomp and lme were not close friends (even if that is some time ago now) Dieter -- View this message in context: http://r.789695.n4.nabble.com/Problems-with-glht-function-for-lme-object-tp3179128p3205045.html Sent from the R help mailing list archive at Nabble.com.
Dear all, Thanks for your input. It works fine now. All it took was for me to tidy up the workspace. Simple solution, that I should have considered earlier. Regards, Andreas -- View this message in context: http://r.789695.n4.nabble.com/Problems-with-glht-function-for-lme-object-tp3179128p3206686.html Sent from the R help mailing list archive at Nabble.com.