ベルトランのパラドクスをPythonでシミュレーションしてみる

確率のミステリー ベルトランのパラドクス について。Pythonでシミュレーションした。

問題

次のような確率の問題を考える：

「円に無作為に弦をひく。このとき弦の長さが、この円に内接する正三角形の一辺の長さよりも大きくなる確率は？」

答え

$\frac13$、あるいは$\frac14$、あるいは$\frac12$。

なんと3つの矛盾する答えを導くことができる。具体的な計算方法は、Wikipediaの記事 ベルトランの逆説 を参照(難しくないです)。いずれももっともらしいけど、一体どれが本当の正解なのか？

シミュレーション

「3つとも正解なんてありえない。ごちゃごちゃ考えずシミュレーションして確かめればいいじゃない」ということでシミュレーションしてみる。

以下、円を$xy$平面上の単位円(半径1の円)を考える。このとき、円に内接する正三角形の一辺の長さは$\sqrt{3}$である。

Case 1

「無作為に弦をひく」という行為は、「円周上から無作為に2点を選ぶ」ことと等しい。

円の半径が1ならば、円周上の点は角度$\theta$を用いて$(\cos \theta, \sin \theta)$と表すことができるので、$0 \leq \theta_1, \theta_2 \leq 2 \pi$となる$\theta_1$, $\theta_2$を無作為に選ぶことが、円周上の2点

\vec{x_1} = (\cos \theta_1, \sin \theta_1), \ \vec{x_2} = (\cos \theta_2, \sin \theta_2)

を無作為に選ぶ行為となる。

さて、無作為に選んだ2点間の長さ(つまり弦の長さ)が、内接する正三角形の一辺の長さよりも大きいかどうか、つまり

|| \vec{x_1} - \vec{x_2} ||  > \sqrt{3}

を満たすかを確かめる。

これを何万回か繰り返して、上記が真となる個数を全試行回数で割れば、上記の3つの確率の値に近いはずだ。(ひょっとするとどの値にも近くないかもしれないし)

上の方針で検証。弦を1000000回ひいて確かめる。

import math

import numpy as np

n = 1000000
EDGE = math.sqrt(3)
LONGER = 'Longer'
SHORTER = 'Shorter'
ls = {True: LONGER, False: SHORTER}


def case_1():
    seeds = np.random.rand(n, 2) * 2 * math.pi
    samples = [[
        np.array([math.cos(x[0]), math.sin(x[0])]),
        np.array([math.cos(x[1]), math.sin(x[1])])] for x in seeds
    ]
    result_count = {}
    for x in samples:
        e = ls[np.linalg.norm(x[0] - x[1]) > EDGE]
        result_count[e] = result_count[e] + 1 if e in result_count else 0
    return dict(sorted((x[0], x[1] / n) for x in result_count.items()))


print('case_1: ', case_1())

出力。

case_1:  {'Longer': 0.333066, 'Shorter': 0.666932}

というわけで、答えは$\frac13$である。おわり。$\frac14$やら$\frac12$という答えは、計算過程でどこかおかしなところがあったのだろう。

Case 2

ちょっと待てよ… 弦一本に対して、円の中心から弦に垂線をひいたときの交点は一対一に対応しているな。つまり、「無作為に弦をひく」行為は「円の中から無作為に一点を選ぶ」行為と同じだ。

円の内部の1点を$(a_1, a_2)$とすると、弦の長さは三平方の定理から

2 \sqrt{1 - (a_1^2 + a_2^2)}

となる。従って、円内の点を無作為に1000000個選んで、

2 \sqrt{1 - (a_1^2 + a_2^2)} > \sqrt{3}

が真となる個数の割合を求めよう。

import math

import numpy as np

n = 1000000
EDGE = math.sqrt(3)
LONGER = 'Longer'
SHORTER = 'Shorter'
ls = {True: LONGER, False: SHORTER}


def case_2():
    samples = []
    while len(samples) < n:
        seed = np.random.rand(2)
        if np.linalg.norm(seed) <= 1: samples.append(seed)
    result_count = {}
    for x in samples:
        e = ls[(2 * math.sqrt(1 - np.linalg.norm(x) ** 2)) > EDGE]
        result_count[e] = result_count[e] + 1 if e in result_count else 0
    return dict(sorted((x[0], x[1] / n) for x in result_count.items()))


print('case_2: ', case_2())

出力。

case_2:  {'Longer': 0.249522, 'Shorter': 0.750476}

$\frac14$が出てきた。

Case 3

ダメ押しで、「円の中から無作為に一点を選ぶ」ときに極座標$(r, \theta)$で考えてみる。

判定式は、やはり三平方の定理から

2 \sqrt{1 - r^2} > \sqrt{3}

となる。

この条件でサンプル数1000000回のシミュレーションをしてみる。

尚、この判定式に$\theta$が現れないので、コードの中に$\theta$に該当するものはない。

import math

import numpy as np

n = 1000000
EDGE = math.sqrt(3)
LONGER = 'Longer'
SHORTER = 'Shorter'
ls = {True: LONGER, False: SHORTER}


def case_3():
    samples = np.random.rand(n)
    result_count = {}
    for x in samples:
        e = ls[2 * math.sqrt(1 - x ** 2) > EDGE]
        result_count[e] = result_count[e] + 1 if e in result_count else 0
    return dict(sorted((x[0], x[1] / n) for x in result_count.items()))

print('case_3: ', case_3())

出力。

case_3:  {'Longer': 0.50085, 'Shorter': 0.499148}

$\frac12$が出てしまった。

考察

シミュレーションをしてみても$\frac13$, $\frac14$, $\frac12$という3つの答えが出てしまった。

ここまで試した3つのCaseは、シミュレーションに用いる1000000個のサンプルを作るときに考えた「無作為に弦をひく」という行為の解釈に違いがある：

• Case 1: 「円周上の2点を無作為に選ぶ」
• Case 2: 「円内の1点を無作為に選ぶ」

Case 3はわかりにくいが、極座標で捉えなおしただけに見えるが、実際行っているのは

• Case 3: 「固定した半径上の1点を無作為に選ぶ」

である。

実際、それぞれのCaseで使った判定式をもとに確率を計算しても、Wikipediaと同じ結論がでてくる。

Case 1 の場合

加法定理やらで注意深く不等式を評価すると

\cos (\theta_1 - \theta_2) < - \frac12

が出てくる。$0 \leq |\theta_1 - \theta_2| < 2 \pi$となる角度だけを考えることにすれば、

\frac23 \pi < |\theta_1 - \theta_2| < \frac43 \pi 

となる。さらに 円周角 = $\frac12$ × 中心角 となることを思い出すと、円周角に関する評価は

\frac13 \pi < \frac12 |\theta_1 - \theta_2| < \frac23 \pi 

となって、これはWikipediaの中の「無作為な端点」方式の解答に相当することがわかる。

Case 2 の場合

これも判定式を変形すれば

0 \leq a_1^2 + a_2^2 < \frac14

となって、これはWikipediaの中の「無作為な中点」方式の解答に相当することがわかる。

Case 3の場合

判定式が

0 \leq r < \frac12

となって、これはWikipediaの中の「無作為な半径」方式の解答に相当することがわかる。

まとめ

3つとも正解でした。つまり、無作為に弦を選ぶ方法毎に、それぞれ正しい答えであった、ということだった。

これを物理的に本当に実験してみたら、どのような結果が出てくるのか。実験を始めようと思うときに、どのように「無作為に弦をひくか」を考えるだろう。地面に円を描いてそこにシャープペンシルの芯を投げ込んで弦を描くか、あるいはボールを円の中に投げ込んでそれを中点とするような弦を描くか、あるいは円周上で回るコマが倒れた方向に弦を描くのか ... と考えると、今回のシミュレーションでサンプルを作るときに考えたことと似てくる。まじりっけなしの純粋な「無作為さ」ってなんだ？

この問題がパラドクス足り得るのは、「無作為に選ぶ」という行為の曖昧さに依っているのでは。普段"無作為"や"ランダムさ"という概念は、文脈に依存せず、唯一つのことを意味しているように思う。しかし今回のシミュレーションでわかったことは、この"ランダムさ"という概念が、文脈に大きく依存して、都度異なる様相を成すことであった。

あと標本空間が無限であることも起因しているのだろう。似たようなパラドクスに水とワインの問題(p.76あたり)がある。

ちなみに、この問題への解答はこの3つの他にもある。この本によると、他にも$\frac{14}{23}, \frac23 ,\frac34, 1 - \frac{\sqrt{3}}{2}$という解答があることが順番に示されて、さらには無限に多くの解答があることまでわかっているらしい。

ランダムに～という言葉が出てきたら一旦深呼吸しましょう、っていうことですかね。

おわり

参考

• 書籍: 確率は迷う: 道標となった古典的な33の問題 面白かった