# 2016 Q4(4)  2024.11.26

# n を同じ値に設定して三つの推定量を比較する

# 必要なライブラリをインポート
import numpy as np
import matplotlib.pyplot as plt

# パラメータ設定
n_common = 10000  # n を全ての推定量で共通に設定
n_trials = 1000  # シミュレーション試行回数

# 推定量のリストを作成
theta_hat_1_values = []
theta_hat_2_values = []
theta_hat_3_values = []

# シミュレーション実行
for _ in range(n_trials):
    # 推定量 θ̂1
    samples_1 = np.random.normal(0, 1, n_common)
    X = np.sum((samples_1 >= 0) & (samples_1 <= 1))
    theta_hat_1_values.append(X / n_common)

    # 推定量 θ̂2
    samples_2 = np.random.normal(0, 1, n_common)
    Y = np.sum(np.abs(samples_2) <= 1)
    theta_hat_2_values.append(Y / (2 * n_common))

    # 推定量 θ̂3
    U = np.random.uniform(0, 1, n_common)
    theta_hat_3 = np.mean((1 / np.sqrt(2 * np.pi)) * np.exp(-U**2 / 2))
    theta_hat_3_values.append(theta_hat_3)

# 各推定量の分散を計算
var_theta_hat_1 = np.var(theta_hat_1_values)
var_theta_hat_2 = np.var(theta_hat_2_values)
var_theta_hat_3 = np.var(theta_hat_3_values)

# グラフ描画
fig = plt.figure(figsize=(12, 8))  # 全体の図を作成
gs = fig.add_gridspec(3, 2, width_ratios=[3, 1])  # ヒストグラムとバーグラフの比率を指定

# 共通の横軸レンジ
x_min, x_max = 0.325, 0.355

# ヒストグラム: θ̂1
ax1 = fig.add_subplot(gs[0, 0])
ax1.hist(theta_hat_1_values, bins=20, density=True, color='blue', alpha=0.7, label='$\\hat{\\theta}_1$')
ax1.axvline(theta, color='red', linestyle='--', label='理論値 $\\theta$')
ax1.set_xlim(x_min, x_max)
ax1.set_title('$\\hat{\\theta}_1$ の分布 (n = 1,000)', fontsize=12)
ax1.set_ylabel('確率密度', fontsize=10)
ax1.legend(fontsize=10)
ax1.grid(alpha=0.3)

# ヒストグラム: θ̂2
ax2 = fig.add_subplot(gs[1, 0])
ax2.hist(theta_hat_2_values, bins=20, density=True, color='green', alpha=0.7, label='$\\hat{\\theta}_2$')
ax2.axvline(theta, color='red', linestyle='--', label='理論値 $\\theta$')
ax2.set_xlim(x_min, x_max)
ax2.set_title('$\\hat{\\theta}_2$ の分布 (n = 1,000)', fontsize=12)
ax2.set_ylabel('確率密度', fontsize=10)
ax2.legend(fontsize=10)
ax2.grid(alpha=0.3)

# ヒストグラム: θ̂3
ax3 = fig.add_subplot(gs[2, 0])
ax3.hist(theta_hat_3_values, bins=20, density=True, color='orange', alpha=0.7, label='$\\hat{\\theta}_3$')
ax3.axvline(theta, color='red', linestyle='--', label='理論値 $\\theta$')
ax3.set_xlim(x_min, x_max)
ax3.set_title('$\\hat{\\theta}_3$ の分布 (n = 1,000)', fontsize=12)
ax3.set_xlabel('推定値', fontsize=10)
ax3.set_ylabel('確率密度', fontsize=10)
ax3.legend(fontsize=10)
ax3.grid(alpha=0.3)

# バーグラフ
ax_bar = fig.add_subplot(gs[:, 1])  # 全行を使う
labels = ['$\\hat{\\theta}_1$', '$\\hat{\\theta}_2$', '$\\hat{\\theta}_3$']
variances = [var_theta_hat_1, var_theta_hat_2, var_theta_hat_3]
colors = ['blue', 'green', 'orange']

ax_bar.bar(labels, variances, color=colors, alpha=0.7)
ax_bar.set_title('各推定量の分散 (n = 1,000)', fontsize=12)
ax_bar.set_ylabel('分散', fontsize=10)
ax_bar.grid(axis='y', alpha=0.3)

# 全体の調整
plt.tight_layout()
plt.show()

# 2016 Q4(4)  2024.11.26

# 必要なライブラリをインポート
import numpy as np
import matplotlib.pyplot as plt

# パラメータ設定
n1, n2, n3 = 10000, 2411, 102  # 各推定量のサンプルサイズ
theta = 0.3413  # 理論値
n_trials = 1000  # シミュレーション試行回数

# 推定量のリストを作成
theta_hat_1_values = []
theta_hat_2_values = []
theta_hat_3_values = []

# シミュレーション実行
for _ in range(n_trials):
    # 推定量 θ̂1
    samples_1 = np.random.normal(0, 1, n1)
    X = np.sum((samples_1 >= 0) & (samples_1 <= 1))
    theta_hat_1_values.append(X / n1)

    # 推定量 θ̂2
    samples_2 = np.random.normal(0, 1, n2)
    Y = np.sum(np.abs(samples_2) <= 1)
    theta_hat_2_values.append(Y / (2 * n2))

    # 推定量 θ̂3
    U = np.random.uniform(0, 1, n3)
    theta_hat_3 = np.mean((1 / np.sqrt(2 * np.pi)) * np.exp(-U**2 / 2))
    theta_hat_3_values.append(theta_hat_3)

# 各推定量の分散を計算
var_theta_hat_1 = np.var(theta_hat_1_values)
var_theta_hat_2 = np.var(theta_hat_2_values)
var_theta_hat_3 = np.var(theta_hat_3_values)

# グラフ描画
fig = plt.figure(figsize=(12, 8))  # 全体の図を作成
gs = fig.add_gridspec(3, 2, width_ratios=[3, 1])  # ヒストグラムとバーグラフの比率を指定

# 共通の横軸レンジ
x_min, x_max = 0.325, 0.355

# ヒストグラム: θ̂1
ax1 = fig.add_subplot(gs[0, 0])
ax1.hist(theta_hat_1_values, bins=20, density=True, color='blue', alpha=0.7, label='$\\hat{\\theta}_1$')
ax1.axvline(theta, color='red', linestyle='--', label='理論値 $\\theta$')
ax1.set_xlim(x_min, x_max)
ax1.set_title('$\\hat{\\theta}_1$ の分布 (n = 10,000)', fontsize=12)
ax1.set_ylabel('確率密度', fontsize=10)
ax1.legend(fontsize=10)
ax1.grid(alpha=0.3)

# ヒストグラム: θ̂2
ax2 = fig.add_subplot(gs[1, 0])
ax2.hist(theta_hat_2_values, bins=20, density=True, color='green', alpha=0.7, label='$\\hat{\\theta}_2$')
ax2.axvline(theta, color='red', linestyle='--', label='理論値 $\\theta$')
ax2.set_xlim(x_min, x_max)
ax2.set_title('$\\hat{\\theta}_2$ の分布 (n = 2,411)', fontsize=12)
ax2.set_ylabel('確率密度', fontsize=10)
ax2.legend(fontsize=10)
ax2.grid(alpha=0.3)

# ヒストグラム: θ̂3
ax3 = fig.add_subplot(gs[2, 0])
ax3.hist(theta_hat_3_values, bins=20, density=True, color='orange', alpha=0.7, label='$\\hat{\\theta}_3$')
ax3.axvline(theta, color='red', linestyle='--', label='理論値 $\\theta$')
ax3.set_xlim(x_min, x_max)
ax3.set_title('$\\hat{\\theta}_3$ の分布 (n = 102)', fontsize=12)
ax3.set_xlabel('推定値', fontsize=10)
ax3.set_ylabel('確率密度', fontsize=10)
ax3.legend(fontsize=10)
ax3.grid(alpha=0.3)

# バーグラフ
ax_bar = fig.add_subplot(gs[:, 1])  # 全行を使う
labels = ['$\\hat{\\theta}_1$', '$\\hat{\\theta}_2$', '$\\hat{\\theta}_3$']
variances = [var_theta_hat_1, var_theta_hat_2, var_theta_hat_3]
colors = ['blue', 'green', 'orange']

ax_bar.bar(labels, variances, color=colors, alpha=0.7)
ax_bar.set_title('各推定量の分散', fontsize=12)
ax_bar.set_ylabel('分散', fontsize=10)
ax_bar.grid(axis='y', alpha=0.3)

# 全体の調整
plt.tight_layout()
plt.show()