zoukankan      html  css  js  c++  java
  • Alias采样算法

    一、离散分布

    离散分布:给你一个概率分布,是离散的,比如[1/2, 1/3, 1/12, 1/12],代表某个变量属于事件A的概率为1/2, 属于事件B的概率为1/3,属于事件C的概率为1/12,属于事件D的概率为1/12。

    离散分布的随机变量的取样问题:

    一个随机事件包含四种情况,每种情况发生的概率分别为: 1/2, 1/3, 1/12, 1/12,问怎么用产生符合这个概率的采样方法。

    二、时间复杂度为o(N)的采样算法

    首先将其概率分布按其概率对应到线段上,像下图这样:

     

    接着产生0~1之间的一个随机数,然后看起对应到线段的的哪一段,就产生一个采样事件。比如落在0~ 1/2之间就是事件A,落在1/2~5/6之间就是事件B,落在5/6~11/12之间就是事件C,落在11/12~1之间就是事件D。 

    构建线段的时间复杂度为o(N),

    如果顺序查找线段的话,查找时间复杂度为O(N),

    如果二分查找的话,查找的时间复杂度为O(logN)。

    三、时间复杂度O(1)的采样算法---Alias

    (1)alias分为两步

    • 做表
    • 根据表采样

    (2)做表

    将概率分布的每个概率乘上N,画出柱状图。N为事件数量。

    其总面积为N,可以看出某些位置面积大于1某些位置的面积小于1,

    将面积大于1的事件多出的面积补充到面积小于1对应的事件中,以确保每一个小方格的面积为1,同时,保证每一方格至多存储两个事件,这样我们就能看到一个1*N的矩形啦。

    表里面有两个数组,一个数组存储的是事件i占第i列矩形的面积的比例,即

    另一个数组存储第i列中不是事件i的另一个事件的编号,即

    做表的时间复杂度是o(N)。

    (3)采样

    产生两个随机数,

    第一个产生一个1~N 之间的整数 i,决定落在哪一列。

    第二个产生一个0~1之间的随机数,判断其与Prab[i]大小,

    如果小于Prab[i],则采样 i,(表示若其小于事件i占第i列矩形的面积的比例,则表示接受事件i

    如果大于Prab[i],则采样Alias[i](否则,接收第i列中不是事件i的另一个事件。)

    这种方式采样的时间复杂度为 o(1) ,且对应事件i的概率,完全对应原来的概率分布。

    计算:选择第一列的概率为1/4,则选择第一列的蓝色部分的概率为1/4*2/3 = 1/6,蓝色部分还有第三、第四列,则蓝色(事件A)的总概率为1/4 * 2/3 * 3 = 1/2。【第一次产生随机数的概率为1/4,第二次产生随机数选择蓝色部分的概率为1/4 *(2/3*3) = 1/2】。--------对应原来的概率分布。

    四、代码

    1、制作表:(create_alias_table函数,O(N))

    (1)概率分布 area_ratio 的每个概率乘上N

    (2)small,large分别存放小于1和大于1的事件下标。

    (3)每次从small,large中各取一个,将大的补充到小的之中,小的出队列,再看大的减去补给之后,如果大于1,继续放入large中,如果等于1,则也出去,如果小于1则放入small中。

     

    (4)获取accept、alias

    accept存放第i列对应的事件i矩形的面积百分比;

    alias存放第i列不是事件i的另外一个事件的标号;

    上例中accept=[2/3,1,1/3,1/3],alias=[2,0,1,1],这里alias[1]的0是默认值,也可默认置为-1避免和事件0冲突;

    2、采样:(alias_sample函数,O(1))

    随机采样1~N 之间的整数i,决定落在哪一列。

    随机采样0~1之间的一个概率值,

    如果小于accept[i],则采样i,

    如果大于accept[i],则采样alias[i];

    import numpy as np
    def create_alias_table(area_ratio):
        """
        area_ratio[i]代表事件i出现的概率
        :param area_ratio: sum(area_ratio)=1
        :return: accept,alias
        """
        N = len(area_ratio)
        accept, alias = [0] * N, [0] * N
        small, large = [], []
      ###------(1)概率 * N ----- area_ratio_
    = np.array(area_ratio) * N

    ###------(2)获取small 、large -----
    for i, prob in enumerate(area_ratio_): if prob < 1.0: small.append(i) else: large.append(i)
      ###------(3)修改柱状图 ----- (4)获取accept和alias -----
    while small and large: small_idx, large_idx = small.pop(), large.pop() accept[small_idx] = area_ratio_[small_idx] alias[small_idx] = large_idx area_ratio_[large_idx] = area_ratio_[large_idx] - (1 - area_ratio_[small_idx]) if area_ratio_[large_idx] < 1.0: small.append(large_idx) else: large.append(large_idx)
    while large: large_idx = large.pop() accept[large_idx] = 1 while small: small_idx = small.pop() accept[small_idx] = 1 return accept, alias def alias_sample(accept, alias): """ :param accept: :param alias: :return: sample index """ N = len(accept) i = int(np.random.random()*N) r = np.random.random() if r < accept[i]: return i else: return alias[i]

    参考文献:

    Alias method:时间复杂度O(1)的离散采样方法

    时间复杂度O(1)的离散采样算法—— Alias method/别名采样方法

    Alias sample(别名采样)

  • 相关阅读:
    高级查询
    简单查询
    CRUD
    T-SQL语句
    数据库规范
    导出含有特定字符串的注册表
    .net、jquery、ajax、wcf实现数据库用户名检测局部刷新
    数据结构实验之二叉树二:遍历二叉树
    传纸条
    数据结构实验之栈与队列五:下一较大值(一)
  • 原文地址:https://www.cnblogs.com/Lee-yl/p/12749070.html
Copyright © 2011-2022 走看看