はじめに
機械学習で登場する確率分布について色々な角度から理解したいシリーズです。
ガンマ関数のグラフをPythonで作成します。
【他の記事一覧】
【この記事の内容】
ガンマ関数の作図
ガンマ関数(Gamma Function)の計算を確認して、グラフを作成します。ガンマ関数については「ガンマ関数の性質の導出 - からっぽのしょこ」を参照してください。
利用するライブラリを読み込みます。
# 利用するライブラリ import numpy as np from scipy.special import gamma, gammaln # ガンマ関数, 対数ガンマ関数 import matplotlib.pyplot as plt
定義式の確認
まずは、ガンマ関数の定義式を確認します。
ガンマ関数は、次の式で定義されます。
ここで、$e$はネイピア数です。
ガンマ関数の計算
続いて、ガンマ関数の計算方法を確認します。
SciPyライブラリのspecialモジュールのガンマ関数gamma()を使って計算します。
# 変数の値を指定 x = 4.0 # ガンマ関数の計算 y = gamma(x) print(y)
6.0
対数ガンマ関数gammaln()を使って計算します。
# 対数ガンマ関数による計算 log_y = gammaln(x) z = np.exp(log_y) print(y)
6.0
計算結果の指数をとるとガンマ関数の計算結果が得られます。
ガンマ関数の計算は非常に大きな値になるため、引数(変数)$x$によっては発散してしまいます。
# ガンマ関数を計算:(発散) print(gamma(171.0)) print(gamma(172.0))
7.257415615308e+306
inf
途中計算などでは、対数ガンマ関数を使うことで計算できます。
# 対数ガンマ関数による計算 print(gammaln(172.0))
711.71472580229
$x$が0または負の整数の場合は計算できません。
# ガンマ関数を計算:(計算できない) print(gamma(0.0)) print(gamma(-3.0))
inf
inf
ガンマ関数の計算を確認できました。次は、ガンマ関数を可視化します。
グラフの作成
MatplotlibライブラリのPyPlotモジュールを利用して、ガンマ関数のグラフを作成します。
非負の値の場合
まずは、0より大きい実数$x$を範囲として、ガンマ関数を計算します。
# xの値を指定:(x > 0) x_vals = np.arange(start=0.01, stop = 6.0, step=0.01) print(x_vals[:5]) # ガンマ関数の計算 y_vals = gamma(x_vals) print(y_vals[:5])
[0.01 0.02 0.03 0.04 0.05]
[99.43258512 49.44221016 32.78499835 24.46095502 19.47008531]
x軸の値を作成してx_valsとします。
作成したx_valsを用いて、gamma()でガンマ関数の計算をします。
ガンマ関数のグラフを作成します。
# ガンマ関数を作図 plt.figure(figsize=(12, 9)) # 図の設定 plt.plot(x_vals, y_vals, color='orange') # 折れ線グラフ plt.xlabel('x') # x軸ラベル plt.ylabel('$\Gamma(x)$') # y軸ラベル plt.suptitle('Gamma Function', fontsize=20) # 全体のタイトル plt.grid() # グリッド線 plt.show() # 描画
$x = 0$の垂線に漸近しているのが分かります。
負の値を含む場合
続いて、負の実数も含めて、ガンマ関数を計算します。
# xの値を指定:(0と負の整数は計算できない) x_vals = np.arange(start=-5.0, stop=5.0, step=0.0005) print(x_vals[:5]) # ガンマ関数の計算 y_vals = gamma(x_vals) print(y_vals[:5])
[-5. -4.9995 -4.999 -4.9985 -4.998 ]
[ inf -16.68089686 -8.34757609 -5.56981089 -4.18093459]
0と負の整数のとき、計算結果がinfになります。
グラフに表示しない値をマスクします。(マスクしておかないと不連続な点を繋いだグラフになります。)
# 閾値を指定 threshold = 10.0 # 閾値外の値をマスク y_mask_vals = np.ma.masked_where((y_vals < -threshold) | (y_vals > threshold), y_vals) print(y_mask_vals[:5])
[-- -- -8.347576090322463 -5.569810889197673 -4.180934591519495]
閾値thresholdを指定して、その範囲外の値の要素をnp.ma.masked_where()でマスクします。第1引数にマスクする条件、第2引数にマスクしない場合に返す値を指定します。この例では、$y < - \mathrm{threshold}$、$\mathrm{threshold} < y$の値をマスクします。複数の条件A, Bの一方でも該当する要素をマスクする場合は、(A) | (B)と書きます。
ガンマ関数のグラフを作成します。
# ガンマ関数を作図 plt.figure(figsize=(12, 9)) # 図の設定 plt.plot(x_vals, y_mask_vals, color='orange') # 折れ線グラフ plt.xlabel('x') # x軸ラベル plt.ylabel('$\Gamma(x)$') # y軸ラベル plt.suptitle('Gamma Function', fontsize=20) # 全体のタイトル plt.xticks(ticks=np.arange(start=np.floor(x_vals.min()), stop=np.ceil(x_vals.max())+1.0)) # x軸目盛 plt.ylim(ymin=-threshold, ymax=threshold) # y軸の表示範囲 plt.grid() # グリッド線 plt.show() # 描画
x軸の値が負の整数と0の点において、不連続なグラフになっているのが分かります。(マスクせずに線が繋がらないようにする引数の設定とかないですか？)
参考文献
- 岩田具治『トピックモデル』(機械学習プロフェッショナルシリーズ)講談社,2015年.
おわりに
もう一段階深く理解するには複素数の知識が要るっぽい？
【次の内容】
