ゼロからわかるモンテカルロ法：第3回 Pythonで乱数を扱う基本操作

イントロダクション

前回は乱数と確率分布の理論を学びました。

今回は、Pythonを使って乱数を自在に扱う方法を詳しく解説します。

モンテカルロ法を実装するためには、乱数を生成し、分布に従ったサンプルを作ることが不可欠です。

ここでは、標準ライブラリとNumPyを使った基本操作をマスターし、次回以降のモンテカルロ積分やMCMCの準備を整えます。

Pythonで乱数を扱うためのライブラリ

Pythonには乱数生成のための標準ライブラリ random と、科学計算向けの numpy.random があります。

random はシンプルで使いやすく、numpy は大量の乱数を効率的に生成できます。

モンテカルロ法では大量のサンプルが必要になるため、NumPyの利用が基本になります。

一様乱数の生成

まずは一様乱数です。区間 [0,1] の一様乱数を生成するには random.random() を使います。

import random
# 0から1の一様乱数を10個生成
uniform_randoms = [random.random() for _ in range(10)]
print("一様乱数:", uniform_randoms)

区間を指定したい場合は random.uniform(a, b) を使います。

# 区間 [5,10] の一様乱数を生成
uniform_randoms_range = [random.uniform(5, 10) for _ in range(10)]
print(uniform_randoms_range)

numpyを使うとさらに簡単です。

import numpy as np
# 区間 [0,1] の一様乱数を100000個生成
samples = np.random.uniform(0, 1, 100000)
print("平均:", np.mean(samples))

大量の乱数を生成して統計量を確認すると、平均は0.5付近、標準偏差は約0.288になるはずです。

正規分布の乱数

モンテカルロ法で頻繁に使われるのは正規分布です。

Pythonでは random.gauss(mu, sigma) や numpy.random.normal(mu, sigma, size) を使います。

# randomモジュールで正規分布乱数
normal_randoms = [random.gauss(0, 1) for _ in range(10)]
print("正規分布乱数:", normal_randoms)
# NumPyで大量の正規分布乱数
samples_normal = np.random.normal(0, 1, 100000)
print("平均:", np.mean(samples_normal))
print("標準偏差:", np.std(samples_normal))

ヒストグラムを描くと、中央に山があるベル型の分布が確認できます。

サイコロや離散分布の乱数

離散分布の例として、サイコロの目をシミュレーションします。

# サイコロの目（1〜6）を100回シミュレーション
dice_rolls = [random.randint(1, 6) for _ in range(100)]
print("サイコロの目:", dice_rolls[:20])

NumPyでは np.random.randint(low, high, size) を使います。

import numpy as np
# high は排他なので 7 にする（1〜6）
np.random.seed(42)
dice = np.random.randint(1, 7, 100000)
print("平均:", dice.mean())  # 3.5 に近い値が得られる

平均は約3.5になるはずです。

乱数の再現性とシード

モンテカルロ法では、結果を再現するために乱数のシードを固定します。

random.seed(42)
print(random.random()) # 毎回同じ値
np.random.seed(42)
print(np.random.rand()) # 毎回同じ値

シードを設定することで、同じ乱数列を再現できます。

これは実験や検証などをする際には不可欠です。

実践：ヒストグラムで分布を確認

乱数の分布を視覚化することで理解が深まります。NumPyとMatplotlibを使ってヒストグラムを描いてみましょう。

import matplotlib.pyplot as plt
# 正規分布の乱数を生成
samples = np.random.normal(0, 1, 100000)
# ヒストグラムを描画
plt.hist(samples, bins=50, density=True, alpha=0.6, color='g')
plt.title("正規分布のヒストグラム")
plt.show()

このグラフを見ると、理論的な正規分布に近い形が確認できます。

まとめ

今回は、Pythonで乱数を扱う基本操作を学びました。

一様乱数、正規分布乱数、離散分布乱数、そしてシードによる再現性の確保まで理解できたと思います。

次回は、これらの乱数を使ってモンテカルロ積分を実際に実装し、理論値との比較を行います。