问题

有这样的一个概率问题：假设我们有 \(M\) 个盒子，\(N\) 个相同的球，可以选择如何把这 \(N\) 个球分配到这 \(M\) 个盒子中。然后我们往外拿球。每次从有球的盒子里面，按照相等的概率随机选一个盒子，然后从这个盒子里拿一个球出来。当只剩一个盒子里面有球的时，停止拿球。问题是：怎样把这 \(N\) 个球放入 \(M\) 个盒子中，使得我们拿出来的球的期望数量最大。

这个问题可以用Markov链来进行建模。首先引入1-范数，对于 \(x=(x_1,\cdots,x_M)\in\mathbb{R}^M\)，定义

\[ \Vert x\Vert_1=\sum_{i=1}^M\vert x_i\vert, \]

1-范数会诱导Manhattan距离，对于 \(x,y\in\mathbb{R}^M\)，定义 \[ d(x,y)=\Vert x-y\Vert_1. \] 考虑一个离散时间Markov链 \((V_t)\). 在格点 \(\mathbb{Z}^M\)上，从满足\(x_1+\cdots+x_M=N,\,x_i\ge0\) 的点 \(x=(x_1,\cdots,x_M)\) 出发，每次随机选一个非0分量进行减1，到达任意一个坐标轴被吸收。停时 \(\tau\) 可以写成 \[ \tau=\inf\left\{t\ge0:\sum_{i=1}^M \mathbf{1}_{V_t^{(i)}\ge 1}=1, V_t=\left(V_t^{(1)},\cdots,V_t^{(M)}\right)\right\}. \] 目标是选择最优的初始位置 \(x\)，来最大化游走的期望时间 \[ \max_{\Vert x\Vert_1=M}\mathrm{E}[\tau]. \] 假设被吸收的位置是 \(y=V_\tau\)，那目标就是最大化 \(\tau=d(x,y)\) 的期望，即 \[ \max_{\Vert x\Vert_1=M}\mathrm{E}[d(x,y)]. \]

这也等价于最小化\(d(0,y)\)的期望，因为\(d(0,y)+d(x,y)=d(0,x)=N\).

猜想

最符合直觉的一个答案是尽量把\(N\)个球平均分成\(M\)堆，让\(x\)的各分量等于\(\lfloor\frac{N}{M}\rfloor\)或\(\lfloor\frac{N}{M}\rfloor+1\)。从x点出发，它的概率密度波沿\((-1,-1,\cdots,-1)\)方向扩散，只是波前的形状不是欧式距离下的圆形而是曼哈顿距离下的菱形。让\((-1,-1,\cdots,-1)\)在\(x\)处对准原点，能让它更晚撞上吸收壁。

想法虽然直观，但是我并不能严格证明这个猜想。下面我会给出一些具体的计算。

计算

动态规划

\((V_t)\)在任何一点被吸收的概率可以直接计算。因为用组合的方法，原点到\(x\)的合法路径数应该是 \[ \frac{N !}{x_{1} ! \cdots x_{k} !} \] 也即Multinomial distribution pmf的系数，但这个计算会比较复杂。更常见的解决思路，是利用Markov性给出一些期望的递推关系，然后用动态规划来解决。

比如设\(\mathrm{E}[d(x,y)\mid V_0=x]=f(x)\)。又设\(e_k\in\mathbb{Z}^M\)的第\(k\)个分量为1，其余分量为0。规定边界条件 \[ f(c\,e_k)=0,\quad(c=0,1,2,\cdots). \] 由\(\{v_t\}\)的Markov性可知 \[ \begin{align*} f(x)&=\mathrm{E}[d(x,y)\mid V_0=x]\\ &=\sum_{k=1}^M\mathrm{P}\left(V_1=\hat{x}_k\mid V_0=x\right)\mathrm{E}[d(V_0,y)\mid V_0=x,V_1=\hat{x}_k]\mathbf{1}_{x_k\ge 1}\\ &=\sum_{k=1}^M\frac{1}{\sum_{k=1}^M\mathbf{1}_{x_k\ge 1}}\,\mathrm{E}[d(V_0,V_1)+d(V_1,y)\mid V_1=\hat{x}_k]\mathbf{1}_{x_k\ge 1}\\ &=\frac{1}{\sum_{k=1}^M\mathbf{1}_{x_k\ge 1}}\sum_{k=1}^M\left(1+\mathrm{E}[d(V_1,y)\mid V_1=\hat{x}_k]\right)\mathbf{1}_{x_k\ge 1}. \end{align*} \] 由此得到递推关系 \[ \begin{align*} f(x_1,\cdots,x_M)&=\frac{1}{\sum\limits_{k=1}^M\mathbf{1}_{x_k\ge 1}}\sum_{k=1}^M \left(1+f(x_1,\cdots,x_k-1,\cdots,x_M)\right)\mathbf{1}_{x_k\ge 1}\\ &=1+\frac{1}{\sum\limits_{k=1}^M\mathbf{1}_{x_k\ge 1}}\sum_{k=1}^M f\left(x_1,\cdots,x_k-1,\cdots,x_M\right)\mathbf{1}_{x_k\ge 1}. \end{align*} \] 进一步规定边界条件 \[ f(x_1,\cdots,x_M)=0 \text{ if }x_k< 0, \] 得到 \[ \begin{align*} f(x_1,\cdots,x_M)&=1+\frac{1}{\sum\limits_{k=1}^M\mathbf{1}_{x_k\ge 1}}\sum_{k=1}^M f\left(x_1,\cdots,x_k-1,\cdots,x_M\right)\\ &=1+\frac{1}{\sum\limits_{k=1}^M\mathbf{1}_{x_k\ge 1}}\sum_{k=1}^M f\left(x-e_k\right). \end{align*} \] 这是一个 \(M\) 维的动态规划.

3维情况

在 \(M=3\) 的情况下，下图展示了 \(f(x_1,x_2,x_3)\) 的一些计算结果。