ch4_practice.R

# 環境の初期化 ----------------------------------------------------------

rm(list = ls())

## 練習問題1 -------------------------------------------------------------------------------------------
#　まず、標本分散の正の平方根を返す関数を定義する
sd_p <- function(x){
  n <- length(x)
  return(
    (var(x) * (n - 1) / n) |> sqrt()
  )
}

## 設定と準備（ここでは一例として、母集団分布を標準正規分布とする）
n <- 20 # サンプルサイズ
iter <- 20000
# 結果を格納するオブジェクト
sample_sd <- rep(0, each = iter)
## シミュレーション
set.seed(123)
for (i in 1:iter) {
  sample_sd[i] <- rnorm(n) |> sd_p() # 自作関数を使用
}
## 結果（標準正規分布に従う確率変数の標準偏差は1だが、それと一致しない）
mean(sample_sd)


## 練習問題2 -------------------------------------------------------------------------------------------
## 設定と準備
n <- 20 # サンプルサイズ
iter <- 100000

k <- 5
theta <- 0.4
mu <- k * theta # 母平均

# 結果を格納するオブジェクト
# 母集団分布がk = 5, θ = 0.4の二項分布のとき、95％信頼区間が母平均を含んだ回数
true_num_binom <- 0
## シミュレーション
set.seed(123)
for (i in 1:iter) {
  rnd <- rbinom(n, size = k, prob = theta)
  t <- (mean(rnd) - mu) / (sd(rnd) / sqrt(n))
  if (qt(p = 0.025, df = n - 1) <= t & t <= qt(p = 0.975, df = n - 1)) {
    true_num_binom <- true_num_binom + 1
  }
}
## 結果
# 母集団分布が = 5, θ = 0.4の二項分布のとき、95％信頼区間が母平均を含んだ割合
true_num_binom / iter


## 練習問題3 -------------------------------------------------------------------------------------------

## 設定と準備
n <- 50 # サンプルサイズ
dat_obs <- rnorm(n) # 母集団分布が標準正規分布なので、標本平均は正規分布に従う
B <- 2000 # ブートストラップ標本数

# 結果を格納するオブジェクト
boot_mean <- rep(0, each = B)

## シミュレーション
set.seed(123)
for (i in 1:B) {
  # 1, 2, ..., n（dat_obsの要素数）の数列から、n個を復元抽出
  # 書籍では相関係数を例に説明しているので、dat_obsはデータフレームだったためnrow()で行数を取得したが、
  # この練習問題ではdat_obsはベクトルなので、length()で要素数を取得している
  row_num <- sample(1:length(dat_obs),
                    size = length(dat_obs),
                    replace = TRUE)
  # row_numに合致する行番号を抽出
  resampled_data <- dat_obs[row_num]
  # ブートストラップ標本から、標本平均の実現値を求める
  boot_mean[i] <- mean(resampled_data)
}

## 信頼区間
quantile(boot_mean, prob = c(0.025, 0.975)) # パーセンタイル信頼区間

t.test(dat_obs)$conf.int[1:2] # t.test()関数を用いて計算した95%信頼区間