Python 3 で数学を。

Python 3 とライブラリで数学の問題を解いていきます。統計学や機械学習はときどき。

Python (Python 3) で数学をやるにはどうしたらいいか。その21. 確率。"俺だってマーク・トウェインの相続人さ!"

このシリーズの過去記事は以下にまとめてある

py3math.hatenablog.com

当記事について

当ブログ筆者 (以下、筆者) も、Python で数学をやるにはどうしたらいいか悩んでいた時期があるから、昔の自分に向けて書いたような記事。

対象とする読者

Python の入門書を一冊か二冊、一通りやった人で、Python の基本的な構文や基本的な用語がだいたいわかっている人。

そして、公式サイトを読む努力をちゃんとする人。

確率をやってみよう

六面体のサイコロを例にして確率の基本的なことをやってみよう。

サイコロの目は、1 2 3 4 5 6 だ。

これは言うまでもないことだろう。

では、この期待値は何だろうか。

そもそも期待値とは何だろうか。

先に結論を書けば、期待値とは平均のことだ。

これを求めてみよう。それぞれの目が出る確率は、1/6 だ。

Python の対話モードから:

>>> 1*(1/6) + 2*(1/6) + 3*(1/6) + 4*(1/6) + 5*(1/6) + 6*(1/6)
3.5

サイコロの期待値 (平均) は、3.5 だということがわかったが、もう少し簡潔に書こう。

>>> (1 + 2 + 3 + 4 + 5 + 6)*(1/6)
3.5

上でもよいが、リスト内包表記で書こう。

>>> res = sum([i*(1/6) for i in range(1, 7)])
>>> res
3.5

SymPy でサイコロの確率分布と期待値 (平均) を求めてみよう。

以下のようにインポートする。

from sympy.stats import Die, density, E

まず六面体のサイコロを用意しよう。

>>> D6 = Die('D6', 6)

そして、確率分布を求める。

>>> density(D6).dict
{1: 1/6, 2: 1/6, 3: 1/6, 4: 1/6, 5: 1/6, 6: 1/6}

結果は各々 1/6 だ。間違いない。

次に、期待値 (平均) を求めてみよう。

>>> E(D6)
7/2

結果は、7/2 だ。これも間違いない。

だが、小数にしたい人もいるかもしれない。

単純にこうすればよい。

>>> 7/2
3.5

最後に大数の法則を確認しておこう (以下では、とくに (擬似) 乱数のタネは使用しない)。

NumPy を使用する (マシンに NumPy がすでに入っているものとして話を進める)。

まず NumPy をインポートしよう。

>>> import numpy as np

以下は、1 から 6 までの整数をランダムに100 万個生成し、その平均を求めている。

>>> x1 = [np.random.randint(1, 7) for _ in np.arange(1000000)]
>>> np.mean(x1)
3.497124

上で書いたコードについて少し細かく書いておこう。

np.random.randint() で 1 から 6 までの整数を指定。

np.arange() で 100 万と指定。

np.mean() で平均を求めている。

それから、for _ としているのは、例えば ループカウンタ変数 i 等は使用しないからこうしている。

さて、結果の 3.4971243.5 はかなり近いといえる。

だが、1 回だけではちょっと不安かもしれない。

上で書いたリスト内包表記を 10 回、for 文で回してみよう (なお、マシンスペックによっては多少時間が掛かるかもしれない)。

>>> for _ in range(10):
...     x2 = [np.random.randint(1, 7) for _ in np.arange(1000000)]
...     np.mean(x2)
...
3.498765
3.50086
3.501293
3.500276
3.499952
3.502131
3.503294
3.497259
3.500632
3.501244

上の結果から、大数の法則を確認することができただろう。

(つづく)。

参考文献 (数式を参考)

高校数学でわかる統計学―本格的に理解するために (ブルーバックス)

高校数学でわかる統計学―本格的に理解するために (ブルーバックス)