Python で
結果を
指数分布とは
そもそも、
個人的には、指数分布 Web 解析
なので
書籍
統計学入門 - 東京大学出版会 には
指数分布は
exponetial distribution
は、(5.10) で 述べたように 確率密度関数
$$ f(x) = \lambda e^{-\lambda x} (x >= 0), 0 (x < 0) $$
で定義される 連続型の 分布である。 幾何分布に
おいて $ q = e^{-\lambda} $ すなわち $ \lambda = -log q $ と おくと、 形式上この 確率分布の 主要な 部分 $ e^{-\lambda x} $ が 得られるから、 この 確率分布は 連続的な 待ち時間分布の 性質を もつことが わかる。
自然科学の
指数分布は、
ポアソン分布に 従って 起きる 事象の 生起感覚を 表現する 分布である。 指数分布は、 さまざまな サービスの 待ち時間を 表現する 分布と して 用いられる。
実世界での 例
実世界で、
- サービスの
待ち行列の 分布 - システムの
寿命 - 災害が
発生するまでの 年数 - 雨滴の
大きさ - Wikipedia への
ページリクエスト間の 時間
Web
指数分布で
* 指数分布 - Wikipedia
指数分布の 計算方法
Python では、
numpy.random.exponential を 使う
パラメータmu(分布の
numpy.random.exponential
は
関数に
引数 scale のみ
第一引数scale
のみ
スカラー値とは
スカラー (数学) - Wikipedia
import numpy as np import pylab as plt target = 500 beta = 1.0/target Y = np.random.exponential(beta) print(beta, Y)
0.002 0.00026848474693625724
引数 scale には リストも 指定可能
引数 scale には、
戻り
import numpy as np import pylab as plt target = 500 beta = [1.0/target, 2.0/target, 3.0/target] Y = np.random.exponential(beta) print(Y) beta = (1.0/target, 2.0/target, 3.0/target) Y = np.random.exponential(beta) print(Y)
[ 0.00035789 0.00079873 0.0059433 ] [ 0.00070439 0.00116696 0.01094726]
size を 指定する
第二引数size
を
scale
は
これが
import numpy as np import pylab as plt # target が ラムダ ある期間に平均して10回起こる target = 10 beta = 1.0/target # 要素数10000 ある期間に平均して10回 指数分布 Y = np.random.exponential(beta, 10000) plt.hist(Y, normed=True, bins=200,lw=0,alpha=.8) # Y の最大値をプロットしておく plt.plot([0,max(Y)],[target,target],'r--') plt.ylim(0,target*1.1) plt.show()
### size に3要素のタプルを指定する size には 3要素の タブル を指定することができます。 戻り値が3次元配列で返ります。
import numpy as np import pylab as plt # target が ラムダ ある期間に平均して10回起こる target = 10 beta = 1.0/target # (m=3, n=5, k=2) を指定する Y = np.random.exponential(beta, (3, 5, 2)) print(Y)
[[[ 0.03303321 0.31016278] [ 0.2045412 0.24228506] [ 0.12344505 0.15427801] [ 0.21190096 0.02285032] [ 0.00356998 0.08352543]] [[ 0.05553856 0.06446413] [ 0.11237701 0.12704291] [ 0.00912467 0.1277612 ] [ 0.01687101 0.19313137] [ 0.0858672 0.08564713]] [[ 0.0119443 0.0315812 ] [ 0.00289961 0.07208328] [ 0.04683225 0.04912672] [ 0.08552624 0.04571101] [ 0.01842625 0.0141593 ]]]
scipy.stats.expon を 使う
続いてscipy.stats.expon
を
モーメントを 求める
expon.stats(moments='mvsk')
で
これは
歪度、
from scipy.stats import expon mean, var, skew, kurt = expon.stats(moments='mvsk') #平均1、分散1、歪度 2、尖度 6 print(mean, var, skew, kurt)
1.0 1.0 2.0 6.0
指数分布の グラフを 描く
scipy.stats.expon — SciPy v1.1.0 Reference Guide を
from scipy.stats import expon import matplotlib.pyplot as plt fig, ax = plt.subplots(1, 1) # 指定した区間をN等分した数列を生成する # N 100 を増減させると、曲線の滑らかさが変わる # ppf は percent point function の略 x = np.linspace(expon.ppf(0.01),expon.ppf(0.99), 100) # expon は、Exponential distribution の expon だと思う ax.plot(x, expon.pdf(x),'r-', lw=5, alpha=0.6, label='expon pdf') # frozen > 固定した rv = expon() ax.plot(x, rv.pdf(x), 'k-', lw=2, label='frozen pdf') # ランダムな要素数1000 の指数分布に従う配列を取得 # loc (指数分布のx軸の位置)、scale (1 / lambda) の rvs に指定可能。 というか他の関数にも指定できる r = expon.rvs(size=1000) # ヒストグラム描画 ax.hist(r, density=True, histtype='stepfilled', alpha=0.2) ax.legend(loc='best', frameon=False) plt.show()
参考
以下参考に
- Manipulating the numpy.random.exponential distribution in Python - Stack Overflow
- scipy.stats.expon — SciPy v1.1.0 Reference Guide
正直、
以上です。
コメント