ベルトランのパラドックスを考えてみよう

ベルトランのパラドックスという確率の問題を知っていますか。考え方により数パターンの答えが現れる面白い問題です。この現象をPythonコードにより確認し、また私なりに考えてみたので共有します。

投稿者は数学が専門ではなく、用語の使い方が正確ではないかもしれません。ご了承の上、お付き合いください。

ベルトランのパラドックスとは

以下の問題を考えます。

「円に内接する正三角形を考える。その円の弦を1本無作為に選ぶ。その弦が正三角形の辺よりも長くなる確率はどれだけか？」¹

この問題の例を図示します。

円上で無作為に選んだ弦(赤線)が、正三角形の各辺(青線)より長くなる場合の確率を考えます。左の図のような弦を選ぶと弦は正三角形の1辺より長く、右のような弦では正三角形の一変より短いことが分かります。

この問題ではパラドックスというだけあり、3つの考え方がよく知られています。それぞれの考え方を案①〜③と表現します。以下にその3種類の案を記載しますが、この問題を初めて知った方はここで先に読むのを止めて、ぜひ一度考えてみてください。

案①　無作為な端点

この案では次のようなアプローチで、弦を引きます。
1. 円周上に2個の点をランダムに選び、配置する。
2. その2点を結び弦を作り、その長さと円に内接する正三角形の1辺の長さと比較する。

この選び方では、弦が1辺より長くなる確率は$1/3$になることが知られています。説明のために1個目の点のみを配置した次のような図を考えます。

円の下部に設定した点を基準にして円の内接正三角形の内部に入る$60^{\circ}$の範囲内では弦は1辺より長く(図の緑線)なります。対して、正三角形の外である円の接線方向から$60^{\circ}$の範囲では弦のほうが短く(図の紫線)なります。選択できる$180^{\circ}$の中で、正三角形の中に入るのは$60^{\circ}$ですので、確率は$1/3$になります。

ちなみに筆者がまず思いついたのも、この考え方でした。

案②　無作為な半径

この案では次のようなアプローチで、弦を引きます。
1. 円の中心からランダムな角度$\theta$を選ぶ。
2. 手順１で選んだ角度上の線の中から、ランダムな点を選ぶ。
　　(ここまでの操作が無作為な半径を選ぶ)
3. 手順２で選んだ点を中点とする弦を引き、三角形の1辺の長さと比較する。

この選び方が案②であり、この場合三角形の1辺より長くなる確率は$1/2$となります。次のように作図するとわかりやすいです。

補助線として、考える弦に平行な辺を持つ内接正三角形を二つ用意しました。円の半径を$r$とした場合に、緑線2本で示した領域($r=\frac{r}{2} + \frac{r}{2}$)内に入れば、正三角形の1辺より長くなります。逆に緑線の外の領域では1辺よりも短くなります。このエリアも($r=\frac{r}{2} + \frac{r}{2}$)の長さとなります。
よって、案②の考え方をすると三角形1辺よりも長くなる確率は$1/2$となることがわかります。

案③　無作為な中点

この案では次のようなアプローチで、弦を引きます。
1. 円の中でランダムな点を選びます。
2. その点を中点とする弦を引き、三角形の1辺の長さと比較します。

この選び方を行うと、正三角形の1辺より長くなる確率は$1/4$となります。これも、次のように作図するとわかりやすいです。

円の内接正三角形にさらに内接する円(緑)を作図します。この円は外側の円の半分の半径になり、この緑の領域内に中点を持つ弦は、三角形の1辺より長くなります。
この例では外側の円の半径が$r$、内側の緑の円の半径が$\frac{2}{r}$です。そのため面積を計算すると、外側の円が

S = \pi r^2

となり、内側の円が

S_{in} = \pi \left(\frac{r}{2}\right)^2 = \frac{\pi r^2}{4}

となります。つまり内側の円の面積は外側の円に比べて$1/4$となります。円の中でランダムに中点とする座標を1か所選び、その領域が緑の領域となる確率は$1/4$というのが、案③になります。

確認

各案の試行を指定した回数行うスクリプトを作成しました。いずれのスクリプトも単位円を検討しており、単位円に内接する3角形の1辺は$\sqrt3$となります。繰り返し弦を引き、弦の長さが$\sqrt3$を超えた割合を求めます。

(※、以下の"▲"マークをクリックすると実装したスクリプトと解説が開きます。)

案①を実行するスクリプト

s1_RandomEndpoints.py

import numpy as np

def main():

    print("Please input N:")
    N = int(input())

    ta = 0
    fa = 0
    for i in range(N):
        if(e2eLength()):
            ta+=1
        else:
            fa+=1
    
    print("N:",N, "True:", ta, "False", fa)
    print("T/N:",ta/N, "& F/N:",fa/N, "%" )

def e2eLength():
    rd1 = np.random.rand()*2*np.pi
    x1 = np.cos(rd1)
    y1 = np.sin(rd1)

    rd2 = np.random.rand()*2*np.pi
    x2 = np.cos(rd2)
    y2 = np.sin(rd2)

    e2eLength = np.sqrt((x2-x1)**2 + (y2-y1)**2)
    
    if e2eLength > (np.sqrt(3)):
        return True
    else:
        return False

if __name__ == "__main__":
    main()

このスクリプトでは、２個の角度をランダムで指定し、その角度から弦を結ぶ2点の座標を取得します。スクリプト実行後にNを入力することで、試行回数を指定できます。

案②を実行するスクリプト

s2_RandomRadius.py

import numpy as np

def main():

    print("Please input N:")
    N = int(input())

    ta = 0
    fa = 0
    for i in range(N):
        if(e2eRadius()):
            ta+=1
        else:
            fa+=1
    
    print("N:",N, "True:", ta, "False", fa)
    print("T/N:",ta/N, "& F/N:",fa/N, "%" )

def e2eRadius():
    #中心からの角度をランダムに
    rd = np.random.rand()*2*np.pi #角度
    xrd = np.cos(rd)
    yrd = np.sin(rd)

    diffLength = np.random.rand() #中心からの長さ

    rx0 = xrd/diffLength
    ry0 = yrd/diffLength

    phi = np.arccos(diffLength)
    alpha = rd - phi
    beta = rd + phi

    x1 = np.cos(alpha)
    y1 = np.sin(alpha)

    x2 = np.cos(beta)
    y2 = np.sin(beta)

    e2eLength = np.sqrt((x2-x1)**2 + (y2-y1)**2)
    
    if e2eLength > (np.sqrt(3)):
        return True
    else:
        return False

if __name__ == "__main__":
    main()

このスクリプトでは1つの角度$\theta$と、中心からの距離$d$(0〜1)をランダムに指定します。単位円の場合には、中心からの距離$d$が決まると、

\phi=\arccos(d)

により計算できる$\phi$を利用して、弦の2点の角度をそれぞれ

\displaylines{
\alpha &=& \theta - \phi \\
\beta &=& \theta + \phi
}

で計算できるので、この角度から端点の位置を計算し、弦の長さを計算しました。

案③を実行するスクリプト

s3_RandomCenter.py

import numpy as np

def main():

    print("Please input N:")
    N = int(input())

    ta = 0
    fa = 0
    for i in range(N):
        if(e2eCenter()):
            ta+=1
        else:
            fa+=1
    
    print("N:",N, "True:", ta, "False", fa)
    print("T/N:",ta/N, "& F/N:",fa/N, "%" )

def e2eCenter():
    while True:
        #ランダム座標を取得
        x0 = np.random.rand()*2-1 #中点　x座標候補
        y0 = np.random.rand()*2-1 #中点　y座標候補
        r0 = np.sqrt(x0**2 + y0**2)
        if r0 < 1:
            break # 採用

    rd = np.arctan2(y0,x0)
    if rd < 0:
        rd += 2*np.pi

    diffLength = r0 #中心からの長さ

    phi = np.arccos(diffLength)
    alpha = rd - phi
    beta = rd + phi

    x1 = np.cos(alpha)
    y1 = np.sin(alpha)

    x2 = np.cos(beta)
    y2 = np.sin(beta)

    e2eLength = np.sqrt((x2-x1)**2 + (y2-y1)**2)
    
    if e2eLength > (np.sqrt(3)):
        return True
    else:
        return False

if __name__ == "__main__":
    main()

このスクリプトでは、座標$(x,y)$をランダムに指定します$x$と$y$はそれぞれ-1〜1の範囲でランダムに選びますが、円の外の座標を選んだ場合には選び直しています。座標が決まれば$\arctan()$を利用して角度を計算できます。また座標位置から中心との距離$d$も計算できますので、案②と同様に弦の長さを計算できます。

それぞれ1,000,000回ずつ試行した結果を表に示します。

パターン	TRUE	FALSE	確率
案①	333896	666104	0.3339%
案②	499702	500298	0.4997%
案③	250623	749377	0.2506%

案①が概ね$1/3$、案②が$1/2$、案③が$1/4$であり、それぞれ期待した確率にほぼ一致していることがわかります。よって考え方の違いから、導き出される確率が異なるという、ベルトランのパラドックスを確認することができました。

考察

ベルトランのパラドックスは確認できましたが、三つの案で導かれる答えの違いはどこから来るのでしょうか。

そこでプロットとヒストグラムから各案の傾向を調べました。

確認方法	試行回数	説明
プロット	5,000	作成した弦の中点をプロットする。
ヒストグラム	1,000,000	作成した弦の中点と原点(円の中心)間の距離のヒストグラム。

なおヒストグラムの$r_0 < 0.5$の領域が、正三角形の1辺より長い領域を示します。

案①　無作為な端点の場合

案①　プロット	案①　ヒストグラム

案①プロット　プログラム

s1Plot_RandomEndpoints.py

import numpy as np
import matplotlib.pyplot as plt

centerX_list = []
centerY_list = []

def main():

    print("Please input N:")
    N = int(input())

    ta = 0
    fa = 0
    for i in range(N):
        if(e2eLength()):
            ta+=1
        else:
            fa+=1
    
    print("N:",N, "True:", ta, "False", fa)
    print("T/N:",ta/N, "& F/N:",fa/N, "%" )

    plt.figure()
    plt.plot(centerX_list,centerY_list, "bo", markersize=1)
    plt.xlim(-1,1)
    plt.ylim(-1,1)
    plt.axis('equal')
    plt.show()
    input()

def e2eLength():
    rd1 = np.random.rand()*2*np.pi
    x1 = np.cos(rd1)
    y1 = np.sin(rd1)

    rd2 = np.random.rand()*2*np.pi
    x2 = np.cos(rd2)
    y2 = np.sin(rd2)

    e2eLength = np.sqrt((x2-x1)**2 + (y2-y1)**2)
    
    centerX = (x2+x1)/2
    centerY = (y2+y1)/2
    r0 = np.sqrt(centerX**2+centerY**2)
    centerX_list.append(centerX)
    centerY_list.append(centerY)

    if e2eLength > (np.sqrt(3)):
        return True
    else:
        return False

if __name__ == "__main__":
    main()

案①ヒストグラム　プログラム

s1Count_RandomEndpoints.py

import numpy as np
import matplotlib.pyplot as plt

centerPointLength_list = []

def main():

    print("Please input N:")
    N = int(input())

    ta = 0
    fa = 0
    for i in range(N):
        if(e2eLength()):
            ta+=1
        else:
            fa+=1
    
    print("N:",N, "True:", ta, "False", fa)
    print("T/N:",ta/N, "& F/N:",fa/N, "%" )

    nplist = np.array(centerPointLength_list)
    
    fig = plt.figure()
    ax = fig.add_subplot(1,1,1)
    ax.hist(nplist, bins=50)
    ax.set_title('pettern 1 histgram')
    ax.set_xlabel('r0')
    ax.set_ylabel('count')
    fig.show()
    input()

def e2eLength():
    rd1 = np.random.rand()*2*np.pi
    x1 = np.cos(rd1)
    y1 = np.sin(rd1)

    rd2 = np.random.rand()*2*np.pi
    x2 = np.cos(rd2)
    y2 = np.sin(rd2)

    e2eLength = np.sqrt((x2-x1)**2 + (y2-y1)**2)
    
    centerX = (x2+x1)/2
    centerY = (y2+y1)/2
    r0 = np.sqrt(centerX**2+centerY**2)
    centerPointLength_list.append(r0)

    if e2eLength > (np.sqrt(3)):
        return True
    else:
        return False

if __name__ == "__main__":
    main()

まずプロットに着目すると円の中心付近に多くのプロットがあります。また円周付近にもプロットが多数確認できます。プロット点の原点からの長さを$r_0$とした場合に、$|r_0|<0.5$の範囲に入る円の中心に近いプロット（具体的には）は弦のほうが長くなります。

次に中心から弦の中点までの距離$r_0$のヒストグラムを確認します。これを見ると1が圧倒的に多いことがわかります。ヒストグラムの横軸$r_0$が0付近は中心付近のプロットを、1付近は円周近くのプロットを示します。ヒストグラムをみると0～0.6程度までは差がないように見えますが、プロットする領域は$r_0$が大きくなるにつれて増加し特に1付近で最大となります。つまり円周付近のプロットが一番多くなることが分かります。

逆にヒストグラムからは$r_0$が0.5以下ではほぼ一定に見えますが、プロットでは円の中心付近に密集して見えるのは、同じ半径の変化量$dr$の場合には、$r_0$小さい領域のほうがより面積は小さくなるので、密集してるように見えます。

案①方式でランダムに選んだ場合には、このような傾向となることが分かりました。

• 中心付近と円周付近にプロットが密になる
• ヒストグラムは$r_0$が大きくなるに連れて指数関数的に増加する

案②　無作為な半径の場合

案②　プロット	案②　ヒストグラム

案②プロット　プログラム

s2Plot_RandomRadius.py

import numpy as np
import matplotlib.pyplot as plt

centerX_list = []
centerY_list = []

def main():

    print("Please input N:")
    N = int(input())

    ta = 0
    fa = 0
    for i in range(N):
        if(e2eRadius()):
            ta+=1
        else:
            fa+=1
    
    print("N:",N, "True:", ta, "False", fa)
    print("T/N:",ta/N, "& F/N:",fa/N, "%" )

    
    plt.figure()
    plt.plot(centerX_list,centerY_list, "bo", markersize=1)
    plt.xlim(-1,1)
    plt.ylim(-1,1)
    plt.axis('equal')
    plt.show()
    input()


def e2eRadius():
    #中心からの角度をランダムに
    rd = np.random.rand()*2*np.pi #角度
    xrd = np.cos(rd)
    yrd = np.sin(rd)

    diffLength = np.random.rand() #中心からの長さ

    phi = np.arccos(diffLength)
    alpha = rd - phi
    beta = rd + phi

    x1 = np.cos(alpha)
    y1 = np.sin(alpha)

    x2 = np.cos(beta)
    y2 = np.sin(beta)

    e2eLength = np.sqrt((x2-x1)**2 + (y2-y1)**2)
    
    rx0 = xrd*diffLength
    ry0 = yrd*diffLength
    r0 = np.sqrt(rx0**2 + ry0**2)
    centerX_list.append(rx0)
    centerY_list.append(ry0)


    if e2eLength > (np.sqrt(3)):
        return True
    else:
        return False

if __name__ == "__main__":
    main()

案②ヒストグラム　プログラム

s2Count_RandomRadius.py

import numpy as np
import matplotlib.pyplot as plt

centerPointLength_list = []

def main():

    print("Please input N:")
    N = int(input())

    ta = 0
    fa = 0
    for i in range(N):
        if(e2eRadius()):
            ta+=1
        else:
            fa+=1
    
    print("N:",N, "True:", ta, "False", fa)
    print("T/N:",ta/N, "& F/N:",fa/N, "%" )

    nplist = np.array(centerPointLength_list)
    
    fig = plt.figure()
    ax = fig.add_subplot(1,1,1)
    ax.hist(nplist, bins=50)
    ax.set_title('pettern 2 histgram')
    ax.set_xlabel('r0')
    ax.set_ylabel('count')
    fig.show()
    input()


def e2eRadius():
    #中心からの角度をランダムに
    rd = np.random.rand()*2*np.pi #角度
    xrd = np.cos(rd)
    yrd = np.sin(rd)

    diffLength = np.random.rand() #中心からの長さ

    phi = np.arccos(diffLength)
    alpha = rd - phi
    beta = rd + phi

    x1 = np.cos(alpha)
    y1 = np.sin(alpha)

    x2 = np.cos(beta)
    y2 = np.sin(beta)

    e2eLength = np.sqrt((x2-x1)**2 + (y2-y1)**2)
    
    rx0 = xrd*diffLength
    ry0 = yrd*diffLength
    r0 = np.sqrt(rx0**2 + ry0**2)
    centerPointLength_list.append(r0)

    if e2eLength > (np.sqrt(3)):
        return True
    else:
        return False

if __name__ == "__main__":
    main()
    

このプロットでも円の中心付近に多くのプロットがあるのは案①と同じですが、より中心付近に密集しています。また案①と異なり円周付近にはプロットが密集しないという特徴があります。中心から外に向かうにつれて中点のプロット密度が下がっていることが確認できます。

ヒストグラムを見ると$r_0$が0から1.0までのすべての領域で、同程度の値となりました。案②では、角度と半径について自由度がありますが、半径は0～1までのいずれの値になる確率も同じなので、円の中心はプロットが濃く、円周近くが薄くなる結果と合致します。またヒストグラムから案①(確率：$1/3$)に比べ$r_0 < 0.5$の領域が多いことから、正三角形の1辺より長い弦となる確率が高い($1/2$である)ことも分かります。

まとめると、案②方式でランダムに選んだ場合には、以下の傾向となります。

• 中心付近のプロットが密になる
• ヒストグラムは$r_0$によらず、いずれの領域も等しい確率で発生する

案③　無作為な中点の場合

案③　プロット	案③　ヒストグラム

案③プロット　プログラム

s3Count_RandomCenter.py

import numpy as np
import matplotlib.pyplot as plt

centerX_list = []
centerY_list = []

def main():

    print("Please input N:")
    N = int(input())

    ta = 0
    fa = 0
    for i in range(N):
        if(e2eCenter()):
            ta+=1
        else:
            fa+=1
    
    print("N:",N, "True:", ta, "False", fa)
    print("T/N:",ta/N, "& F/N:",fa/N, "%" )

    plt.figure()
    plt.plot(centerX_list,centerY_list, "bo", markersize=1)
    plt.xlim(-1,1)
    plt.ylim(-1,1)
    plt.axis('equal')
    plt.show()
    input()

def e2eCenter():
    while True:
        #ランダム座標を取得
        x0 = np.random.rand()*2-1 #中点　x座標候補
        y0 = np.random.rand()*2-1 #中点　y座標候補
        r0 = np.sqrt(x0**2 + y0**2)
        if r0 < 1:
            break # 採用

    rd = np.arctan2(y0,x0)
    if rd < 0:
        rd += 2*np.pi

    diffLength = r0 #中心からの長さ

    phi = np.arccos(diffLength)
    alpha = rd - phi
    beta = rd + phi

    x1 = np.cos(alpha)
    y1 = np.sin(alpha)

    x2 = np.cos(beta)
    y2 = np.sin(beta)

    e2eLength = np.sqrt((x2-x1)**2 + (y2-y1)**2)

    centerX_list.append(x0)
    centerY_list.append(y0)
    
    if e2eLength > (np.sqrt(3)):
        return True
    else:
        return False

if __name__ == "__main__":
    main()

案③ヒストグラム　プログラム

s3Plot_RandomCenter.py

import numpy as np
import matplotlib.pyplot as plt

centerPointLength_list = []

def main():

    print("Please input N:")
    N = int(input())

    ta = 0
    fa = 0
    for i in range(N):
        if(e2eCenter()):
            ta+=1
        else:
            fa+=1
    
    print("N:",N, "True:", ta, "False", fa)
    print("T/N:",ta/N, "& F/N:",fa/N, "%" )

    nplist = np.array(centerPointLength_list)
    
    fig = plt.figure()
    ax = fig.add_subplot(1,1,1)
    ax.hist(nplist, bins=50)
    ax.set_title('pettern 3 histgram')
    ax.set_xlabel('r0')
    ax.set_ylabel('count')
    fig.show()
    input()

def e2eCenter():
    while True:
        #ランダム座標を取得
        x0 = np.random.rand()*2-1 #中点　x座標候補
        y0 = np.random.rand()*2-1 #中点　y座標候補
        r0 = np.sqrt(x0**2 + y0**2)
        if r0 < 1:
            break # 採用

    rd = np.arctan2(y0,x0)
    if rd < 0:
        rd += 2*np.pi

    diffLength = r0 #中心からの長さ

    phi = np.arccos(diffLength)
    alpha = rd - phi
    beta = rd + phi

    x1 = np.cos(alpha)
    y1 = np.sin(alpha)

    x2 = np.cos(beta)
    y2 = np.sin(beta)

    e2eLength = np.sqrt((x2-x1)**2 + (y2-y1)**2)

    centerPointLength_list.append(r0)
    
    if e2eLength > (np.sqrt(3)):
        return True
    else:
        return False

if __name__ == "__main__":
    main()

案３は円周内に対してランダムで中点を決定するので、プロットでは密な領域がありません。またヒストグラムで確認する$r_0 = 1.0$付近が最大であり、$r_0$に比例して単純増加する傾向がわかります。

当初私は円の面積なのでヒストグラムは$r_0$の2乗に比例するのかと思ったのですが、$r_0$に比例しているようです。例として、ある領域に着目し$r_2 > r_1$の面積$S$に着目すると

S = \pi r_2^2 - \pi r_1^2

と表現できます。ここで、$r_2 = r$、$r_1 = r - dr$という仮定を置くと

\displaylines{
S &=& \pi r^2 - \pi (r-dr)^2 \\
 &=& 2\pi rdr - dr^2 
}

となり、$r$の2乗項が消え、また$dr^2$は微小領域だと無視できるため、$r$に比例する項$2\pi rdr$により$r$に比例することが分かります。

まとめると、案③方式でランダムに選んだ場合には、以下の傾向となります。

• プロットは一様となり、密な部分はない
• ヒストグラムは、$r_0$の大きさに比例した確率で発生する

最後に

この記事ではベルトランのパラドックスとして知られる問題について、3種の解法を比較しました。各案それぞれで無作為に弦を選ぶ操作を行っていますが、プロットやヒストグラムで比較すると分布が異なる、すなわち確率密度分布が異なることが確認できました。

ちなみに、このような定義により、異なる解が複数存在する数学の問題を"well-defined"な問題と呼ぶそうです。

ここまで読んでくださって、ありがとうございました!

1. Wikipedia:ベルトランの逆説 ↩