Hi Walmes, The problem goes back to the singularity caused by the incomplete design. The summary() and coef() methods for aov class objects hide it rather than printing NA as with lm(), but look at just:
m1 it says at the end "1 out of 6 effects not estimable". Now when you use glht() you pass a 1 x 5 matrix, but there are really *6* effects. glht drops inestimable parameters, and it knows the 6th should be dropped, so it is trying to subset your linear hypothesis matrix to only those effects that are estimable. The code basically does: > c(rep(TRUE, 5), FALSE) # TRUE being estimable effects [1] TRUE TRUE TRUE TRUE TRUE FALSE > matrix(c(1,1,0,10,0), nrow = 1) # your hypothesis matrix [,1] [,2] [,3] [,4] [,5] [1,] 1 1 0 10 0 > matrix(c(1,1,0,10,0), nrow = 1)[, c(rep(TRUE, 5), FALSE)] # how glht tries to > subset it Error: (subscript) logical subscript too long # the error you get The solution is relatively straight forward, just include a value in your H matrix for the effect that cannot be estimated. For example: glht(m1, linfct=matrix(c(1,1,0,10,0,0), nrow=1)) Cheers, Josh On Sat, Aug 6, 2011 at 7:48 AM, Walmes Zeviani <walmeszevi...@gmail.com> wrote: > Hi R users, > > I sent a message yesterday about NA in model estimates ( > http://r.789695.n4.nabble.com/How-set-lm-to-don-t-return-NA-in-summary-td3722587.html). > If I use aov() instead of lm() I get no NA in model estimates and I use > gmodels::estimable() without problems. Ok! > Now I'm performing a lot of contrasts and I need correcting for > multiplicity. So, I can use multcomp::glht() for this. However, glht() > return an error message that is not compatible with my expectations. Someone > know I or has a suggestion for? Below some reproducible code. > > # toy data > adi <- expand.grid(cult=gl(1,3,la=LETTERS[1]), fert=101) > fat <- expand.grid(cult=gl(2,3,la=LETTERS[2:3]), fert=seq(50,150,50)) > da <- rbind(adi, fat) > da$y <- da$fert+rnorm(nrow(da),0,10) > > # plot > require(lattice) > xyplot(y~fert|cult, da) > > # adjust complete factorial > m0 <- aov(y~cult*fert, subset(da, cult!="A")) > summary(m0) > coef(m0) > > # adjust incomplete factorial > m1 <- aov(y~cult*fert, da) > summary(m1) > coef(m1) > > require(multcomp) > glht(m0, linfct=matrix(c(1,1,10,0), nrow=1)) # work > glht(m1, linfct=matrix(c(1,1,0,10,0), nrow=1)) # don't work > > Thank you. > Walmes. > > ========================================================================== > Walmes Marques Zeviani > LEG (Laboratório de Estatística e Geoinformação, 25.450418 S, 49.231759 W) > Departamento de Estatística - Universidade Federal do Paraná > fone: (+55) 41 3361 3573 > VoIP: (3361 3600) 1053 1173 > e-mail: wal...@ufpr.br > twitter: @walmeszeviani > homepage: http://www.leg.ufpr.br/~walmes > linux user number: 531218 > ========================================================================== > > [[alternative HTML version deleted]] > ______________________________________________ > R-help@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. -- Joshua Wiley Ph.D. Student, Health Psychology Programmer Analyst II, ATS Statistical Consulting Group University of California, Los Angeles https://joshuawiley.com/ ______________________________________________ R-help@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.