ホーム » Qchan の投稿 (ページ 8)

作者アーカイブ: Qchan

投稿一覧

2016 Q3(2)

線形モデルの未知の母数βの最小二乗推定量を求め、それが不偏推定量であることを確認しました。

2016 Q3(1)

線形モデルの未知の母数βの2つの推定量が不偏推定量であることを確認しました。

2016 Q2(4)

指数分布の上側確率の推定量が、ガンマ分布に従う確率変数の式で表され、それが不偏推定量であることを示しました。

2016 Q2(4)

指数分布の和がガンマ分布になることを示しました。

2016 Q2(3)

指数分布のλと、上側α点を返す関数の最尤推定量を求め、不偏であるか確認しました。

2016 Q2(2)

指数分布の上側α点を返す関数を求めました。

2016 Q2(1)

指数分布の期待値を求めました。

独立でない確率変数の積の期待値について考えました。

独立でない確率変数の積の期待値

2016 Q1(4)

正規分布の母平均μに対してθ=exp(μ)をパラメータとする尤度関数からフィッシャー情報量を求め、その逆数が不偏推定量の分散と一致するか否かを確かめた。

 

コード

V[\tilde{\theta}] > \frac{1}{I(\theta)}となることを確認するために、両辺をグラフで比較してみます。

# 2016 Q1(4)  2024.11.14

import numpy as np
import matplotlib.pyplot as plt

# θ の範囲設定
theta_values = np.linspace(0.1, 5.0, 10)  # 0.1から5.0までの範囲で計算
n_samples = 50  # サンプルサイズ

# 分散 V[θ~] と 1/I(θ) の計算
V_theta_tilde = theta_values**2 * (np.exp(1 / n_samples) - 1)
Fisher_info_inverse = theta_values**2 / n_samples  # 1 / I(θ)

# 折れ線グラフの描画
plt.figure(figsize=(10, 6))
plt.plot(theta_values, V_theta_tilde, label=r'$V[\tilde{\theta}]$ (不偏推定量の分散)', color='blue', linestyle='-', marker='o')
plt.plot(theta_values, Fisher_info_inverse, label=r'$1/I(\theta)$ (フィッシャー情報量の逆数)', color='red', linestyle='--', marker='s')
plt.title('不偏推定量の分散とフィッシャー情報量の逆数の比較', fontsize=16)
plt.xlabel(r'$\theta$', fontsize=14)
plt.ylabel('分散', fontsize=14)
plt.legend(fontsize=12)
plt.grid(True, linestyle='--', alpha=0.7)
plt.show()

わずかにV[\tilde{\theta}] > \frac{1}{I(\theta)}となっているようです。

次に、2式の差V[\tilde{\theta}] - \frac{1}{I(\theta)}をプロットし、その差を確認してみます。

# 2016 Q1(4)  2024.11.14

import numpy as np
import matplotlib.pyplot as plt

# θ の範囲設定
theta_values = np.linspace(0.1, 5.0, 10)  # 0.1から5.0までの範囲で計算
n_samples = 50  # サンプルサイズ

# 分散 V[θ~] と 1/I(θ) の計算
V_theta_tilde = theta_values**2 * (np.exp(1 / n_samples) - 1)
Fisher_info_inverse = theta_values**2 / n_samples  # 1 / I(θ)

# V[θ~] - 1/I(θ) を計算
diff_V_Fisher = V_theta_tilde - Fisher_info_inverse

# 差分をプロット
plt.figure(figsize=(10, 6))
plt.plot(theta_values, diff_V_Fisher, color='purple', linestyle='-', marker='o', label=r'$V[\tilde{\theta}] - \frac{1}{I(\theta)}$')
plt.title('不偏推定量の分散とフィッシャー情報量逆数の差', fontsize=16)
plt.xlabel(r'$\theta$', fontsize=14)
plt.ylabel(r'$V[\tilde{\theta}] - \frac{1}{I(\theta)}$', fontsize=14)
plt.axhline(0, color='black', linestyle='--', linewidth=1)  # y=0 の線を追加
plt.legend(fontsize=12)
plt.grid(True, linestyle='--', alpha=0.7)
plt.show()

2式の差V[\tilde{\theta}] - \frac{1}{I(\theta)}が存在することが確認できました。この差は、θが増加するにつれて大きくなる傾向が見られます。

次に、乱数を用いたシミュレーションを行い、V[\tilde{\theta}]を重ねてプロットしてみます。

# 2016 Q1(4)  2024.11.14

import numpy as np
import matplotlib.pyplot as plt

# 再定義:thetaの範囲
theta_values = np.linspace(0.1, 5.0, 100)
n_samples = 50
n_trials = 10000  # シミュレーション試行回数

# 理論分散 V[θ~] と 1/I(θ) の計算
V_theta_tilde = theta_values**2 * (np.exp(1 / n_samples) - 1)
Fisher_info_inverse = theta_values**2 / n_samples

# シミュレーションによる分散の再計算
simulated_V_theta_tilde = []

for theta in theta_values:
    mse_trials = []

    for _ in range(n_trials):
        sample = np.random.normal(loc=np.log(theta), scale=1, size=n_samples)  # 平均 log(theta), 分散 1
        X_bar = np.mean(sample)
        theta_tilde = np.exp(X_bar - 1 / (2 * n_samples))  # 不偏推定量
        mse_trials.append((theta_tilde - theta)**2)  # 二乗誤差

    simulated_var = np.mean(mse_trials)  # シミュレーションによる分散 (MSE)
    simulated_V_theta_tilde.append(simulated_var)

# グラフの描画
plt.figure(figsize=(10, 6))
plt.plot(theta_values, V_theta_tilde, label='理論分散 $V[\\tilde{\\theta}]$', linestyle='-', color='blue', alpha=0.7)
plt.scatter(theta_values, simulated_V_theta_tilde, label='シミュレーション分散 $V[\\tilde{\\theta}]$', marker='s', color='red', alpha=0.6)
plt.plot(theta_values, Fisher_info_inverse, label='$1/I(\\theta)$ (フィッシャー情報量の逆数)', linestyle='--', color='green', alpha=0.8)
plt.title('不偏推定量の分散とフィッシャー情報量の逆数の比較', fontsize=16)
plt.xlabel('$\\theta$', fontsize=14)
plt.ylabel('分散', fontsize=14)
plt.legend(fontsize=12)
plt.grid(True, linestyle='--', alpha=0.7)
plt.show()

V[\tilde{\theta}]は、わずかに\frac{1}{I(\theta)}を下回る場合があるようです。これは、シミュレーションにおける乱数の誤差によるものだと考えられます。

2016 Q1(3)

正規分布の母平均μに対するθ=exp(μ)の最尤推定量の平均二乗誤差が0になることを示しました。

 

コード

シミュレーションを行います。\hat{\theta} = e^{\bar{X}}に対して、MSE[\hat{\theta}] = E\left[(\hat{\theta} - \theta)^2\right]を計算し、理論値とともにプロットします。サンプルサイズnを増加させて、MSEの変化を観察します。

# 2016 Q1(3)  2024.11.13

import numpy as np
import matplotlib.pyplot as plt

# 真の値とシミュレーション設定
mu_true = 1.0  # 真の平均 (mu)
theta_true = np.exp(mu_true)  # 真の θ = e^mu
n_trials = 10000  # シミュレーション試行回数
sample_sizes = np.arange(5, 101, 5)  # サンプルサイズの範囲

# MSEを格納するリスト
mse_simulated = []

# 理論的なMSEの計算関数
def theoretical_mse(n, mu):
    return np.exp(2 * mu) * (np.exp(2 / n) - 2 * np.exp(1 / (2 * n)) + 1)

mse_theoretical = [theoretical_mse(n, mu_true) for n in sample_sizes]

# シミュレーションでのMSE計算
for n_samples in sample_sizes:
    mse_trials = []

    for _ in range(n_trials):

        sample = np.random.normal(loc=mu_true, scale=1, size=n_samples)
        X_bar = np.mean(sample)
        
        theta_hat = np.exp(X_bar)  # 最尤推定量
        mse_trials.append((theta_hat - theta_true) ** 2)  # MSEの計算

    mse_simulated.append(np.mean(mse_trials))  # 試行平均を格納

# グラフの描画
plt.figure(figsize=(10, 6))
plt.plot(sample_sizes, mse_simulated, marker='o', linestyle='-', color='blue', label='シミュレーション MSE')
plt.plot(sample_sizes, mse_theoretical, marker='s', linestyle='--', color='red', label='理論的 MSE')
plt.title('MSEのシミュレーションと理論値', fontsize=16)
plt.xlabel('サンプルサイズ $n$', fontsize=14)
plt.ylabel('MSE', fontsize=14)
plt.legend(fontsize=12)
plt.grid(True)
plt.show()

サンプルサイズnが増加するにつれて、MSE[\hat{\theta}]は0に近づきました。これにより、\lim_{n \to \infty} MSE[\hat{\theta}] = 0になることが理論とシミュレーションの両方で確認できました。また、この結果から最尤推定量 \hat{\theta} が漸近的一致性を持つことが分かりました。