# 2016 Q5(2)  2024.11.28

import numpy as np
import matplotlib.pyplot as plt

# シミュレーションのパラメータ設定
n = 100  # 全体のサンプル数
repeats = 10000  # シミュレーションの繰り返し回数

# シミュレーション結果を格納するリスト
d2_direct_simulations = []
d2_simplified_simulations = []

# シミュレーションを実行
np.random.seed(42)  # 再現性のため乱数の種を固定
for _ in range(repeats):
    # 群間平方和 (SSB), 群内平方和 (SSW) を正の値でランダム生成
    SSB = np.random.uniform(0.1, 10)
    SSW = np.random.uniform(0.1, 10)
    
    # SST = SSB + SSW の関係を利用
    SST = SSB + SSW
    
    # F値を計算
    F = SSB / (SSW / (n - 2))
    
    # d^2 の直接計算
    d2_direct = SSB / (SST / (n - 1))
    d2_direct_simulations.append(d2_direct)
    
    # d^2 の簡略化された式による計算
    d2_simplified = (n - 1) * F / (n - 2 + F)
    d2_simplified_simulations.append(d2_simplified)

# 平均絶対差を計算
mean_absolute_difference = np.mean(np.abs(np.array(d2_direct_simulations) - np.array(d2_simplified_simulations)))

# 結果を出力
print(f"シミュレーションによる d^2 の平均絶対差: {mean_absolute_difference:.2e}")

# 数式を含む散布図をプロット
plt.figure(figsize=(8, 6))
plt.scatter(d2_direct_simulations, d2_simplified_simulations, alpha=0.5, s=10, label="シミュレーション結果")
plt.plot([min(d2_direct_simulations), max(d2_direct_simulations)],
         [min(d2_direct_simulations), max(d2_direct_simulations)], 'r--', label="$y=x$")  # y=xの基準線
plt.title("$d^2$ の二式の比較", fontsize=14)
plt.xlabel(r"$d^2 = \frac{\mathrm{SSB}}{\mathrm{SST} / (n-1)}$", fontsize=12)
plt.ylabel(r"$d^2 = \frac{(n-1)F}{n-2+F}$", fontsize=12)
plt.legend(fontsize=12)  # 凡例を追加
plt.grid(alpha=0.3)  # グリッドを薄く表示
plt.show()

シミュレーションによる d^2 の平均絶対差: 3.97e-15