1 超幾何分布

非復元無作為抽出

1.1 赤6個、白4個から5個抽出

A=c(0,0,0,0,1,1,1,1,1,1)#母集団
#X <- sample(A, size = 5, replace = FALSE)#標本
par(mfrow = c(1, 1))#グラフを1つ並べる
X = c()
for (i in 1:10000) {
  Sample_5 <- sample(A, size = 5, replace = FALSE)#非復元
  #print(Sample_5)
  X <-c(X,sum(Sample_5))
}
table(X)/10000
## X
##      1      2      3      4      5 
## 0.0230 0.2373 0.4756 0.2435 0.0206
hist(X,prob=T,breaks=c(0,1,2,3,4,5,6))#breaksは値の区切り

「ランダムなイベント」とは大雑把に言うと「起こる確率が常に一定である」ようなイベントのことです。

1.2 乱数 dhyper （超幾何分布）

par(mfrow = c(1, 2))
x=c(1:5)
plot(dhyper(x,6,10-6,5),type="l",ylim=c(0,0.5),xlim=c(0,6),ylab="")
hist(X,prob=T,breaks=c(0,1,2,3,4,5,6))

#print(dhyper(x,6,10-6,5))

乱数と抽出が一致していることがわかる。

2 二項分布に従う乱数

2.1 オッズ比の分布

2.1.1 オッズ比の算出

A <- rbinom(6, 1, 1/3)
odsA <- sum(A)/(6-sum(A))
B <- rbinom(4, 1, 3/4)
odsB <- sum(B)/(4-sum(B))
OR <- odsA/odsB
OR
## [1] 3

偏りががあるのでオッズ比は１から離れている。

2.1.2 偏りのないとき

OR=c()
for (i in 1:1000){
  A <- rbinom(6, 1, 1/3)
  sum(A)+(6-sum(A))
  odsA <- sum(A)/(6-sum(A))
  B <- rbinom(4, 1, 1/3)
  sum(B)+(4-sum(B))
  odsB <- sum(B)/(4-sum(B))
  or <- odsA/odsB
  OR <- c(OR,or)
}
head(OR)
## [1] Inf 1.5 1.0 1.5 0.6 2.0
mean(OR<Inf,na.rm = TRUE)#非数値を除外
## [1] 0.8082051
par(mfrow = c(1, 2))
hist(OR,prob=T)
hist(log(OR),prob=T)

偏りがないので１に近くなる。

2.1.3 サンプル数を１０倍

OR=c()
for (i in 1:1000){
  A <- rbinom(600, 1, 1/3)
  odsA <- sum(A)/(600-sum(A))
  B <- rbinom(400, 1, 1/3)
  odsB <- sum(B)/(400-sum(B))
  or <- odsA/odsB
  OR <- c(OR,or)
}
head(OR)
## [1] 1.1311866 1.2450566 0.8507463 1.0798358 0.9462712 1.2333420
mean(OR)
## [1] 1.013274
par(mfrow = c(1, 2))
hist(OR,prob=T)
hist(log(OR),prob=T)

OR=c()
for (i in 1:1000){
  A <- rbinom(60, 1, 1/3)
  odsA <- sum(A)/(60-sum(A))
  B <- rbinom(40, 1, 3/4)
  odsB <- sum(B)/(40-sum(B))
  or <- odsA/odsB
  OR <- c(OR,or)
}
head(OR)
## [1] 0.19298246 0.15632754 0.08571429 0.08270677 0.20720721 0.13185913
str(OR)
##  num [1:1000] 0.193 0.1563 0.0857 0.0827 0.2072 ...
mean(OR)
## [1] 0.173435
par(mfrow = c(1, 2))
hist(OR,prob=T)
hist(log(OR),prob=T)

オッズ比はデータ数に依存しない。 ことがわかる。

2.2 カイ二乗値

\(T=\dfrac{(ad-bc)^2N}{(a+b)(c+d)(a+c)(b+d)}\).

K=c()
for (i in 1:100) {
  A <- rbinom(6, 1, 1/3)
  a <- sum(A)
  b <- (6-sum(A))
  B <- rbinom(4, 1, 1/3)
  c <- sum(B)
  d <- (4-sum(B))
  k <- (abs(a*d-b*c))^2*10/((a+b)*(c+d)*(a+c)*(b+d))
  K <- c(K,k)
}
hist(K,prob=T)
curve(dchisq(x, 1),add=T)

L=c()
for (i in 1:100) {
  A <- rbinom(6, 1, 1/3)
  a <- sum(A)
  b <- (6-sum(A))
  B <- rbinom(4, 1, 3/4)
  c <- sum(B)
  d <- (4-sum(B))
  k <- (a*d-b*c)^2*10/((a+b)*(c+d)*(a+c)*(b+d))
  L <- c(L,k)
}
par(mfrow = c(1, 2))
hist(K,prob=T,ylim=c(0,0.5),xlim=c(0,12))
hist(L,prob=T,ylim=c(0,0.5),xlim=c(0,12))
curve(dchisq(x, 1),add=T)

sum(L>=3.84)/100#検出力
## [1] NA
#curve(dchisq(x, 1),add=T)