Генерация случайных чисел с заданным (числовым) распределением
у меня есть файл с некоторыми вероятностями для разных значений, например:
1 0.1
2 0.05
3 0.05
4 0.2
5 0.4
6 0.2
Я хотел бы генерировать случайные числа, используя это распределение. Существует ли существующий модуль, который обрабатывает это? Это довольно просто кодировать самостоятельно (построить функцию кумулятивной плотности, сгенерировать случайное значение [0,1] и выбрать соответствующее значение), но похоже, что это должно быть общей проблемой, и, вероятно, кто-то создал функцию/модуль для нее.
мне это нужно потому что я хочу создать список дней рождения (которые не следуют за любым распределением в стандарте random
модуль).
12 ответов:
scipy.stats.rv_discrete
может быть то, что вы хотите. Вы можете предоставить свои вероятности через
начиная с 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})
преимуществом создания списка с помощью CDF является то, что вы можете использовать двоичный поиск. Хотя вам нужно O(n) время и пространство для предварительной обработки, вы можете получить k чисел в O(K log n). Поскольку обычные списки Python неэффективны, вы можете использовать
array
модуль.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
может быть, уже поздно. Но вы можете использовать
numpy.random.choice()
, передает
(хорошо, я знаю, что вы просите термоусадочную пленку, но, возможно, эти домашние решения просто не были достаточно краткими по своему вкусу. : -)
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))
возможно, вы захотите взглянуть на 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
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
основываясь на других решениях, вы генерируете накопительное распределение (как целое или плавающее, что вам нравится), затем вы можете использовать bisect, чтобы сделать его быстрым
Это простой пример (я использовал здесь целых)
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),
the
get_cdf
функция преобразует его из 20, 60, 10, 10 в 20, 20+60, 20+60+10, 20+60+10+10Теперь мы выбираем случайное число до 20+60+10+10 используя
random.randint
затем мы используем bisect, чтобы быстро получить фактическое значение
ни один из этих ответов не является особенно ясным или простым.
здесь ясный, простой метод, который гарантированно работает.
accumulate_normalize_probabilities берет словарь
p
что сопоставляет символы с вероятностями или частот. Он выводит полезный список кортежей, из которых нужно сделать выбор.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 в случае первого символа). Эти интервалы можно использовать для выбора из (и, таким образом, выборки предоставленного распределения), просто шагая по списку, пока случайное число в интервале 0.0 -> 1.0 (подготовленное ранее) не станет меньше или равно конечной точке интервала текущего символа.
The нормализация освобождает нас от необходимость убедиться, что все складывается в некоторую ценность. После нормализации "вектор" вероятностей суммируется до 1,0.
The остальная часть кода для выбора и генерации произвольно длинной выборки из распределения ниже:
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
здесь более эффективным способом для этого:
просто вызовите следующую функцию с массивом "Весов" (предполагая, что индексы являются соответствующими элементами) и нет. необходимых образцов. Эта функция может быть легко изменена для обработки упорядоченной пары.
возвращает индексы (или элементы), отобранные/выбранные (с заменой), используя их соответствующие вероятности:
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
короткая заметка о концепции, используемой в то время петля. Мы уменьшаем вес текущего элемента от кумулятивного бета, который является кумулятивным значением, построенным равномерно случайным образом, и увеличиваем текущий индекс, чтобы найти элемент, вес которого соответствует значению бета.