Hi,
I'm having difficulty to replicate in R a nested (hierarchical) anova example
found in p. 308 of
Zar, J.H. 1996. Bostatistical Analysis. Prentice Hall. 3rd ed.
The example (15.1) is as follows:
The variable is blood cholesterol concentration in women (in mg/100 ml of
plasma). This
variable was measured after the administration of one of three different drugs,
each drug
having been obtained from two sources
DRUG 1 DRUG 2 DRUG 3
A Q D B L S
102 103 108 1091 104 105
104 104 110 1082 106 107
Zar results are these:
----------------------------------------------------------------------------
Source of variation SS DF MS F P
----------------------------------------------------------------------------
*Total 71.67 11
*Among all subgroups 62.67 5
*Groups 61.17 2 30.58 61.16 0.0037
*Subgroups 1.50 3 0.50 0.33 0.80
*Error 9.00 6 1.50
----------------------------------------------------------------------------
It is straightforward to follow the equations from the book step by step to
obtain those results.
Now, when trying the 'aov()' function in R:
## creating the data frame
zar <- data.frame(
cChol=c(102, 104, 103, 104, 108, 110, 109, 108, 104, 106, 105, 107),
drug=gl(3, 4, 12),
dsource=gl(6, 2, 12, labels=c("A", "Q", "D", "B", "L", "S")) )
## This gives the relevant SS, DF and MS values, but only the F and P
corresponding to drug
## effect (Zar's Groups)
summary( aov(cChol ~ drug + Error(dsource/drug), data=zar) )
## The F and P values corresponding to subgroups can be obtained by:
0.5/1.5 # F
pf(0.5/1.5, 3, 6, lower.tail=FALSE) # P
## or from the second row in the two way anova:
summary( aov(cChol ~ drug + dsource, data=zar) )
In order to avoid confussion in our students (and myself!) I wonder if we can
obtain the
desired results in one step?
Thank you in advance for your kind help.
Héctor
--
Héctor Villalobos <hvill...@ipn.mx>
CICIMAR - IPN
A.P. 592. Col. Centro
La Paz, Baja California Sur, MÉXICO. 23000
Tels. (+52 612) 122 53 44; 123 46 58; 123 47 34 ext. 82425
Fax. (+52 612) 122 53 22
[[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.
