はじめに
機械学習で登場する確率分布について色々な角度から理解したいシリーズです。
この記事では、1次元(一変量)スチューデントのt分布のグラフを作成します。
【前の内容】
【他の記事一覧】
【この記事の内容】
1次元スチューデントのt分布の作図
1次元スチューデントのt分布(Student's t-Distribution)のグラフを作成します。t分布については「1次元スチューデントのt分布の定義式の確認 - からっぽのしょこ」を参照してください。
利用するライブラリを読み込みます。
# 利用ライブラリ import numpy as np from scipy.stats import t import matplotlib.pyplot as plt from matplotlib.animation import FuncAnimation
分布の変化をアニメーション(gif画像)で確認するのにmatplotlibライブラリのanimationモジュールのFuncAnimation関数を利用します。不要であれば省略してください。
定義式の確認
まずは、標準化されたスチューデントのt分布と一般化された(標準化されていない)スチューデントのt分布の定義式を確認します。(標準化・一般化という表現は雰囲気です。正しい言い方があれば教えてください。)
標準化t分布は、次の式で定義されます。
$$
\mathrm{St}(x | \nu)
= \frac{
\Gamma(\frac{\nu + 1}{2})
}{
\sqrt{\pi \nu}
\Gamma(\frac{\nu}{2})
}
\left(
1 + \frac{x^2}{\nu}
\right)^{-\frac{(\nu+1)}{2}}
$$
ここで、$\nu$は自由度、$\pi$は円周率、$\Gamma(x)$はガンマ関数です。自由度は形状パラメータとも呼ばれます。
確率変数の実現値$x$は、実数となります。$\nu$は、正の整数を満たす必要があります。
また、一般化t分布は、次の2つの式で定義されます。
$$
\begin{aligned}
\mathrm{St}(x | \nu, \mu, \sigma)
&= \frac{
\Gamma(\frac{\nu + 1}{2})
}{
\Gamma(\frac{\nu}{2})
}
\frac{1}{\sqrt{\pi \nu} \sigma}
\left\{
1 + \frac{1}{\nu} \left(
\frac{x - \mu}{\sigma}
\right)^2
\right\}^{-\frac{(\nu+1)}{2}}
\\
\mathrm{St}(x | \nu, \mu, \lambda)
&= \frac{
\Gamma(\frac{\nu + 1}{2})
}{
\Gamma(\frac{\nu}{2})
}
\left(
\frac{\lambda}{\pi \nu}
\right)^{\frac{1}{2}}
\left\{
1
+ \frac{\lambda (x - \mu)^2}{\nu}
\right\}^{-\frac{(\nu+1)}{2}}
\end{aligned}
$$
ここで、$\mu$は位置パラメータ、$\sigma$はスケールパラメータ、$\lambda$は逆スケールパラメータです。
$\mu$は実数、$\sigma$は正の実数、$\lambda$は正の実数を満たす必要があります。また、$\sigma = \frac{1}{\sqrt{\lambda}}$、$\lambda = \frac{1}{\sigma^2}$の関係が成り立ちます。
$\mu = 0$、$\sigma = \lambda = 1$のとき、標準化t分布と一般化t分布は一致します。
$$
\mathrm{St}(x | \nu)
= \mathrm{St}(x | \nu, \mu = 0, \sigma = 1)
= \mathrm{St}(x | \nu, \mu = 0, \lambda = 1)
$$
t分布の期待値・分散・最頻値は、次の式で計算できます。詳しくはいつか書きます。
$$
\begin{aligned}
\mathbb{E}[x]
&= \mu
\quad (\nu > 1)
\\
\mathbb{V}[x]
&= \sigma^2
\frac{\nu}{\nu - 2}
\quad (\nu > 2)
\\
&= \frac{1}{\lambda}
\frac{\nu}{\nu - 2}
\quad (\nu > 2)
\\
\mathrm{mode}[x]
&= \mu
\end{aligned}
$$
これらの計算を行いグラフを作成します。
グラフの作成
Matplotlibライブラリを利用して、スチューデントのt分布のグラフを作成します。t分布の確率や統計量の計算については「【Python】1次元スチューデントのt分布の計算 - からっぽのしょこ」を参照してください。
t分布のパラメータ$\nu, \mu, \sigma$を設定します。
# 形状パラメータ(自由度)を指定 nu = 5 # 位置パラメータを指定 mu = 2.0 # 尺度パラメータを指定 sigma = 0.5
形状パラメータ(自由度)(正の整数)$\nu$、位置パラメータ(実数)$\mu$、スケールパラメータ(正の実数)$\sigma$を指定します。
t分布の確率変数がとり得る値$x$を作成します。
# xの値を作成 x_vals = np.linspace(start=mu-sigma*5, stop=mu+sigma*5, num=251) print(x_vals[:10])
[-0.5 -0.48 -0.46 -0.44 -0.42 -0.4 -0.38 -0.36 -0.34 -0.32]
$x$の値を作成してx_valsとします。この例では、muを中心にスケールパラメータの5倍を範囲とします。
$x$の値ごとに確率密度を計算します。
# スチューデントのt分布を計算 density = t.pdf(x=x_vals, df=nu, loc=mu, scale=sigma) print(density[:10])
[0.00351488 0.00365871 0.00380925 0.00396684 0.00413185 0.00430467
0.0044857 0.00467538 0.00487417 0.00508255]
t分布の確率密度は、SciPyライブラリのt分布のクラスtの確率密度メソッドpdf()で計算できます。確率変数の引数x
にx_vals、形状パラメータ(自由度)の引数dfにnu、位置パラメータの引数locにmu、スケールパラメータの引数scaleにsigmaを指定します。逆スケールパラメータ$\lambda$を使う場合は、scale引数に$\lambda$の逆数の平方根$\frac{1}{\sqrt{\lambda}}$を指定します。
t分布のグラフを作成します。
# スチューデントのt分布のグラフを作成 plt.figure(figsize=(12, 9)) # 図の設定 plt.plot(x_vals, density, color='#00A968', linewidth=2.5) # 折れ線グラフ plt.xlabel('x') # x軸ラベル plt.ylabel('density') # y軸ラベル plt.suptitle("Student's t Distribution", fontsize=20) # 全体のタイトル plt.title('$\\nu='+str(nu) + ', \mu='+str(mu) + ', \sigma='+str(sigma)+'$', loc='left') # タイトル plt.grid() # グリッド線 plt.show() # 描画
この分布に統計量の情報を重ねて表示します。
# 補助線用の統計量を計算 E_x = mu s_x = np.sqrt(sigma**2 * nu / (nu - 2)) # 統計量を重ねたt分布のグラフを作成 plt.figure(figsize=(12, 9)) # 図の設定 plt.plot(x_vals, density, color='#00A968', linewidth=2.5, label='$p(x | \\nu, \mu, \sigma)$') # 分布 plt.axvline(x=E_x, color='blue', linewidth=2.5, linestyle='--', label='$E[x]$') # 期待値 plt.axvline(x=E_x-s_x, color='orange', linewidth=2.5, linestyle=':', label='$E[x]\pm \sqrt{V[x]}$') # 期待値 - 標準偏差 plt.axvline(x=E_x+s_x, color='orange', linewidth=2.5, linestyle=':') # 期待値 + 標準偏差 plt.xlabel('x') # x軸ラベル plt.ylabel('density') # y軸ラベル plt.suptitle("Student's t Distribution", fontsize=20) # 全体のタイトル plt.title('$\\nu='+str(nu) + ', \mu='+str(mu) + ', \lambda='+str(lmd)+'$', loc='left') # タイトル plt.legend() # 凡例 plt.grid() # グリッド線 plt.show() # 描画
$\nu > 1$のとき期待値と最頻値が一致します。$\nu \leq 0$のときは期待値を計算できません。
ガンマ分布のグラフを描画できました。以降は、ここまでの作図処理を用いて、パラメータの影響を確認していきます。
パラメータと分布の関係を並べて比較
複数のパラメータのグラフを比較することで、パラメータの値と分布の形状の関係を確認します。
複数の$\nu$を指定し、$\mu, \sigma$を固定して、それぞれ分布を作図します。
# 自由度として利用する値を指定 nu_vals = np.array([1, 2, 3, 4, 5, 10, 100]) # 固定するパラメータを指定 mu = 0.0 sigma = 1.0 # xの値を作成 x_vals = np.linspace(start=mu-sigma*5, stop=mu+sigma*5, num=251) # スチューデントのt分布のグラフを作成 plt.figure(figsize=(12, 9)) # 図の設定 for nu in nu_vals: density = t.pdf(x=x_vals, df=nu, loc=mu, scale=sigma) # 確率密度を計算 plt.plot(x_vals, density, label='$\\nu='+str(nu) + ', \mu='+str(mu) + ', \sigma='+str(sigma)+'$') # 折れ線グラフ plt.xlabel('x') # x軸ラベル plt.ylabel('density') # y軸ラベル plt.suptitle("Student's t Distribution", fontsize=20) # 全体のタイトル plt.legend() # 凡例 plt.grid() # グリッド線 plt.show() # 描画
$\nu$が小さいほど裾が厚くなり、$\nu$が大きいほどガウス分布に近付きます。
複数の$\mu$を指定し、$\nu, \sigma$を固定して、それぞれ分布を作図します。
# 位置パラメータとして利用する値を指定 mu_vals = np.array([-3.5, -2.0, -0.5, 1.0, 2.5]) # 固定するパラメータを指定 nu = 1 sigma = 1.0 # xの値を作成 x_vals = np.linspace(start=mu_vals.min()-sigma, stop=mu_vals.max()+sigma, num=251) # スチューデントのt分布のグラフを作成 plt.figure(figsize=(12, 9)) # 図の設定 for mu in mu_vals: density = t.pdf(x=x_vals, df=nu, loc=mu, scale=sigma) # 確率密度を計算 plt.plot(x_vals, density, label='$\\nu='+str(nu) + ', \mu='+str(mu) + ', \sigma='+str(sigma)+'$') # 折れ線グラフ plt.xlabel('x') # x軸ラベル plt.ylabel('density') # y軸ラベル plt.suptitle("Student's t Distribution", fontsize=20) # 全体のタイトル plt.legend() # 凡例 plt.grid() # グリッド線 plt.show() # 描画
$\mu$が大きいほど、山が右に移動します。
複数の$\sigma$を指定し、$\nu, \mu$を固定して、それぞれ分布を作図します。
# スケールパラメータとして利用する値を指定 sigma_vals = np.array([0.5, 1.0, 2.0, 4.0, 8.0]) # 固定するパラメータを指定 nu = 1 mu = 0.0 # xの値を作成 x_vals = np.linspace(start=mu-sigma_vals.max(), stop=mu+sigma_vals.max(), num=251) # スチューデントのt分布のグラフを作成 plt.figure(figsize=(12, 9)) # 図の設定 for sigma in sigma_vals: density = t.pdf(x=x_vals, df=nu, loc=mu, scale=sigma) # 確率密度を計算 plt.plot(x_vals, density, label='$\\nu='+str(nu) + ', \mu='+str(mu) + ', \sigma='+str(sigma)+'$') # 折れ線グラフ plt.xlabel('x') # x軸ラベル plt.ylabel('density') # y軸ラベル plt.suptitle("Student's t Distribution", fontsize=20) # 全体のタイトル plt.legend() # 凡例 plt.grid() # グリッド線 plt.show() # 描画
$\sigma$が小さいほど、山が細く高くなります。
複数の$\lambda$を指定し、$\nu, \mu$を固定して、それぞれ分布を作図します。
# 逆スケールパラメータとして利用する値を指定 lambda_vals = np.array([0.25, 0.5, 1.0, 2.0, 4.0]) # 固定するパラメータを指定 nu = 1 mu = 0.0 # xの値を作成 sigma = 1.0 / np.sqrt(lambda_vals.min()) x_vals = np.linspace(start=mu-sigma*4, stop=mu+sigma*4, num=251) # スチューデントのt分布のグラフを作成 plt.figure(figsize=(12, 9)) # 図の設定 for lmd in lambda_vals: density = t.pdf(x=x_vals, df=nu, loc=mu, scale=1.0/np.sqrt(lmd)) # 確率密度を計算 plt.plot(x_vals, density, label='$\\nu='+str(nu) + ', \mu='+str(mu) + ', \lambda='+str(lmd)+'$') # 折れ線グラフ plt.xlabel('x') # x軸ラベル plt.ylabel('density') # y軸ラベル plt.suptitle("Student's t Distribution", fontsize=20) # 全体のタイトル plt.legend() # 凡例 plt.grid() # グリッド線 plt.show() # 描画
$\lambda$が大きいほど、山が細く高くなります。
パラメータと分布の関係をアニメーションで可視化
前節では、複数のパラメータのグラフを並べて比較しました。次は、パラメータの値を少しずつ変化させて、分布の形状の変化をアニメーションで確認します。
$\nu$の値を変化させ、$\mu, \sigma$を固定して、それぞれ分布を作図します。
・作図コード(クリックで展開)
# 自由度として利用する値を指定 nu_vals = np.arange(1, 31) print(len(nu_vals)) # フレーム数 # 固定するパラメータを指定 mu = 0.0 sigma = 1.0 # xの値を作成 x_vals = np.linspace(start=mu-sigma*5, stop=mu+sigma*5, num=251) # 確率密度(y軸)の最大値を計算 dens_max = t.pdf(x=x_vals, df=nu_vals.max(), loc=mu, scale=sigma).max() # 図を初期化 fig = plt.figure(figsize=(12, 9), facecolor='white') # 図の設定 fig.suptitle("Student's t Distribution", fontsize=20) # 全体のタイトル # 作図処理を関数として定義 def update(i): # 前フレームのグラフを初期化 plt.cla() # i番目の値を取得 nu = nu_vals[i] # スチューデントのt分布を計算 density = t.pdf(x=x_vals, df=nu, loc=mu, scale=sigma) # スチューデントのt分布を作図 plt.plot(x_vals, density, color='#00A968', linewidth=2.5) # 折れ線グラフ plt.xlabel('x') # x軸ラベル plt.ylabel('density') # y軸ラベル plt.title('$\\nu='+str(np.round(nu, 2)) + ', \mu='+str(mu) + ', \sigma='+str(sigma)+'$', loc='left') # タイトル plt.grid() # グリッド線 plt.ylim(-0.01, dens_max*1.1) # y軸の表示範囲 # gif画像を作成 dens_anime = FuncAnimation(fig, func=update, frames=len(nu_vals), interval=100) # gif画像を保存 dens_anime.save('t_dens_nu.gif')
30
$\mu$の値を変化させ、$\nu, \sigma$を固定して、それぞれ分布を作図します。
・作図コード(クリックで展開)
# 位置パラメータとして利用する値を指定 mu_vals = np.arange(-3.0, 3.1, 0.1) print(len(mu_vals)) # フレーム数 # 固定するパラメータを指定 nu = 1 sigma = 1.0 # xの値を作成 x_vals = np.linspace(start=mu_vals.min()-sigma, stop=mu_vals.max()+sigma, num=251) # 確率密度(y軸)の最大値を計算 dens_max = t.pdf(x=x_vals, df=nu, loc=np.median(mu_vals), scale=sigma).max() # 図を初期化 fig = plt.figure(figsize=(12, 9), facecolor='white') # 図の設定 fig.suptitle("Student's t Distribution", fontsize=20) # 全体のタイトル # 作図処理を関数として定義 def update(i): # 前フレームのグラフを初期化 plt.cla() # i番目の値を取得 mu = mu_vals[i] # スチューデントのt分布を計算 density = t.pdf(x=x_vals, df=nu, loc=mu, scale=sigma) # スチューデントのt分布を作図 plt.plot(x_vals, density, color='#00A968', linewidth=2.5) # 折れ線グラフ plt.xlabel('x') # x軸ラベル plt.ylabel('density') # y軸ラベル plt.title('$\\nu='+str(nu) + ', \mu='+str(np.round(mu, 2)) + ', \sigma='+str(sigma)+'$', loc='left') # タイトル plt.grid() # グリッド線 plt.ylim(-0.01, dens_max*1.1) # y軸の表示範囲 # gif画像を作成 dens_anime = FuncAnimation(fig, func=update, frames=len(mu_vals), interval=100) # gif画像を保存 dens_anime.save('t_dens_mu.gif')
61
$\sigma$の値を変化させ、$\nu, \mu$を固定して、それぞれ分布を作図します。
・作図コード(クリックで展開)
# スケールパラメータとして利用する値を指定 sigma_vals = np.arange(0.5, 3.1, 0.1) print(len(sigma_vals)) # フレーム数 # 固定するパラメータを指定 nu = 1 mu = 0.0 # xの値を作成 x_vals = np.linspace(start=mu-sigma_vals.max()*2, stop=mu+sigma_vals.max()*2, num=251) # 確率密度(y軸)の最大値を計算 dens_max = t.pdf(x=x_vals, df=nu, loc=mu, scale=sigma_vals.min()).max() # 図を初期化 fig = plt.figure(figsize=(12, 9), facecolor='white') # 図の設定 fig.suptitle("Student's t Distribution", fontsize=20) # 全体のタイトル # 作図処理を関数として定義 def update(i): # 前フレームのグラフを初期化 plt.cla() # i番目の値を取得 sigma = sigma_vals[i] # スチューデントのt分布を計算 density = t.pdf(x=x_vals, df=nu, loc=mu, scale=sigma) # スチューデントのt分布を作図 plt.plot(x_vals, density, color='#00A968', linewidth=2.5) # 折れ線グラフ plt.xlabel('x') # x軸ラベル plt.ylabel('density') # y軸ラベル plt.title('$\\nu='+str(nu) + ', \mu='+str(mu) + ', \sigma='+str(np.round(sigma, 2))+'$', loc='left') # タイトル plt.grid() # グリッド線 plt.ylim(-0.01, dens_max*1.1) # y軸の表示範囲 # gif画像を作成 dens_anime = FuncAnimation(fig, func=update, frames=len(sigma_vals), interval=100) # gif画像を保存 dens_anime.save('t_dens_sigma.gif')
26
$\lambda$の値を変化させ、$\nu, \mu$を固定して、それぞれ分布を作図します。
・作図コード(クリックで展開)
# 逆スケールパラメータとして利用する値を指定 lambda_vals = np.arange(0.1, 3.1, 0.1) print(len(lambda_vals)) # フレーム数 # 固定するパラメータを指定 nu = 1 mu = 0.0 # xの値を作成 sigma = 1.0 / np.sqrt(lambda_vals.min()) x_vals = np.linspace(start=mu-sigma*2, stop=mu+sigma*2, num=251) # 確率密度(y軸)の最大値を計算 dens_max = t.pdf(x=x_vals, df=nu, loc=mu, scale=1.0/np.sqrt(lambda_vals.max())).max() # 図を初期化 fig = plt.figure(figsize=(12, 9), facecolor='white') # 図の設定 fig.suptitle("Student's t Distribution", fontsize=20) # 全体のタイトル # 作図処理を関数として定義 def update(i): # 前フレームのグラフを初期化 plt.cla() # i番目の値を取得 lmd = lambda_vals[i] # スチューデントのt分布を計算 density = t.pdf(x=x_vals, df=nu, loc=mu, scale=1.0/np.sqrt(lmd)) # スチューデントのt分布を作図 plt.plot(x_vals, density, color='#00A968', linewidth=2.5) # 折れ線グラフ plt.xlabel('x') # x軸ラベル plt.ylabel('density') # y軸ラベル plt.title('$\\nu='+str(nu) + ', \mu='+str(mu) + ', \lambda='+str(np.round(lmd, 2))+'$', loc='left') # タイトル plt.grid() # グリッド線 plt.ylim(-0.01, dens_max*1.1) # y軸の表示範囲 # gif画像を作成 dens_anime = FuncAnimation(fig, func=update, frames=len(lambda_vals), interval=100) # gif画像を保存 dens_anime.save('t_dens_lambda.gif')
30
パラメータと統計量の関係をアニメーションで可視化
ここまでは、パラメータと分布の関係を確認しました。次は、パラメータと統計量の関係をアニメーションで確認します。
$\nu$の値を変化させ、$\mu, \sigma$を固定して、それぞれ統計量を計算し、分布を作図します。
・作図コード(クリックで展開)
# 自由度として利用する値を指定 nu_vals = np.arange(1, 31) print(len(nu_vals)) # フレーム数 # 固定するパラメータを指定 mu = 0.0 sigma = 1.0 # xの値を作成 x_vals = np.linspace(start=mu-sigma*5, stop=mu+sigma*5, num=251) # 確率密度(y軸)の最大値を計算 dens_max = t.pdf(x=x_vals, df=nu_vals.max(), loc=mu, scale=sigma).max() # 図を初期化 fig = plt.figure(figsize=(12, 9), facecolor='white') # 図の設定 fig.suptitle("Student's t Distribution", fontsize=20) # 全体のタイトル # 作図処理を関数として定義 def update(i): # 前フレームのグラフを初期化 plt.cla() # i番目の値を取得 nu = nu_vals[i] # スチューデントのt分布を計算 density = t.pdf(x=x_vals, df=nu, loc=mu, scale=sigma) # スチューデントのt分布を作図 plt.plot(x_vals, density, color='#00A968', linewidth=2.5, label='$p(x | \\nu, \mu, \sigma)$') # 分布 if nu > 1: plt.axvline(x=mu, color='blue', linewidth=2.5, linestyle='--', label='$E[x]$') # 期待値 if nu > 2: s_x = np.sqrt(sigma**2 * nu / (nu - 2)) # 標準偏差を計算 plt.axvline(x=mu-s_x, color='orange', linewidth=2.5, linestyle=':', label='$E[x]\pm \sqrt{V[x]}$') # 期待値 - 標準偏差 plt.axvline(x=mu+s_x, color='orange', linewidth=2.5, linestyle=':') # 期待値 + 標準偏差 plt.xlabel('x') # x軸ラベル plt.ylabel('density') # y軸ラベル plt.title('$\\nu='+str(nu) + ', \mu='+str(np.round(mu, 2)) + ', \sigma='+str(sigma)+'$', loc='left') # タイトル plt.legend() # 凡例 plt.grid() # グリッド線 plt.ylim(-0.01, dens_max*1.1) # y軸の表示範囲 # gif画像を作成 dens_anime = FuncAnimation(fig, func=update, frames=len(nu_vals), interval=100) # gif画像を保存 dens_anime.save('t_stat_nu.gif')
30
$\mu$の値を変化させ、$\nu, \sigma$を固定して、それぞれ統計量を計算し、分布を作図します。
・作図コード(クリックで展開)
# 位置パラメータとして利用する値を指定 mu_vals = np.arange(-3.0, 3.1, 0.1) print(len(mu_vals)) # フレーム数 # 固定するパラメータを指定 nu = 5 sigma = 1.0 # xの値を作成 x_vals = np.linspace(start=mu_vals.min()-sigma, stop=mu_vals.max()+sigma, num=251) # 確率密度(y軸)の最大値を計算 dens_max = t.pdf(x=x_vals, df=nu, loc=np.median(mu_vals), scale=sigma).max() # 図を初期化 fig = plt.figure(figsize=(12, 9), facecolor='white') # 図の設定 fig.suptitle("Student's t Distribution", fontsize=20) # 全体のタイトル # 作図処理を関数として定義 def update(i): # 前フレームのグラフを初期化 plt.cla() # i番目の値を取得 mu = mu_vals[i] # スチューデントのt分布を計算 density = t.pdf(x=x_vals, df=nu, loc=mu, scale=sigma) # スチューデントのt分布を作図 plt.plot(x_vals, density, color='#00A968', linewidth=2.5, label='$p(x | \\nu, \mu, \sigma)$') # 分布 if nu > 1: plt.axvline(x=mu, color='blue', linewidth=2.5, linestyle='--', label='$E[x]$') # 期待値 if nu > 2: s_x = np.sqrt(sigma**2 * nu / (nu - 2)) # 標準偏差を計算 plt.axvline(x=mu-s_x, color='orange', linewidth=2.5, linestyle=':', label='$E[x]\pm \sqrt{V[x]}$') # 期待値 - 標準偏差 plt.axvline(x=mu+s_x, color='orange', linewidth=2.5, linestyle=':') # 期待値 + 標準偏差 plt.xlabel('x') # x軸ラベル plt.ylabel('density') # y軸ラベル plt.title('$\\nu='+str(nu) + ', \mu='+str(np.round(mu, 2)) + ', \sigma='+str(sigma)+'$', loc='left') # タイトル plt.legend() # 凡例 plt.grid() # グリッド線 plt.ylim(-0.01, dens_max*1.1) # y軸の表示範囲 # gif画像を作成 dens_anime = FuncAnimation(fig, func=update, frames=len(mu_vals), interval=100) # gif画像を保存 dens_anime.save('t_stat_mu.gif')
61
$\sigma$の値を変化させ、$\nu, \mu$を固定して、それぞれ統計量を計算し、分布を作図します。
・作図コード(クリックで展開)
# スケールパラメータとして利用する値を指定 sigma_vals = np.arange(0.5, 3.1, 0.1) print(len(sigma_vals)) # フレーム数 # 固定するパラメータを指定 nu = 5 mu = 0.0 # xの値を作成 s_x = np.sqrt(sigma_vals.max()**2 * nu / (nu - 2)) x_vals = np.linspace(start=mu-s_x*2, stop=mu+s_x*2, num=251) # 確率密度(y軸)の最大値を計算 dens_max = t.pdf(x=x_vals, df=nu, loc=mu, scale=sigma_vals.min()).max() # 図を初期化 fig = plt.figure(figsize=(12, 9), facecolor='white') # 図の設定 fig.suptitle("Student's t Distribution", fontsize=20) # 全体のタイトル # 作図処理を関数として定義 def update(i): # 前フレームのグラフを初期化 plt.cla() # i番目の値を取得 sigma = sigma_vals[i] # スチューデントのt分布を計算 density = t.pdf(x=x_vals, df=nu, loc=mu, scale=sigma) # スチューデントのt分布を作図 plt.plot(x_vals, density, color='#00A968', linewidth=2.5, label='$p(x | \\nu, \mu, \sigma)$') # 分布 if nu > 1: plt.axvline(x=mu, color='blue', linewidth=2.5, linestyle='--', label='$E[x]$') # 期待値 if nu > 2: s_x = np.sqrt(sigma**2 * nu / (nu - 2)) # 標準偏差を計算 plt.axvline(x=mu-s_x, color='orange', linewidth=2.5, linestyle=':', label='$E[x]\pm \sqrt{V[x]}$') # 期待値 - 標準偏差 plt.axvline(x=mu+s_x, color='orange', linewidth=2.5, linestyle=':') # 期待値 + 標準偏差 plt.xlabel('x') # x軸ラベル plt.ylabel('density') # y軸ラベル plt.title('$\\nu='+str(nu) + ', \mu='+str(np.round(mu, 2)) + ', \sigma='+str(sigma)+'$', loc='left') # タイトル plt.legend() # 凡例 plt.xlabel('x') # x軸ラベル plt.ylabel('density') # y軸ラベル plt.title('$\\nu='+str(nu) + ', \mu='+str(mu) + ', \sigma='+str(np.round(sigma, 2))+'$', loc='left') # タイトル plt.grid() # グリッド線 plt.ylim(-0.01, dens_max*1.1) # y軸の表示範囲 # gif画像を作成 dens_anime = FuncAnimation(fig, func=update, frames=len(sigma_vals), interval=100) # gif画像を保存 dens_anime.save('t_stat_sigma.gif')
26
$\lambda$の値を変化させ、$\nu, \mu$を固定して、それぞれ統計量を計算し、分布を作図します。
・作図コード(クリックで展開)
# 逆スケールパラメータとして利用する値を指定 lambda_vals = np.arange(0.1, 3.1, 0.1) print(len(lambda_vals)) # フレーム数 # 固定するパラメータを指定 nu = 5 mu = 0.0 # xの値を作成 sigma = 1.0 / np.sqrt(lambda_vals.min()) x_vals = np.linspace(start=mu-sigma*2, stop=mu+sigma*2, num=251) # 確率密度(y軸)の最大値を計算 dens_max = t.pdf(x=x_vals, df=nu, loc=mu, scale=1.0/np.sqrt(lambda_vals.max())).max() # 図を初期化 fig = plt.figure(figsize=(12, 9), facecolor='white') # 図の設定 fig.suptitle("Student's t Distribution", fontsize=20) # 全体のタイトル # 作図処理を関数として定義 def update(i): # 前フレームのグラフを初期化 plt.cla() # i番目の値を取得 lmd = lambda_vals[i] # スチューデントのt分布を計算 density = t.pdf(x=x_vals, df=nu, loc=mu, scale=1.0/np.sqrt(lmd)) # スチューデントのt分布を作図 plt.plot(x_vals, density, color='#00A968', linewidth=2.5, label='$p(x | \\nu, \mu, \lambda)$') # 分布 if nu > 1: plt.axvline(x=mu, color='blue', linewidth=2.5, linestyle='--', label='$E[x]$') # 期待値 if nu > 2: s_x = np.sqrt(nu / (nu - 2) / lmd) # 標準偏差を計算 plt.axvline(x=mu-s_x, color='orange', linewidth=2.5, linestyle=':', label='$E[x]\pm \sqrt{V[x]}$') # 期待値 - 標準偏差 plt.axvline(x=mu+s_x, color='orange', linewidth=2.5, linestyle=':') # 期待値 + 標準偏差 plt.xlabel('x') # x軸ラベル plt.ylabel('density') # y軸ラベル plt.title('$\\nu='+str(nu) + ', \mu='+str(mu) + ', \lambda='+str(np.round(lmd, 2))+'$', loc='left') # タイトル plt.legend() # 凡例 plt.grid() # グリッド線 plt.ylim(-0.01, dens_max*1.1) # y軸の表示範囲 # gif画像を作成 dens_anime = FuncAnimation(fig, func=update, frames=len(lambda_vals), interval=100) # gif画像を保存 dens_anime.save('t_stat_lambda.gif')
30
この記事では、1次元スチューデントのt分布のグラフを作成しました。
参考文献
- C.M.ビショップ著,元田 浩・他訳『パターン認識と機械学習 上』,丸善出版,2012年.
おわりに
Pythonの方もやっていかないとなぁ。
【次の内容】
つづく
