Генерация нормально распределённой случайной величины в Python

#python #python_3x

Я хочу для генерации нормально распределённой случайной величины на языке Python
воспользоваться методом random.normalvariate(mu, sigma) из модуля random. Однако я
не сумел найти как этот метод реализован. Я надеюсь что там используется метод Бокса-Мюллера
в одной из двух вариаций, расписанных на Википедии (потому что ещё есть методы, основанные
на ЦПТ, которые мне не нравятся). Основной вопрос: это так или нет? Где можно посмотреть
реализацию этого метода?
    

Ответы

Ответ 1

В исходниках 2.7 и в 3.5 (ссылка есть в другом ответе) написано следующее:

def normalvariate(self, mu, sigma):
    """Normal distribution.
    mu is the mean, and sigma is the standard deviation.
    """
    # mu = mean, sigma = standard deviation

    # Uses Kinderman and Monahan method. Reference: Kinderman,
    # A.J. and Monahan, J.F., "Computer generation of random
    # variables using the ratio of uniform deviates", ACM Trans
    # Math Software, 3, (1977), pp257-260.


Ответ 2

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

Поскольку мне лень копаться там, я просто скину пример с гитхаба. А там сами смотрите,
что это за метод используется.

#include 
#include 

#define NV_MAGICCONST 1.7155277699214135 /* (4 * exp(-0.5) / sqrt(2.0)); */
#define MAX_RANDOM 2147483647 /* (2 ** 31) - 1; */

/*
 * normalvariate ported from python's stdlib (random.normalvariate).
 * released under the same licence as that (python license).
 * 
 * remember to initialize the random engine before calling this function.
 */
double normalvariate(double mu, double sigma) {
    double u1, u2, z, zz;
    for (;;) {
        u1 = ((float)random()) / MAX_RANDOM;
        u2 = 1.0 - (((float)random()) / MAX_RANDOM);

        z = NV_MAGICCONST * (u1 - 0.5) / u2;
        zz = z * z / 4.0;

        if (zz <= -(log(u2))) {
            break;
        }
    }

    return mu + z * sigma;
}