水塘抽樣(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;