I've tried without success to calculate interaction and main effects
contrasts using R. I've found the functions C(), contrasts(),
se.contrasts() and fit.contrasts() in package gmodels. Given the url
for a small dataset and the two-way anova model below, I'd like to
reproduce the results from appended SAS code. Thanks. --Dale.
## the dataset (from Montgomery)
twoway <- read.table("http://dsteele.veryspeedy.net/sta501/twoway.txt",
col.names=c('material', 'temp','voltage'),colClasses=c('factor',
'factor', 'numeric'))
## the model
fit <- aov(voltage ~ material*temp, data=twoway)
/* SAS code */
proc glm data=twoway;
class material temp;
model voltage = material temp material*temp;
contrast '21-22-31+32' material*temp 0 0 0 1 -1 0 -1 1 0;
estimate '21-22-31+32' material*temp 0 0 0 1 -1 0 -1 1 0;
contrast 'material1-material2' material 1 -1 0;
estimate 'material1-material2' material 1 -1 0;
contrast 'temp50 - temp80' temp 1 0 -1;
estimate 'temp50 - temp80' temp 1 0 -1;
run;
Finding Interaction and main effects contrasts for two-way ANOVA
5 messages · Dale Steele, Thompson, David (MNR), Gregory R. Warnes +1 more
Dale, Other than the first SAS contrast, does the following demonstrate what your asking for?
summary(twoway)
material temp voltage
1:12 50:12 Min. : 20
2:12 65:12 1st Qu.: 70
3:12 80:12 Median :108
Mean :106
3rd Qu.:142
Max. :188
contrasts(twoway$material)
2 3 1 0 0 2 1 0 3 0 1
contrasts(twoway$temp)
65 80 50 0 0 65 1 0 80 0 1
fit <- aov(voltage ~ material*temp, data=twoway) summary.aov(fit)
Df Sum Sq Mean Sq F value Pr(>F) material 2 10684 5342 7.91 0.0020 ** temp 2 39119 19559 28.97 1.9e-07 *** material:temp 4 9614 2403 3.56 0.0186 * Residuals 27 18231 675 --- Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1 # setting (partial) contrasts
contrasts(twoway$material) <- c(1,-1,0) # ignoring the second
available df
contrasts(twoway$temp) <- c(0,1,-1) # ignoring the second available df contrasts(twoway$material)
[,1] [,2] 1 1 -0.41 2 -1 -0.41 3 0 0.82
contrasts(twoway$temp)
[,1] [,2] 50 0 -0.82 65 1 0.41 80 -1 0.41
summary.aov(fit, split=list(material=list('m1-m2'=1), temp=list('t50 -
t80'=1)))
Df Sum Sq Mean Sq F value Pr(>F)
material 2 10684 5342 7.91 0.00198 **
material: m1-m2 1 3800 3800 5.63 0.02506 *
temp 2 39119 19559 28.97 1.9e-07 ***
temp: t50 - t80 1 11310 11310 16.75 0.00035 ***
material:temp 4 9614 2403 3.56 0.01861 *
material:temp: m1-m2.t50 - t80 1 4970 4970 7.36 0.01146 *
Residuals 27 18231 675
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
# other examples of setting contrasts
# compare m1 vs m2 and m2 vs m3
contrasts(twoway$material) <- matrix(c(1,-1,0,1,1,-2), nrow=3) contrasts(twoway$material)
[,1] [,2] 1 1 0 2 -1 1 3 0 -1 # compare m1 vs m2 and m1+m2 vs m3
contrasts(twoway$material) <- matrix(c(1,-1,0,1,1,-2), nrow=3) contrasts(twoway$material)
[,1] [,2] 1 1 1 2 -1 1 3 0 -2 I'm not sure if 'summary.aov' is the only lm-family summary method with the split argument. DaveT. ************************************* Silviculture Data Analyst Ontario Forest Research Institute Ontario Ministry of Natural Resources david.john.thompson at ontario.ca http://ofri.mnr.gov.on.ca *************************************
-----Original Message-----
From: Steele [mailto:dale.w.steele at gmail.com]
Sent: March 6, 2008 09:08 PM
To: r-help at stat.math.ethz.ch
Subject: [R] Finding Interaction and main effects contrasts
for two-way ANOVA
I've tried without success to calculate interaction and main effects
contrasts using R. I've found the functions C(), contrasts(),
se.contrasts() and fit.contrasts() in package gmodels. Given the url
for a small dataset and the two-way anova model below, I'd like to
reproduce the results from appended SAS code. Thanks. --Dale.
## the dataset (from Montgomery)
twoway <- read.table("http://dsteele.veryspeedy.net/sta501/twoway.txt",
col.names=c('material', 'temp','voltage'),colClasses=c('factor',
'factor', 'numeric'))
## the model
fit <- aov(voltage ~ material*temp, data=twoway)
/* SAS code */
proc glm data=twoway;
class material temp;
model voltage = material temp material*temp;
contrast '21-22-31+32' material*temp 0 0 0 1 -1 0 -1 1 0;
estimate '21-22-31+32' material*temp 0 0 0 1 -1 0 -1 1 0;
contrast 'material1-material2' material 1 -1 0;
estimate 'material1-material2' material 1 -1 0;
contrast 'temp50 - temp80' temp 1 0 -1;
estimate 'temp50 - temp80' temp 1 0 -1;
run;
Dale, You might find it fruitful to look at the help pages for fit.contrast () and estimble() functions in the gmodels package, and the contrast () functions in the Hmisc package. -Greg
On Mar 7, 2008, at 4:20PM , Thompson, David ((MNR)) wrote:
Dale, Other than the first SAS contrast, does the following demonstrate what your asking for?
summary(twoway)
material temp voltage
1:12 50:12 Min. : 20
2:12 65:12 1st Qu.: 70
3:12 80:12 Median :108
Mean :106
3rd Qu.:142
Max. :188
contrasts(twoway$material)
2 3 1 0 0 2 1 0 3 0 1
contrasts(twoway$temp)
65 80 50 0 0 65 1 0 80 0 1
fit <- aov(voltage ~ material*temp, data=twoway) summary.aov(fit)
Df Sum Sq Mean Sq F value Pr(>F) material 2 10684 5342 7.91 0.0020 ** temp 2 39119 19559 28.97 1.9e-07 *** material:temp 4 9614 2403 3.56 0.0186 * Residuals 27 18231 675 --- Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1 # setting (partial) contrasts
contrasts(twoway$material) <- c(1,-1,0) # ignoring the second
available df
contrasts(twoway$temp) <- c(0,1,-1) # ignoring the second available df contrasts(twoway$material)
[,1] [,2] 1 1 -0.41 2 -1 -0.41 3 0 0.82
contrasts(twoway$temp)
[,1] [,2] 50 0 -0.82 65 1 0.41 80 -1 0.41
summary.aov(fit, split=list(material=list('m1-m2'=1), temp=list
('t50 -
t80'=1)))
Df Sum Sq Mean Sq F value Pr(>F)
material 2 10684 5342 7.91 0.00198 **
material: m1-m2 1 3800 3800 5.63 0.02506 *
temp 2 39119 19559 28.97 1.9e-07 ***
temp: t50 - t80 1 11310 11310 16.75 0.00035 ***
material:temp 4 9614 2403 3.56 0.01861 *
material:temp: m1-m2.t50 - t80 1 4970 4970 7.36 0.01146 *
Residuals 27 18231 675
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
# other examples of setting contrasts
# compare m1 vs m2 and m2 vs m3
contrasts(twoway$material) <- matrix(c(1,-1,0,1,1,-2), nrow=3) contrasts(twoway$material)
[,1] [,2] 1 1 0 2 -1 1 3 0 -1 # compare m1 vs m2 and m1+m2 vs m3
contrasts(twoway$material) <- matrix(c(1,-1,0,1,1,-2), nrow=3) contrasts(twoway$material)
[,1] [,2] 1 1 1 2 -1 1 3 0 -2 I'm not sure if 'summary.aov' is the only lm-family summary method with the split argument. DaveT. ************************************* Silviculture Data Analyst Ontario Forest Research Institute Ontario Ministry of Natural Resources david.john.thompson at ontario.ca http://ofri.mnr.gov.on.ca *************************************
-----Original Message-----
From: Steele [mailto:dale.w.steele at gmail.com]
Sent: March 6, 2008 09:08 PM
To: r-help at stat.math.ethz.ch
Subject: [R] Finding Interaction and main effects contrasts
for two-way ANOVA
I've tried without success to calculate interaction and main effects
contrasts using R. I've found the functions C(), contrasts(),
se.contrasts() and fit.contrasts() in package gmodels. Given the url
for a small dataset and the two-way anova model below, I'd like to
reproduce the results from appended SAS code. Thanks. --Dale.
## the dataset (from Montgomery)
twoway <- read.table("http://dsteele.veryspeedy.net/sta501/
twoway.txt",
col.names=c('material', 'temp','voltage'),colClasses=c('factor',
'factor', 'numeric'))
## the model
fit <- aov(voltage ~ material*temp, data=twoway)
/* SAS code */
proc glm data=twoway;
class material temp;
model voltage = material temp material*temp;
contrast '21-22-31+32' material*temp 0 0 0 1 -1 0 -1 1 0;
estimate '21-22-31+32' material*temp 0 0 0 1 -1 0 -1 1 0;
contrast 'material1-material2' material 1 -1 0;
estimate 'material1-material2' material 1 -1 0;
contrast 'temp50 - temp80' temp 1 0 -1;
estimate 'temp50 - temp80' temp 1 0 -1;
run;
______________________________________________ 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.
Thanks to those who have replied to my original query. However, I'm
still confused on how obtain estimates, standard error and F-tests for
main effect and interaction contrasts which agree with the SAS code
with output appended below.
for example,
## Given the dataset (from Montgomery)
twoway <- read.table("http://dsteele.veryspeedy.net/sta501/twoway.txt",
col.names=c('material', 'temp','voltage'),colClasses=c('factor',
'factor', 'numeric'))
## the model
fit <- aov(voltage ~ material*temp, data=twoway)
material.means <- tapply(twoway$voltage, twoway$material, mean)
temp.means <- tapply(twoway$voltage, twoway$temp, mean)
cell.means <- tapply(twoway$voltage, twoway[,1:2], mean)
Contrasts of Interest ....
cell.means[2,1] - cell.means[2,2] - cell.means[3,1] + cell.means[3,2]
[1] 37.75
material.means[1] - material.means[2]
1 -25.16667
temp.means[1] - temp.means[3]
50 80.66667 I expected the following code to provide the estimates above for (material 1 - material 2) and (temp1 - temp3), but get unexpected results...
library(gmodels)
fit.contrast(fit, "material", rbind("(1 - 2)" =c(1, -1, 0) ))
Estimate Std. Error t value Pr(>|t|) material(1 - 2) -21 18.37407 -1.142915 0.2631074
fit.contrast(fit, "temp", rbind("50 - 80" =c(1, 0, -1) ))
Estimate Std. Error t value Pr(>|t|)
temp50 - 80 77.25 18.37407 4.204294 0.0002572756
Thanks. --Dale
/* SAS code */
proc glm data=twoway;
class material temp;
model voltage = material temp material*temp;
contrast '21-22-31+32' material*temp 0 0 0 1 -1 0 -1 1 0;
estimate '21-22-31+32' material*temp 0 0 0 1 -1 0 -1 1 0;
contrast 'material1-material2' material 1 -1 0;
estimate 'material1-material2' material 1 -1 0;
contrast 'temp50 - temp80' temp 1 0 -1;
estimate 'temp50 - temp80' temp 1 0 -1;
run;
SAS output
Contrast DF Contrast SS Mean Square F Value Pr > F
21-22-31+32 1 1425.06250 1425.06250 2.11 0.1578
material1-material2 1 3800.16667 3800.16667 5.63 0.0251
temp50 - temp80 1 39042.66667 39042.66667 57.82 <.0001
Standard
Parameter Estimate Error t Value Pr > |t|
21-22-31+32 37.7500000 25.9848603 1.45 0.1578
material1-material2 -25.1666667 10.6082748 -2.37 0.0251
temp50 - temp80 80.6666667 10.6082748 7.60 <.0001
On Sat, Mar 8, 2008 at 11:02 AM, Gregory Warnes <gregory.warnes at mac.com> wrote:
Dale, You might find it fruitful to look at the help pages for fit.contrast () and estimble() functions in the gmodels package, and the contrast () functions in the Hmisc package. -Greg On Mar 7, 2008, at 4:20PM , Thompson, David ((MNR)) wrote:
> Dale, > > Other than the first SAS contrast, does the following demonstrate what > your asking for?
>> summary(twoway)
> material temp voltage > 1:12 50:12 Min. : 20 > 2:12 65:12 1st Qu.: 70 > 3:12 80:12 Median :108 > Mean :106 > 3rd Qu.:142 > Max. :188
>> contrasts(twoway$material)
> 2 3 > 1 0 0 > 2 1 0 > 3 0 1
>> contrasts(twoway$temp)
> 65 80 > 50 0 0 > 65 1 0 > 80 0 1
>> fit <- aov(voltage ~ material*temp, data=twoway) >> summary.aov(fit)
> Df Sum Sq Mean Sq F value Pr(>F) > material 2 10684 5342 7.91 0.0020 ** > temp 2 39119 19559 28.97 1.9e-07 *** > material:temp 4 9614 2403 3.56 0.0186 * > Residuals 27 18231 675 > --- > Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1 > > > # setting (partial) contrasts
>> contrasts(twoway$material) <- c(1,-1,0) # ignoring the second
> available df
>> contrasts(twoway$temp) <- c(0,1,-1) # ignoring the second >> available df >> contrasts(twoway$material)
> [,1] [,2] > 1 1 -0.41 > 2 -1 -0.41 > 3 0 0.82
>> contrasts(twoway$temp)
> [,1] [,2] > 50 0 -0.82 > 65 1 0.41 > 80 -1 0.41 >
>> summary.aov(fit, split=list(material=list('m1-m2'=1), temp=list
>> ('t50 -
> t80'=1))) > Df Sum Sq Mean Sq F value Pr(>F) > material 2 10684 5342 7.91 0.00198 ** > material: m1-m2 1 3800 3800 5.63 0.02506 * > temp 2 39119 19559 28.97 1.9e-07 *** > temp: t50 - t80 1 11310 11310 16.75 0.00035 *** > material:temp 4 9614 2403 3.56 0.01861 * > material:temp: m1-m2.t50 - t80 1 4970 4970 7.36 0.01146 * > Residuals 27 18231 675 > --- > Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1 > > # other examples of setting contrasts > # compare m1 vs m2 and m2 vs m3
>> contrasts(twoway$material) <- matrix(c(1,-1,0,1,1,-2), nrow=3) >> contrasts(twoway$material)
> [,1] [,2] > 1 1 0 > 2 -1 1 > 3 0 -1 > # compare m1 vs m2 and m1+m2 vs m3
>> contrasts(twoway$material) <- matrix(c(1,-1,0,1,1,-2), nrow=3) >> contrasts(twoway$material)
> [,1] [,2] > 1 1 1 > 2 -1 1 > 3 0 -2 > > I'm not sure if 'summary.aov' is the only lm-family summary method > with > the split argument. > > DaveT. > ************************************* > Silviculture Data Analyst > Ontario Forest Research Institute > Ontario Ministry of Natural Resources > david.john.thompson at ontario.ca > http://ofri.mnr.gov.on.ca > *************************************
>> -----Original Message-----
>> From: Steele [mailto:dale.w.steele at gmail.com]
>> Sent: March 6, 2008 09:08 PM
>> To: r-help at stat.math.ethz.ch
>> Subject: [R] Finding Interaction and main effects contrasts
>> for two-way ANOVA
>>
>> I've tried without success to calculate interaction and main effects
>> contrasts using R. I've found the functions C(), contrasts(),
>> se.contrasts() and fit.contrasts() in package gmodels. Given the url
>> for a small dataset and the two-way anova model below, I'd like to
>> reproduce the results from appended SAS code. Thanks. --Dale.
>>
>> ## the dataset (from Montgomery)
>> twoway <- read.table("http://dsteele.veryspeedy.net/sta501/
>> twoway.txt",
>> col.names=c('material', 'temp','voltage'),colClasses=c('factor',
>> 'factor', 'numeric'))
>>
>> ## the model
>> fit <- aov(voltage ~ material*temp, data=twoway)
>>
>> /* SAS code */
>> proc glm data=twoway;
>> class material temp;
>> model voltage = material temp material*temp;
>> contrast '21-22-31+32' material*temp 0 0 0 1 -1 0 -1 1 0;
>> estimate '21-22-31+32' material*temp 0 0 0 1 -1 0 -1 1 0;
>> contrast 'material1-material2' material 1 -1 0;
>> estimate 'material1-material2' material 1 -1 0;
>> contrast 'temp50 - temp80' temp 1 0 -1;
>> estimate 'temp50 - temp80' temp 1 0 -1;
>> run;
>
______________________________________________ > 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.
On 3/8/2008 3:05 PM, Dale Steele wrote:
Thanks to those who have replied to my original query. However, I'm
still confused on how obtain estimates, standard error and F-tests for
main effect and interaction contrasts which agree with the SAS code
with output appended below.
for example,
## Given the dataset (from Montgomery)
twoway <- read.table("http://dsteele.veryspeedy.net/sta501/twoway.txt",
col.names=c('material', 'temp','voltage'),colClasses=c('factor',
'factor', 'numeric'))
## the model
fit <- aov(voltage ~ material*temp, data=twoway)
material.means <- tapply(twoway$voltage, twoway$material, mean)
temp.means <- tapply(twoway$voltage, twoway$temp, mean)
cell.means <- tapply(twoway$voltage, twoway[,1:2], mean)
Contrasts of Interest ....
cell.means[2,1] - cell.means[2,2] - cell.means[3,1] + cell.means[3,2]
[1] 37.75
material.means[1] - material.means[2]
1 -25.16667
temp.means[1] - temp.means[3]
50 80.66667 I expected the following code to provide the estimates above for (material 1 - material 2) and (temp1 - temp3), but get unexpected results...
library(gmodels)
fit.contrast(fit, "material", rbind("(1 - 2)" =c(1, -1, 0) ))
Estimate Std. Error t value Pr(>|t|) material(1 - 2) -21 18.37407 -1.142915 0.2631074
fit.contrast(fit, "temp", rbind("50 - 80" =c(1, 0, -1) ))
Estimate Std. Error t value Pr(>|t|) temp50 - 80 77.25 18.37407 4.204294 0.0002572756 Thanks. --Dale
Here is one way to reproduce the SAS contrasts:
## the dataset (from Montgomery)
twoway <- read.table("http://dsteele.veryspeedy.net/sta501/twoway.txt",
col.names=c('material', 'temp','voltage'),colClasses=c('factor',
'factor', 'numeric'))
library(gmodels)
fm <- lm(voltage ~ material:temp + 0, data = twoway)
cm <- rbind(
'21-22-31+32 ' = c( 0, 1, -1, 0, -1,1, 0, 0, 0),
'material1-material2' = c(1/3,-1/3, 0,1/3,-1/3,0, 1/3,-1/3, 0),
'temp50-temp80 ' = c(1/3, 1/3,1/3, 0, 0,0,-1/3,-1/3,-1/3))
estimable(fm, cm)
Estimate Std. Error t value DF Pr(>|t|)
21-22-31+32 37.75000 25.98486 1.452769 27 1.578118e-01
material1-material2 -25.16667 10.60827 -2.372362 27 2.505884e-02
temp50-temp80 80.66667 10.60827 7.604127 27 3.525248e-08
This formulates the model so that each coefficient corresponds to one
of the 9 cell means. For me, that makes specifying the contrasts much
easier.
/* SAS code */
proc glm data=twoway;
class material temp;
model voltage = material temp material*temp;
contrast '21-22-31+32' material*temp 0 0 0 1 -1 0 -1 1 0;
estimate '21-22-31+32' material*temp 0 0 0 1 -1 0 -1 1 0;
contrast 'material1-material2' material 1 -1 0;
estimate 'material1-material2' material 1 -1 0;
contrast 'temp50 - temp80' temp 1 0 -1;
estimate 'temp50 - temp80' temp 1 0 -1;
run;
SAS output
Contrast DF Contrast SS Mean Square F Value Pr > F
21-22-31+32 1 1425.06250 1425.06250 2.11 0.1578
material1-material2 1 3800.16667 3800.16667 5.63 0.0251
temp50 - temp80 1 39042.66667 39042.66667 57.82 <.0001
Standard
Parameter Estimate Error t Value Pr > |t|
21-22-31+32 37.7500000 25.9848603 1.45 0.1578
material1-material2 -25.1666667 10.6082748 -2.37 0.0251
temp50 - temp80 80.6666667 10.6082748 7.60 <.0001
On Sat, Mar 8, 2008 at 11:02 AM, Gregory Warnes <gregory.warnes at mac.com> wrote:
Dale, You might find it fruitful to look at the help pages for fit.contrast () and estimble() functions in the gmodels package, and the contrast () functions in the Hmisc package. -Greg On Mar 7, 2008, at 4:20PM , Thompson, David ((MNR)) wrote:
> Dale, > > Other than the first SAS contrast, does the following demonstrate what > your asking for?
>> summary(twoway)
> material temp voltage > 1:12 50:12 Min. : 20 > 2:12 65:12 1st Qu.: 70 > 3:12 80:12 Median :108 > Mean :106 > 3rd Qu.:142 > Max. :188
>> contrasts(twoway$material)
> 2 3 > 1 0 0 > 2 1 0 > 3 0 1
>> contrasts(twoway$temp)
> 65 80 > 50 0 0 > 65 1 0 > 80 0 1
>> fit <- aov(voltage ~ material*temp, data=twoway) >> summary.aov(fit)
> Df Sum Sq Mean Sq F value Pr(>F) > material 2 10684 5342 7.91 0.0020 ** > temp 2 39119 19559 28.97 1.9e-07 *** > material:temp 4 9614 2403 3.56 0.0186 * > Residuals 27 18231 675 > --- > Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1 > > > # setting (partial) contrasts
>> contrasts(twoway$material) <- c(1,-1,0) # ignoring the second
> available df
>> contrasts(twoway$temp) <- c(0,1,-1) # ignoring the second >> available df >> contrasts(twoway$material)
> [,1] [,2] > 1 1 -0.41 > 2 -1 -0.41 > 3 0 0.82
>> contrasts(twoway$temp)
> [,1] [,2] > 50 0 -0.82 > 65 1 0.41 > 80 -1 0.41 >
>> summary.aov(fit, split=list(material=list('m1-m2'=1), temp=list
>> ('t50 -
> t80'=1))) > Df Sum Sq Mean Sq F value Pr(>F) > material 2 10684 5342 7.91 0.00198 ** > material: m1-m2 1 3800 3800 5.63 0.02506 * > temp 2 39119 19559 28.97 1.9e-07 *** > temp: t50 - t80 1 11310 11310 16.75 0.00035 *** > material:temp 4 9614 2403 3.56 0.01861 * > material:temp: m1-m2.t50 - t80 1 4970 4970 7.36 0.01146 * > Residuals 27 18231 675 > --- > Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1 > > # other examples of setting contrasts > # compare m1 vs m2 and m2 vs m3
>> contrasts(twoway$material) <- matrix(c(1,-1,0,1,1,-2), nrow=3) >> contrasts(twoway$material)
> [,1] [,2] > 1 1 0 > 2 -1 1 > 3 0 -1 > # compare m1 vs m2 and m1+m2 vs m3
>> contrasts(twoway$material) <- matrix(c(1,-1,0,1,1,-2), nrow=3) >> contrasts(twoway$material)
> [,1] [,2] > 1 1 1 > 2 -1 1 > 3 0 -2 > > I'm not sure if 'summary.aov' is the only lm-family summary method > with > the split argument. > > DaveT. > ************************************* > Silviculture Data Analyst > Ontario Forest Research Institute > Ontario Ministry of Natural Resources > david.john.thompson at ontario.ca > http://ofri.mnr.gov.on.ca > *************************************
>> -----Original Message-----
>> From: Steele [mailto:dale.w.steele at gmail.com]
>> Sent: March 6, 2008 09:08 PM
>> To: r-help at stat.math.ethz.ch
>> Subject: [R] Finding Interaction and main effects contrasts
>> for two-way ANOVA
>>
>> I've tried without success to calculate interaction and main effects
>> contrasts using R. I've found the functions C(), contrasts(),
>> se.contrasts() and fit.contrasts() in package gmodels. Given the url
>> for a small dataset and the two-way anova model below, I'd like to
>> reproduce the results from appended SAS code. Thanks. --Dale.
>>
>> ## the dataset (from Montgomery)
>> twoway <- read.table("http://dsteele.veryspeedy.net/sta501/
>> twoway.txt",
>> col.names=c('material', 'temp','voltage'),colClasses=c('factor',
>> 'factor', 'numeric'))
>>
>> ## the model
>> fit <- aov(voltage ~ material*temp, data=twoway)
>>
>> /* SAS code */
>> proc glm data=twoway;
>> class material temp;
>> model voltage = material temp material*temp;
>> contrast '21-22-31+32' material*temp 0 0 0 1 -1 0 -1 1 0;
>> estimate '21-22-31+32' material*temp 0 0 0 1 -1 0 -1 1 0;
>> contrast 'material1-material2' material 1 -1 0;
>> estimate 'material1-material2' material 1 -1 0;
>> contrast 'temp50 - temp80' temp 1 0 -1;
>> estimate 'temp50 - temp80' temp 1 0 -1;
>> run;
______________________________________________ > 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.
______________________________________________ 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.
Chuck Cleland, Ph.D. NDRI, Inc. 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