﻿rm(list = ls())  # ワークスペースをクリア
pb <- 0.05  # 男の子がこの球技に取り組む確率
pg <- 0.00  # 女の子がこの球技に取り組む確率
n <- 100000  # 試行回数
denominator <- 0 ; numerator <- 0  # 分母，分子
set.seed(1)  # 乱数のシード値を固定
for(i in 1:n) {
  n1 <- 0  # 上の子がこの球技に取り組む場合は1，それ以外は0
  n2 <- 0  # 下の子がこの球技に取り組む場合は1，それ以外は0
  x <- sample(c("b","g"), 1)  # 上の子の性別 (bは男, gは女)
  if(x == "b") {  # 上の子が男の場合
    # 確率pbでn1に1を代入，確率(1-pb)でn1に0を代入
    n1 <- sample(c(1,0), 1, TRUE, c(pb,1-pb)) 
  } else {  # 上の子が女の場合
    # 確率pgでn1に1を代入，確率(1-pg)でn1に0を代入
    n1 <- sample(c(1,0), 1, TRUE, c(pg,1-pg))
  }
  y <- sample(c("b","g"), 1)  # 下の子の性別 (bは男, gは女)
  if(y == "b") {  # 下の子が男の場合
    # 確率pbでn2に1を代入，確率(1-pb)でn2に0を代入
    n2 <- sample(c(1,0), 1, TRUE, c(pb,1-pb)) 
  } else {  # 下の子が女の場合
    # 確率pgでn2に1を代入，確率(1-pg)でn2に0を代入
    n2 <- sample(c(1,0), 1, TRUE, c(pg,1-pg))
  }
  if(n1+n2 >= 1) {  # 少なくとも1人はこの球技に取り組む場合
    denominator <- denominator + 1
    if(x=="b" && y=="b") {  # 2人とも男の場合
      numerator <- numerator + 1
    }
  }
}
# 2人の子供の性別がどちらも男である確率(モンテカルロ法)
numerator / denominator