linear contrasts with anova

classic Classic list List threaded Threaded
4 messages Options
Reply | Threaded
Open this post in threaded view
|

linear contrasts with anova

Posta Univ. Cagliari
I have some doubts about the validity of my procedure to estimeate linear contrasts ina a factorial design.
For sake of semplicity, let's imagine a one way ANOVA with three levels. I am interested to test the significance of the difference between the first and third level (called here contrast C1) and between the first and the seconda level (called here contrast C2). I used the following procedure:


------------------- reading data from a text file -----------------------

> ar <-read.table("C:/Programmi/R/myworks/contrasti/cont1.txt",header=TRUE)

> ar

     CC GROUP

1   3.0     0

2   3.0     0

3   4.0     0

4   5.0     0

5   6.0     0

6   7.0     0

7   3.0     0

8   2.0     0

9   1.0     1

10  6.0     1

11  5.0     1

12  7.0     1

13  2.0     1

14  3.0     1

15  1.5     1

16  1.7     1

17 17.0     2

18 12.0     2

19 15.0     2

20 16.0     2

21 12.0     2

22 23.0     2

23 19.0     2

24 21.0     2

 

------------------- creating a new array of data-----------------------

> ar<-data.frame(GROUP=factor(ar$GROUP),DIP=ar$CC)

> ar

   GROUP  DIP

1      0  3.0

2      0  3.0

3      0  4.0

4      0  5.0

5      0  6.0

6      0  7.0

7      0  3.0

8      0  2.0

9      1  1.0

10     1  6.0

11     1  5.0

12     1  7.0

13     1  2.0

14     1  3.0

15     1  1.5

16     1  1.7

17     2 17.0

18     2 12.0

19     2 15.0

20     2 16.0

21     2 12.0

22     2 23.0

23     2 19.0

24     2 21.0

 

------------------- creating two dummy variables (C1 and C2) for linear contrasts-----------------------

> ar<-data.frame(GROUP=factor(ar$GROUP),C1=factor(ar$GROUP),C2=factor(ar$GROUP),DIP=ar$DIP)

> ar

   GROUP C1 C2  DIP

1      0  0  0  3.0

2      0  0  0  3.0

3      0  0  0  4.0

4      0  0  0  5.0

5      0  0  0  6.0

6      0  0  0  7.0

7      0  0  0  3.0

8      0  0  0  2.0

9      1  1  1  1.0

10     1  1  1  6.0

11     1  1  1  5.0

12     1  1  1  7.0

13     1  1  1  2.0

14     1  1  1  3.0

15     1  1  1  1.5

16     1  1  1  1.7

17     2  2  2 17.0

18     2  2  2 12.0

19     2  2  2 15.0

20     2  2  2 16.0

21     2  2  2 12.0

22     2  2  2 23.0

23     2  2  2 19.0

24     2  2  2 21.0

 

------------------- selecting the contrast levels-----------------------

> ar$C1 <- C(ar$C1, c(1,0,-1), how.many = 1)

> ar$C2 <- C(ar$C2, c(1,-1,0), how.many = 1)

 

 

------------------- contrast analysis of C2 -----------------------

> r.aov8 <-aov(DIP ~  C2 + GROUP , data = ar)

> anova(r.aov8)

Analysis of Variance Table

 

Response: DIP

          Df Sum Sq Mean Sq  F value    Pr(>F)    

C2         1   2.10    2.10   0.2622     0.614    

GROUP      1 917.00  917.00 114.3460 5.915e-10 ***

Residuals 21 168.41    8.02                      

---

Signif. codes:  0 `***' 0.001 `**' 0.01 `*' 0.05 `.' 0.1 ` ' 1

 

------------------- contrast analysis of C1 -----------------------

> r.aov9 <-aov(DIP ~  C1 + GROUP , data = ar)

> anova(r.aov9)

Analysis of Variance Table

 

Response: DIP

          Df Sum Sq Mean Sq F value    Pr(>F)    

C1         1 650.25  650.25  81.083 1.175e-08 ***

GROUP      1 268.85  268.85  33.525 9.532e-06 ***

Residuals 21 168.41    8.02                      

---

Signif. codes:  0 `***' 0.001 `**' 0.01 `*' 0.05 `.' 0.1 ` ' 1

 

------------------- anova of the global design -----------------------

> r.aov10 <-aov(DIP ~  GROUP , data = ar)

> anova(r.aov10)

Analysis of Variance Table

 

Response: DIP

          Df Sum Sq Mean Sq F value    Pr(>F)    

GROUP      2 919.10  459.55  57.304 3.121e-09 ***

Residuals 21 168.41    8.02                      

---

Signif. codes:  0 `***' 0.001 `**' 0.01 `*' 0.05 `.' 0.1 ` ' 1









I would like to know if there is a more economic procedure with R to do linear contrasts.

Every comments will be well accepted.



Thank you very much and best regards



Marco Tommasi

        [[alternative HTML version deleted]]

______________________________________________
[hidden email] mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Reply | Threaded
Open this post in threaded view
|

Re: linear contrasts with anova

Steffen Katzner
group=factor(rep(c(0:2), each = 8))
ar = data.frame(group, dip)

con = matrix(c(1, -1, 0, 1, 0, -1), nrow=3, ncol=2, byrow=F)
contrasts(ar$group)=con

aovRes = aov(dip~group, ar)

 > summary.aov(aovRes, split=list(group = list("0 vs 1" = 1, "0 vs 3" = 2)))

                 Df Sum Sq Mean Sq  F value    Pr(>F)
group            2 919.10  459.55  57.3041 3.121e-09 ***
   group: 0 vs 1  1   2.10    2.10   0.2622     0.614
   group: 0 vs 3  1 917.00  917.00 114.3460 5.915e-10 ***
Residuals       21 168.41    8.02
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1


I only don't know why it does not work with within-subject designs.
I have posted this question before, does anybody know?

-steffen






> I have some doubts about the validity of my procedure to estimeate linear contrasts ina a factorial design.
> For sake of semplicity, let's imagine a one way ANOVA with three levels. I am interested to test the significance of the difference between the first and third level (called here contrast C1) and between the first and the seconda level (called here contrast C2). I used the following procedure:
>
>
> ------------------- reading data from a text file -----------------------
>
>> ar <-read.table("C:/Programmi/R/myworks/contrasti/cont1.txt",header=TRUE)
>
>> ar
>
>      CC GROUP
>
> 1   3.0     0
>
> 2   3.0     0
>
> 3   4.0     0
>
> 4   5.0     0
>
> 5   6.0     0
>
> 6   7.0     0
>
> 7   3.0     0
>
> 8   2.0     0
>
> 9   1.0     1
>
> 10  6.0     1
>
> 11  5.0     1
>
> 12  7.0     1
>
> 13  2.0     1
>
> 14  3.0     1
>
> 15  1.5     1
>
> 16  1.7     1
>
> 17 17.0     2
>
> 18 12.0     2
>
> 19 15.0     2
>
> 20 16.0     2
>
> 21 12.0     2
>
> 22 23.0     2
>
> 23 19.0     2
>
> 24 21.0     2
>
>  
>
> ------------------- creating a new array of data-----------------------
>
>> ar<-data.frame(GROUP=factor(ar$GROUP),DIP=ar$CC)
>
>> ar
>
>    GROUP  DIP
>
> 1      0  3.0
>
> 2      0  3.0
>
> 3      0  4.0
>
> 4      0  5.0
>
> 5      0  6.0
>
> 6      0  7.0
>
> 7      0  3.0
>
> 8      0  2.0
>
> 9      1  1.0
>
> 10     1  6.0
>
> 11     1  5.0
>
> 12     1  7.0
>
> 13     1  2.0
>
> 14     1  3.0
>
> 15     1  1.5
>
> 16     1  1.7
>
> 17     2 17.0
>
> 18     2 12.0
>
> 19     2 15.0
>
> 20     2 16.0
>
> 21     2 12.0
>
> 22     2 23.0
>
> 23     2 19.0
>
> 24     2 21.0
>
>  
>
> ------------------- creating two dummy variables (C1 and C2) for linear contrasts-----------------------
>
>> ar<-data.frame(GROUP=factor(ar$GROUP),C1=factor(ar$GROUP),C2=factor(ar$GROUP),DIP=ar$DIP)
>
>> ar
>
>    GROUP C1 C2  DIP
>
> 1      0  0  0  3.0
>
> 2      0  0  0  3.0
>
> 3      0  0  0  4.0
>
> 4      0  0  0  5.0
>
> 5      0  0  0  6.0
>
> 6      0  0  0  7.0
>
> 7      0  0  0  3.0
>
> 8      0  0  0  2.0
>
> 9      1  1  1  1.0
>
> 10     1  1  1  6.0
>
> 11     1  1  1  5.0
>
> 12     1  1  1  7.0
>
> 13     1  1  1  2.0
>
> 14     1  1  1  3.0
>
> 15     1  1  1  1.5
>
> 16     1  1  1  1.7
>
> 17     2  2  2 17.0
>
> 18     2  2  2 12.0
>
> 19     2  2  2 15.0
>
> 20     2  2  2 16.0
>
> 21     2  2  2 12.0
>
> 22     2  2  2 23.0
>
> 23     2  2  2 19.0
>
> 24     2  2  2 21.0
>
>  
>
> ------------------- selecting the contrast levels-----------------------
>
>> ar$C1 <- C(ar$C1, c(1,0,-1), how.many = 1)
>
>> ar$C2 <- C(ar$C2, c(1,-1,0), how.many = 1)
>
>  
>
>  
>
> ------------------- contrast analysis of C2 -----------------------
>
>> r.aov8 <-aov(DIP ~  C2 + GROUP , data = ar)
>
>> anova(r.aov8)
>
> Analysis of Variance Table
>
>  
>
> Response: DIP
>
>           Df Sum Sq Mean Sq  F value    Pr(>F)    
>
> C2         1   2.10    2.10   0.2622     0.614    
>
> GROUP      1 917.00  917.00 114.3460 5.915e-10 ***
>
> Residuals 21 168.41    8.02                      
>
> ---
>
> Signif. codes:  0 `***' 0.001 `**' 0.01 `*' 0.05 `.' 0.1 ` ' 1
>
>  
>
> ------------------- contrast analysis of C1 -----------------------
>
>> r.aov9 <-aov(DIP ~  C1 + GROUP , data = ar)
>
>> anova(r.aov9)
>
> Analysis of Variance Table
>
>  
>
> Response: DIP
>
>           Df Sum Sq Mean Sq F value    Pr(>F)    
>
> C1         1 650.25  650.25  81.083 1.175e-08 ***
>
> GROUP      1 268.85  268.85  33.525 9.532e-06 ***
>
> Residuals 21 168.41    8.02                      
>
> ---
>
> Signif. codes:  0 `***' 0.001 `**' 0.01 `*' 0.05 `.' 0.1 ` ' 1
>
>  
>
> ------------------- anova of the global design -----------------------
>
>> r.aov10 <-aov(DIP ~  GROUP , data = ar)
>
>> anova(r.aov10)
>
> Analysis of Variance Table
>
>  
>
> Response: DIP
>
>           Df Sum Sq Mean Sq F value    Pr(>F)    
>
> GROUP      2 919.10  459.55  57.304 3.121e-09 ***
>
> Residuals 21 168.41    8.02                      
>
> ---
>
> Signif. codes:  0 `***' 0.001 `**' 0.01 `*' 0.05 `.' 0.1 ` ' 1
>
>
>
>
>
>
>
>
>
> I would like to know if there is a more economic procedure with R to do linear contrasts.
>
> Every comments will be well accepted.
>
>
>
> Thank you very much and best regards
>
>
>
> Marco Tommasi
>
> [[alternative HTML version deleted]]
>

______________________________________________
[hidden email] mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Reply | Threaded
Open this post in threaded view
|

Re: linear contrasts with anova

Christoph Buser
In reply to this post by Posta Univ. Cagliari
Dear Marco

If you are interested in a comparison of a reference level
against each other level (in your case level 1 against level 2
and level 1 against level 3), you can use the summary.lm(),
since this is the default contrast (see ?contr.treatment)

ar <- data.frame(GROUP = factor(rep(1:3, each = 8)),
                 DIP = c(3.0, 3.0, 4.0, 5.0, 6.0, 7.0, 3.0, 2.0, 1.0, 6.0, 5.0,
                   7.0, 2.0, 3.0, 1.5, 1.7, 17.0, 12.0, 15.0, 16.0, 12.0, 23.0,
                   19.0, 21.0))


r.aov10 <- aov(DIP ~  GROUP, data = ar)
anova(r.aov10)
summary.lm(r.aov10)

As result you will get the comparison GROUP 2 against GROUP 1,
denoted by GROUP2 and the comparison GROUP 3 against GROUP 1,
denoted by GROUP3.

Be careful. In your example you include both GROUP and C1 or C2,
respectively in your model. This results in a over parameterized
model and you get a warning that not all coefficients have been
estimated, due to singularities.

It is possible to use other contrasts than contr.treatment,
contr.sum, contr.helmert, contr.poly, but then you have to
specify the correctly in the model.

Regards,

Christoph Buser

--------------------------------------------------------------
Christoph Buser <[hidden email]>
Seminar fuer Statistik, LEO C13
ETH (Federal Inst. Technology) 8092 Zurich SWITZERLAND
phone: x-41-44-632-4673 fax: 632-1228
http://stat.ethz.ch/~buser/
--------------------------------------------------------------

Posta Univ. Cagliari writes:
 > I have some doubts about the validity of my procedure to estimeate linear contrasts ina a factorial design.
 > For sake of semplicity, let's imagine a one way ANOVA with three levels. I am interested to test the significance of the difference between the first and third level (called here contrast C1) and between the first and the seconda level (called here contrast C2). I used the following procedure:
 >
 >
 > ------------------- reading data from a text file -----------------------
 >
 > > ar <-read.table("C:/Programmi/R/myworks/contrasti/cont1.txt",header=TRUE)
 >
 > > ar
 >
 >      CC GROUP
 >
 > 1   3.0     0
 >
 > 2   3.0     0
 >
 > 3   4.0     0
 >
 > 4   5.0     0
 >
 > 5   6.0     0
 >
 > 6   7.0     0
 >
 > 7   3.0     0
 >
 > 8   2.0     0
 >
 > 9   1.0     1
 >
 > 10  6.0     1
 >
 > 11  5.0     1
 >
 > 12  7.0     1
 >
 > 13  2.0     1
 >
 > 14  3.0     1
 >
 > 15  1.5     1
 >
 > 16  1.7     1
 >
 > 17 17.0     2
 >
 > 18 12.0     2
 >
 > 19 15.0     2
 >
 > 20 16.0     2
 >
 > 21 12.0     2
 >
 > 22 23.0     2
 >
 > 23 19.0     2
 >
 > 24 21.0     2
 >
 >  
 >
 > ------------------- creating a new array of data-----------------------
 >
 > > ar<-data.frame(GROUP=factor(ar$GROUP),DIP=ar$CC)
 >
 > > ar
 >
 >    GROUP  DIP
 >
 > 1      0  3.0
 >
 > 2      0  3.0
 >
 > 3      0  4.0
 >
 > 4      0  5.0
 >
 > 5      0  6.0
 >
 > 6      0  7.0
 >
 > 7      0  3.0
 >
 > 8      0  2.0
 >
 > 9      1  1.0
 >
 > 10     1  6.0
 >
 > 11     1  5.0
 >
 > 12     1  7.0
 >
 > 13     1  2.0
 >
 > 14     1  3.0
 >
 > 15     1  1.5
 >
 > 16     1  1.7
 >
 > 17     2 17.0
 >
 > 18     2 12.0
 >
 > 19     2 15.0
 >
 > 20     2 16.0
 >
 > 21     2 12.0
 >
 > 22     2 23.0
 >
 > 23     2 19.0
 >
 > 24     2 21.0
 >
 >  
 >
 > ------------------- creating two dummy variables (C1 and C2) for linear contrasts-----------------------
 >
 > > ar<-data.frame(GROUP=factor(ar$GROUP),C1=factor(ar$GROUP),C2=factor(ar$GROUP),DIP=ar$DIP)
 >
 > > ar
 >
 >    GROUP C1 C2  DIP
 >
 > 1      0  0  0  3.0
 >
 > 2      0  0  0  3.0
 >
 > 3      0  0  0  4.0
 >
 > 4      0  0  0  5.0
 >
 > 5      0  0  0  6.0
 >
 > 6      0  0  0  7.0
 >
 > 7      0  0  0  3.0
 >
 > 8      0  0  0  2.0
 >
 > 9      1  1  1  1.0
 >
 > 10     1  1  1  6.0
 >
 > 11     1  1  1  5.0
 >
 > 12     1  1  1  7.0
 >
 > 13     1  1  1  2.0
 >
 > 14     1  1  1  3.0
 >
 > 15     1  1  1  1.5
 >
 > 16     1  1  1  1.7
 >
 > 17     2  2  2 17.0
 >
 > 18     2  2  2 12.0
 >
 > 19     2  2  2 15.0
 >
 > 20     2  2  2 16.0
 >
 > 21     2  2  2 12.0
 >
 > 22     2  2  2 23.0
 >
 > 23     2  2  2 19.0
 >
 > 24     2  2  2 21.0
 >
 >  
 >
 > ------------------- selecting the contrast levels-----------------------
 >
 > > ar$C1 <- C(ar$C1, c(1,0,-1), how.many = 1)
 >
 > > ar$C2 <- C(ar$C2, c(1,-1,0), how.many = 1)
 >
 >  
 >
 >  
 >
 > ------------------- contrast analysis of C2 -----------------------
 >
 > > r.aov8 <-aov(DIP ~  C2 + GROUP , data = ar)
 >
 > > anova(r.aov8)
 >
 > Analysis of Variance Table
 >
 >  
 >
 > Response: DIP
 >
 >           Df Sum Sq Mean Sq  F value    Pr(>F)    
 >
 > C2         1   2.10    2.10   0.2622     0.614    
 >
 > GROUP      1 917.00  917.00 114.3460 5.915e-10 ***
 >
 > Residuals 21 168.41    8.02                      
 >
 > ---
 >
 > Signif. codes:  0 `***' 0.001 `**' 0.01 `*' 0.05 `.' 0.1 ` ' 1
 >
 >  
 >
 > ------------------- contrast analysis of C1 -----------------------
 >
 > > r.aov9 <-aov(DIP ~  C1 + GROUP , data = ar)
 >
 > > anova(r.aov9)
 >
 > Analysis of Variance Table
 >
 >  
 >
 > Response: DIP
 >
 >           Df Sum Sq Mean Sq F value    Pr(>F)    
 >
 > C1         1 650.25  650.25  81.083 1.175e-08 ***
 >
 > GROUP      1 268.85  268.85  33.525 9.532e-06 ***
 >
 > Residuals 21 168.41    8.02                      
 >
 > ---
 >
 > Signif. codes:  0 `***' 0.001 `**' 0.01 `*' 0.05 `.' 0.1 ` ' 1
 >
 >  
 >
 > ------------------- anova of the global design -----------------------
 >
 > > r.aov10 <-aov(DIP ~  GROUP , data = ar)
 >
 > > anova(r.aov10)
 >
 > Analysis of Variance Table
 >
 >  
 >
 > Response: DIP
 >
 >           Df Sum Sq Mean Sq F value    Pr(>F)    
 >
 > GROUP      2 919.10  459.55  57.304 3.121e-09 ***
 >
 > Residuals 21 168.41    8.02                      
 >
 > ---
 >
 > Signif. codes:  0 `***' 0.001 `**' 0.01 `*' 0.05 `.' 0.1 ` ' 1
 >
 >
 >
 >
 >
 >
 >
 >
 >
 > I would like to know if there is a more economic procedure with R to do linear contrasts.
 >
 > Every comments will be well accepted.
 >
 >
 >
 > Thank you very much and best regards
 >
 >
 >
 > Marco Tommasi
 >
 > [[alternative HTML version deleted]]
 >
 > ______________________________________________
 > [hidden email] mailing list
 > https://stat.ethz.ch/mailman/listinfo/r-help
 > PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html

______________________________________________
[hidden email] mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Reply | Threaded
Open this post in threaded view
|

Re: linear contrasts with anova

joseclaudio.faria
In reply to this post by Posta Univ. Cagliari
> Date: Mon, 23 Jan 2006 13:25:33 +0100
> From: Christoph Buser <[hidden email]>
> Subject: Re: [R] linear contrasts with anova
> To: "Posta Univ. Cagliari" <[hidden email]>
> Cc: [hidden email]
> Message-ID: <[hidden email]>
> Content-Type: text/plain; charset=us-ascii
>
> Dear Marco
>
> If you are interested in a comparison of a reference level
> against each other level (in your case level 1 against level 2
> and level 1 against level 3), you can use the summary.lm(),
> since this is the default contrast (see ?contr.treatment)
>
> ar <- data.frame(GROUP = factor(rep(1:3, each = 8)),
>                  DIP = c(3.0, 3.0, 4.0, 5.0, 6.0, 7.0, 3.0, 2.0, 1.0, 6.0, 5.0,
>                    7.0, 2.0, 3.0, 1.5, 1.7, 17.0, 12.0, 15.0, 16.0, 12.0, 23.0,
>                    19.0, 21.0))
>
>
> r.aov10 <- aov(DIP ~  GROUP, data = ar)
> anova(r.aov10)
> summary.lm(r.aov10)
>
> As result you will get the comparison GROUP 2 against GROUP 1,
> denoted by GROUP2 and the comparison GROUP 3 against GROUP 1,
> denoted by GROUP3.
>
> Be careful. In your example you include both GROUP and C1 or C2,
> respectively in your model. This results in a over parameterized
> model and you get a warning that not all coefficients have been
> estimated, due to singularities.
>
> It is possible to use other contrasts than contr.treatment,
> contr.sum, contr.helmert, contr.poly, but then you have to
> specify the correctly in the model.
>
> Regards,
>
> Christoph Buser
>
> --------------------------------------------------------------
> Christoph Buser <[hidden email]>
> Seminar fuer Statistik, LEO C13
> ETH (Federal Inst. Technology) 8092 Zurich SWITZERLAND
> phone: x-41-44-632-4673 fax: 632-1228
> http://stat.ethz.ch/~buser/
> --------------------------------------------------------------

Dear Marco

Try also the the below:

# Loading packages
library(gplots)
library(gmodels)

# Testing the nature of dF
is.data.frame(dF)
is.factor(dF$GROUP)
is.numeric(dF$DIP)

#Plotting GROUPS
win.graph(w = 4, h = 5)
plotmeans(DIP ~ GROUP, data = dF, mean.labels = TRUE,
           digits = 3, col = 'blue', connect = FALSE,
           ylab = 'DIP', xlab = 'GROUP', pch='')

# Contrasts
attach(dF)
                              #1   2   3 GROUP
        cmat = rbind('1 vs. 3' = c( 1,  0, -1,),
                     '1 vs. 2' = c( 1, -1,  0,))

   av  = aov(DIP ~ GROUP, data = dF, contrasts = list(GROUP =
make.contrasts(cmat)))
   sav = summary(av1, split = list(GROUP=list('1 vs. 3'=1,
                                             '1 vs. 2'=2)))
   sav
detach(dF)

# Another option

attach(dF)
                             #A   B   C
   fit.contrast(av, GROUP, c( 1,  0, -1)) # from gmodels
   fit.contrast(av, GROUP, c( 1, -1,  0))
detach(dF)

HTH, []s,
--
Jose Claudio Faria
Brasil/Bahia/Ilhéus/UESC/DCET
Estatística Experimental/Prof. Adjunto
mails:
    [hidden email]
    [hidden email]
    [hidden email]

______________________________________________
[hidden email] mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html