# 2018 Q5(1)  2024.10.21

import numpy as np
import matplotlib.pyplot as plt

# 定義域
y_values = np.linspace(0, 1, 500)

# 確率密度関数の定義
f1_y = 3 * (1 - y_values) ** 2
f3_y = 3 * y_values ** 2

# シミュレーションの設定
n_samples = 10000   # サンプル数

# 一様分布から3つの変数をサンプリング
samples = np.random.uniform(0, 1, (n_samples, 3))

# 順序統計量を求める
Y1_samples = np.min(samples, axis=1)  # 最小値
Y3_samples = np.max(samples, axis=1)  # 最大値

# 期待値の計算
expected_Y1 = np.mean(Y1_samples)
expected_Y3 = np.mean(Y3_samples)

# 理論値との比較
theoretical_Y1 = 1/4
theoretical_Y3 = 3/4

# 期待値を表示
print(f"シミュレーションによる E[Y1]（最小値の期待値）: {expected_Y1:.4f}, 理論値: {theoretical_Y1:.4f}")
print(f"シミュレーションによる E[Y3]（最大値の期待値）: {expected_Y3:.4f}, 理論値: {theoretical_Y3:.4f}")

# ヒストグラムの作成
plt.figure(figsize=(8, 6))

# 理論密度関数のプロット
plt.plot(y_values, f1_y, label='$f_1(y) = 3(1-y)^2$', color='blue')
plt.plot(y_values, f3_y, label='$f_3(y) = 3y^2$', color='green')

# サンプルのヒストグラムを追加 (確率密度として正規化)
plt.hist(Y1_samples, bins=50, density=True, alpha=0.5, color='blue', label='Y1 samples (min)')
plt.hist(Y3_samples, bins=50, density=True, alpha=0.5, color='green', label='Y3 samples (max)')

# グラフの設定
plt.title('確率密度関数 $f_1(y)$ と $f_3(y)$ (シミュレーションと理論)', fontsize=14)
plt.xlabel('y', fontsize=12)
plt.ylabel('密度', fontsize=12)
plt.legend(loc='upper right', fontsize=10)
plt.grid(True)

# グラフの表示
plt.show()

シミュレーションによる E[Y1]（最小値の期待値）: 0.2499, 理論値: 0.2500
シミュレーションによる E[Y3]（最大値の期待値）: 0.7487, 理論値: 0.7500