• Alias sample(别名采样)


    应用场景:加权采样,即按照随机事件出现的概率抽样

    具体算法:

    as

    举例如上,随机事件出现的概率依次是1/2,1/3,1/12,1/12;记随机事件的个数为N,则所有事件概率乘以N后概率为2,4/3,1/3,1/3;

    记队列small,large分别存放小于1和大于1的事件下标(例子中small=[0,1],large=[2,3]);

    记accept存放第i列对应的事件i矩形的面积百分比;alias存放第i列不是事件i的另外一个事件的标号;

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

    上例中accept=[2/3,1,1/3,1/3],alias=[1,0,1,1],这里alias[1]的0是默认值,也可默认置为-1避免和事件0冲突;以上部分在源码create_alias_table函数中,时间复杂度是O(N);

    随机采样1~N 之间的整数i,决定落在哪一列。随机采样0~1之间的一个概率值,如果小于accept[i],则采样i,如果大于accept[i],则采样alias[i];这一部分源码在alias_sample;

    下面附上具体的源码:

    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 = [], []
        area_ratio_ = np.array(area_ratio) * N
        for i, prob in enumerate(area_ratio_):
            if prob < 1.0:
                small.append(i)
            else:
                large.append(i)
    
        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]
    
  • 相关阅读:
    HDOJ1004
    HDOJ1001
    HDOJ1000
    HDOJ1003
    HDOJ1005
    新手如何正确使用CLion之输出hello world
    hihoCoder#1032 : 最长回文子串
    P3805 【模版】manacher算法(manacher)
    P1198 [JSOI2008]最大数(单调栈)
    P1351 联合权值
  • 原文地址:https://www.cnblogs.com/arachis/p/alias_sample.html
Copyright © 2020-2023  润新知