demo1 <- read.csv("http://www.ats.ucla.edu/stat/data/demo1.csv")
demo1 <- within(demo1, {
group <- factor(group)
time <- factor(time)
id <- factor(id)
})
par(cex = .6)
with(demo1, interaction.plot(time, group, pulse,
ylim = c(5, 20), lty= c(1, 12), lwd = 3,
ylab = "mean of pulse", xlab = "time", trace.label = "group"))
demo1.aov <- aov(pulse ~ group * time + Error(id), data = demo1)
summary(demo1.aov)
-
- データはどんなもの?
- 全部で24測定。それは8人について3時刻(8*3=24)の測定結果であり、8人は2グループに分かれている
> demo1
id group pulse time
1 1 1 10 1
2 1 1 10 2
3 1 1 10 3
4 2 1 10 1
5 2 1 10 2
6 2 1 10 3
7 3 1 10 1
8 3 1 10 2
9 3 1 10 3
10 4 1 10 1
11 4 1 10 2
12 4 1 10 3
13 5 2 15 1
14 5 2 15 2
15 5 2 15 3
16 6 2 15 1
17 6 2 15 2
18 6 2 15 3
19 7 2 16 1
20 7 2 15 2
21 7 2 15 3
22 8 2 15 1
23 8 2 15 2
24 8 2 15 3
-
- どんなモデル?
- pulseという量的な従属変数が、時刻によって変わるか、しかも、その変わり方はgroupに影響されるか、が知りたい。ただし、8人はばらばらなので、その8人のバラバラ加減は織り込みたい
pulse ~ group * time + Error(id)
-
-
- ただし、説明変数(とバラバラ織り込み変数)はすべてファクター扱い
demo1 <- within(demo1, {
group <- factor(group)
time <- factor(time)
id <- factor(id)
})
> summary(demo1.aov)
Error: id
Df Sum Sq Mean Sq F value Pr(>F)
group 1 155.04 155.04 3721 1.3e-09 ***
Residuals 6 0.25 0.04
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Error: Within
Df Sum Sq Mean Sq F value Pr(>F)
time 2 0.0833 0.04167 1 0.397
group:time 2 0.0833 0.04167 1 0.397
Residuals 12 0.5000 0.04167
-
- Sum Sqが155.02,0.25,0.0833,0.0833,0.5000に分かれているのがわかる
- これらを足すと全平方和になる
> summary(demo1.aov)[[1]][[1]][,2]
[1] 155.0417 0.2500
> summary(demo1.aov)[[2]][[1]][,2]
[1] 0.08333333 0.08333333 0.50000000
> sum(summary(demo1.aov)[[2]][[1]][,2])
[1] 0.6666667
> sum(summary(demo1.aov)[[1]][[1]][,2])
[1] 155.2917
> sum(summary(demo1.aov)[[1]][[1]][,2])+sum(summary(demo1.aov)[[2]][[1]][,2])
[1] 155.9583
> sum((demo1$pulse-mean(demo1$pulse))^2)
[1] 155.9583