У меня есть функция плотности вероятности, как это:Как имитировать из (произвольного) непрерывного распределения вероятностей?
def p1(x):
return (sin(x) ** (-0.75))/(4.32141 * (x ** (1/5)))
Я хочу denerate случайного значения на [0; 1] с этим pdf. Как я могу сделать случайное значение?
Как уже упоминалось Фрэнсисом, вам лучше знать cdf вашего дистрибутива. В любом случае scipy предоставляет удобный способ определения пользовательских дистрибутивов. Это выглядит очень похоже, что
from scipy import stats
class your_distribution(stats.rv_continuous):
def _pdf(self, x):
return (sin(x) ** (-0.75))/(4.32141 * (x ** (1/5)))
distribution = your_distribution()
distribution.rvs()
(+1), хороший ответ. И если это окажется слишком медленным (так как оно проходит через общие коды кода), то только тогда стоит интегрировать/интерполировать cdf, и он обратный для случайной выборки. –
@Gioelelm, спасибо. – Denis
Без использования SciPy и дал численную выборку вашей PDF, вы можете попробовать с помощью кумулятивного распределения и линейной интерполяции. В приведенном ниже коде предполагается равное расстояние в x. Он может быть изменен для интеграции для произвольно отбираемого PDF-файла. Обратите внимание, что он перенормирует PDF в 1 в диапазоне от x.
import numpy as np
def randdist(x, pdf, nvals):
"""Produce nvals random samples from pdf(x), assuming constant spacing in x."""
# get cumulative distribution from 0 to 1
cumpdf = np.cumsum(pdf)
cumpdf *= 1/cumpdf[-1]
# input random values
randv = np.random.uniform(size=nvals)
# find where random values would go
idx1 = np.searchsorted(cumpdf, randv)
# get previous value, avoiding division by zero below
idx0 = np.where(idx1==0, 0, idx1-1)
idx1[idx0==0] = 1
# do linear interpolation in x
frac1 = (randv - cumpdf[idx0])/(cumpdf[idx1] - cumpdf[idx0])
randdist = x[idx0]*(1-frac1) + x[idx1]*frac1
return randdist
Сэмплирование из произвольного распределения может быть выполнено путем выборочного измерения равномерно в [0, 1], а затем с использованием обратной функции кумулятивной плотности. Если вы не можете вычислить этот обратный-cdf аналитически, вы можете использовать численное интегрирование своего pdf (например, используя правило трапеции) и сохранить значения в списке; то вы можете выполнить двоичный поиск вашего образца [0, 1], чтобы найти образец вашего pdf. –
@FrancisColas, это типичный способ решить мою задачу. Но действительно ли это, что я могу решить свою проблему, используя классический способ, когда есть такой мощный инструмент, как scipy? – Denis
Ну, я не знаю какого-либо конкретного способа использования numpy/scipy. Но обратите внимание, что это просто комментарий, чтобы дать вам подсказку для решения, другие могут отправлять более удовлетворительные ответы. –