Repeated Measures ANOVA

  • リンクを集めているページ(こちら)
  • まずはここ(UCLA)のサイトに行ってみよう
  • 例1
demo1 <- read.csv("http://www.ats.ucla.edu/stat/data/demo1.csv")
## Convert variables to factor
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