水塘抽样(Reservoir Sampling)问题

水塘抽样是一系列的随机算法，其目的在于从包含n个项目的集合S中选取k个样本，其中n为一很大或未知的数量，尤为适用于不能把全部n个项目都存放到主内存的状况。

在高德纳的计算机程序设计艺术中，有以下问题：能否在一未知大小的集合中，随机取出一元素？。或者是Google面试题： I have a linked list of numbers of length N. N is very large and I don’t know in advance the exact value of N. How can I most efficiently write a function that will return k completely random numbers from the list（中文简化的意思就是：在不知道文件总行数的状况下，如何从文件中随机的抽取一行？）。两题的核心意思都是在总数不知道的状况下如何等几率地从中抽取一行？便是说若是最后发现文字档共有N行，则每一行被抽取的几率均为1/N？

咱们能够：定义取出的行号为choice，第一次直接以第一行做为取出行 choice ，然后第二次以二分之一律率决定是否用第二行替换 choice ，第三次以三分之一的几率决定是否以第三行替换 choice ……，以此类推。由上面的分析咱们能够得出结论，在取第n个数据的时候，咱们生成一个0到1的随机数p，若是p小于1/n，保留第n个数。大于1/n，继续保留前面的数。直到数据流结束，返回此数，算法结束。

 

问题一

首先考虑k为1的状况，即：给定一个长度很大或者长度未知数据流，限定对每一个元素只能访问一次，写出一个随机选择算法，使得全部元素被选中的几率相等。

设当前读取的是第n个元素，采用概括法分析以下：

1. n = 1 时，只有一个元素，直接返回便可，几率为1。
2. n = 2 时，须要等几率返回前两个元素，显然几率为1/2。能够生成一个0～1之间的随机数p，p < 0.5 时返回第一个，不然返回第二个。
3. n = 3 时，要求每一个元素返回的几率为1/3。注意此时前两个元素留下来的几率均为1/2。作法是：生成一个0～1之间的随机数，若<1/3，则返回第三个，不然返回上一步留下的那个。元素1和2留下的几率均为：1/2 * (1 - 1/3) = 1/3，即上一步留下的几率乘以这一步留下（即元素3不留下）的几率。
4. 假设 n = m 时，前n个元素留下的几率均为：1/n = 1/m；
5. 那么 n = m+1 时，生成0～1之间的随机数并判断是否<1/(m+1)，如果则留下元素m+1，不然留下上一步留下的元素。这样一来，元素m+1留下的几率为1/(m+1)，前m个元素留下来的几率均为：1/m * (1 - 1/(m+1)) = 1/(m+1)，也就是1/n。
6. 综上可知，算法成立。

 

问题二

将问题一中的条件变为，k为任意整数的状况，即要求最终返回的元素有k个，这就是水塘抽样（Reservoir Sampling）问题。要求是：取到第n个元素时，前n个元素被留下的概率相等，即k/n。

算法同上面思路相似，将1/n换乘k/n便可。在取第n个数据的时候，咱们生成一个0到1的随机数p，若是p小于k/n，替换池中任意一个为第n个数。大于k/n，继续保留前面的数。直到数据流结束，返回此k个数。可是为了保证计算机计算分数额准确性，通常是生成一个0到n的随机数，跟k相比，道理是同样的。

一样采用概括法来分析：

1. 初始状况 n <= k：此时每一个元素留下的几率均为1。
2. 当 n = k+1 时，第k+1个元素留下的几率为k/(k+1)，前k个元素留下的几率均为：k/k * (1 - k/(k+1) * 1/k) = k/(k+1)，即上一步留下的几率乘以这一步留下的几率。
3. 假设 n = m 时，每一个元素留下的几率均为 k/n = k/m。
4. 那么，当 n = m+1 时，第m+1个元素留下的几率为1/(m+1)，前m个元素留下的几率均为：k/m * (1 - k/(m+1) * 1/k) = k/(m+1)，其中：k/m为上一步留下来的几率，k/(m+1) * 1/k 为这一步不能留下来的几率（第m+1个留下来，同时池中一个元素被踢出的几率）。
5. 综上可知，算法成立。

 

伪代码以下：

 1 //stream表明数据流
 2 //reservoir表明返回长度为k的池塘
 3 
 4 //从stream中取前k个放入reservoir；
 5 for ( int i = 1; i < k; i++)
 6     reservoir[i] = stream[i];
 7 for (i = k; stream != null; i++) {
 8     p = random(0, i);
 9     if (p < k) reservoir[p] = stream[i];
10 return reservoir;