Генерация случайных чисел вручную

30

Как я могу вручную сгенерировать случайное число из данного распределения, например, 10 реализаций из стандартного нормального распределения?

Дэн иш
источник
10
Можете ли вы объяснить, почему вы хотите это сделать и какие другие ограничения существуют на вас?
mdewey
2
Вы можете получить таблицу случайных чисел (ранд используется для публикации их среди других).
Бэтмен
4
@Batman: да, и книга RAND из 10 of случайных чисел получила 655 комментариев на Amazon. Все предсказуемо конечно.
Сиань
10
(Я удивлен, что никто не прокомментировал это ранее) Если это не абсолютно по существу, не следует пытаться использовать пользовательские реализации для генерации случайных чисел. Да, это здорово знать, как это сделать, и, вероятно, это первое, что вам покажут в классе « Методы МС » (и это здорово!), Но не делайте этого ни в каких реальных проектах. По разным причинам существуют специализированные процедуры генерации случайных чисел, реализация базового алгоритма с нуля - это потеря времени, источник ошибок и плохая работа по повышению осведомленности на местах.
usεr11852 говорит восстановить Monic
2
@ Сиань: Да, я согласен, что это безопасная презумпция. Как уже упоминалось, это всего лишь комментарий, чтобы предостеречь людей от использования собственного ГСБ, не осознавая, что проектирование ГСЧ является очень серьезным бизнесом. Чтобы процитировать фон Нейман : «Любой, кто считает арифметические методы получения случайных чисел является, конечно же , в состоянии греха.»
usεr11852 говорит восстановить Monic

Ответы:

46

Если «вручную» включает в себя «механический», то у вас есть много вариантов, доступных для вас. Чтобы смоделировать переменную Бернулли с половиной вероятности, мы можем бросить монету: для хвостов, 1 для голов. Чтобы смоделировать геометрическое распределение, мы можем посчитать, сколько бросков монеты необходимо, прежде чем мы получим головы. Чтобы смоделировать биномиальное распределение, мы можем бросить нашу монету n раз (или просто бросить n монет) и сосчитать головы. « В шахматном порядке» или «боб машина» или «Гальтон окно» является более кинетическую альтернативу - почему бы не установить один из них в действие и увидеть для себя ? Кажется, нет такой вещи, как «взвешенная монета»01nnно если мы хотим изменить параметр вероятности нашей Бернулли или биномиальной переменной на значения, отличные от , игла Жоржа-Луи Леклерка, графа де Буффона , позволит нам сделать это. Чтобы смоделировать дискретное равномерное распределение на { 1 , 2 , 3 , 4 , 5 , 6 }, мы бросаем шестигранный кристалл. Любители ролевых игр столкнутся с более экзотическими играми в кости , например, с тетраэдрическими игральными костями, которые можно получить из { 1 , 2 , 3 , 4 }p=0.5{1,2,3,4,5,6}{1,2,3,4}в то время как со спиннером или колесом рулетки можно пойти еще дальше. ( Изображение предоставлено )

Разнообразие игральных костей

Должны ли мы быть сумасшедшими, чтобы генерировать случайные числа таким способом сегодня, когда всего одна команда находится на компьютерной консоли, или, если у нас есть подходящая таблица случайных чисел, один набег на более темные углы книжной полки? Ну, возможно, хотя в физическом эксперименте есть что-то приятное, тактильное. Но для людей, работавших до компьютерного века, даже до широко доступных крупномасштабных таблиц случайных чисел (о которых позже), имитация случайных величин вручную имела большее практическое значение. Когда Буффон расследовал петербургский парадокс- знаменитая игра с подбрасыванием монет, в которой сумма выигрыша игрока удваивается при каждом подбрасывании головы, игрок проигрывает на первых хвостах, и чья ожидаемая выплата нелогично бесконечна - ему нужно было смоделировать геометрическое распределение с помощью . Для этого, кажется, он нанял ребенка, чтобы бросить монетку, чтобы имитировать 2048 пьес петербургской игры, записывая, сколько бросков до окончания игры. Это смоделированное геометрическое распределение воспроизведено в Стиглер (1991) :пзнак равно0,5

Tosses Frequency
1      1061
2      494
3      232
4      137
5      56
6      29
7      25
8      8
9      6

В том же эссе, где он опубликовал это эмпирическое исследование петербургского парадокса, Буффон также представил знаменитую « иглу Буффона ». Если плоскость разделена на полосы параллельными линиями на расстоянии друг от друга, и на нее опущена игла длиной l d , вероятность того, что стрелка пересекает одну из линий, равна 2 л.dLd .2Lπd

Игольчатый эксперимент Буффона

Поэтому иглу Буффона можно использовать для имитации случайной величины илиXбином(n,2лИкс~Бернулли(2Lπd), и мы можем отрегулировать вероятность успеха, изменив длину наших игл или (возможно, более удобно) расстояние, на котором мы правим линиями. Альтернативное использование игл Буффона - ужасно неэффективный способ найти вероятностное приближение дляπ. Изображение (вкредит) показывает 17 спичек, 11 из которых пересекают линию. Когда расстояние между линиями разметки установлено равным длине спички, как здесь, ожидаемая доля пересекающихся спичек равна2Икс~бином(N,2Lπd)π иследовательномы можем оценить π в два раза обратной наблюдаемой дроби: здесь мы получаем π =2172ππ^. В 1901 году Марио Лаццарини утверждал, что проводили экспериментиспользованием 2,5 см иглы с линиями 3 см другдруга, а после3408 бросков получается П =355π^знак равно217113,1 . Это хорошо известное рациональноеπ сточностью до шести десятичных знаков. Badger (1994) предоставляет убедительные доказательства того, что это было мошенничеством, и не в последнюю очередь, чтобы быть на 95% уверенным в шести десятичных разрядах точности с использованием аппарата Лаццарини, необходимо выбросить терпеливые 134 триллиона игл! Конечно, игла Буффона более полезна в качестве генератора случайных чисел, чем в качестве метода оценкиπ.π^знак равно355113ππ


До сих пор наши генераторы были удручающе дискретными. Что если мы хотим смоделировать нормальное распределение? Один из вариантов - получить случайные цифры и использовать их для формирования хороших дискретных приближений к равномерному распределению на , а затем выполнить некоторые вычисления, чтобы преобразовать их в случайные нормальные отклонения. Вращающееся колесо или колесо рулетки может давать десятичные цифры от нуля до девяти; игральные кости могут генерировать двоичные цифры; если наши арифметические навыки могут справиться с более интересной основой, подойдет даже стандартный набор кубиков. Другие ответы более подробно рассмотрели этот вид подхода, основанного на преобразованиях; Я откладываю дальнейшее обсуждение этого до конца.[0,1]

К концу девятнадцатого века полезность нормального распределения была хорошо известна, и поэтому были статистики, заинтересованные в моделировании случайных нормальных отклонений. Излишне говорить, что длительные ручные вычисления были бы неуместны, если бы сначала не был настроен процесс симуляции. Как только это было установлено, генерация случайных чисел должна была быть относительно быстрой и легкой. Стиглер (1991) перечисляет методы, использованные тремя статистиками этой эпохи. Все они исследовали методы сглаживания: случайные нормальные отклонения представляли очевидный интерес, например, для имитации ошибки измерения, которую необходимо сгладить.

Замечательный американский статистик Эрастус Лиман Де Форест интересовался сглаживанием таблиц жизни и столкнулся с проблемой, которая требовала моделирования абсолютных значений нормальных отклонений. В том, что докажет текущую тему, De Forest действительно выбрал половину нормального распределения . Кроме того, вместо того чтобы использовать стандартное отклонение один ( , мы привыкли называть «стандарт»), Де Форест хочет «вероятную ошибку» (медианное отклонение) одного. Это была форма, приведенная в таблице «Вероятность ошибок» в приложениях к «Руководству по сферической и практической астрономии, том II»Z~N(0,12)Уильям Шовен . Из этой таблицы Де Форест интерполировал квантили полунормального распределения от до p = 0,995 , которые он считал «ошибками одинаковой частоты».пзнак равно0,005пзнак равно0,995

De Forest таблица ошибок одинаковой частоты

Если вы хотите смоделировать нормальное распределение, следуя Де Форесту, вы можете распечатать эту таблицу и разрезать ее. Де Форест (1876) писал, что ошибки «были записаны на 100 битах картона одинакового размера, которые встряхивали в коробке и все вытягивали одну за другой».

Астроном и метеоролог сэр Джордж Ховард Дарвин (сын натуралиста Чарльза) изменил положение вещей, разработав то, что он назвал «рулеткой» для генерации случайных нормальных отклонений. Дарвин (1877) описывает , как:

Икс720π0Иксе-Икс2dИкс+-+-

«Индекс» следует понимать здесь как «указатель» или «индикатор» (см. «Указательный палец»). Стиглер отмечает, что Дарвин, как и Де Форест, использовал половинное нормальное кумулятивное распределение вокруг диска. Впоследствии использование монеты для произвольного прикрепления знака приводит к полному нормальному распределению. Стиглер отмечает, что неясно, насколько точно шкала была градуирована, но предполагает, что инструкция для ручного останова середины вращения диска заключалась в том, чтобы «уменьшить потенциальное смещение к одной части диска и ускорить процедуру».

Сэр Фрэнсис Гальтон , кстати, двоюродный брат Чарльза Дарвина, уже упоминался в связи с его квинкунксом. Хотя это механически имитирует биномиальное распределение, которое по теореме де Моавра-Лапласа имеет поразительное сходство с нормальным распределением (и иногда используется в качестве учебного пособия по этой теме), Гальтон фактически создал гораздо более сложную схему, когда он хотел образец из нормального распределения. Даже более необычный, чем нетрадиционные примеры в верхней части этого ответа, Гальтон разработал нормально распределенные костиИли, точнее, набор игральных костей, которые дают превосходное дискретное приближение к нормальному распределению со средним отклонением. Эти кости, датируемые 1890 годом, хранятся в коллекции Гальтона в Университетском колледже Лондона.

Гальтон нормальные кости

В 1890 году в журнале « Nature» Гальтон писал, что:

Как инструмент для случайного выбора, я не нашел ничего лучше кости. Наиболее утомительно перемешивать карты между каждым последующим розыгрышем, а метод перемешивания и перемешивания помеченных шариков в сумке все еще более утомителен. Для них предпочтительнее тройник или какая-то разновидность рулетки, но кости лучше, чем все. Когда их трясут и бросают в корзину, они так по-разному бьются друг о друга и о ребра плетеной корзины, что они дико кувыркаются, и их позиции с самого начала не дают ощутимого ключа к тому, чем они станут после даже один хорошо встряхнуть и бросить. Шансы, предоставляемые штампом, более разнообразны, чем принято считать; есть 24 равных возможности, а не только 6, потому что каждая грань имеет четыре ребра, которые можно использовать, как я покажу.

+-114

Гальтон нормальный дизайн игры в кости

Лаборатория математических статистических экспериментов Раазеша Сайнудина включает студенческий проект из Кентерберийского университета, Новая Зеландия, по воспроизведению кости Гальтона . Проект включает в себя эмпирическое исследование по бросанию костей много раз (включая эмпирический CDF, который выглядит обнадеживающе «нормальным») и адаптацию оценок костей, чтобы они следовали стандартному нормальному распределению. Используя исходные оценки Гальтона, есть также график дискретного нормального распределения, которому фактически следуют оценки костей.

Дискретное распределение костей Galton


В целом, если вы готовы растянуть «механическое» на электрическое, обратите внимание, что эпическое «Миллион случайных цифр» RAND с 100 000 нормальных девиатов RAND было основано на своего рода электронном моделировании колеса рулетки. Из технического отчета (Джордж У. Браун, первоначально июнь 1949 г.) мы находим:

Таким образом, мотивированные сотрудники RAND при содействии инженерного персонала Douglas Aircraft Company разработали колесо электрической рулетки на основе варианта, предложенного Сесилом Хастингсом. Для целей этого разговора будет достаточно краткого описания. Источник импульсов случайной частоты управлялся импульсом постоянной частоты, примерно раз в секунду, обеспечивая в среднем около 100 000 импульсов в секунду. Схемы стандартизации импульсов передавали импульсы в пятизначный двоичный счетчик, так что в принципе машина похожа на колесо рулетки с 32 положениями, совершая в среднем около 3000 оборотов за каждый оборот. Использовалось преобразование двоичной системы в десятичную, отбрасывая 12 из 32 позиций, и полученная в результате случайная цифра подавалась в перфорацию IBM, давая таблицы перфокарт случайных цифр.

χ2Проверка частот нечетных и четных цифр показала, что некоторые партии имели небольшой дисбаланс. В некоторых партиях это было хуже, чем в других, что говорит о том, что «машина не работала в течение месяца с момента ее настройки ... Признаки того, что эта машина требовала чрезмерного обслуживания, чтобы поддерживать ее в идеальном состоянии». Однако был найден статистический способ решения этих проблем:

В этот момент у нас были наши исходные миллионные цифры, 20 000 карт IBM с 50 цифрами на карту, с небольшим, но ощутимым нечетно-четным смещением, выявленным статистическим анализом. Теперь было решено переупорядочить таблицу или, по крайней мере, изменить ее, сыграв с ней небольшую рулетку, чтобы убрать нечетно-четное смещение. Мы добавили (мод 10) цифры в каждой карте, цифра за цифрой, к соответствующим цифрам предыдущей карты. Полученная таблица из миллиона цифр была затем подвергнута различным стандартным тестам, частотным тестам, серийным тестам, покерным тестам и т. Д. Эти миллионы цифр имеют чистый счет здоровья и были приняты в качестве современной таблицы случайных цифр RAND.

Конечно, имелись веские основания полагать, что процесс добавления принесет пользу. В общем случае основным механизмом является предельный подход сумм случайных величин по модулю единичного интервала в прямоугольном распределении таким же образом, как неограниченные суммы случайных величин приближаются к нормальности. Этот метод использовался Хортоном и Смитом из Межгосударственной торговой комиссии для получения некоторых хороших партий явно случайных чисел из больших партий плохо неслучайных чисел.

[0,1]U[0,1]FF-1(U)

Зиккурат для полунормальных

Ссылки


(*)«Любой, кто рассматривает арифметические методы получения случайных цифр, конечно же, находится в состоянии греха. Поскольку, как уже несколько раз указывалось, не существует такого понятия, как случайное число - существуют только методы для получения случайных чисел. и строгая арифметическая процедура, конечно, не такой метод ".

тарпон
источник
5
Вы, должно быть, шутите ... (но, да, это действительно вручную +1)
nalzok
8
Сколько работы вложено в это!
Ричард Харди
3
Обратите внимание, что нестандартные кости могут быть несправедливыми, поэтому сначала хорошо их протестировать, например, youtube.com/watch?v=VI3N4Qg-JZM
Тим
2
@RichardHardy Как это ни парадоксально, на самом деле мне проще записать все эти вещи, пока они передо мной и хотя бы временно остались в моей памяти, чем пытаться нести все это в моей голове!
Серебряная рыба
2
В любом случае, я нахожу это впечатляющим!
Ричард Харди
44

Иксзнак равно-2журналU1соз(2πU2)
Yзнак равно-2журналU1грех(2πU2)Икс

Например, на моей ОС Linux я могу проверить

$ date +%s.%N
1479733744.077762986
$ date +%s.%N
1479733980.615056616

U1знак равно0,077762986, U2знак равно0,615056616
Икс
> sqrt(-2*log(.077762986))*cos(2*pi*.615056616)
[1] -1.694815

Приложение: поскольку вычисление логарифмов и косинусов может считаться недостаточно ручным , существует вариант Бокса-Мюллера, в котором не используются эти трансцендентные функции (см. Упражнение 2.9 в нашей книге « Статистические методы Монте-Карло» ):

вариант Бокса-Мюллера

Теперь можно поспорить с этой версией из-за экспоненциальных вариаций. Но также существует очень умный способ моделирования этих вариаций без вызова трансцендентных функций , благодаря фон Нейману, как показано в этом алгоритме, воспроизведенном из неоднородных случайных вариантов Люка Деврой :

введите описание изображения здесь

По общему признанию, это требует вычисления 1 / e, но только один раз.

(0,1)2(0,1)(0,2π)

Использование CLT для аппроксимации нормальности, безусловно, не тот метод, который я бы посоветовал, поскольку (1) вам все еще нужны другие переменные, чтобы подавать среднее значение, поэтому можно также использовать униформу в алгоритме Бокса-Мюллера, и (2) точность растет довольно медленно с количеством симуляций. Особенно если использовать дискретную случайную величину, например, результат игры в кости, даже с более чем шестью гранями . Цитировать из Thomas et al. (2007), обзор плюсов и минусов гауссовских генераторов случайных чисел:

Центральная предельная теорема, конечно, является примером «приближенного» метода - даже если используется идеальная арифметика, для конечного K результат не будет гауссовым.

Вот быстрый эксперимент, чтобы проиллюстрировать проблему: я получил в 100 раз среднее значение 30 результатов:

dies=apply(matrix(sample(1:6,30*100,rep=TRUE),ncol=30),1,mean)

затем нормализовали эти средние значения до среднего нуля - дисперсия один меняется

stdies=(dies-3.5)/sqrt(35/12/30)

и посмотрел на нормальное соответствие [или его отсутствие] этого образца:

qqnorm (stdies, Col = "Gold2", PCH = 19); abline (а = 0, Ь = 1, Col = "SteelBlue", LWD = 2, LTY = 2)

dies76/30122/30Dя

Uзнак равноΣязнак равно1КDя-16я
Кзнак равно15
dies=matrix(apply(matrix(sample(0:5,15*200,rep=TRUE),nrow=15)/6^(1:15),2,sum),ncol=2) 
norma=sqrt(-2*log(dies[,1]))*c(cos(2*pi*dies[,2]),sin(2*pi*dies[,2]))

аппроксимация настолько хороша, насколько можно ожидать для нормальной выборки размером 200 (просто нанесите еще одну для истинной нормальной выборки norma=rnorm(100)):

введите описание изображения здесь

как далее показано с помощью теста Колмогорова-Смирнова:

> ks.test(norma,pnorm)

        One-sample Kolmogorov-Smirnov test

data:  norma
D = 0.06439, p-value = 0.3783
alternative hypothesis: two-sided
Сиань
источник
3

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

Используйте свой телефон для настройки хронометра. После хороших 10 секунд остановите его (чем дольше вы ждете, тем больше вы приближаетесь к действительно «случайному» результату, но 10 секунд - это нормально). Возьмите последние цифры (например, 10,67 сек даст вам 67). Примените таблицу процентилей для нормального распределения. В этом примере вам просто нужно найти 0,67, и вы найдете номер. В этом случае ваше значение составляет около 0,45. Это не совсем точно, но это даст вам надежную оценку.

Если вы опуститесь ниже 50, просто сделайте 100- [Ваш результат] и используйте таблицу. Ваш результат будет таким же, со знаком минус, из-за симметрии N (0,1).

Thalantas
источник
3

N+1 если головы, -1если хвосты. ПослеN подбрасывание монет, мы делим счетчик на N, Используя центральную предельную теорему , еслиN достаточно велик, то мы должны иметь "приближенную реализацию" нормализованного гауссова N(0,1),


Зачем? Позволять

ИксКзнак равно{+1 если Кмонета броска головы-1 если Кбросок монеты хвостами

быть случайными переменными Бернулли с п(ИксКзнак равно±1)знак равно12, Следовательно,

Е(ИксК)знак равно0Var(ИксК)знак равно1

Позволять Yзнак равноИкс1+Икс2++ИксN, Следовательно,

Е(Y)знак равно0Var(Y)знак равноN

Нормализация,

Zзнак равноYN

мы получаем случайную величину с единичной дисперсией

Е(Z)знак равно0Var(Z)знак равно1
Родриго де Азеведо
источник
1
Как обсуждалось в моем ответе, я боюсь, что CLT не является эффективным использованием случайности в бросках монеты: их лучше использовать, так как двоичные цифры псевдослучайного U (0,1) изменяются с простой или двойной точностью.
Сиань
@ Сиань Я прочитал твой ответ, прежде чем опубликовать свой, и твои возражения против CLT, похоже, основывались на медленной конвергенции. Поскольку это мысленный эксперимент, подбрасывание монеты 10 миллиардов раз ничего не стоит. И это действительно ручная процедура, не требующая компьютеров, вычислений логарифмов, квадратных корней или косинусов. Конечно, можно использовать правило слайдов, но это может зайти слишком далеко.
Родриго де Азеведо
1
:}}: 10000000000 монета переворачивается не звучит очень ручной для меня ...!
Сиань
Ручной = вручную. Вычисление логарифмов и косинусов вручную также требует времени.
Родриго де Азеведо
0

Стоит отметить, что как только вы можете сгенерировать униформу (0,1), вы можете сгенерировать любую случайную переменную, для которой обратный cdf рассчитывается, просто подключив равномерную случайную переменную к обратному CDF.

Итак, как можно вычислить форму (0,1) вручную? Ну, как упомянул @Silverfish, есть множество костей, используемых традиционными игроками RPG. Один из них - десятигранный кубик. Предполагая, что это справедливый кубик, теперь мы можем сгенерировать дискретную униформу (0, 9).

Мы также можем использовать эту форму (0,9), чтобы представить одну цифру случайной величины. Поэтому, если мы используем две кости, мы получим равномерную случайную величину, которая может принимать значения0,01,0.02,,,,,0,99,1,00, С тремя кубиками мы можем получить равномерное распределение0,001,0,002,,,,,0,999,+1,000,

Таким образом, мы можем очень близко подойти к непрерывной униформе (0,1), аппроксимировав ее точным равномерно распределенным дискретным распределением с несколькими 10-гранными кубиками. Затем его можно подключить к обратному CDF для получения интересующей случайной величины.

Клифф AB
источник