t検定

  • こちらでt.test()の中のハンドリングのことが話題になっていた
    • 2次データ(各標本の値がわかっているわけでなくて、検体数・平均・(不偏)分散など)だけを用いて検定とも比較していた。ちなみに1次データは「各標本の値が全部わかっているデータ」
  • 等分散を仮定した対応のないt検定のt統計量はどうなっているのか、確認
# t.test(x,y,var.equal=TRUE)とanova(lm(x,y))とが一致すること、そのときのt.test()の統計量を手計算できることを確認
# 平方和(ss/nxが標本分散)
ss <- function(x){
	sum((x-mean(x))^2)
}
# データセットを作って一致することをn.iter回繰り返して確認しよう
n.iter <- 100
ps <- matrix(0,n.iter,3)
for(i in 1:n.iter){
x <- rnorm(10)
y <- rnorm(5)

nx <- length(x)
ny <- length(y)
df <- nx+ny-2
mx <- mean(x)
my <- mean(y)
ssx <- ss(x)
ssy <- ss(y)
# 2群の群内平方和の和
vxy <- (ssx+ssy)
# それを全標本数・群数・群ごとの標本数で調整する
z <- vxy/df*(1/nx+1/ny)
diff.mean <- abs(mx-my)
# 統計量は平均値の差の絶対値とzとの関係
tstat <- diff.mean/sqrt(z)
p <- pt(tstat,df,lower.tail=FALSE)*2
#tstat
ps[i,1] <- p
t.test.out <- t.test(x,y,var.equal=TRUE)
ps[i,2] <- t.test.out[[3]]
anova.out <- anova(lm(c(x,y)~c(rep(1,nx),rep(2,ny))))
ps[i,3] <- anova.out[[5]][1]
}