import sympy as sp

# 定義
lambda_var = sp.symbols('lambda')
alpha, beta = sp.symbols('alpha beta', positive=True)

# ガンマ分布の確率密度関数
gamma_pdf = (beta**alpha / sp.gamma(alpha)) * lambda_var**(alpha-1) * sp.exp(-beta * lambda_var)

# 期待値 E[Λ] の計算
expected_value = sp.integrate(lambda_var * gamma_pdf, (lambda_var, 0, sp.oo)).simplify()

# E[Λ^2] の計算
expected_value_Lambda2 = sp.integrate(lambda_var**2 * gamma_pdf, (lambda_var, 0, sp.oo)).simplify()

# 分散 V[Λ] の計算
variance = (expected_value_Lambda2 - expected_value**2).simplify()

# 結果を辞書形式で表示
{
    "期待値 E[Λ]": expected_value,
    "分散 V[Λ]": variance
}

{'期待値 E[Λ]': alpha/beta, '分散 V[Λ]': alpha/beta**2}

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

# パラメータ α と β の設定
alpha = 3.0
beta = 2.0
# サンプルサイズ
sample_size = 10000

# ガンマ分布に従う乱数を生成
samples = np.random.gamma(alpha, 1/beta, sample_size)

# 期待値と分散の計算
expected_value = np.mean(samples)
variance = np.var(samples)

# 結果の表示
print(f"期待値のシミュレーション結果: {expected_value}")
print(f"分散のシミュレーション結果: {variance}")

# ヒストグラムの描画
plt.hist(samples, bins=50, density=True, alpha=0.75, color='blue', edgecolor='black')

# 理論的なガンマ分布の確率密度関数をプロット
x = np.linspace(0, max(samples), 100)
gamma_pdf = gamma.pdf(x, alpha, scale=1/beta)
plt.plot(x, gamma_pdf, 'r', linestyle='-', label='理論的なガンマ分布')

# グラフのタイトルとラベル
plt.title('ガンマ分布のシミュレーション結果')
plt.xlabel('値')
plt.ylabel('確率密度')
plt.legend()

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

期待値のシミュレーション結果: 1.4947527551017439
分散のシミュレーション結果: 0.7452822182213243

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

def custom_exponential(beta):
    # [0, 1) の一様乱数を生成
    u = np.random.rand()
    # 指数分布に従う乱数を生成
    return -np.log(u) / beta

def custom_gamma(alpha, beta, sample_size):
    samples = []
    for _ in range(sample_size):
        sample = sum(custom_exponential(beta) for _ in range(int(alpha)))
        samples.append(sample)
    return samples

# パラメータ α と β の設定
alpha = 3.0
beta = 2.0
# サンプルサイズ
sample_size = 10000

# カスタム関数を使ってガンマ分布に従う乱数を生成
samples = custom_gamma(alpha, beta, sample_size)

# 期待値と分散の計算
expected_value = np.mean(samples)
variance = np.var(samples)

# 理論値
theoretical_expected_value = alpha / beta
theoretical_variance = alpha / beta**2

# 結果の表示
print(f"シミュレーション結果 - 期待値: {expected_value}, 理論値: {theoretical_expected_value}")
print(f"シミュレーション結果 - 分散: {variance}, 理論値: {theoretical_variance}")

# ヒストグラムの描画
plt.hist(samples, bins=50, density=True, alpha=0.75, color='blue', edgecolor='black')

# 理論的なガンマ分布の確率密度関数をプロット
x = np.linspace(0, max(samples), 100)
gamma_pdf = gamma.pdf(x, alpha, scale=1/beta)
plt.plot(x, gamma_pdf, 'r', linestyle='-', label='理論的なガンマ分布')

# グラフのタイトルとラベル
plt.title('ガンマ分布のシミュレーション結果')
plt.xlabel('値')
plt.ylabel('確率密度')
plt.legend()

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

シミュレーション結果 - 期待値: 1.4966126702823181, 理論値: 1.5
シミュレーション結果 - 分散: 0.7242038142114708, 理論値: 0.75