Генерировать случайные числа с заданным (численным) распределением

Ответ 1

scipy.stats.rv_discrete может быть тем, что вы хотите. Вы можете указать свои вероятности с помощью параметра values. Затем вы можете использовать метод rvs() объекта распределения для генерации случайных чисел.

Как отметил Евгений Пахомов в комментариях, вы также можете передать параметр ключевого слова 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() из стандартной библиотеки - см. Ответ Марка Дикинсона.

Ответ 2

Начиная с 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})

Ответ 3

Преимущество создания списка с помощью CDF заключается в том, что вы можете использовать двоичный поиск. Хотя вам нужно O (n) время и пространство для предварительной обработки, вы можете получить k чисел в O (k log n). Поскольку обычные списки 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

Ответ 4

Может быть, уже поздно. Но вы можете использовать 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])

Ответ 5

(Хорошо, я знаю, что вы просите термоусадочную пленку, но, возможно, эти домашние решения просто не были достаточно краткими для вашей симпатии.: -)

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

Ответ 6

вы можете посмотреть на NumPy Случайные распределения выборки

Ответ 7

Составьте список элементов на основе их 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

Оптимизация может заключаться в нормализации сумм наибольшим общим делителем, чтобы сделать список целей меньшим.

Кроме того, этот может быть интересным.

Ответ 8

Другой ответ, возможно, быстрее:)

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  

Ответ 9

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

Ответ 10

основываясь на других решениях, вы генерируете накопительное распределение (как целое или плавающее, как вам нравится), тогда вы можете использовать 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),

функция get_cdf преобразует ее из 20, 60, 10, 10 в 20, 20 + 60, 20 + 60 + 10, 20 + 60 + 10 + 10

теперь мы выбираем случайное число до 20 + 60 + 10 + 10, используя random.randint, затем мы используем bisect для быстрого получения фактического значения

Ответ 11

Ни один из этих ответов не является особенно ясным или простым.

Вот простой и понятный метод, который гарантированно работает.

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 (подготовленное ранее) не будет меньше или равно текущей конечной точке интервала символов.

нормализация освобождает нас от необходимости убедиться, что все суммируется до некоторого значения. После нормализации "вектор" вероятностей суммируется до 1,0.

Остальная часть кода для выбора и создания произвольно длинного образца из дистрибутива ниже:

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

Ответ 12

Я написал решение для рисования случайных выборок из пользовательского непрерывного распределения.

Мне это нужно для аналогичного случая использования (например, генерация случайных дат с заданным распределением вероятности).

Вам просто требуется функция random_custDist и строка samples=random_custDist(x0,x1,custDist=custDist,size=1000). Остальное украшение ^^.

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

Производительность этого решения наверняка улучшена, но я предпочитаю удобочитаемость.

Ответ 13

Вот более эффективный способ:

Просто вызовите следующую функцию с массивом 'weightights' (при условии, что индексы соответствуют соответствующим элементам) и no. необходимых образцов. Эта функция может быть легко модифицирована для обработки упорядоченной пары.

Возвращает индексы (или элементы), отобранные/выбранные (с заменой), используя их соответствующие вероятности:

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. Мы уменьшаем текущий вес элемента от кумулятивной бета-версии, которая является совокупным значением, построенным равномерно случайным образом, и увеличиваем текущий индекс, чтобы найти элемент, вес которого соответствует значению бета-версии.