SciPy потребовалось очень много времени для генерации гамма-распределения в Python 3.2.

Мне нужно создать усеченную кривую pdf гамма-распределения и гистограмму в Python 3.2 на win7.

import numpy as np
import matplotlib.pyplot as plt
import scipy.special as sps

shape, scale = 2., 2. # mean and dispersion
counter =0
s = []
upper_bound = 4
lower_bound  = 0.5
while (counter <= 1000):
    t = np.random.gamma(shape, scale, 1)
    if (lower_bound <= t <= upper_bound) :   
       s.append(t)
       counter += 1 

 count, bins, ignored = plt.hist(s, 50, normed=True)

 // this part take s very long time
 y = bins**(shape-1)*(np.exp(-bins/scale) /
                  (sps.gamma(shape)*scale**shape))

 plt.plot(bins, y, linewidth=2, color='r')
 plt.show()

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

"y = bins**(shape-1)*(np.exp(-bins/scale)/(sps.gamma(shape)*scale**shape))"

Если я уберу нижнюю и верхнюю границы гамма-распределения, оно будет работать очень быстро.

Любая помощь будет оценена по достоинству.


person user3356689    schedule 13.03.2014    source источник


Ответы (1)


В строке с np.random.gamma вам не нужен аргумент size=1, вам просто нужен float, так что просто сделайте:

t = np.random.gamma(shape, scale)

На данный момент это занимает много времени, потому что каждый генерируемый вами t представляет собой массив с 1 элементом, а ваш s представляет собой большой вложенный массив.

То, что вы сделали в цикле while, уже усекает ваш дистрибутив! Хотя гораздо быстрее было бы получить желаемый диапазон впоследствии, т.е. заменить весь цикл while на:

t = np.random.gamma(shape, scale, 1000) # 1000 entries in a gamma distribution
s = filter( lambda x: upper_bound>x>lower_bound, t ) # get entries within bounds
person Higgs    schedule 14.03.2014
comment
оно работает. Почему ? И мне нужно сгенерировать усеченное гамма-распределение, не могли бы вы рассказать мне, как это сделать? Спасибо ! - person user3356689; 14.03.2014
comment
как я уже сказал, если вы включите аргумент size=1, это даст вам сказать [0.5] вместо 0.5, и ваш s станет [[0.5],[0.4],...] вместо [0.5,0.4], что вы и хотите. вы всегда можете попробовать распечатать/построить s с аргументом size=1 или без него, чтобы увидеть разницу. документацию для более подробной информации: вы скопировать код из примера? - person Higgs; 14.03.2014
comment
Я внес некоторые изменения в код из примера. Кстати, как создать усеченное гамма-распределение? Спасибо ! - person user3356689; 14.03.2014
comment
посмотрите мое редактирование для усечения, убедитесь, что вы понимаете код =) наконец, добро пожаловать в ТАК, вам просто нужно принять мой ответ, чтобы сказать спасибо - person Higgs; 14.03.2014