码迷,mamicode.com
首页 > 其他好文 > 详细

【编程珠玑】【第一章】生成随机数、随机取样的问题

时间:2018-01-17 23:24:17      阅读:255      评论:0      收藏:0      [点我收藏+]

标签:turn   put   gen   insert   问题:   奇数   答案   整数   res   

一、利用随机数函数生成随机数

问题1(《编程珠玑》习题12.1后半段):

给定一个rand(),可以产生从0RAND_MAX的随机数,其中RAND_MAX很大(常见值:16int能表示的最大整数32767),写出利用rand()生成[a,b]中任意整数的函数,其中a>=0, b<=RAND_MAX,且b-a<<RAND_MAX。

分析:

这是在编程工作最常见的随机函数的应用,在这里做一个起点再合适不过。把随机数区间的起点从0变为a,同时把一共RAND_MAX+1个数的区间缩小至只含有b-a+1个数的区间,写为 a + rand()%(b-a+1),此时显然最大值是a+(b-a)= b进一步地,这个b-a<<RAND_MAX的条件虽然看上去不起眼,其实很重要。

附加思考:

如果b-aRAND_MAX很接近会发生什么情况?读者不妨先做思考,问题2的分析会做出解答。这个rand(),其实相当于《编程珠玑》提到的bigrand()

 

问题2(笔试题变形):

给定一个随机数函数rand7(),它能以等概率生成1~77个整数。请根据rand7()写出类似的rand5()

分析:

如果直接像问题1中一样,把1+rand7()%5作为rand5()会有什么情况发生?这时确实能产生1~5的随机数没错,可是各个数的概率相等吗?

对于随机数2,既有可能来自于1+1%5,也有可能来自于1+6%5,显然其概率是2/7,而不是1/5,不满足rand5()等概率产生各随机数这一隐含要求。不同于问题1,问题1中一个很大的区间收缩成较小区间时,各个元素映射后的新元素概率虽然概率可能不完全一样,但却是近似相同的。

为了满足等概率的要求,可以这么做:

int rand5() {
    int res;
    do {
        res = rand7();
    } while(t>=6);
    return res;
}

虽然保证了1~5的概率都变成了1/5,但是有一个无法避免的缺点是,每当产生了6或者7都要抛弃,相当于这一次运行是“空转”,浪费了时间。如果对1/5这个概率不明白,可以有两种理解:每次产生67就被抛弃,剩下数的概率相等,必然为1/5;或者用更严密的推理:产生1~5的随机数,最终得到某一个的概率为:1/7+(1/7)*(2/7)+(1/7)*(2/7)2+...,无限项求和,结果是1/5

问题3(笔试题原题)

给定一个随机数函数rand7(),它能以等概率生成1~77个整数。请根据rand7()写出类似的rand10()

分析:

有了问题2的概率基础,把rand7()变成rand10()仅仅需要一点点思考了。

int rand10() {
    int t1,t2,res;
    do {
        t1 = rand7();      //t1以等概率1/5成为{1,2,3,4,5}中的一个数字
    } while(t1>=6);
   
    do {
        t2 = rand7();    //t2以1/2概率成为奇数或偶数,因为{1,2,3,4,5,6]中3奇3偶
    } while(t2==7); 
    res = t1+5*(t2%2);       //res是1~10中的任意一个数的概率都是1/10,注意到%和*具有相同的优先级,这里去掉括号结果是错的    
return res;
}

解释:

t1的概率是{1~5}中每个元素概率为1/5t21/2概率为奇或者偶,所以(t2%2)是以1/2概率选择0或者1,因此5*(t2%2)是以1/2概率为0或者5。这样,选择1的概率是1/5 * 1/2=1/10...,选择6的概率是1/5 * 1/21/10....由此可见。

问题2和问题3是对《编程珠玑》上使用范围不大的randint(i,j)生成其他范围随机数方法的解答。在掌握了问题2和问题3的解法后,你已经学会随机数区间的收缩和扩张,类似的问题迎刃而解。

问题4(《编程珠玑》习题12.1前半段)

C库函数rand()常返回15个随机位,写出bigrand(),能够返回30个随机位。

分析:

其实问题4和问题3有点像,但是不同之处在于,这次我们的视角是从位出发的,把rand()看做了将15个位每一个位以1/2概率设为01,从而生成0~RAND_MAX。从某种意义上来说,按这种理解方式来解这个问题更轻松一些,但是仅局限于2的幂减1这样的数值的区间,比如从011...11。在此基础之上把两部分的位拼接起来即可。

//《编程珠玑》的答案
//int bigrand() {return RAND_MAX *rand() +rand();}
//最大值不是30个1,在另一个笔记中有叙述,怀疑有错我的答案,把先生成的部分左移15位作为高15位
long bigrand(){return ((long)RAND_MAX+1)*rand() +rand();}

 扩展:拒绝采样

这是一个简单粗暴有效的方法,使用这个方法你可以不用考虑复杂的区间扩展和收缩的问题。以rand7()产生rand10()为例,构造一个二维数组并把它填成:

技术分享图片

然后使用两次rand7()分别生成行号和列号,如果对应元素是0,抛弃重来;否则即是结果。其实这种方法依然用到了区间扩展:把7扩展成7*7,并把不符合的部分抛弃。这样横向生成一个1~7的随机数,纵坐标生成1~7的随机数,这样两次随机数便能定位一个数组元素,若该元素为1~10之间的数,则被选中,否则重新选择。这能保证1~10之间的每个数被选中的概率都是1/10从这里可以看出这种方法其实还是有缺陷的:如果用rand7()生成rand50(),那这50个数可是填满这个表格后还有填不下的。

 

二、利用随机数函数产生随机事件

用随机事件表示随机事件?这个问题具体化为:使用rand()表示以M/N的概率发生的随机事件,其中M<=N,并可用作:if(事件A发生)  ,其中P(A) = M/N,那么表示为:

if(rand()%N< M)

    ...

通过这种方式,我们就可以做出让程序“以M/N的概率执行某个命令”这样的设计了。

 

三、取样问题:n个元素中随机选取m

 

1.从概率角度出发

考虑整数0,1,2,...,n-1,可以用上节的方法以概率m/n选取0(推导方式略)。但是对于1,必须考虑之前0是否被选取而以(m-1)/(n-1)m/(n-1)的概率选取1,后续就更加麻烦。好在迭代是计算机的长项,只需要把这个是否选择当前数的随机事件稍作修改即可,使之变成从r个其余的数中选择s个:

int gen(int m,int n){
    int i, select = m,remaining = n;
    for(i=0;i<n;i++) {
        if(rand() % remaining <select) {
            printf("%d\n",i);
            select--;
        }
        remaining--;
    }
    return 0;
}

其概率证明可见于KnuthThe Art of Computer Programming2卷。进一步地可以优化为

int genknuth(int m,int n){
    int i;
    for(i=0;i<n;i++)
        if(rand()%(n -i) < m) {
            printf("%d\n",i);
            m--;
        }
    return 0;
}

《编程珠玑》提示,这种算法是“所有n的所有m元子集被选中的概率相等”,这个条件强于“所有元素被选中的概率相同”。下面是习题12.2中提到的“以等概率选择有元素,但是有些m元子集被选中的概率比其他子集大”的算法:直接选择1个数,则这个m元集合为它本身即后续的一共m个数,可能包括回绕。

对于这种方法,总要产生n次随机数。进一步可以写为for(i=0;i<n&&m>0;i++),但程序的平均运行时间是否变得更快需要权衡。对应地,习题12.7提供了一种稍微减少随机数产生次数的递归函数:

int randselect(int m,int n) {
    int r;
    //assert(m<=n && m>=0);
    if(m>0) {
        r = rand()%n;
        if(r < m) {
            printf("%d\n",n-1);
            randselect(m-1,n-1);
        }
        else
            randselect(m,n-1);
    }
    return 0;
}

 

2.从集合插入出发

由于集合元素不重复,如果按等概率选择一个随机数,不在集合中就把它插入,反之直接抛弃,直到集合元素个数达到m个,同样可以满足要求,并且用C++STL很容易实现:

void gensets(int m,int n) {
    set<int> S;
    while(S.size() < m)
        S.insert(rand()%n);
    set<int>::iterator i;
    for(i = S.begin();i!=S.end();++i)
        cout<<*i<<"\n";
}

这个算法的主要问题是,如果抛弃已存在的元素的次数过多,相当于多次产生随机数并进行集合操作,性能将明显下降。比如当n=100m=99,取第99个元素时,算法“闭着眼睛乱猜整数,直到偶然碰上正确的那个为止”(《编程珠玑(续)》,13.1节)。虽然这种情况会在“从一般到特殊”提供解决方案,但下面的Floyd算法明显规避了产生随机数超过m次的问题。

习题12.9提供了一种基于STL集合的随机数取样方法,可以在最坏情况下也只产生m个随机数:限定当前从中取值的区间的大小,每当产生重复的随机数,就把这一次迭代时不会产生的第一个随机数拿来替换。

int genfloyd(int m,int n){
    set<int> S;
    set<int>::iterator i;
    for(int j = n-m; j<n;j++) {
        int t = rand()%(j+1);
        if(S.find(t) == S.end())
            S.insert(t);
        else
            S.insert(j);
    }
    for(i=S.begin();i!=S.end();++i)
        cout<<*i<<"\n";
}

不必基于集合而实现,我自己写了一个原理类似的算法,思想是把“用过”的元素“扔”到下一次迭代的随机数取样区间之外:

int gen(int m,int n) {
    int *array;
    int i,j;
    array = (int *)malloc(sizeof(int) * n);
    for(i=0;i<n;i++)
        array[i] = i;
    while(m>=1) {
        j = rand()%n;
        printf("%d\n",array[j]);
        if(j<n-1)
            array[j] = array[n-1];
        m--;
        n--;
    }
    return 0;
} 

3“打乱顺序”出发

这是个来源于实际的想法:将所有n个元素打乱,取出前m个。更快的做法是,打乱前m个即可。对应的C++代码如下:

int genshuf(int m,int n){
    int i,j;
    int *x = new int[n];
    for(i = 0;i<n;i++)
        x[i] = i;
    for(i = 0;i<m;i++) {
        j = randint(i,n-1);  //randint产生i到n-1之间的随机数
        int t = x[i];x[i] = x[j];x[j] = t;
    }
    //sort(x,x+m);            //sort是为了按序输出 
    for(i=0;i<m;i++)
        cout<<x[i]<<"\n";
}       

sort注释掉的这段代码,可以作为随机不重复序列产生器。类似的还有Floyd的算法P。(《编程珠玑(续)》,13.3节)

4从一般到特殊

以上讨论的几种方式都不限定mn的取值,只需m<=n即可。对于特殊的取值,有特殊的解决方案,以下是编程珠玑上的两例:

  1.n = 106m=n-10,这时m很大较为接近n可以先生成10个元素,然后输出其余的元素。进行这种处理的额外代码可以提高算法的平均速度。

  2.n=231m = 107m<<n,这时可以先生成1.1*107个元素,排序后去掉重复的(由于n很大,m中出现重复元素的概率很低),得到 107个元素的样本。

附:

(《编程珠玑(续)》习题13.5Doug McIlory的求N个元素中取M种的第G种情况的组合的算法,原书这个算法我没有理解,也没有看到比较满意的解释,可能用到了某些我不熟悉的组合数性质。原书中没有对其进行进一步解释,并且似乎只是作者的题外话而已。如果看着有困难,这部分代码完全可以跳过。介绍这个算法的原因是用它能在获得随机数G后,直接获得第GN个元素取M个的取法,相当于只产生了一次随机数。

int comb(int n,int m,int g,int* array){
    int d = 1,t;
    while(m>0) {
        t = combination(n-d,m-1);
        //combination(n,m)是计算n中取m个的组合数
        if(g-t<0) {
            m--;
            array[m] = d;
            printf("%d\n",d);
        }
        else
            g-=t;
        d++;
    }
    return 0;
}

三、取样问题:从未知总数的元素中选择一个

从事先未知总数为n的对象中随机选择一个的方法。有两种常见的具体问题:

  1.读取一个未知行数的文件,随机输出其中的的一行,同时最多只能缓冲一行的内容(《编程珠玑(续)》习题15.3利用了这种形式);

  2.对链表进行一趟遍历,随机输出其中的一个结点的元素,只能使用一个临时指针。

解法是,以概率1选择第一个元素存入缓存,以概率1/2用第二个元素替换掉缓存,...,直至遍历完所有元素,最后输出缓存的内容。可以分析出此时所有元素留在缓存的概率均为1/n,比如1,是1*(1/2)*(2/3)*...*((n-1)/n)

int random_select(void){
    int i,num=1;
    for(i=1;i<n;i++)
       if(rand()%i ==0)     //等于0或者任何一个0~i-1之间的数,表示概率1/i。
           num = i;         //以1/i的概率替换num的缓存值。
    return num;
}

对于这个问题的常见应用:马尔科夫文本生成器里选择哈希表中某项对应链表的任意一个结点。

 

【编程珠玑】【第一章】生成随机数、随机取样的问题

标签:turn   put   gen   insert   问题:   奇数   答案   整数   res   

原文地址:https://www.cnblogs.com/blastbao/p/8306501.html

(0)
(0)
   
举报
评论 一句话评论(0
登录后才能评论!
© 2014 mamicode.com 版权所有  联系我们:gaon5@hotmail.com
迷上了代码!