单次遍历，带权随机选取问题（一）

扩展：选取 m 个元素，概率理论值

来看看选取多个元素的问题。当选取多个元素时，可以认为选取过程是逐步进行的，即无放回的多次选取。每一次选取时，任何一个元素被选中的概率都与其权重成正比，但总的权重则又剩余的元素集合决定。

当 m=2 的时候，第 i 个元素被选中可以是两种情况：第一次就被选中；第一次未被选中，第二次被选中。可以得到其概率为这两种情况的概率之和，即：

值得注意的是，即便 w_i 和 w 不变，如果其他元素的概率分布不同，最后得到的结果也不同，因此上面这个式子无法把其中的求和化简掉。

从另一方面来看，第 i 个元素被选中的概率等于 1 减去它不被选中的概率。用 \(\bar p\) 表示不被选中的概率，则有：

显然，\(p_i(2)+\bar p_i(2)=1\)。

当 m>2 时，其概率表达式将会变得异常复杂，因为跟概率分布有关，所以算式无法化简。未被选中的概率计算式要稍微简单些，大概是这个样子的：

其中，\(\forall 1\leq k\leq m,j_k\notin\{i,j_1,j_2,\cdots,j_{k-1}\}\)。

对于给定的一组权重，可以用下面这段程序计算出任意 m、i（程序中的 i 是从 0 开始的）对应的概率数值（请无视其 coding style）：

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18

def Foo(weights, ids, totalweight, m, i, times):
  if times == m: return 1
  p = 0.0
  for j in ids:
    ids.remove(j)
    p += float(weights[j]) / totalweight \
         * Foo(weights, ids, totalweight - weights[j], m, i, times + 1)
    ids.add(j)
  return p

def CalcSampleProbability(weights, m, i):
  n = len(weights)
  assert 0 <= i < n, 'invalid i'
  assert 0 < m <= n, 'invalid m'
  ids = set(xrange(n))
  ids.remove(i)
  p = Foo(weights, ids, sum(weights), m, i, 0)
  return 1 - p

可惜算法的复杂度非常高，CalcSampleProbability 需要 O(n^m) 时间来完成一次计算。期待高手改进。

来看一下等权重、等差数列权重和等比数列权重的 n 选 m 概率分布图（图中 i 依旧采用 1 <= i <= n 的取值范围）：

Mathematica 提供了 RandomSample 函数，支持带权选取，当然它是在遍历之前就已经知道元素个数的。给它一组等差分布的权重，可以看出十万次随机选取后得到的概率分布与上面的理论分布非常接近。

count = 10;
weights = Range[count];
elems = Range[count];
retry = 100000;
map = Table[
    freq = ConstantArray[0, count];
    For [i = 0, i < retry, i++,
        freq += BinCounts[RandomSample[weights -> elems, m], {1, count + 1, 1}]
    ];
    freq, {m, 1, count, 1}
];
ListLinePlot[map / retry, PlotMarkers -> Automatic]

Mathematica RandomSample 随机选取 m 个元素，第 i 个元素被选中的概率

苦苦思考了好几天，但并没有想到一个直观的方法可以给之前的 RandomSample 加上权重处理。因为那概率式子太复杂，实在不知道该怎么去凑。不过在下一篇文章中将会介绍一个神奇的算法（当然不是我想出来的），并且会给出我的证明。

现在发文章的速度越来越慢了，实在因为能力有限，每次为了一两个式子都要演算很久。再接再厉。