Генерация случайных чисел после нормального распределения в C / C++



кто-нибудь знает, как я могу легко генерировать случайные числа после нормального распределения в C/C++ ?



http://www.mathworks.com/access/helpdesk/help/toolbox/stats/normrnd.html



Я не хочу использовать какой-либо из Boost.



Я знаю, что кнут говорит об этом, но у меня нет книги под рукой сейчас.

944   17  

17 ответов:

The Box-Muller преобразование-это то, что обычно используется. Это правильно создает значения с нормальным распределением.

http://en.wikipedia.org/wiki/Normal_distribution#Generating_values_from_normal_distribution

http://en.wikipedia.org/wiki/Box_Muller_transform

математика-это легко. Вы генерируете два однородных числа, и из них вы получаете два нормально распределенных числа. Верните один, сохраните другой для следующего запроса случайного числа.

EDIT: С 12 августа 2011 года C++11, который напрямую предлагает std::normal_distribution, именно так я бы пошел сегодня.

вот ответ:

вот некоторые решения, упорядоченные по возрастанию сложности.

  1. добавить 12 однородных чисел от 0 до 1 и вычесть 6. Это будет соответствовать среднему и стандартному отклонению нормальной переменной. Очевидным недостатком является что диапазон ограничен до +/-6 - в отличие от нормального распределения.

  2. Box - Muller transform - был перечислен выше, и относительно прост в реализации. Однако если вам нужны очень точные образцы, имейте в виду, что преобразование бокса-Мюллера в сочетании с некоторыми однородными генераторами страдает от аномалии, называемой эффектом Нива.

    Нив, " об использовании преобразования бокса-Мюллера с мультипликативным конгруэнтным псевдослучайным числом генераторы, " Прикладная статистика, 22, 92-97, 1973

  3. на наилучшая точность предлагаю чертеж формы и применение обратного кумулятивного нормального распределения для получения нормально распределенных вариаций. Вы можете найти очень хороший алгоритм для обратного кумулятивного нормального распределения в

https://web.archive.org/web/20151030215612/http://home.online.no/~pjacklam/notes/invnorm/

надеюсь, что это поможет

Петр

быстрый и простой метод-это просто суммировать количество равномерно распределенных случайных чисел и взять их среднее значение. Смотрите Центральная Предельная Теорема для полного объяснения, почему это работает.

Я создал C++ open source project для нормально распределенного теста генерации случайных чисел.

Он сравнивает несколько алгоритмов, в том числе

  • метод центральной предельной теоремы
  • Box - Muller transform
  • Марсалья полярный способ
  • зиккурат
  • метод обратного преобразования выборки.
  • cpp11random использует C++11 std::normal_distribution С std::minstd_rand (это на самом деле бокс-Мюллер трансформировать в лязг).

результаты одиночной точности (float) версия на iMac [email protected] , clang 6.1, 64-бит:

normaldistf

для корректности программа проверяет среднее, стандартное отклонение, асимметрию и эксцесс образцов. Было обнаружено, что метод CLT путем суммирования 4, 8 или 16 однородных чисел не имеет хорошего эксцесса, как и другие методы.

алгоритм зиккурата имеет лучшую производительность, чем другие. Однако он не подходит для параллелизма SIMD, поскольку ему нужен поиск таблицы и ветви. Бокса-Мюллера с поддержкой SSE2/с AVX набор инструкций намного быстрее (Х1.79, Х2.99) чем non-Симд версии алгоритм зиггурата.

поэтому я предложу использовать Box-Muller для архитектуры с наборами команд SIMD и может быть зиккуратом в противном случае.


С. П. бенчмарк использует простой ЖКИ ГПСЧ для генерации равномерно распределенной случайной величины. Поэтому этого может быть недостаточно для некоторые приложения. Но сравнение производительности должно быть справедливым, потому что все реализации используют один и тот же PRNG, поэтому тест в основном проверяет производительность преобразования.

вот пример C++, основанный на некоторых ссылках. Это быстро и грязно,вам лучше не изобретать заново и использовать библиотеку boost.

#include "math.h" // for RAND, and rand
double sampleNormal() {
    double u = ((double) rand() / (RAND_MAX)) * 2 - 1;
    double v = ((double) rand() / (RAND_MAX)) * 2 - 1;
    double r = u * u + v * v;
    if (r == 0 || r > 1) return sampleNormal();
    double c = sqrt(-2 * log(r) / r);
    return u * c;
}

вы можете использовать график Q-Q для изучения результатов и посмотреть, насколько хорошо он аппроксимирует реальное нормальное распределение (ранжируйте ваши образцы 1..X, в свою очередь, ряды в пропорции общее число х, т. е.. сколько образцов, получить z-значения и построить их. Желаемый результат-прямая линия вверх).

использовать std::tr1::normal_distribution.

пространство имен std::tr1 не является частью boost. Это пространство имен, которое содержит дополнения библиотеки из технического отчета C++ 1 и доступно в современных компиляторах Microsoft и gcc, независимо от boost.

вот как вы создаете образцы на современном компиляторе C++.

#include <random>
...
std::mt19937 generator;
double mean = 0.0;
double stddev  = 1.0;
std::normal_distribution<double> normal(mean, stddev);
cerr << "Normal: " << normal(generator) << endl;

можно использовать GSL. Некоторые приведены полные примеры чтобы продемонстрировать, как использовать его.

посмотрите на: http://www.cplusplus.com/reference/random/normal_distribution/. это самый простой способ получения нормальных распределений.

Если вы используете C++11, вы можете использовать std::normal_distribution:

#include <random>

std::default_random_engine generator;
std::normal_distribution<double> distribution(/*mean=*/0.0, /*stddev=*/1.0);

double randomNumber = distribution(generator);

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

я следовал определению PDF, приведенному в http://www.mathworks.com/help/stats/normal-distribution.html и придумал такое:

const double DBL_EPS_COMP = 1 - DBL_EPSILON; // DBL_EPSILON is defined in <limits.h>.
inline double RandU() {
    return DBL_EPSILON + ((double) rand()/RAND_MAX);
}
inline double RandN2(double mu, double sigma) {
    return mu + (rand()%2 ? -1.0 : 1.0)*sigma*pow(-log(DBL_EPS_COMP*RandU()), 0.5);
}
inline double RandN() {
    return RandN2(0, 1.0);
}

возможно, это не лучший подход, но это довольно просто.

реализация Box-Muller:

#include <cstdlib>
#include <cmath>
#include <ctime>
#include <iostream>
using namespace std;
 // return a uniformly distributed random number
double RandomGenerator()
{
  return ( (double)(rand()) + 1. )/( (double)(RAND_MAX) + 1. );
}
 // return a normally distributed random number
double normalRandom()
{
  double y1=RandomGenerator();
  double y2=RandomGenerator();
  return cos(2*3.14*y2)*sqrt(-2.*log(y1));
}

int main(){
double sigma = 82.;
double Mi = 40.;
  for(int i=0;i<100;i++){
double x = normalRandom()*sigma+Mi;
    cout << " x = " << x << endl;
  }
  return 0;
}

1) графически интуитивно понятный способ генерации гауссовых случайных чисел заключается в использовании чего-то похожего на метод Монте-Карло. Вы могли бы создать случайную точку в поле вокруг гауссовой кривой, используя генератор псевдослучайных чисел в C. Вы можете вычислить, находится ли эта точка внутри или под гауссовым распределением, используя уравнение распределения. Если эта точка находится внутри гауссовского распределения, то вы получили свое гауссовское случайное число как значение x точка.

этот метод не идеален, потому что технически гауссова кривая идет к бесконечности, и вы не можете создать поле, которое приближается к бесконечности в X-измерении. Но кривая Guassian подходов 0 в y-измерении довольно быстро, так что я не волновался бы об этом. Ограничение размера ваших переменных в C может быть более ограничивающим фактором для вашей точности.

2) Другой способ - использовать центральную предельную теорему, которая гласит, что когда добавляются независимые случайные величины, которые образуют нормальное распределение. Имея в виду эту теорему, вы можете аппроксимировать гауссово случайное число, добавив большое количество независимых случайных величин.

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

взгляните на то, что я нашел.

этой библиотека использует алгоритм Зиггурата.

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

вы можете взглянуть на него:http://c-faq.com/lib/gaussian.html

существуют различные алгоритмы для обратного кумулятивного нормального распределения. Самые популярные в количественном финансировании тестируются на http://chasethedevil.github.io/post/monte-carlo--inverse-cumulative-normal-distribution/

кроме того, он показывает недостаток зиккурат, как подходы.

компьютер является детерминированным устройством. В расчетах нет случайности. Кроме того, арифметическое устройство в ЦП может оценивать суммирование по некоторому конечному набору целых чисел (выполняя оценку в конечном поле) и конечному набору вещественных рациональных чисел. А также выполнял побитовые операции. Математика имеет дело с более большими множествами, такими как [0.0, 1.0] с бесконечным числом точек.

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

существуют алгоритмы, называемые-псевдослучайный генератор. Как я чувствовал, цель псевдослучайного генератора-эмулировать случайность. И критерии goodnes-это: - эмпирическое распределение сходится (в некотором смысле - точечный, равномерный, L2) к теоретическому - значения, которые вы получаете от случайного генератора, кажутся idependent. Конечно, это не так с "реальной точки зрения", но мы предполагаем, что это правда.

один из популярных методов-вы можете суммировать 12 i.r.v с равномерными распределениями....Но, честно говоря, при выводе центральной предельной теоремы с помощью преобразования Фурье, ряда Тейлора, необходимо иметь N->+INF предположения пару раз. так, например, теоретизации - Лично я не понимаю, как люди выполняют суммирование 12 i.r.v. с равномерным распределением.

У меня была теория вероятностей в университете. И особенно для меня это просто математический вопрос. В университете я видел следующую модель:


double generateUniform(double a, double b)
{
  return uniformGen.generateReal(a, b);
}

double generateRelei(double sigma)
{
  return sigma * sqrt(-2 * log(1.0 - uniformGen.generateReal(0.0, 1.0 -kEps)));
}
double generateNorm(double m, double sigma)
{
  double y2 = generateUniform(0.0, 2 * kPi);
  double y1 = generateRelei(1.0);
  double x1 = y1 * cos(y2);
  return sigma*x1 + m;
}

таким образом, как todo это был просто пример, я думаю, что существуют другие способы его реализации.

подтверждение, что это правильно, можно найти в этой книге "Москва, МГТУ, 2004: вероятность XVI в. Теория, пример 6.12, стр. 246-247 " of Крищенко Александр Петрович ISBN 5-7038-2485-0

к сожалению, я не знаю о существовании перевода этой книги на английский язык.

Comments

    Ничего не найдено.