﻿rm(list = ls())  # ワークスペースをクリア
m <- 100000  # 試行回数
n <- 5  # 景品の種類
Yn <- numeric(m)  # 景品を全種類集めるまでの購入回数を格納するベクトル
set.seed(1)  # 乱数のシード値を固定
for (i in 1:m) {
  num_of_purchases <- 0  # 購入回数
  coupon_vec <- integer(n)  # n次元ベクトル(初期値は(0,0,...,0))
  while (sum(coupon_vec) != n) {
    k <- sample.int(n, 1)  # 1からnまでの整数から無作為に1つの整数を取得
    coupon_vec[k] <- 1  # k番目の要素を1とし，k番目の景品を取得済みとみなす
    num_of_purchases <- num_of_purchases + 1
  }
  Yn[i] <- num_of_purchases
}
# 平均と分散の理論値を計算
t <- sum(1 / 1:n)
mean_theoretical <- n * t
u <- sum((1:(n - 1)) / (n - (1:(n - 1)))^2)
var_theoretical <- n * u
# 平均E(Yn)を表示
mean(Yn)  # モンテカルロ法
mean_theoretical  # 参考：理論値
# 分散V(Yn)を表示
var(Yn)  # モンテカルロ法
var_theoretical  # 参考：理論値