Reservoir Sampling 学习笔记

前言: Reservoir Sampling 中文翻译是水塘抽样, 之的是一组随机抽样算法, 而不是某一个具体的算法. 这类算法主要用于解决这样一个问题: 当样本总体很大或者是在数据流上进行采样的时候, 我们往往无法预知总体的样本实例个数N. 那么Reservoir Sampling 就是这样一组算法, 即使不知道N, 也能保证每个样本实例被采样到的概率依然相等. 下面是一个其中一个非常典型的在线(online)线性算法, 可以用于数据无法装入内存或者在流数据上使用:

代码:

const int k = 10;

vector<int> rsampleing(stream sin) {
    int val = 0;
    int cnt = 0;
    vector<int> result(k, -1);
    while(sin >> val) {
        if(cnt < k) {
            result[cnt] = val;
        }
        else {
            // generate a random int within [0, cnt]
            int rnd = rand(0, cnt + 1); // !!! rand(0, cnt) is incorrect when k = 1
            if(rnd < k)
                result[rnd] = val;
        }
        ++cnt;
    }
    return result;
}

分析:
上面这个算法的步骤是这样, 从总体S中抽取前面的k个实例放入预制的数组中, 这个数组就是我们最后要返回的抽样结果, 对于后面的所有样本实例, 从cnt = k个开始, 我们对每一个生成一个[0, cnt]的随机数rnd, 若 rnd < k, 那我们就用当前的val替换掉result[cnd]. 这样做为什么能保证每个实例被抽到的概率当等而且概率为k/(n + 1)呢?

可以这样分析: 我们知道对于第cnt个实例, 当算法遇到他的时候, 他被选中进入result的概率是k / (cnt + 1), 那么他依然出现在最后的result的情况是, cnt后面所有的实例都没有取代掉他, cnt后面任何第t > cnt个实例取代掉他的概率是k / (t + 1) / k = 1 / (t + 1), 即 t被选中(k / (t + 1)), 而且选中取代原来cnt所在的位置(k / (t+1) / k). 所以后面任意一个实例不取代cnt的概率就是 1-1/(t + 1), 那么要所有的的情况都发生, 最后cnt才能留在result中, 这样就是一个连乘的结果:

(k/ (cnt + 1)) * (1 - 1/(cnt+2)) * (1 - 1/(cnt+3)...*(1 - 1/(N + 1))) = k/(N + 1)

特别注意: 和Wiki里面的概率值不一样, 主要是因为实现的时候非常tricky, 注意代码里面的这一句int rnd = rand(0, cnt + 1);我们产生的随机数的的范围必须是一个闭区间[0, cnt], 而不能是开区间[0, cnt), 为什么呢, 考虑k = 1的时候, 这个时候当我们放进第一个实例的时候, 接下来的所有sample都会进到else这个分支, 如果我们是开区间的话, 第二个sample总是会取代掉第一个sample, 也就是说第一个sample永远不会被抽样到, 可以自己分析加深理解, 当然当N很大的时候我们的概率k/(N+1)和k/N没有区别, 并且我们还是保证每一个实例被抽中的概率相等.

参考: http://zh.wikipedia.org/wiki/水塘抽樣

(全文完,转载时请注明作者和出处)


(转载本站文章请注明作者和出处 烟客旅人 sigmainfy — http://www.sigmainfy.com,请勿用于任何商业用途)

Written on September 11, 2013