本記事に関する注意
素人の落書きです。誤り、用語の不適切、証明なしなどが多々見られると思います。
何卒ご了承ください。
本記事の目的
連続確率分布により生成された標本から、元となる確率分布の平均情報量を求めたい。
理論
確率密度関数$,,f,$をもつ連続確率分布の平均情報量$,h(X),$は式 1 で表されるが、標本数$,N,$が十分大きければ、式 2 のように各標本$,,x_i,$の情報量から求められる（はず）。
\begin{align}
&h(X) = \int_{\chi}^{}f(x)\,log\,f(x)\,dx\qquad\qquad・・・式１\\
&h(X) \approx \frac{1}{N}\sum_{i=1}^{N}\,-logP(\,X = x_i\,)\,\qquad・・・式２
\end{align}
式２から$,h(X)$ を求めるためには、各標本$,,x_i,$に関して$,P(,X = x_i,),$を求める必要がある。
ここから、$P(,X = x_i,),$と$,h(X),$の求め方について少し説明していく。
まず、いくつかの量を定義する。
・ 各標本間の距離を$,d(x_i,,x_j),$とおく$^*$。
・ 各標本$,,x_i,$に関して、$x_i,$との距離$,,d,$が$,r,$以下である標本の総数を$,n_i,$とおく（自身含む）。
・ ある$,,x,$からの距離$,,d,$が$,r,$以下である領域の体積を$,V(r),$とおく。
この定義に従い、適切な$,r,$を決めることで$,P(,X = x_i,),$を式３から近似的に求められる。
P(\,X=x_i\,)\, \approx \frac{n_i}{NV(r)}\qquad・・・式３
式２に式３を代入することで、式４が得られる。
\begin{align}
h(X,r) &\approx \frac{1}{N}\sum_{i=1}^{N}\,-log\frac{n_i}{NV(r)}\\
&＝ \,logV(r) + logN - \frac{1}{N}\sum_{i=1}^{N}\,log\,n_i\qquad・・・式４
\end{align}
式３の近似を成立させるために、$r,$はなるべく小さな値をとる方が好ましい。しかし、標本数が有限である以上、$r,$が極端に小さいと大数の法則を満たせず、式３の近似に綻びが出る。そのため、$r,$の適切な決め方は実際のデータを眺めながら考える必要がある。
$^*\lim_{d(x_i,,x_j) \to 0}P(,X=x_i,)=lim_{d(x_i,,x_j) \to 0}P(,X=x_j,),,$さえ満たせば何でもいいはず
応用
適当な二次元ガウス分布から標本を生成して、$h(X),$を求めてみる。
import numpy as np
from matplotlib import pyplot as plt
def calc_d(x):
N = len(x)
x_tiled = np.tile(x, (N, 1, 1))
d = np.linalg.norm(x_tiled - x_tiled.transpose((1, 0, 2)), axis=2)
return d
# 次元数が 2 であるため円の面積の公式を適用
def calc_v(r):
v = np.pi * np.power(r, 2)
return v
def calc_h(d, v, N, r):
n = np.sum(d <= r, axis=0)
h = np.log(v) + np.log(N) - np.sum(np.log(n)) / N
return h
# 適当な 2 次元ガウス分布からデータを生成する
data = np.random.normal(0, 1, (1000, 2))
# r を変化させながら h(X) を計算する
r_list = [(i + 1) * 0.01 for i in range(10000)] # r の範囲は適当に決めた
d = calc_d(data)
N = len(data)
h_list = [calc_h(d, calc_v(r), N, r) for r in r_list]
# グラフを描画する
# 計算値を青い実線でプロットする
plt.figure(0)
plt.plot(r_list, h_list, color='blue', linestyle='solid')
# 標本分散から計算した値を青い点線でプロットする
Z = np.cov(data[:, 0], data[:, 1])
h_s = 0.5 * np.log(np.linalg.det(2 * np.pi * np.e * Z))
plt.plot(r_list, [h_s for _ in range(len(r_list))], color='blue', linestyle='dotted')
# 母集団分散から計算した値をオレンジ色の点線でプロットする
h_u = np.log(2 * np.pi * np.e)
plt.plot(r_list, [h_u for _ in range(len(r_list))], color='orange', linestyle='dotted')
plt.xlim([0, 3])
plt.ylim([0, 5])
plt.show()
実行すると、このようなグラフが得られる。
横軸は$,r,$、縦軸は$,h(X,r),$を表す。
理論で説明したように、$r,$が小さいほど$,h(X,r),$は真の値に近づくが、小さすぎると今度は負の無限大に向かって発散していく。
グラフを見ると、傾きが最も小さくなるように$,r,$を決めるのがなんとなく良さそうに見える。実際に、式３の近似が成り立つ場合$,\frac{\partial}{\partial r}h(X,r)=0,$も成り立つことを考慮すると、この決め方もそこまでおかしくは無いと思うが、何か証明があるわけでも無いので過信は禁物。