# 2019 Q4(2)  2024.9.24

import numpy as np
from scipy.integrate import quad

# コーシー分布の確率密度関数 (theta = 1)
def cauchy_pdf_theta_1(x):
    return 1 / (np.pi * (1 + (x - 1)**2))

# 積分範囲（棄却域R = (1, 3)）
lower_bound = 1
upper_bound = 3

# 積分を実行
power, error = quad(cauchy_pdf_theta_1, lower_bound, upper_bound)

# 結果を表示（小数第3位まで表示）
print(f"検出力 (1 - β): {power:.3f}")

検出力 (1 - β): 0.352

# 2019 Q4(2)  2024.9.24

import numpy as np
import matplotlib.pyplot as plt
from scipy.stats import cauchy

# シミュレーションのパラメータ
np.random.seed(43)
num_trials = 100000

# 棄却域 R = (1, 3)
lower_bound = 1
upper_bound = 3

# コーシー分布 (θ = 1) からのサンプルを生成
samples = cauchy.rvs(loc=1, scale=1, size=num_trials)

# -8 から 8 の範囲にサンプルを制限して外れ値を除外
filtered_samples = samples[(samples > -8) & (samples < 8)]

# 棄却域に入っているかどうかを判定
reject = (samples > lower_bound) & (samples < upper_bound)

# 棄却域に入った割合が検出力 (1 - β)
power_simulated = np.mean(reject)

# 結果を表示
print(f"シミュレーションによる検出力 (1 - β): {power_simulated:.3f}")

# ヒストグラムのプロット（フィルタリングしたサンプルを使う）
plt.figure(figsize=(10, 6))
plt.hist(filtered_samples, bins=100, density=True, alpha=0.5, color='blue', label='シミュレーション結果（フィルタリング済み）')

# 理論的なコーシー分布の確率密度関数 (theta = 1) をプロット
x = np.linspace(-8, 8, 1000)
pdf_theta_1 = cauchy.pdf(x, loc=1)  # 対立仮説の理論曲線 (θ = 1)
plt.plot(x, pdf_theta_1, 'r-', lw=2, label='理論的コーシー分布 (θ = 1)')

# 理論的なコーシー分布の確率密度関数 (theta = 0) をプロット
pdf_theta_0 = cauchy.pdf(x, loc=0)  # 帰無仮説の理論曲線 (θ = 0)
plt.plot(x, pdf_theta_0, 'g--', lw=2, label='理論的コーシー分布 (θ = 0)')

# 棄却域を塗りつぶす
plt.axvspan(lower_bound, upper_bound, color='red', alpha=0.3, label='棄却域 (1 < x < 3)')

# グラフの設定
plt.xlim(-8, 8)
plt.xlabel('x')
plt.ylabel('密度')
plt.title(f'コーシー分布のシミュレーションと棄却域\nシミュレーションによる検出力 (1 - β) = {power_simulated:.3f}')
plt.legend()
plt.grid(True)

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

シミュレーションによる検出力 (1 - β): 0.355