対数正規分布のシミュレーション
中心極限定理に基づいて正規分布は実証的分析における確率分布としての妥当性が説明されることがある。中心極限定理では、要因は加法的に働くが、これを乗法的な形式に変換すると対数正規分布の妥当性が示されることになる。
いま、同一の分布に従う独立な確率変数を、,・・・とおく。中心極限定理は次のことを主張するものである。
すなわち
は、平均
、分散
の正規分布に近づく。
指数関数で変換すると、
は、平均
、分散
の正規分布に対応する対数正規分布に近づく。
となる。
いま、
と置けば、
は、対数正規分布に近づく
となる。
積が対数正規分布に近づくことを、が連続確率変数の場合と、離散確率変数の場合についてシミュレーションで確認する。
連続確率変数の場合
連続確率変数として
の場合について考える。パラメータと
は、
の平均値と分散が、
と
となるものである。
積の個数はとした。シミュレーションのスクリプトは、リスト1に示すものである。
リスト1 確率変数が連続型の場合
import numpy as np
import matplotlib.pyplot as plt
import seaborn as sb
rng = np.random.default_rng()
def lognormal_pdf(x, mu, sgm):
v1 = sgm * x *
((2*np.pi)**0.5)
v2 = np.exp(-((np.log(x) -
mu)**2) / (2 * (sgm**2)))
return v2 / v1
def param_to_mu_sgm(mean, var):
sgm_2 = np.log(var *
(mean**(-2)) + 1.0)
mu = np.log(mean) - sgm_2 /
2.0
return mu, sgm_2**0.5
mean = 1.2
var = 0.1**2
print(mean, var)
mu, sgm = param_to_mu_sgm(mean, var)
print(mu, sgm)
n_simu = 100_000
xs = []
N = 100
for i_simu in range(n_simu):
if i_simu % 10000 == 0:
print(i_simu)
v = 1.0
for i in range(N):
x =
rng.normal(mu, sgm)
v *=
np.exp(x)
xs.append(v)
plt.figure(figsize=(14,5))
plt.hist(xs, density=True,
bins=np.linspace(0, np.max(xs), 1000))
xcoord = np.linspace(0.001,
np.max(xs), 1000)
ycoord = lognormal_pdf(xcoord, N*mu,
sgm * (N**0.5))
plt.plot(xcoord,ycoord, c='r')
plt.show()
リスト1を実行すると、図1の結果を得る。
図1
100000回のシミュレーションの結果のヒストグラムと対応する対数正規分布のグラフである。積による効果が対数正規分布でよく表されていることがわかる。
離散型確率変数の場合
各確率変数の値を
と設定した。
確率変数の個数はとして、リスト2に示すシミュレーションのスクリプトを作成した。
リスト2 離散型確率変数の場合
import numpy as np
import scipy.stats as ss
import matplotlib.pyplot as plt
import seaborn as sb
rng = np.random.default_rng()
def lognormal_pdf(x, mu, sgm):
v1 = sgm * x *
((2*np.pi)**0.5)
v2 = np.exp(-((np.log(x) -
mu)**2) / (2 * (sgm**2)))
return v2 / v1
n_simu = 100_000
xs = []
N = 100
for i_simu in range(n_simu):
if i_simu % 10000 == 0:
print(i_simu)
v = 1.0
ck = rng.binomial(1, 0.25,
size=N)
for i in range(N):
v *=
1.0*(1-ck[i]) + 1.1*ck[i] # Yi = 1.0 or 1.1
xs.append(v)
plt.figure(figsize=(14,5))
plt.hist(xs, density=True)
sb.kdeplot(xs, label='kdeplot')
mu = 0.25 * np.log(1.1)
sgm_2 = 0.75 * (((-0.25 *
np.log(1.1))**2)) + \
0.25 * ((np.log(1.1) - 0.25 * np.log(1.1)) ** 2)
xcoord = np.linspace(0.001,
np.max(xs), 1000)
ycoord = lognormal_pdf(xcoord, N*mu,
(N*sgm_2)**0.5)
plt.plot(xcoord, ycoord,
label='lognormal pdf')
plt.legend(fontsize=16)
plt.show()
リスト2を実行すると、図2の結果を得た。
図2
100,000回のシミュレーションの結果のヒストグラムと、KDE密度関数(赤色の曲線)、および対数正規分布密度関数(緑色の曲線)が描かれている。ヒストグラムおよびKDE密度関数は対数正規分布の密度関数によく合っているといえる。
参考文献
Johnson, N. L., Kotz, S., & Balakrishnan, N. (1994) Continuous Univariate Distributions, volume 1, second edition. John Wiley & Sons, Inc.
國仲寛人・小林奈央樹・松下貢(2011)複雑系にひそむ規則性――対数正規分布を軸にして――.日本物理学会誌、66、658−665.