2要因とも対応がない場合（CRFpqデザイン）Ⅰ [データ数が等しいケース]
ソースコード datはデータ行列,pは要因Aの水準数,qは要因Bの水準数
anova.CRFpq1 <- function(dat,p,q){
n <- length(dat[,1])
ABjk <- apply(dat, 2, sum)
mt <- matrix(ABjk, ncol=q, byrow = TRUE)
Aj <- apply(mt, 1, sum)
Bk <- apply(mt, 2, sum)
ABS <- sum(dat ^ 2)
X <- sum(mt) ^ 2 / (n * p * q)
A <- sum(Aj ^ 2 / (n * q))
B <- sum(Bk ^ 2 / (n * p))
AB <- sum(ABjk ^ 2 / n)
SSa <- A - X
SSb <- B - X
SSab <- AB - A - B + X
SSwc <- ABS - AB
SSt <- ABS - X
dfA <- p - 1
dfB <- q - 1
dfAB <- (p - 1) * (q - 1)
dfWC <- p * q * (n - 1)
dfT <- n * p * q - 1
MSwc <- SSwc / dfWC
MSa <- SSa / dfA
MSb <- SSb / dfB
MSab <- SSab / dfAB
Fa <- MSa / MSwc
Fb <- MSb / MSwc
Fab <- MSab / MSwc
pA <- 1 - pf(Fa, dfA, dfWC)
pB <- 1 - pf(Fb, dfB, dfWC)
pAB <- 1 - pf(Fab, dfAB, dfWC)
result.table <- matrix(c(
SSa, dfA, MSa, Fa, pA,
SSb, dfB, MSb, Fb, pB,
SSab, dfAB, MSab, Fab, pAB,
SSwc, dfWC, MSwc, NA, NA,
SSt, dfT, NA, NA, NA),ncol = 5, byrow = TRUE)
dimnames(result.table) <- list(c(
"行効果:A","列効果:B","交互作用:A*B","誤差:WC","総:T"),
c("SS","df","MS","F-value","p-value"))
result.table
}
～使用例～
> dat
[,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8]
[1,] 3 4 6 5 3 2 3 2
[2,] 3 3 6 7 5 6 2 3
[3,] 1 4 6 8 2 3 3 3
[4,] 3 5 4 7 4 6 6 4
[5,] 5 7 8 9 6 4 5 6
> anova.CRFpq1(dat,2,4) #要因Aの水準数が2,要因Bの水準数が4である
SS df MS F-value p-value
行効果:A 16.9 1 16.900000 7.041667 0.01229668
列効果:B 19.7 3 6.566667 2.736111 0.05972018
交互作用:A*B 30.5 3 10.166667 4.236111 0.01249606
誤差:WC 76.8 32 2.400000 NA NA
総:T 143.9 39 NA NA NA
