SymPyの使い方18 ～統計モジュール2

　前回に続きSymPyの統計計算ツールを利用できるStatsモジュールについて紹介する。

atatat.hatenablog.com

1. 乱数の生成
2. 任意の確率分布から乱数生成
3. まとめ

1. 乱数の生成

　前回はシンボリックな乱数の定義ついて紹介したが、シンボリックな乱数から乱数を生成することもできる。サイコロの目の乱数Dieをインポートし、シンボリックな乱数Xを以下のように生成する。

from sympy.stats import Die, sample, Normal
X=Die('X',6)

シンボリックな乱数Xから乱数(1～6)を生成するにはsampleを用いる。

sample(X)
>3
sample(X)
>1
sample(X)
>3

sampleにより乱数を生成するごとに値が変わっていることがわかる。

正規分布の場合も、同様でsampleにより乱数生成ごとに値が変わる。

Z=Normal('Z',0,1)
sample(Z)
>-0.717809169119500
sample(Z)
>-0.655047199602641
sample(Z)
>0.902397063242097

2. 任意の確率分布から乱数生成

　正規分布のようにあらかじめ定義されている関数以外に任意の確率分布に従ったシンボリックな乱数を生成することができる。例えば確率密度関数exp(-x)で区間[ 0, ∞ ] で定義される確率分布を生成するにはContinuousRV(変数,確率密度関数, 区間)を用いて以下のようにする。

from sympy import exp, Symbol, Interval, oo
from sympy.stats import ContinuousRV, E, P
x=Symbol('x')
pdf=exp(-x)
Z = ContinuousRV(x, pdf, set=Interval(0, oo))

シンボリックな乱数は、DieやNormalのように期待値や確率を以下のように求めることができる。

E(Z)
>1
P(Z>5)
>exp(-5)

引数が自然数に限定するような離散的な場合は、DiscreteRVを用いて任意の確率分布を、ContinuousRVと同様に生成できる。なおDiscreteRVは最新版1.7.1でないとロードできない。

from sympy import S
from sympy.stats import DiscreteRV, P, E
p=S(1)/2
x = Symbol('x', integer=True, positive=True)
pdf = p*(1 - p)**(x - 1)
D = DiscreteRV(x, pdf, set=S.Naturals)

以下のように期待値、確率が計算できる。

E(D)
>2
P(D>3)
>1/8

引数、確率がいずれも離散的な場合はFiniteRVを用いる。以下確率密度関数をDictionaryで定義して、確率分布に従ったシンボリックな乱数を生成する。

from sympy.stats import FiniteRV
from sympy import Rational
pmf = {1: Rational(1, 3), 2: Rational(1, 6), 3: Rational(1, 4), 4: Rational(1, 4)}
X = FiniteRV('X', pmf)
E(X)
>29/12
P(X>3)
>1/4