【モンテカルロ法】ビュフォンの針の問題

はじめに

こんにちは。よっしーです。
前回の記事ではモンテカルロシミュレーションというシミュレーション手法を紹介しました。
今回はそれを用いて確率に  \pi が出現するような面白い事象を紹介したいと思います。
確率に円周率が出てくると言うのは日常生活ではなかなか考えづらいので一つの知識として知っておいていただけますと幸いです。

ビュフォンの針

多数の平行線を引き、そこに針を落とした時、いづれかの線と針が交差する確率はいくつか。

こちらはジョルジュ=ルイ・ルクレール・ド・ビュフォンが提示した「ビュフォンの針の問題」と言われている問題です。
事象としては図を見たほうが分かりやすいので図を載せます。


5 本の線(外枠を除くグレーの横線)と針(赤線と青線)を描画してみました。
赤線は平行線と交差した針、青線は平行線と交差していない針を表現しています。

10 本中 3 本が平行線と交差していることが分かります。

ところで、針と線が交わる確率は針の長さと平行線の間隔で決まります。
なぜなら、とてつもなく長い針を用いたり、極端に狭い幅の平行線だとすぐに針が平行線に触れてしまうからです。(短い針、間隔が広い場合も同様です)

針が長く幅が狭い場合

針が短く幅が広い場合

上図の 2 つでは確率に差が出ますので、それも考慮して確率を考えていく必要があります。

ここで、以前行ったモンテカルロシミュレーションを使えば、確率計算せずとも色々なシチュエーションで確率を求めることができます。

実際にシミュレーションしてみる

実際にシミュレーションしてみましょう。
今後、針の長さを l , 平行線の幅を d と置きます。

 l = 10 ,  d = 20 で試してみます。
10000 本の針を落としてみました。

10000 本中 3186 本平行線と交差したので確率は 0.3186 でした。
逆に、何本に 1 本の割合で平行線にヒットしたかを考えると 3.1387…でした。
これは確率の逆数を取ったと言えます。
どこか見覚えのある数字の気がする。
そうです。これこそが円周率  \pi なのです。

回数を増やして見てみますと以下のようになります。
結果は以下となりました。

試行回数 交差した本数 確率 確率の逆数
10000 3186 0.3186 3.1387
1000000 317980 0.3180 3.145
100000000 3184038 0.3184 3.1407

このことをまとめると l = 10 ,  d = 20 の時、確率  P
 P \sim \dfrac{1}{\pi}
となることが分かります。

ここで、 l  d を変更してしてみて、1000000 回試行してみます。
両方一度に変えるとややこしいので  l のみ変えてみます。

l d 交差した本数 確率 確率の逆数
2.5 20 79649 0.0796 12.5551
5 20 158908 0.1589 6.2929
7.5 20 238626 0.2386 4.1907
10 20 317980 0.3180 3.1449
12.5 20 398002 0.398 2.5126
15 20 477211 0.4772 2.0955
17.5 20 556991 0.557 1.7954
20 20 636654 0.6367 1.5707
22.5 20 691036 0.691 1.4471
25 20 727357 0.7274 1.3748
27.5 20 756181 0.7562 1.3224
30 20 778936 0.7789 1.2838
32.5 20 796710 0.7967 1.2552
35 20 812240 0.8122 1.2312
37.5 20 825567 0.8256 1.2113
40 20 837685 0.8377 1.1938

 d = 20 で固定した時の  l と確率 P の関係は以下のようです。

 l が 20 になるまではグラフが比例しているように見えます。

 l < d の場合と  d \le l の場合で関係性が違いそうです。
 d \le l の場合は複雑なので、 l < d の場合を考えると、
 P \propto l
と置いておきます。

上記の関係性を裏付けるために  d = 30 の時を見てみると、

l d 交差した本数 確率 確率の逆数
2.5 30 53343 0.0533 18.7466
5 30 106089 0.1061 9.4260
7.5 30 159141 0.1591 6.2837
10 30 212292 0.2123 4.7105
12.5 30 264526 0.2645 3.7803
15 30 317836 0.3178 3.1463
17.5 30 371694 0.3717 2.6904
20 30 423372 0.4234 2.3620
22.5 30 476238 0.4762 2.0998
25 30 530611 0.5306 1.8846
27.5 30 584066 0.5841 1.7121
30 30 636691 0.6367 1.5706
32.5 30 675306 0.6753 1.4808
35 30 704432 0.7044 1.4196
37.5 30 728430 0.7284 1.3728
40 30 747569 0.7476 1.3377
42.5 30 763607 0.7636 1.3096
45 30 778750 0.7788 1.2841
47.5 30 791145 0.7911 1.2640
50 30 803157 0.8032 1.2451
52.5 30 812986 0.8130 1.2300
55 30 821634 0.8216 1.2171
57.5 30 830030 0.8300 1.2048
60 30 837049 0.8370 1.1947


 d = 20 の時の関係と一致していそうです。

今度は l を固定して、 d を変化させてみましょう。

l d 交差した本数 確率 確率の逆数
10 2 936056 0.9361 1.0683
10 4 870416 0.8704 1.1489 dfrac
10 6 802577 0.8026 1.2460
10 8 727014 0.7270 1.3755
10 10 636667 0.6367 1.5707
10 12 530250 0.5303 1.8859
10 14 448384 0.4484 2.2302
10 16 389057 0.3891 2.5703
10 18 345457 0.3455 2.8947
10 20 317980 0.3180 3.1449

 l = 10 で固定した時の  d と確率  P の関係は以下のようです。

こちらも  l < d の場合と  d \le l の場合で関係性が違いそうです。
 l < d の場合は反比例してるように見えますが図では分かりづらいため、確率の逆数との関係を見てみると、

ということで、以下の関係性があると考えられます。
 P \propto \dfrac{1}{d}

上記のことをまとめると、  l < d の場合は以下の関係性となりそうです。
 P \propto \dfrac{l}{d}
つまり、
 P \sim \alpha\dfrac{l}{d}
と書くことができます。

揃っているデータを代入すると係数  \alpha が分かります。

l d P α
10 12 0.5303 0.6364
10 14 0.4484 0.6278
10 16 0.3891 0.6226
10 18 0.3455 0.6219
10 20 0.3178 0.6356
2.5 20 0.0796 0.6368
5 20 0.1589 0.6356
7.5 20 0.2386 0.6363
10 20 0.3181 0.6362
12.5 20 0.3980 0.6368
15 20 0.4772 0.6363
17.5 20 0.5570 0.6366
20 20 0.6367 0.6367
2.5 30 0.0533 0.6396
5 30 0.1061 0.6366
7.5 30 0.1591 0.6364
10 30 0.2123 0.6369
12.5 30 0.2645 0.6348
15 30 0.3178 0.6356
17.5 30 0.3717 0.6372
20 30 0.4234 0.6351
22.5 30 0.4762 0.6349
25 30 0.5306 0.6367
27.5 30 0.5841 0.6372
30 30 0.6367 0.6367

上記データから  \alpha の中央値を取ると
 \alpha = 0.6364 \sim \dfrac{2}{\pi} です。

結論を総合すると、
 l < d の場合、
 P \sim \dfrac{2l}{\pi d}
となります。

今回利用したプログラムは以下です。

考察

確率に  \pi が出現した理由

今回の確率に  \pi が出てきたことは角度が関係あります。
理系の方はよく知っていると思いますが、弧度法における角度はラジアンで表されて、 360^{\circ} = 2\pi です。
コードにも記載しておりますが、ランダムで角度を設定する際に上記ラジアンを用いています。

そのためランダムな角度で針を落下させるということを考えると、確率に  \pi が出てくることも納得ができると思います。

針の長さと平行線の幅の長さによる場合分け

 l < d の場合と  d \le l の場合で振る舞いが変わった理由に関しては凄く単純化して考えてみましょう。

 l < d の場合
針の長さがほぼ平行線の幅と同じ時を考えます。

上図にあるように、以下の 2 パターンが存在しています。

  1. 平行線と平行線のちょうど間の点 O1 に針の中心を落とした場合はどの角度で落としても平行線と交わることは無い。
    そのように、どの角度で針を落としても平行線と交わることの無い点が存在する領域を領域 A と呼ぶ
  2. 点 O2 に張りの中心を落とした場合は特定の角度の場合は平行線と交わる。
    そのように針を落とした時、特定の角度で平行線と交わる点が存在する領域を領域 B と呼ぶ

 l \le d の場合

上図にあるように、領域 A が存在せずすべての空間が領域 B となり、性質が異なるために差が生じることになります。

おわりに

今回はモンテカルロシミュレーションで日常ではあまり馴染みの無い値が出る確率に関して紹介しました。
難しい数式をできるだけ省いて、シミュレーションでできる範囲を述べているので、より詳しく知りたい人は調べてみると数学的に面白いことが分かったりするかもしれません。
ではまた次回も見ていただけますと幸いです。

参考図書

直感を裏切る数学 「思い込み」にだまされない数学的思考法
https://www.amazon.co.jp/dp/B00RXCOQDK/

最近の記事

  • 関連記事
  • おすすめ記事
  • 特集記事

アーカイブ

カテゴリー

PAGE TOP