为什么同一个概率问题会有3种答案？用Python模拟伯特兰悖论

在数学的世界里，概率问题往往看似简单却暗藏玄机。1889年，法国数学家约瑟夫·伯特兰提出了一个令人费解的问题：在一个单位圆内随机画一条弦，这条弦的长度超过内接等边三角形边长的概率是多少？令人惊讶的是，这个问题竟然有三种看似合理的解答——1/3、1/2和1/4。本文将用Python代码模拟这三种解法，带你直观理解这个经典概率悖论背后的数学原理。

1. 伯特兰悖论的三重解答

伯特兰悖论的核心在于"随机"选择弦的方式不同会导致不同的概率结果。让我们先理解这三种经典解法：

1.1 圆周均匀取点法（概率：1/3）

这种方法假设弦的两个端点在圆周上均匀随机分布。想象固定一个端点，让另一个端点在圆周上随机移动：

import numpy as np
import matplotlib.pyplot as plt

def method1(num_samples):
    angles = np.random.uniform(0, 2*np.pi, num_samples)
    return np.where(angles < 2*np.pi/3, 1, 0)

samples = method1(10000)
print(f"方法1概率估计: {np.mean(samples):.4f}")

运行结果示例：方法1概率估计: 0.3327

1.2 半径均匀取点法（概率：1/2）

这种方法先在半径上随机取一点，然后画垂直于半径的弦：

def method2(num_samples):
    positions = np.random.uniform(0, 1, num_samples)
    return np.where(positions < 0.5, 1, 0)

samples = method2(10000)
print(f"方法2概率估计: {np.mean(samples):.4f}")

运行结果示例：方法2概率估计: 0.5013

1.3 圆内均匀取点法（概率：1/4）

这种方法在圆内随机取一点作为弦的中点：

def method3(num_samples):
    # 使用极坐标确保均匀分布
    r = np.sqrt(np.random.uniform(0, 1, num_samples))
    return np.where(r < 0.5, 1, 0)

samples = method3(10000)
print(f"方法3概率估计: {np.mean(samples):.4f}")

运行结果示例：方法3概率估计: 0.2489

2. 可视化三种方法的差异

理解这三种方法的关键在于观察它们生成弦的分布差异。我们可以用Matplotlib创建可视化：

def plot_chords(method, num_samples=100):
    fig, ax = plt.subplots(figsize=(6,6))
    circle = plt.Circle((0,0), 1, fill=False)
    ax.add_patch(circle)
    
    if method == 1:
        angles = np.random.uniform(0, 2*np.pi, num_samples)
        for angle in angles:
            x1, y1 = np.cos(angle), np.sin(angle)
            x2, y2 = np.cos(angle + np.pi/1.5), np.sin(angle + np.pi/1.5)
            ax.plot([x1,x2], [y1,y2], 'b-', alpha=0.3)
    
    elif method == 2:
        positions = np.random.uniform(-1, 1, num_samples)
        for p in positions:
            ax.plot([p,p], [-np.sqrt(1-p**2), np.sqrt(1-p**2)], 'g-', alpha=0.3)
    
    elif method == 3:
        r = np.sqrt(np.random.uniform(0, 1, num_samples))
        theta = np.random.uniform(0, 2*np.pi, num_samples)
        for i in range(num_samples):
            x, y = r[i]*np.cos(theta[i]), r[i]*np.sin(theta[i])
            # 计算弦的端点
            if abs(x) < 1e-10 and abs(y) < 1e-10:
                continue  # 避免除以零
            slope = -x/y if abs(y) > 1e-10 else 1e10
            dx = np.sqrt(1 - r[i]**2) / np.sqrt(1 + slope**2)
            x1, y1 = x + dx, y + slope*dx
            x2, y2 = x - dx, y - slope*dx
            ax.plot([x1,x2], [y1,y2], 'r-', alpha=0.3)
    
    ax.set_aspect('equal')
    ax.set_title(f'方法{method}生成的弦分布')
    plt.xlim(-1.1,1.1)
    plt.ylim(-1.1,1.1)
    plt.show()

plot_chords(1)
plot_chords(2)
plot_chords(3)

这三种可视化会展示完全不同的弦分布模式，解释了为什么概率结果不同。

3. 深入理解悖论本质

伯特兰悖论之所以成为悖论，是因为它揭示了"随机"这个概念在无限样本空间中的模糊性。关键在于：

• 样本空间的选择：每种方法对应不同的样本空间
• 均匀分布的定义：在不同几何维度上"均匀"的含义不同
• 不变性要求：Jaynes提出的解决方案强调尺度和平移不变性

让我们用Python验证Jaynes的观点：

def verify_invariance(method):
    # 测试尺度不变性
    scale_prob = []
    for scale in [0.5, 1, 2]:
        if method == 1:
            angles = np.random.uniform(0, 2*np.pi, 10000)
            prob = np.mean(angles < 2*np.pi/3)
        elif method == 2:
            positions = np.random.uniform(0, scale, 10000)
            prob = np.mean(positions < scale/2)
        elif method == 3:
            r = scale * np.sqrt(np.random.uniform(0, 1, 10000))
            prob = np.mean(r < scale/2)
        scale_prob.append(prob)
    
    print(f"方法{method}尺度不变性测试:", scale_prob)
    
    # 测试平移不变性(仅对方法2有意义)
    if method == 2:
        trans_prob = []
        for dx in [-0.3, 0, 0.3]:
            positions = np.random.uniform(0, 1, 10000) + dx
            positions = positions % 1  # 保持在线段内
            prob = np.mean(positions < 0.5)
            trans_prob.append(prob)
        print(f"方法2平移不变性测试:", trans_prob)

verify_invariance(1)
verify_invariance(2)
verify_invariance(3)

典型输出：

方法1尺度不变性测试: [0.3332, 0.3347, 0.3341]
方法2尺度不变性测试: [0.501, 0.4983, 0.4997]
方法2平移不变性测试: [0.4989, 0.4983, 0.5018]
方法3尺度不变性测试: [0.2516, 0.2493, 0.2502]

4. 实际应用与教学启示

伯特兰悖论不仅是一个有趣的数学谜题，它对我们理解概率和统计有重要启示：

1. 明确定义随机过程：在解决任何概率问题时，必须明确说明随机机制
2. 警惕无限样本空间：无限情况下，"均匀分布"需要特别小心处理
3. 物理系统的建模：在物理和工程中，选择适当的概率模型至关重要

以下是一个教学演示代码，可以让学生交互式探索悖论：

from ipywidgets import interact

def interactive_bertrand(method, num_samples=1000):
    if method == 1:
        angles = np.random.uniform(0, 2*np.pi, num_samples)
        prob = np.mean(angles < 2*np.pi/3)
    elif method == 2:
        positions = np.random.uniform(0, 1, num_samples)
        prob = np.mean(positions < 0.5)
    elif method == 3:
        r = np.sqrt(np.random.uniform(0, 1, num_samples))
        prob = np.mean(r < 0.5)
    
    fig, ax = plt.subplots(figsize=(6,6))
    circle = plt.Circle((0,0), 1, fill=False)
    ax.add_patch(circle)
    
    # 绘制内接三角形
    triangle_x = [np.cos(2*np.pi*i/3) for i in range(4)]
    triangle_y = [np.sin(2*np.pi*i/3) for i in range(4)]
    ax.plot(triangle_x, triangle_y, 'k-')
    
    # 绘制弦
    if method == 1:
        angles = np.random.uniform(0, 2*np.pi, num_samples)
        for angle in angles[:100]:  # 只画100条避免太密集
            x1, y1 = np.cos(angle), np.sin(angle)
            x2, y2 = np.cos(angle + np.pi/1.5), np.sin(angle + np.pi/1.5)
            color = 'red' if angle < 2*np.pi/3 else 'blue'
            ax.plot([x1,x2], [y1,y2], color=color, alpha=0.5)
    
    elif method == 2:
        positions = np.random.uniform(-1, 1, num_samples)
        for p in positions[:100]:
            y_max = np.sqrt(1 - p**2)
            color = 'red' if abs(p) < 0.5 else 'blue'
            ax.plot([p,p], [-y_max, y_max], color=color, alpha=0.5)
    
    elif method == 3:
        r = np.sqrt(np.random.uniform(0, 1, num_samples))
        theta = np.random.uniform(0, 2*np.pi, num_samples)
        for i in range(100):
            x, y = r[i]*np.cos(theta[i]), r[i]*np.sin(theta[i])
            if abs(x) < 1e-10 and abs(y) < 1e-10:
                continue
            slope = -x/y if abs(y) > 1e-10 else 1e10
            dx = np.sqrt(1 - r[i]**2) / np.sqrt(1 + slope**2)
            x1, y1 = x + dx, y + slope*dx
            x2, y2 = x - dx, y - slope*dx
            color = 'red' if r[i] < 0.5 else 'blue'
            ax.plot([x1,x2], [y1,y2], color=color, alpha=0.5)
    
    ax.set_aspect('equal')
    ax.set_title(f'方法{method}: 估计概率={prob:.4f} (理论值={1/3 if method==1 else 0.5 if method==2 else 0.25})')
    plt.xlim(-1.1,1.1)
    plt.ylim(-1.1,1.1)
    plt.show()

interact(interactive_bertrand, method=[1,2,3], num_samples=(100,10000,100))

这个交互式工具可以让学生直观看到不同方法产生的弦分布差异，以及对应的概率估计值如何收敛到不同的理论值。