Генератор случайных чисел, который производит степенное распределение?

Я пишу несколько тестов для приложения Linux для командной строки C++. Я хотел бы сгенерировать кучу целых чисел с степенным распределением / длинным хвостом. То есть я получаю некоторые цифры очень часто, но большинство из них относительно редко.

В идеале это были бы магические уравнения, которые я мог бы использовать с rand() или одной из случайных функций stdlib. Если нет, то простой в использовании кусок C/C++ был бы великолепен.

Спасибо!

Источник

4 ответа

Решение

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

Краткий ответ (вывод по вышеуказанной ссылке):

x = [(x1^(n+1) - x0^(n+1))*y + x0^(n+1)]^(1/(n+1))

где y - равномерная переменная, n - степень распределения, x0 и x1 определяют диапазон распределения, а x - ваша степенная распределенная переменная.

Источник

Если вы знаете, какое распределение вы хотите (называемое функцией распределения вероятностей (PDF)) и правильно ли оно нормализовано, вы можете интегрировать его для получения функции накопительного распределения (CDF), а затем инвертировать CDF (если это возможно), чтобы получить преобразование, которое вы нужно из униформы [0,1] раздача по вашему желанию.

Итак, вы начинаете с определения дистрибутива, который вы хотите.

P = F(x)

(для х в [0,1]) затем интегрируется, чтобы дать

C(y) = \int_0^y F(x) dx

Если это можно перевернуть, вы получите

y = F^{-1}(C)

Так что звоните rand() и включите результат в виде C в последней строке и используйте y.

Этот результат называется фундаментальной теоремой выборки. Это хлопотно из-за требования нормализации и необходимости аналитически инвертировать функцию.

В качестве альтернативы вы можете использовать технику отклонения: бросить число равномерно в желаемом диапазоне, затем бросить другое число и сравнить с PDF в месте, не определенном вашим первым броском. Отклонить, если второй бросок превышает PDF. Как правило, неэффективно для PDF-файлов с большим количеством областей с низкой вероятностью, например, с длинными хвостами...

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

Настоящий вонючий здесь так просто x^-n распределения не нормируются в диапазоне [0,1], поэтому вы не можете использовать теорему выборки. Попробуйте (x+1)^-n вместо...

Источник

Я просто хотел провести реальное моделирование как дополнение к (справедливо) принятому ответу. Хотя в R этот код настолько прост, что представляет собой (псевдо)-псевдокод.

Одно крошечное различие между формулой Wolfram MathWorld в принятом ответе и другими, возможно, более общими уравнениями заключается в том, что показатель степенного закона n (который обычно обозначается как альфа) не имеет явного отрицательного знака. Таким образом, выбранное альфа-значение должно быть отрицательным, обычно от 2 до 3.

x0 а также x1 обозначают нижний и верхний пределы распределения.

Итак, вот оно:

x1 = 5           # Maximum value
x0 = 0.1         # It can't be zero; otherwise X^0^(neg) is 1/0.
alpha = -2.5     # It has to be negative.
y = runif(1e5)   # Number of samples
x = ((x1^(alpha+1) - x0^(alpha+1))*y + x0^(alpha+1))^(1/(alpha+1))
hist(x, prob = T, breaks=40, ylim=c(0,10), xlim=c(0,1.2), border=F, 
col="yellowgreen", main="Power law density")
lines(density(x), col="chocolate", lwd=1)
lines(density(x, adjust=2), lty="dotted", col="darkblue", lwd=2)

или построены в логарифмическом масштабе:

h = hist(x, prob=T, breaks=40, plot=F)
     plot(h$count, log="xy", type='l', lwd=1, lend=2, 
     xlab="", ylab="", main="Density in logarithmic scale")

Вот краткое изложение данных:

> summary(x)
   Min.   1st Qu.  Median    Mean   3rd Qu.    Max. 
  0.1000  0.1208  0.1584    0.2590  0.2511   4.9388
Источник

Я не могу комментировать математику, необходимую для получения степенного распределения (другие посты содержат предложения), но я бы посоветовал вам ознакомиться со средствами случайных чисел в стандартной библиотеке TR1 C++ в <random>, Они обеспечивают больше функциональности, чем std::rand а также std::srand, Новая система задает модульный API для генераторов, движков и дистрибутивов и предоставляет набор пресетов.

Включены предустановки распространения:

  • uniform_int
  • bernoulli_distribution
  • geometric_distribution
  • poisson_distribution
  • binomial_distribution
  • uniform_real
  • exponential_distribution
  • normal_distribution
  • gamma_distribution

Когда вы определите свое распределение по степенным законам, вы сможете подключить его к существующим генераторам и двигателям. Книга Питера Беккера " Расширения стандартной библиотеки C++" имеет большую главу <random>,

Вот статья о том, как создавать другие дистрибутивы (с примерами для Коши, Chi-squared, Student t и Snedecor F)

Источник
Другие вопросы по тегам