Hi Dylan, Chuck, >> contrast(l, a=list(f=levels(d$f)[1:3], x=0), b=list(f=levels(d$f)[4], >> x=0))
There is a subtlety here that needs to be emphasized. Setting the interacting variable (x) to zero is reasonable in this case, because the mean value of rnorm(n) is zero. However, in the real world of real data it is somewhat unusual for zero to be a reasonable value for an interacting variable. In fact, it is this setting-to-zero that causes multcomp's "bells to ring" with a word of warning. ## glht(l, linfct=mcp(f="Seq"))$linfct Regards, Mark. Chuck Cleland wrote: > > On 2/16/2009 10:18 PM, Dylan Beaudette wrote: >> On Mon, Feb 16, 2009 at 5:28 PM, Patrick Giraudoux >> <patrick.giraud...@univ-fcomte.fr> wrote: >>> Greg Snow a écrit : >>>> One approach is to create your own contrasts matrix: >>>> >>>> >>>>> mycmat <- diag(8) >>>>> mycmat[ row(mycmat) == col(mycmat) + 1 ] <- -1 >>>>> mycmati <- solve(mycmat) >>>>> contrasts(agefactor) <- mycmati[,-1] >>>>> >>>> Now when you use agefactor, the intercept will be the first age group >>>> and the slopes will be the differences between the pairs of groups >>>> (make sure that the order of the levels of agefactor is correct). >>>> >>>> The difference between this method and the contr.sdif function in MASS >>>> is how the intercept will end up being interpreted (and the dimnames). >>>> >>>> Hope this helps, >>>> >>>> >>> Actually, exactly what I needed including the reference to contr.sdif in >>> MASS I did not spot before (although I am a faithful reader of the >>> yellow book... but so many things still escape to me). Again thanks a >>> lot. >>> >>> Patrick >>> >> >> Glad you were able to solve your problem. Frank replied earlier with >> the suggestion to use contrast.Design() to perform the tests after the >> initial model has been fit. Perhaps I am a little denser than normal, >> but I cannot figure out how to apply contrast.Design() to a similar >> model with several levels of a factor. Example data: >> >> # need these >> library(lattice) >> library(Design) >> >> # replicate an important experimental dataset >> set.seed(10101010) >> x <- rnorm(100) >> y1 <- x[1:25] * 2 + rnorm(25, mean=1) >> y2 <- x[26:50] * 2.6 + rnorm(25, mean=1.5) >> y3 <- x[51:75] * 2.9 + rnorm(25, mean=5) >> y4 <- x[76:100] * 3.5 + rnorm(25, mean=5.5) >> >> d <- data.frame(x=x, y=c(y1,y2,y3,y4), f=factor(rep(letters[1:4], >> each=25))) >> >> # plot >> xyplot(y ~ x, groups=f, data=d, auto.key=list(columns=4, title='Beard >> Type', lines=TRUE, points=FALSE, cex=0.75), type=c('p','r'), >> ylab='Number of Pirates', xlab='Distance from Land') >> >> # standard comparison to base level of f >> summary(lm(y ~ x * f, data=d)) >> >> # compare to level 4 of f >> summary(lm(y ~ x * C(f, base=4), data=d)) >> >> # now try with contrast.Design(): >> dd <- datadist(d) >> options(datadist='dd') >> l <- ols(y ~ x * f, data=d) >> >> # probably the wrong approach, as the results do not look too familiar: >> contrast(l, list(f=levels(d$f))) >> >> x f Contrast S.E. Lower Upper t Pr(>|t|) >> -0.3449623 a 0.3966682 0.1961059 0.007184856 0.7861515 2.02 0.0460 >> -0.3449623 b 0.5587395 0.1889173 0.183533422 0.9339456 2.96 0.0039 >> -0.3449623 c 4.1511677 0.1889170 3.775962254 4.5263730 21.97 0.0000 >> -0.3449623 d 4.3510407 0.1888820 3.975904726 4.7261766 23.04 0.0000 >> >> This is adjusting the output to a given value of 'x'... Am I trying to >> use contrast.Design() for something that it was not intended for? Any >> tips Frank? > > Does this help? > > # Compare with summary(lm(y ~ x * f, data=d)) > contrast(l, a=list(f=levels(d$f)[2:4], x=0), > b=list(f=levels(d$f)[1], x=0)) > > f Contrast S.E. Lower Upper t Pr(>|t|) > b 0.3673455 0.2724247 -0.1737135 0.9084046 1.35 0.1808 > c 4.1310015 0.2714011 3.5919754 4.6700275 15.22 0.0000 > d 4.4308653 0.2731223 3.8884207 4.9733098 16.22 0.0000 > > Error d.f.= 92 > > # Compare with summary(lm(y ~ x * C(f, base=4), data=d)) > contrast(l, a=list(f=levels(d$f)[1:3], x=0), > b=list(f=levels(d$f)[4], x=0)) > > f Contrast S.E. Lower Upper t Pr(>|t|) > a -4.4308653 0.2731223 -4.9733098 -3.8884207 -16.22 0.0000 > b -4.0635198 0.2782704 -4.6161888 -3.5108507 -14.60 0.0000 > c -0.2998638 0.2772684 -0.8505427 0.2508151 -1.08 0.2823 > > Error d.f.= 92 > >> Cheers, >> Dylan >> >> ______________________________________________ >> 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. > > -- > Chuck Cleland, Ph.D. > NDRI, Inc. (www.ndri.org) > 71 West 23rd Street, 8th floor > New York, NY 10010 > tel: (212) 845-4495 (Tu, Th) > tel: (732) 512-0171 (M, W, F) > fax: (917) 438-0894 > > ______________________________________________ > 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. > > -- View this message in context: http://www.nabble.com/Comparison-of-age-categories-using-contrasts-tp22031426p22057400.html Sent from the R help mailing list archive at Nabble.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.