高维概率问题

这个问题是视频《Darts in Higher Dimensions》后的思考题。

视频是由3Blue1Brown的Manim (Mathematical Animation Engine) 引擎制作

问题

等概率地从区间 $[0,1]$ 中不停地取数，直到这些数的和超过 $1$ 为止。记 $x_i$ 为第 $i$ 次取出的数，$n$ 为这些数的个数，即有 $\sum\limits_{i=1}^n x_i>1$ ，且 $\sum\limits_{i=1}^{n-1} x_i \le 1$ 。求个数 $n$ 的数学期望 $E(n)$ 。

解答

显然个数 $n$ 的数学期望

$$\tag{1} E(n)=\sum_{k=1}^\infty k \cdot P(n=k)$$

问题转变为求 $n=k$ 时的概率 $P(n=k)$ 。

记 $P_k=P(x_1+x_2+\cdots+x_k \le 1)$ ，则 $P_k$ 表示个数 $n$ 不小于 $k+1$ 的概率，即 $P_k=P(n>k)$ 。这样有

$$\tag{2} P(n=k) = P(n>k-1) - P(n>k)$$

将式$(2)$代入式$(1)$可得

$$\tag{3} \begin{aligned} E(n)&=\sum_{k=1}^\infty k \cdot \big[P(n>k-1) - P(n>k)\big]\\[5pt] &=\sum_{k=1}^\infty P(n>k-1)=\sum_{k=1}^\infty P_{k-1} \end{aligned}$$

问题又转变为求概率 $P_k$ 。

$P_k=P(x_1+x_2+\cdots+x_k \le 1)$ ，显然，$P_0=1$ 、$P_1=1$

对于二维空间，$P_2=\frac{\text{三角形面积}}{\text{正方形面积}}=0.5$ ；对于三维空间，$P_3=\frac{\text{三棱锥体积}}{\text{正方体体积}}=0.1\dot{6}$ 。实际上概率 $P_k$ 就对应 $k$ 维空间分割。

将 $k$ 维切片得到 $k-1$ 维，将 $k-1$ 维的测度积分得到 $k$ 维的测度（属实抽象）

下面采用纯概率的思想来计算 $P_k$ 。

$$\overset{0}{|} \overbrace{\cdots\cdots}^{\displaystyle x_1} | \underbrace{\overbrace{\cdots\cdots}^{\displaystyle x_2} | \underbrace{\overbrace{\cdots\cdots}^{\displaystyle x_3}}_{1-x_1-x_2}}_{1-x_1} \overset{1}{|}$$

如上图所示，以三维为例，计算 $P_3$ 的值：

• 首先确定取出 $x_1$，则从长度为 $1$ 的区间中取到 $x_1$ 的概率为 $\mathrm{d}x_1\big/1$ ；
• 接着取出 $x_2$，为满足要求，则只能从长度为 $1-x_1$ 的区间中取值，同样的概率为 $\mathrm{d}x_2\big/1$ ；
• 最后取出 $x_3$，而 $x_3$ 的可选区间由 $x_1$、$x_2$ 确定，其长度为 $1-x_1-x_2$ 。

于是，$P_3$ 的计算式为

$$\tag{4} P_3=\Big(\int_0^1 \frac{\mathrm{d}x_1}{1}\Big)\Big(\int_0^{1-x_1} \frac{\mathrm{d}x_2}{1}\Big)\big(1-x_1-x_2\big)=\frac{1}{6}$$

应用上述计算思想来求 $P_k$ ：

$$\tag{5} \begin{aligned} P_k&=\Big(\int_0^1 \frac{\mathrm{d}x_1}{1}\Big)\cdots\Big(\int_0^\textcolor{coral}{1-x_1-\cdots-x_{n-2}} \frac{\mathrm{d}x_{n-1}}{1}\Big)\big(\textcolor{coral}{1-x_1-\cdots-x_{n-2}}-x_{n-1}\big) \\[5pt] &=\Big(\int_0^1 \frac{\mathrm{d}x_1}{1}\Big)\cdots\Big(\int_0^\textcolor{coral}{1-x_1-\cdots-x_{n-3}} \frac{\mathrm{d}x_{n-2}}{1}\Big)\frac{\big(\textcolor{coral}{1-x_1-\cdots-x_{n-3}}-x_{n-2}\big)^2}{2!} \\[5pt] &=\Big(\int_0^1 \frac{\mathrm{d}x_1}{1}\Big)\cdots\Big(\int_0^\textcolor{coral}{1-x_1-\cdots-x_{n-4}} \frac{\mathrm{d}x_{n-3}}{1}\Big)\frac{\big(\textcolor{coral}{1-x_1-\cdots-x_{n-4}}-x_{n-3}\big)^3}{3!} \\[5pt] &=\cdots\cdots\cdots \\[5pt] &=\int_0^1\frac{(1-x_1)^{k-1}}{(k-1)!} \mathrm{d}x_1 \\[7pt] &=\frac{1}{k!} \end{aligned}$$

式$(5)$的推导用到了下式

$$\int_0^t \frac{(t-x)^k}{k!}\mathrm{d}x=\frac{t^{k+1}}{(k+1)!}$$ 式中 $t$ 代表式$(5)$中被标记为橘色的部分

将式$(5)$代入式$(3)$可得

$$\tag{6} \begin{aligned} E(n)&=\sum_{k=1}^\infty \frac{1}{(k-1)!}\\[5pt] &=1+\frac{1}{1!}+\frac{1}{2!}+\frac{1}{3!}+\frac{1}{4!}+\cdots\\[5pt] &=\mathrm{e} \end{aligned}$$

个数 $n$ 的数学期望为自然常数 $\mathrm{e}$ ，约为2.71828182845904

实验验证

MATLAB程序：

% 由实验计算个数n的数学期望
clc;clear;

num=10^6;% 随机实验次数
record1=zeros(num,1);% 记录每次实验得到的个数

for i=1:num
    s=0;k=0;
    while s<=1
        temp=rand(1);% 返回一个在区间(0,1)内均匀分布的随机数
        s=s+temp;
        k=k+1;
    end
    record1(i)=k;
end

record2=tabulate(record1);% 统计个数k的频数
EN=sum(record2(:,1).*record2(:,3)/100);% 个数n的数学期望
disp(['期望 = ',sprintf('%1.5f',EN)]);

随机实验次数与 $E(n)$ 的关系：

实验次数	数学期望 $E(n)$
$10^2$	2.78000
$10^3$	2.70000
$10^4$	2.71710
$10^5$	2.71647
$10^6$	2.71826
$10^7$	2.71881
$10^8$	2.71842

---