1 分割表の検定

1.1 四分表の検定

最も基本的な場合として、四分表の独立検定を考える。
X=1 であることが Yの判断に影響しないこと、即ち、条件付きの分布が周辺分布と一致しているときに独立となる。
独立である場合の四分表のセルの比率と実現値の比率を比べ、確率的に有意であるほどに乖離がみられるかどうかを考えれば良いので、四分表のセルについて適合度検定を行えば良い。
クラス数が２の場合の適合度検定は次の通りである。
関数chisq.testの出力からカイ二乗値を取り出してみよう。

X = c()
s <- rmultinom(n=1000, size = 10, prob = c(0.5,0.5)) 
for (i in 1:1000) {
  res <- chisq.test(s[,i], p = c(0.5,0.5))
  rs <- res$ statistic
  X <-c(X,rs)
}
hist(X,prob=T)
curve(dchisq(x,1),add=T)

市原先生の例題について考えよう。

1.2 例4-5

定義通りカイ二乗値を計算する。

X = c()
for (i in 1:1000) {
  a1 <- rmultinom(n=1, size = 72, prob = c(80,60)/140)
  a2 <- rmultinom(n=1, size = 68, prob = c(80,60)/140)
  x <-c(a1,a2)
  y <- c(c(80,60)*(72/140),c(80,60)*68/140)
  s <- sum((x-y)^2/y)
  X <-c(X,s)
}
hist(X,prob=T)
curve(dchisq(x,1),add=T)

シミュレーションによる出力は、理論上の関数より少し大きいように見える。
は一致しているようである。
### イエーツの補正

X1 = c()
X2 = c()
for (i in 1:1000) {
  a1 <- rmultinom(n=1, size = 72, prob = c(80,60)/140)
  a2 <- rmultinom(n=1, size = 68, prob = c(80,60)/140)
  x <-c(a1,a2)
  y <- c(c(80,60)*(72/140),c(80,60)*68/140)
  s1 <- sum((x-y)^2/y)
  s2 <- sum((abs(x-y)-0.5)^2/y)
  X1 <-c(X1,s1)
  X2 <-c(X2,s2)
}
hist(X1,prob=T)
hist(X2,prob=T,add=TRUE,col="lightblue",alpha=0.1)
curve(dchisq(x,1),add=T)

補正により当てはまりが良くなっていることがグラフ上で確認できる。

参考

r <- rmultinom(n=1, size = 140, prob = c(72,68)/140);#r
c <- rmultinom(n=1, size = 140, prob = c(80,60)/140);#c
a1 <- rmultinom(n=1, size = r[1,], prob = c(80,60)/140);#a1
a2 <- rmultinom(n=1, size = r[2,], prob = c(80,60)/140);#a2
x <-c(a1,a2)
y <- c(c(80,60)*(72/140),c(80,60)*68/140);y
## [1] 41.14286 30.85714 38.85714 29.14286
sum((x-y)^2/y)
## [1] 1.871051
a <-x[1];b <-x[2];c <-x[3];d <-x[4]
(a*d-b*c)^2*140/((a+b)*(c+d)*(a+c)*(b+d))
## [1] 0.141065
x1 <-c(48,24,32,36)
a1 <-x1[1];b1 <-x1[2];c1 <-x1[3];d1 <-x1[4]
(a1*d1-b1*c1)^2*140/((a1+b1)*(c1+d1)*(a1+c1)*(b1+d1))
## [1] 5.490196