Kia ora list members:
I'm having a little difficulty getting the correct standard errors from a glm.object
(R 1.9.0 under Windows XP 5.1). predict() will gives standard errors of the predicted
values, but I am wanting the standard errors of the mean.
To clarify:
Assume I have a 4x3x2 factorial with 2 complete replications (i.e. 48 observations,
I've appended a dummy set of data at the end of this message). Call the treatments
trt1 (4 levels), trt2 (3 levels) and trt3 (2 levels) and the replications rep - all
are factors. The observed data is S. Then:
temp.aov <- aov(S~rep+trt1*trt2*trt3, data=dummy.data)
model.tables(temp.aov, type='mean', se=T)
Returns the means, but states "Design is unbalanced - use se.contrasts for se's" which
is a little surprising since the design is balanced. Nevertheless, se.contrast gives
what I'd expect:
se.contrast(temp.aov, list(trt1==0, trt1==1), data=dummy.data)
[1] 5.960012
i.e. standard error of mean is 5.960012/sqrt(2) = 4.214, which is the
sqrt(anova(temp.aov)[9,3]/12) as expected. Similarly for interactions, e.g.:
se.contrast(temp.aov, list(trt1==0 & trt2==0, trt1==1 & trt2==1),
data=dummy.data)/sqrt(2)
[1] 7.299494
How do I get the equivalent of these standard errors if I have used lm(), and by
extension glm()? I think I should be able to get these using predict(...,
type='terms', se=T) or coef(summary()) but can't quite see how.
predict(lm(S~rep+trt1*trt2*trt3, data=dummy.data), type='terms', se=T)
or
predict(glm(cbind(S, 100-S)~rep+trt1*trt2*trt3, data=dummy.data, family='binomial'),
type='terms', se=T)
or, as in my case,
predict(glm(cbind(S, 100-S)~rep+trt1*trt2*trt3, data=dummy.data,
family='quasibinomial'), type='terms', se=T)
Thanks ........
Peter Alspach
HortResearch
dummy.data
trt1 trt2 trt3 rep S
0 0 0 1 33
0 0 0 2 55
0 0 1 1 18
0 0 1 2 12
0 1 0 1 47
0 1 0 2 16
0 1 1 1 22
0 1 1 2 33
0 2 0 1 22
0 2 0 2 18
0 2 1 1 60
0 2 1 2 40
1 0 0 1 38
1 0 0 2 24
1 0 1 1 8
1 0 1 2 14
1 1 0 1 69
1 1 0 2 42
1 1 1 1 42
1 1 1 2 44
1 2 0 1 48
1 2 0 2 26
1 2 1 1 46
1 2 1 2 33
2 0 0 1 48
2 0 0 2 46
2 0 1 1 24
2 0 1 2 8
2 1 0 1 69
2 1 0 2 33
2 1 1 1 22
2 1 1 2 33
2 2 0 1 33
2 2 0 2 18
2 2 1 1 26
2 2 1 2 42
3 0 0 1 12
3 0 0 2 42
3 0 1 1 16
3 0 1 2 22
3 1 0 1 14
3 1 0 2 60
3 1 1 1 40
3 1 1 2 55
3 2 0 1 47
3 2 0 2 38
3 2 1 1 18
3 2 1 2 44
______________________________________________________
The contents of this e-mail are privileged and/or confidenti...{{dropped}}
______________________________________________
[EMAIL PROTECTED] mailing list
https://www.stat.math.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
