# 2016 Q2(2)  2024.11.16

import numpy as np
import matplotlib.pyplot as plt

# パラメータ設定
lambda_value = 2  # λの値
sample_size = 10000  # シミュレーションに使用するサンプル数
c = 1.0  # 上側確率を計算するための定数 c
alpha_values = [0.1, 0.25, 0.5, 0.75, 0.9]  # 上側 100α% 点を確認する α の値

# 1. 上側確率 Q(c) のシミュレーション
random_samples = np.random.exponential(scale=1/lambda_value, size=sample_size)
empirical_prob_c = np.mean(random_samples > c)  # シミュレーションによる上側確率
theoretical_prob_c = np.exp(-lambda_value * c)  # 理論的な上側確率

# 結果表示
print("上側確率 Q(c) のシミュレーション結果")
print(f"理論的な Q(c) = {theoretical_prob_c:.4f}")
print(f"シミュレーションによる Q(c) = {empirical_prob_c:.4f}\n")

# 2. 上側 100α% 点 u(α) の確認
alpha_results = []
for alpha in alpha_values:
    u_alpha = -np.log(alpha) / lambda_value  # 理論的な上側 100α% 点
    empirical_prob_alpha = np.mean(random_samples > u_alpha)  # u(α) を超える割合
    alpha_results.append((alpha, u_alpha, empirical_prob_alpha))

# 結果表示
print("上側 100α% 点 u(α) のシミュレーション結果")
for alpha, u_alpha, empirical_prob_alpha in alpha_results:
    print(f"α = {alpha:.2f} | 理論的な u(α) = {u_alpha:.4f} | 理論的確率P(X>u(α)) = {alpha:.2f} | シミュレーション確率 = {empirical_prob_alpha:.4f}")

# ヒストグラムのプロット
plt.figure(figsize=(10, 6))
plt.hist(random_samples, bins=50, density=True, alpha=0.6, color='skyblue', label='サンプルのヒストグラム')

# 理論的な確率密度関数
x = np.linspace(0, max(random_samples), 1000)
theoretical_pdf = lambda_value * np.exp(-lambda_value * x)
plt.plot(x, theoretical_pdf, 'r-', label='理論的な確率密度関数', linewidth=2)

# 上側 100α% 点 u(α) のプロット
for alpha, u_alpha, _ in alpha_results:
    plt.axvline(u_alpha, linestyle='dashed', label=f'u({alpha}) = {u_alpha:.2f}')

# グラフ設定
plt.title(f'指数分布 Exp({lambda_value}) のシミュレーションと理論分布')
plt.xlabel('x')
plt.ylabel('密度')
plt.legend(fontsize=10)
plt.grid(True)

# プロット表示
plt.show()

上側確率 Q(c) のシミュレーション結果
理論的な Q(c) = 0.1353
シミュレーションによる Q(c) = 0.1365

上側 100α% 点 u(α) のシミュレーション結果
α = 0.10 | 理論的な u(α) = 1.1513 | 理論的確率P(X>u(α)) = 0.10 | シミュレーション確率 = 0.1035
α = 0.25 | 理論的な u(α) = 0.6931 | 理論的確率P(X>u(α)) = 0.25 | シミュレーション確率 = 0.2480
α = 0.50 | 理論的な u(α) = 0.3466 | 理論的確率P(X>u(α)) = 0.50 | シミュレーション確率 = 0.4962
α = 0.75 | 理論的な u(α) = 0.1438 | 理論的確率P(X>u(α)) = 0.75 | シミュレーション確率 = 0.7513
α = 0.90 | 理論的な u(α) = 0.0527 | 理論的確率P(X>u(α)) = 0.90 | シミュレーション確率 = 0.9011

# 2016 Q2(1)  2024.11.15

import numpy as np
import matplotlib.pyplot as plt

# シミュレーションとプロット
def simulate_and_plot_exponential(lambda_value=2, sample_size=50000, bins=100):
    # 指数分布に従う乱数を生成
    random_samples = np.random.exponential(scale=1/lambda_value, size=sample_size)

    # 理論的な確率密度関数を計算
    x = np.linspace(0, 5, 1000)
    theoretical_pdf = lambda_value * np.exp(-lambda_value * x)

    # サンプルのヒストグラムを描画
    plt.figure(figsize=(10, 6))
    plt.hist(random_samples, bins=bins, density=True, alpha=0.6, color='skyblue', label='サンプルのヒストグラム')

    # 理論的な確率密度関数を重ねる
    plt.plot(x, theoretical_pdf, 'r-', label='理論的な確率密度関数', linewidth=2)

    # 期待値を計算
    simulated_mean = np.mean(random_samples)
    theoretical_mean = 1 / lambda_value

    # 期待値を直線でプロット
    plt.axvline(simulated_mean, color='green', linestyle='dashed', linewidth=2, label=f'シミュレーション期待値: {simulated_mean:.2f}')
    plt.axvline(theoretical_mean, color='blue', linestyle='dotted', linewidth=2, label=f'理論期待値: {theoretical_mean:.2f}')

    # グラフの設定
    plt.title(f'指数分布のシミュレーションと理論分布 (λ={lambda_value})', fontsize=16)
    plt.xlabel('x', fontsize=14)
    plt.ylabel('密度', fontsize=14)
    plt.legend(fontsize=12)
    plt.grid(True)

    # プロットを表示
    plt.show()

# 実行
simulate_and_plot_exponential()

# 2016 Q2(1)  2024.11.15

import numpy as np
import matplotlib.pyplot as plt

# シミュレーションとプロット
def simulate_and_plot_adjusted_lambdas(lambdas=[0.5, 1, 2, 3], sample_size=50000, bins=100):
    fig, axes = plt.subplots(2, 2, figsize=(12, 10), sharex=True, sharey=True)
    x = np.linspace(0, 6, 1000)  # 横軸の範囲を0～6に限定

    for i, lambda_value in enumerate(lambdas):
        row, col = divmod(i, 2)
        ax = axes[row, col]
        
        # 乱数生成と理論値計算
        random_samples = np.random.exponential(scale=1/lambda_value, size=sample_size)
        theoretical_pdf = lambda_value * np.exp(-lambda_value * x)
        
        # ヒストグラムと理論分布
        ax.hist(random_samples, bins=bins, range=(0, 6), density=True, alpha=0.6, color='skyblue', label='サンプルのヒストグラム')
        ax.plot(x, theoretical_pdf, 'r-', label='理論的な確率密度関数', linewidth=2)
        
        # 期待値を計算してプロット
        simulated_mean = np.mean(random_samples)
        theoretical_mean = 1 / lambda_value
        ax.axvline(simulated_mean, color='green', linestyle='dashed', linewidth=2, label=f'シミュレーション期待値: {simulated_mean:.2f}')
        ax.axvline(theoretical_mean, color='blue', linestyle='dotted', linewidth=2, label=f'理論期待値: {theoretical_mean:.2f}')
        
        # グラフ設定
        ax.set_xlim(0, 6)  # 横軸を0～6に限定
        ax.set_title(f'λ={lambda_value}', fontsize=14)
        ax.set_xlabel('x', fontsize=12)
        ax.set_ylabel('密度', fontsize=12)
        ax.legend(fontsize=10)
        ax.grid(True)
    
    # 全体のタイトルと調整
    plt.suptitle('異なるλにおける指数分布のシミュレーションと理論分布', fontsize=16)
    plt.tight_layout(rect=[0, 0, 1, 0.96])  # タイトルとの重なり防止
    plt.show()

# 実行
simulate_and_plot_adjusted_lambdas()