モンテカルロシミュレーション

モンテカルロシミュレーションとは

きちんとした定義は色々な文献で散見されるので、ここでは自分なりの理解を述べます。

モンテカルロ法とは、①論理的に数値を直接導きにくい場合に、②確率的な問題に置き換え、③大量の試行回数を用いて近似的な数値を求める手法です。

以下、この定義をもとに例を考えてみます。

①論理的に数値を直接導きにくい場合

例えば、円周率$π$の値がわからないとします。
数学者ではないので、論理的に$π$の値を導くのは厳しいとしましょう。

そこで、確率的に求める方法を考えます。

②確率的な問題に置き換える

面積の関係を利用して、確率の問題に置き換えてみます。
2x2の正方形の中に半径1の正円があるとして、この正方形内にランダムな点を打ち、それが円の中にある確率を$p$とします。
ここで以下のように面積の関係を利用すると：

• 正方形の面積：$4(=2×2)$
• 円の面積：$𝜋⋅𝑟^2=𝜋⋅1^2=π$
• 点が円の中に入る確率：$p=\frac{円の面積}{正方形の面積}=\frac{𝜋}{4}$

$p$についての式を導けます。

この確率$p$を大量の試行によって近似的に求めることで、円周率$π$を計算できます。

2x2の正方形と半径1の正円

③大量の試行回数を用いて近似的な数値を求める

最後にコンピューターによるシミュレーションで$p$の値を求めます。
人の手ではできないので、シュミレーションであることが重要です。

さて、当然プログラムを書くわけですが、書きやすいように条件を次のように読み替えます：

\text{正方形内のランダムな点}\Leftrightarrow 範囲が[-1,1]のランダムな二つの値

円の内側\Leftrightarrow　(原点からの距離) \leq 1

コードは以下のようになります↓
（※プログラムは$p$を求めてから、πを求めるところまで書いてます。）

import numpy as np
num_samples = 1000000  # 試行回数
radius = 1            # 円の半径

# ランダムに点を生成 [-1, 1) の範囲
x = np.random.uniform(-1, 1, num_samples)
y = np.random.uniform(-1, 1, num_samples)

# 丸め込み誤差で1.0を含むようにする
x_fl32 = np.float32(x)
y_fl32 = np.float32(y)

# 円の内側にある点を判定
inside_circle = (x_fl32**2 + y_fl32**2) <= radius**2

# 円周率を近似
probability = np.sum(inside_circle) / num_samples  # 円内の点の割合
approx_pi = probability * 4  # 正方形の面積比からπを計算
print(approx_pi)

$\text{approx_pi}\approx3.14...$に近似ができました。

終わり