developer tip

주어진 (숫자) 분포로 난수 생성

optionbox 2020. 8. 11. 08:22
반응형

주어진 (숫자) 분포로 난수 생성


다른 값에 대한 확률이있는 파일이 있습니다. 예 :

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

이 분포를 사용하여 난수를 생성하고 싶습니다. 이를 처리하는 기존 모듈이 있습니까? 직접 코딩하는 것은 매우 간단합니다 (누적 밀도 함수를 빌드하고 임의의 값 [0,1]을 생성하고 해당 값을 선택). 이것은 일반적인 문제인 것처럼 보이며 아마도 누군가에 대한 함수 / 모듈을 만들었을 것입니다. 그것.

생일 목록 (표준 random모듈의 배포를 따르지 않음)을 생성하고 싶기 때문에 이것이 필요합니다 .


scipy.stats.rv_discrete당신이 원하는 것일 수 있습니다. values매개 변수 를 통해 확률을 제공 할 수 있습니다 . 그런 다음 rvs()분포 객체 방법 을 사용하여 난수를 생성 할 수 있습니다.

코멘트에 유진 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()하는 경우 표준 라이브러리에서 사용할 수 있습니다 . Mark Dickinson답변을 참조하십시오 .


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) 시간과 공간이 필요하지만 O (k log n)에서 k 개의 숫자를 얻을 수 있습니다. 일반 Python 목록은 비효율적이므로 arraymodule 을 사용할 수 있습니다 .

일정한 공간을 고집하면 다음을 수행 할 수 있습니다. 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

아마 늦었을 것입니다. 그러나 매개 변수를 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])

(좋아, 나는 당신이 수축 포장을 요구하고 있다는 것을 알고 있지만 아마도 그 자체 개발 솔루션은 당신의 취향에 충분히 간결하지 않았을 것입니다. :-)

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

사용자 지정 연속 분포에서 무작위 샘플그리는 솔루션을 작성했습니다 .

나는 당신과 비슷한 사용 사례를 위해 이것이 필요했습니다 (즉 주어진 확률 분포로 무작위 날짜 생성).

당신은 단지 기능 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()

연속 사용자 지정 분포 및 이산 표본 분포

이 솔루션의 성능은 확실히 즉흥적이지만 가독성을 선호합니다.


NumPy Random 샘플링 분포를 보고 싶을 수도 있습니다.


다음을 기반으로 항목 목록을 작성하십시오 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),

get_cdf함수 (20)에 20, 60, 10, 10에서 변환 할 + 60 20 20 + 60 + 10 20 + 60 + 10 + 10

이제 우리는 20 + 60 + 10 + 10까지 임의의 숫자를 random.randint고른 다음 bisect를 사용하여 빠른 방법으로 실제 값을 얻습니다.


이 답변 중 특별히 명확하거나 간단한 것은 없습니다.

다음은 작동이 보장되는 명확하고 간단한 방법입니다.

accumulate_normalize_probabilities takes a dictionary p that maps symbols to probabilities OR frequencies. It outputs usable list of tuples from which to do selection.

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

Yields:

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

Why it works

The accumulation step turns each symbol into an interval between itself and the previous symbols probability or frequency (or 0 in the case of the first symbol). These intervals can be used to select from (and thus sample the provided distribution) by simply stepping through the list until the random number in interval 0.0 -> 1.0 (prepared earlier) is less or equal to the current symbol's interval end-point.

The normalization releases us from the need to make sure everything sums to some value. After normalization the "vector" of probabilities sums to 1.0.

The rest of the code for selection and generating a arbitrarily long sample from the distribution is below :

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

Usage :

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

Here is a more effective way of doing this:

Just call the following function with your 'weights' array (assuming the indices as the corresponding items) and the no. of samples needed. This function can be easily modified to handle ordered pair.

Returns indexes (or items) sampled/picked (with replacement) using their respective probabilities:

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 루프에서 사용되는 개념에 대한 짧은 메모입니다. 무작위로 균일하게 구성된 누적 값 인 누적 베타에서 현재 항목의 가중치를 줄이고, 베타 값과 일치하는 항목을 찾기 위해 현재 인덱스를 증가시킵니다.

참고 URL : https://stackoverflow.com/questions/4265988/generate-random-numbers-with-a-given-numerical-distribution

반응형