CC BY-NC 4.0 - 非営利目的のみ許可

モンテカルロ法

概要

モンテカルロ法は、乱数を用いてシミュレーションを行う数値計算手法です。

• 計算が複雑だったり解析的に解が得られない問題に強みがある
• 確率的な現象を扱う問題に適している
• サンプリング数を増やすことで精度を向上させることができる
• 1940年代にスタニスワフ・ウラムとジョン・フォン・ノイマンによって開発された

適用例

• 数値積分（高次元積分の計算）
• 物理シミュレーション（粒子の挙動、熱伝導など）
• 金融工学（オプション価格の評価、リスク分析）
• 統計推定（ベイズ推定、マルコフ連鎖モンテカルロ法）

モンテカルロ法による円周率の求め方

アルゴリズム

1. 設定: 一辺の長さが\(2\)の正方形の中に、半径\(1\)の円を内接させる
2. ランダムサンプリング: 正方形内にランダムに点を打ち、その点が円の内部か外部かを調べる
3. 比率計算: 円の内部に入った点の数と全点の数の比率から円周率を推定する
4. 精度向上: 試行回数（打つ点の数）を増やすほど精度が高くなる

数学的原理

• 点が均一にばらまかれていれば、円の面積と正方形の面積の比率に比例する
• 正方形の面積: \(4\)（一辺が\(2\)なので \(2 \times 2 = 4\)）
• 円の面積: \(\pi \times 1^2 = \pi\)
• 円の内部に入った点の比率: \(\frac{\text{円内の点数}}{\text{全点数}} \approx \frac{\pi}{4}\)
• したがって、\(\pi \approx 4 \times \frac{\text{円内の点数}}{\text{全点数}}\)

精度の評価

試行回数を\(N\)とすると、推定誤差は約\(\frac{1}{\sqrt{N}}\)に比例します。

• 試行回数を100倍にすると、誤差は約10分の1になる
• より正確な結果を得るには、大量のサンプリングが必要

プログラム実装での判定方法

ある点\((x, y)\)が円の内部にあるかどうかは、原点からの距離で判定します。

• 座標系: 原点を正方形の中心に設定し、\(-1 \leq x \leq 1, -1 \leq y \leq 1\)の範囲でランダムに点を生成
• 原点からの距離: \(\sqrt{x^2 + y^2}\)
• 円の内部判定: \(\sqrt{x^2 + y^2} \leq 1\)
• 計算効率化: 通常は二乗した値で判定: \(x^2 + y^2 \leq 1\)（平方根の計算を回避）

乱数生成方法

1. 物理的な乱数生成

サイコロ

• 20面体など特殊なサイコロも利用可能
• 10面体では製造が難しいため、0〜9が複数記載された20面体を用いることもある
• 確率が均一になるよう設計されている
• 主に教育目的や小規模な実験で使用される

乱数表

• 0〜9までの数字が不規則に並んでおり、出現確率が等しくなるよう工夫されている
• 使い方: 任意の場所からスタートし、順に乱数を拾う（縦、横、斜めなど、どの方向でも可）
• 現在は主に疑似乱数が使われるため、一般的には利用されていない

物理的過程

• 原子核の崩壊、ダイオードの電気的ノイズ等の確率的現象を利用する
• 真の乱数（真性乱数）を生成できるが、生成速度が遅い
• 暗号学的用途など、高度なセキュリティが必要な場合に使用される
• ハードウェア乱数生成器（HRNG）として実装されることが多い

2. コンピュータによる疑似乱数生成

コンピュータによる疑似乱数生成には、平方根中法、線形合同法、Mersenne Twisterなど様々なアルゴリズムが存在します。
疑似乱数は完全な乱数ではないが、統計的性質が良く、高速に生成できるため、モンテカルロ法で広く利用されています。

疑似乱数の評価基準

• 周期性: 乱数列が繰り返し始めるまでの長さ
• 統計的独立性: 前後の値の間に相関がないこと
• 均一性: 各値が等しい確率で出現すること
• 計算速度: 生成にかかる時間

平方根中法（中位平方根法）

アルゴリズム

適当な\(n\)桁の数字を2乗し、中位の\(n\)桁を取り出すことで乱数を生成する方法です。
1949年にジョン・フォン・ノイマンによって提案された、一様分布を生成する古典的な方法です。

手順例（n=4の場合）

1. 初期値として4桁の数値を選ぶ（例: 1234）
2. 2乗する: \(1234^2 = 1522756\)（7桁）
3. 中位4桁を取り出す:
   • 全体が7桁の場合、中央から4桁を取る（2275）
   • 具体的には、両端の桁を除いた中央部分を取る
   • 桁数が足りない場合は0を前に補う
4. 次の乱数として使用し、同様の操作を繰り返す

計算例の続き

• 初期値: 1234
• 1回目: \(1234^2 = 1522756\)（7桁）→ 中位4桁: 2275
• 2回目: \(2275^2 = 5175625\)（7桁）→ 中位4桁: 1756
• 3回目: \(1756^2 = 3083536\)（7桁）→ 中位4桁: 0835 → 835（先頭の0は通常省略、または835として扱う）
• （以下同様に続く）

注意: 平方根中法では、実装によって「中位」の取り方が若干異なる場合があります。また、0が出現すると以後すべて0になる問題があります。

長所

• 簡単で分かりやすい
• 計算が速い
• 実装が容易

短所

• 0が出ると、それ以降の乱数がすべて0になってしまう（収束問題）
• 周期が解析的にわかっていない
• 統計的性質が悪く、相関が生じやすい
• 現在ではほとんど使用されていない（歴史的な手法）

線形合同法（Linear Congruential Generator, LCG）

アルゴリズム

線形合同法は、1951年にD.H. Lehmerによって提案された、以下の漸化式で疑似乱数を生成する方法です：

\[x_{n+1} = (a \times x_n + b) \bmod M\]

ここで：

• \(a\): 乗数（multiplier）
• \(b\): 加数（increment）
• \(M\): 法（modulus）
• \(a, b, M\)は整数でなければならない
• 初期値\(x_0\)も整数

具体例

例1: 小規模な例（教育用）

• \(a = 5, b = 3, M = 16, x_0 = 1\)
• \(x_1 = (5 \times 1 + 3) \bmod 16 = 8\)
• \(x_2 = (5 \times 8 + 3) \bmod 16 = 43 \bmod 16 = 11\)
• \(x_3 = (5 \times 11 + 3) \bmod 16 = 58 \bmod 16 = 10\)
• （以下同様に続く）

例2: 実用的なパラメータ（Glibcの古いバージョンで使用）

• \(a = 1103515245, b = 12345, M = 2^{31}\)
• 周期が長く、統計的性質が比較的良好
• 注意: 現在のGlibcでは、より高度な乱数生成器を使用している

特徴

• 乱数の周期を最大にするための値の選び方が理論的に研究されている
• 計算が単純で高速（加算、乗算、剰余演算のみ）
• 広く使用されている基本的な疑似乱数生成法
• メモリ使用量が少ない（前の値1つだけ保持）

欠点

• 係数の与え方によって規則的だったりパターン化することがある
• 初期値が極端に小さいと、生成される値が偏る可能性がある
• 高次元での相関が問題になることがある（超平面に点が集中する）
• 下位ビットの周期性が短い（特に\(M\)が2のべき乗の場合）

パラメータの選び方（Hull-Dobellの定理）

最長周期\(M\)を得るための条件：

1. \(b\)と\(M\)は互いに素（最大公約数が1）
2. \(a-1\)は\(M\)のすべての素因数で割り切れる
3. \(M\)が4の倍数の場合、\(a-1\)も4の倍数

推奨されるパラメータの選び方:

• Hull-Dobellの定理の条件を満たす範囲で\(a\)を選ぶ
• 初期値は適切な値を選ぶ（極端に小さすぎると、最初のいくつかの値が小さな値に偏る可能性がある）
• \(M\)は通常、2のべき乗または素数を用いる
• \(b\)と\(M\)は互いに素であることが必要（Hull-Dobellの定理の条件1）
• 実用的には、既知の良好なパラメータセットを使用することが推奨される

線形合同法の変種

• 混合合同法: \(b \neq 0\)の場合（上記の標準的な形式）
• 乗合同法: \(b = 0\)の場合（周期が短くなるが計算が速い）

その他の疑似乱数生成法

Mersenne Twister

• 1997年に開発された現代的な疑似乱数生成法
• 周期が非常に長い（\(2^{19937}-1\)）
• 統計的性質が優れている
• 多くのプログラミング言語で標準的に使用されている（Python、Rなど）
• メモリ使用量が比較的多い（約2.5KB）

その他の手法

• Xorshift: ビット演算のみで高速な乱数生成
• WELL (Well Equidistributed Long-period Linear): Mersenne Twisterの改良版
• PCG (Permuted Congruential Generator): 線形合同法を改良した手法

まとめ

モンテカルロ法は、乱数を利用した強力な数値計算手法です。
円周率の計算を例として、乱数の生成方法から実際の応用まで理解することが重要です。
適切な乱数生成法を選択し、十分な試行回数を確保することで、高い精度の結果を得ることができます。

実装時の注意点

1. 乱数生成器の選択: 用途に応じて適切な手法を選ぶ
2. 初期化: シード値（初期値）を適切に設定する
3. 周期の考慮: 十分に長い周期を持つ生成器を使う
4. 統計的検定: 生成された乱数の品質を評価する
5. 再現性: 同じシード値で同じ乱数列を生成できるようにする（デバッグ時に有用）