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

我有一个文件,有不同的值的概率,例如:

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

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

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

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]) 

使用CDF生成列表的优点是可以使用二分search。 当需要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 

(好吧,我知道你们正在寻求收缩包装,但是也许这些自制的解决scheme并不足以满足你的喜好。:-)

 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) 

我通过观察这个expression式的输出来伪证实了这一点:

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

也许这是晚了。 但是你可以使用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]) 

自Python 3.6以来,在Python的标准库中有一个解决scheme,即random.choices

示例用法:让我们build立一个匹配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}) 

你可能想看看NumPy 随机抽样分布

根据weights制定项目清单:

 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 

优化可以是通过最大公约数来标准化量,以使目标列表更小。

另外, 这可能是有趣的。

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

 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 

基于其他解决scheme,你可以产生累积分布(如你喜欢的整数或浮点数),那么你可以使用平分来快速

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

 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), 

get_cdf函数将它从get_cdf转换为20,20 + get_cdf + 60 + 10,20 + 60 + 10 + 10

现在我们用random.randintselect一个20 + 60 + 10 + 10的随机数,然后我们用对分来快速得到实际值

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

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

accumulate_normalize_probabilities需要一个将符号映射到概率频率的字典p 。 它输出可用于select的元组列表的可用列表。

 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)] 

为什么它的作品

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

规范化使我们摆脱了确保一切都达到一定价值的需要。 归一化之后,概率的“向量”总和为1.0。

下面的代码供select和生成任意长度的样本:

 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 

用法:

 >>> gen_random (['a','b','c','d'],10,[100,300,400,200]) ['d', 'b', 'b', 'a', 'c', 'c', 'b', 'c', 'c', 'c'] #<--- some of the time 
 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] 

validation:

 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 

这是一个更有效的方法

只需要调用你的'权重'数组(假设索引作为相应的项目)和没有。 需要的样品。 这个function可以很容易地修改来处理有序对。

返回索引(或项目)采样/采摘(与replace)使用各自的概率:

 def resample(weights, n): beta = 0 # Caveat: Assign max weight to max*2 for best results max_w = max(weights)*2 # Pick an item uniformly at random, to start with current_item = random.randint(0,n-1) result = [] for i in range(n): beta += random.uniform(0,max_w) while weights[current_item] < beta: beta -= weights[current_item] current_item = (current_item + 1) % n # cyclic else: result.append(current_item) return result 

关于while循环中使用的概念的简短说明。 我们从累积贝塔(这是一个随机均匀构造的累积值)中减去当前物品的重量,并增加当前索引以find与贝塔值相匹配的物品。