生成具有给定(数字)分布的随机数

问题描述 投票:94回答:13

我有一个文件,其中包含一些不同值的概率,例如:

1 0.1
2 0.05
3 0.05
4 0.2
5 0.4
6 0.2

我想用这个发行版生成随机数。是否存在处理此问题的现有模块?自己编写代码相当简单(构建累积密度函数,生成随机值[0,1]并选择相应的值)但似乎这应该是一个常见问题,可能有人创建了一个函数/模块它。

我需要这个,因为我想生成一个生日列表(不遵循标准random模块中的任何分布)。

python module random
13个回答
88
投票

scipy.stats.rv_discrete可能是你想要的。您可以通过values参数提供概率。然后,您可以使用分布对象的rvs()方法生成随机数。

正如Eugene Pakhomov在评论中指出的那样,您也可以将p关键字参数传递给numpy.random.choice(),例如:

numpy.random.choice(numpy.arange(1, 7), p=[0.1, 0.05, 0.05, 0.2, 0.4, 0.2])

如果您使用的是Python 3.6或更高版本,则可以使用标准库中的random.choices() - 请参阅answer by Mark Dickinson


1
投票
this

验证:

distribution = [(1, 0.2), (2, 0.3), (3, 0.5)]  
# init distribution  
dlist = []  
sumchance = 0  
for value, chance in distribution:  
    sumchance += chance  
    dlist.append((value, sumchance))  
assert sumchance == 1.0 # not good assert because of float equality  

# get random value  
r = random.random()  
# for small distributions use lineair search  
if len(distribution) < 64: # don't know exact speed limit  
    for value, sumchance in dlist:  
        if r < sumchance:  
            return value  
else:  
    # else (not implemented) binary search algorithm  

1
投票

基于其他解决方案,您可以生成累积分布(整数或浮动任何您喜欢的),然后您可以使用bisect使其快速

这是一个简单的例子(我在这里使用了整数)

from __future__ import division
import random
from collections import Counter


def num_gen(num_probs):
    # calculate minimum probability to normalize
    min_prob = min(prob for num, prob in num_probs)
    lst = []
    for num, prob in num_probs:
        # keep appending num to lst, proportional to its probability in the distribution
        for _ in range(int(prob/min_prob)):
            lst.append(num)
    # all elems in lst occur proportional to their distribution probablities
    while True:
        # pick a random index from lst
        ind = random.randint(0, len(lst)-1)
        yield lst[ind]

gen = num_gen([(1, 0.1), (2, 0.05), (3, 0.05), (4, 0.2), (5, 0.4), (6, 0.2)]) lst = [] times = 10000 for _ in range(times): lst.append(next(gen)) # Verify the created distribution: for item, count in Counter(lst).iteritems(): print '%d has %f probability' % (item, count/times) 1 has 0.099737 probability 2 has 0.050022 probability 3 has 0.049996 probability 4 has 0.200154 probability 5 has 0.399791 probability 6 has 0.200300 probability 函数将它从20,60,10,10转换为20,20 + 60,20 + 60 + 10,20 + 60 + 10 + 10

现在我们使用l=[(20, 'foo'), (60, 'banana'), (10, 'monkey'), (10, 'monkey2')] def get_cdf(l): ret=[] c=0 for i in l: c+=i[0]; ret.append((c, i[1])) return ret def get_random_item(cdf): return cdf[bisect.bisect_left(cdf, (random.randint(0, cdf[-1][0]),))][1] cdf=get_cdf(l) for i in range(100): print get_random_item(cdf), 选择一个高达20 + 60 + 10 + 10的随机数,然后我们使用bisect以快速的方式获得实际值


0
投票

这些答案都不是特别清楚或简单。

这是一个明确,简单的方法,保证工作。

accumulate_normalize_probabilities采用字典get_cdf,将符号映射到概率OR频率。它输出可用于进行选择的元组列表。

random.randint

产量:

p

为什么会这样

累积步骤将每个符号变为其自身与先前符号概率或频率之间的间隔(或者在第一符号的情况下为0)。这些间隔可用于通过简单地逐步遍历列表来选择(并因此对所提供的分布进行采样),直到间隔0.0-> 1.0(先前准备的)中的随机数小于或等于当前符号的间隔终点。

规范化使我们无需确保一切都达到某种价值。归一化后,概率的“向量”总和为1.0。

用于从分布中选择和生成任意长样本的其余代码如下:

def accumulate_normalize_values(p):
        pi = p.items() if isinstance(p,dict) else p
        accum_pi = []
        accum = 0
        for i in pi:
                accum_pi.append((i[0],i[1]+accum))
                accum += i[1]
        if accum == 0:
                raise Exception( "You are about to explode the universe. Continue ? Y/N " )
        normed_a = []
        for a in accum_pi:
                normed_a.append((a[0],a[1]*1.0/accum))
        return normed_a

用法:

>>> accumulate_normalize_values( { 'a': 100, 'b' : 300, 'c' : 400, 'd' : 200  } )
[('a', 0.1), ('c', 0.5), ('b', 0.8), ('d', 1.0)]

-1
投票

这是一种更有效的方法:

只需使用'weights'数组调用以下函数(假设索引为相应的项)和no。需要的样品。可以轻松修改此功能以处理有序对。

使用各自的概率返回采样/拾取(替换)的索引(或项):

def select(symbol_intervals,random):
        print symbol_intervals,random
        i = 0
        while random > symbol_intervals[i][1]:
                i += 1
                if i >= len(symbol_intervals):
                        raise Exception( "What did you DO to that poor list?" )
        return symbol_intervals[i][0]


def gen_random(alphabet,length,probabilities=None):
        from random import random
        from itertools import repeat
        if probabilities is None:
                probabilities = dict(zip(alphabet,repeat(1.0)))
        elif len(probabilities) > 0 and isinstance(probabilities[0],(int,long,float)):
                probabilities = dict(zip(alphabet,probabilities)) #ordered
        usable_probabilities = accumulate_normalize_values(probabilities)
        gen = []
        while len(gen) < length:
                gen.append(select(usable_probabilities,random()))
        return gen

关于while循环中使用的概念的简短说明。我们从累积beta减去当前项目的权重,累积beta是随机统一构造的累积值,并且增加当前索引以便找到项目,其权重与beta的值匹配。


75
投票

从Python 3.6开始,在Python的标准库中有一个解决方案,即random.choices

用法示例:让我们设置一个与OP问题匹配的总体和权重:

>>> from random import choices
>>> population = [1, 2, 3, 4, 5, 6]
>>> weights = [0.1, 0.05, 0.05, 0.2, 0.4, 0.2]

现在choices(population, weights)生成一个样本:

>>> choices(population, weights)
4

可选的仅关键字参数k允许一次请求多个样本。这很有价值,因为在生成任何样本之前,random.choices每次调用时都必须做一些准备工作;通过一次生成许多样本,我们只需要做一次准备工作。在这里,我们生成了一百万个样本,并使用collections.Counter来检查我们得到的分布大致与我们给出的权重相匹配。

>>> million_samples = choices(population, weights, k=10**6)
>>> from collections import Counter
>>> Counter(million_samples)
Counter({5: 399616, 6: 200387, 4: 200117, 1: 99636, 3: 50219, 2: 50025})

25
投票

使用CDF生成列表的一个优点是您可以使用二进制搜索。虽然您需要O(n)时间和空间进行预处理,但您可以在O(k log n)中获得k个数字。由于普通的Python列表效率低下,因此可以使用array模块。

如果你坚持不变的空间,你可以做以下事情; O(n)时间,O(1)空间。

def random_distr(l):
    r = random.uniform(0, 1)
    s = 0
    for item, prob in l:
        s += prob
        if s >= r:
            return item
    return item  # Might occur because of floating point inaccuracies

14
投票

也许有点晚了。但你可以使用numpy.random.choice(),传递p参数:

val = numpy.random.choice(numpy.arange(1, 7), p=[0.1, 0.05, 0.05, 0.2, 0.4, 0.2])

12
投票

(好吧,我知道你要求收缩包装,但也许那些本土解决方案根本不够简洁,不符合你的喜好。:-)

pdf = [(1, 0.1), (2, 0.05), (3, 0.05), (4, 0.2), (5, 0.4), (6, 0.2)]
cdf = [(i, sum(p for j,p in pdf if j < i)) for i,_ in pdf]
R = max(i for r in [random.random()] for i,c in cdf if c <= r)

我伪确认这是通过观察这个表达式的输出来起作用的:

sorted(max(i for r in [random.random()] for i,c in cdf if c <= r)
       for _ in range(1000))

2
投票

我写了一个从自定义连续分布中抽取随机样本的解决方案。

我需要这个用于你的类似用例(即生成具有给定概率分布的随机日期)。

你只需要功能qazxsw poi和行qazxsw poi。其余的是装饰^^。

random_custDist

samples=random_custDist(x0,x1,custDist=custDist,size=1000)

这个解决方案的性能肯定是可以改进的,但我更喜欢可读性。


1
投票

你可能想看看NumPy import numpy as np #funtion def random_custDist(x0,x1,custDist,size=None, nControl=10**6): #genearte a list of size random samples, obeying the distribution custDist #suggests random samples between x0 and x1 and accepts the suggestion with probability custDist(x) #custDist noes not need to be normalized. Add this condition to increase performance. #Best performance for max_{x in [x0,x1]} custDist(x) = 1 samples=[] nLoop=0 while len(samples)<size and nLoop<nControl: x=np.random.uniform(low=x0,high=x1) prop=custDist(x) assert prop>=0 and prop<=1 if np.random.uniform(low=0,high=1) <=prop: samples += [x] nLoop+=1 return samples #call x0=2007 x1=2019 def custDist(x): if x<2010: return .3 else: return (np.exp(x-2008)-1)/(np.exp(2019-2007)-1) samples=random_custDist(x0,x1,custDist=custDist,size=1000) print(samples) #plot import matplotlib.pyplot as plt #hist bins=np.linspace(x0,x1,int(x1-x0+1)) hist=np.histogram(samples, bins )[0] hist=hist/np.sum(hist) plt.bar( (bins[:-1]+bins[1:])/2, hist, width=.96, label='sample distribution') #dist grid=np.linspace(x0,x1,100) discCustDist=np.array([custDist(x) for x in grid]) #distrete version discCustDist*=1/(grid[1]-grid[0])/np.sum(discCustDist) plt.plot(grid,discCustDist,label='custom distribustion (custDist)', color='C1', linewidth=4) #decoration plt.legend(loc=3,bbox_to_anchor=(1,0)) plt.show()


1
投票

根据他们的Continuous custom distribution and discrete sample distribution制作一个项目列表:

Random sampling distributions

优化可以是通过最大公约数对量进行归一化,以使目标列表更小。

此外,weights可能很有趣。


1
投票

另一个答案,可能更快:)

items = [1, 2, 3, 4, 5, 6]
probabilities= [0.1, 0.05, 0.05, 0.2, 0.4, 0.2]
# if the list of probs is normalized (sum(probs) == 1), omit this part
prob = sum(probabilities) # find sum of probs, to normalize them
c = (1.0)/prob # a multiplier to make a list of normalized probs
probabilities = map(lambda x: c*x, probabilities)
print probabilities

ml = max(probabilities, key=lambda x: len(str(x)) - str(x).find('.'))
ml = len(str(ml)) - str(ml).find('.') -1
amounts = [ int(x*(10**ml)) for x in probabilities]
itemsList = list()
for i in range(0, len(items)): # iterate through original items
  itemsList += items[i:i+1]*amounts[i]

# choose from itemsList randomly
print itemsList
© www.soinside.com 2019 - 2024. All rights reserved.