Модуль random. Часть 1 | Программирование на Питоне

Модуль random. Часть 1

Все курсы > Программирование на Питоне > Занятие 10 (часть 1)

Библиотека Numpy предоставляет широкие возможности по работе со случайными числами. Для этого в ней предусмотрен модуль random.

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

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

Содержание занятия

Случайная величина

Случайная величина (random variable) — это величина, которая в зависимости от случая может принимать определенное значение с заданной вероятностью. Разберем это определение более подробно и параллельно введем несколько новых терминов.

Предположим, что у нас есть некоторое событие (event), например, выпадение двойки или тройки на игральной кости. Это событие случайно (random process), если его невозможно предугадать.

Тогда испытанием (trial) мы назовем процесс бросания игральной кости, а исходами (outcomes) испытания — возможные результаты, то есть выпадение чисел от 1 до 6.

Вероятность случайной величины

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

Классическая (теоретическая) вероятность

Пусть $k$ — количество благоприятствующих событию A исходов (то есть исходов, которые влекут за собой наступление события $A$), а $n$ — общее количество равновозможных исходов, тогда теоретической вероятностью события $A$ будет

$$ P(A) = \frac{k}{n} $$

Разберем это определение на том же примере бросания игральной кости.

Выпадение двойки или тройки благоприятствует (приводит к наступлению) случайного события ($A$). Таких исходов два (либо двойка, либо тройка) и соответственно $k = 2$. Всего исходов шесть и вероятность их наступления одинакова (они равновозможны). Значит $n = 6$. Рассчитаем теоретическую вероятность такого события.

$$ P(A) = \frac{2}{6} = \frac{1}{3} $$

Аналогично, вероятность выпадения орла или решки на «честной» монете (fair coin) равна 1/2, а вероятность вытащить червы из колоды, в которой 36 карт, равна 9/36 или 1/4.

Вероятность любого события находится в диапазоне от 0 до 1. При этом событие с нулевой вероятностью называется невозможным, а с вероятностью равной единице — достоверным. В промежутке находятся случайные события.

Сумма вероятностей событий, образующих полную группу, всегда равна единице. Выпадение орла и выпадение решки образуют полную группу (других событий при подбрасывании монеты быть не может). Вероятность каждого из этих событий равна 1/2, сумма этих вероятностей — единице.

Относительная (эмпирическая) вероятность

Эмпирическая вероятность события — это отношение частоты наступления события А к общему числу испытаний. По большому счету, мы бросаем кости несколько раз и считаем, сколько раз выпала двойка или тройка по отношению к общему количеству бросков.

Предположим, мы бросали кость 10 раз и двойка или тройка выпали 4 раза. Рассчитаем эмпирическую вероятность такого события.

$$ P(A) = \frac{4}{10} = \frac{2}{5} $$

Испытаний (бросков) было проведено очень мало и эмпирическая вероятность довольно сильно отличается от теоретической. Что же будет, если кратно увеличить количество бросков?

Закон больших чисел

Закон больших чисел (Law of large numbers) утверждает, что при проведении множества испытаний эмпирическая вероятность начнет приближаться к теоретической. Другими словами, если бросать кость достаточное количество раз, то эмпирическая вероятность выпадения двойки или тройки приблизится к 1/3.

Откроем ноутбук к этому занятию

Модуль random библиотеки Numpy

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

Подбрасывание монеты

Подбрасывание монеты можно имитировать, обозначив орел и решку через 0 и 1 и случайно генерируя соответствующие числа. Для этого подойдет функция np.random.randint().

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

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

Бросание игральной кости

Бросание игральной кости можно имиритировать, изменив диапазон на [1, 7).

Метод Монте-Карло

Теперь давайте вернемся к Закону больших чисел. Для его проверки воспользуемся методом Монте-Карло.

Метод Монте-Карло позволяет исследовать какой-либо случайный процесс, многократно имитируя его с помощью компьютера.

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

Будем создавать нашу программу поэтапно.

Шаг 1. Метод Монте-Карло

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

Создать такой список, вернее массив Numpy, можно с помощью функции np.logspace(). Она похожа на уже изученную нами функцию np.linspace() с той лишь разницей, что возвращает заданное количество чисел, которые равномерно распределены по логарифмической шкале в пределах заданного диапазона.

Пусть наша последовательность будет состоять из чисел 10, возведенных в степень от 1 до 6. То есть первым числом будет $10^1 = 10$, а последним $10^6 = 1 000 000$. Всего таких чисел будет 50 (значение функции np.logspace() по умолчанию).

Промежуточные числа Numpy рассчитает самостоятельно с помощью десятичного логарифма. Например, для получения числа 790604 мы возвели 10 в степень 5.897959007291967.

Убедимся в верности вычислений.

На графике такая последовательность демонстрирует экспоненциальный рост.

функция np.linspace()

Ровно то, что нам нужно.

Шаг 2. Напишем код для проведения испытаний

Используем два цикла:

  • С помощью первого (внешнего) цикла for мы можем имитировать серии бросков (всего 50 испытаний).
  • Второй (внутренний) цикл отвечает за броски каждой отдельной серии. В первый раз мы бросим кость 10 раз, во второй — 12 и так далее.

Пока ограничимся одной серией. Для этого в конце внешнего цикла поставим оператор break.

Шаг 3. Рассчитаем долю успешных испытаний

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

Обратите внимание, эмпирическая вероятность первой серии, в которой было только 10 бросков (1/2), далека от теоретической (1/3).

Шаг 4. Проведем полноценные испытания алгоритма

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

Помимо этого будем замерять время исполнения кода с помощью магической команды %%time.

Общее время выполнения кода составило 26 секунд.

Шаг 5. Оценим результат визуально

метод Монте-Карло

Как мы видим, если в первых сериях эмпирическая вероятность довольно сильно отличалась от теоретической, с увеличением количества бросков она вплотную приблизилась к 1/3.

Мы убедились в верности Закона больших чисел с помощью метода Монте-Карло.

Пример векторизованного кода

На прошлом занятии мы поговорили про векторизацию кода и выяснили, что векторизованный код исполняется быстрее, чем циклы for. Давайте посмотрим, нет ли возможности оптимизировать код с двумя циклами, написанный ранее.

Векторизованный код исполнился за 87,1 миллисекунды, то есть почти в 321 раз (26/0,081) быстрее, чем без использования векторизации.

Убедимся, что результат идентичен.

Задача о двух конвертах

Задача о двух конвертах (two envelopes problem) — популярная задача из области теории вероятностей. Сегодня мы рассмотрим один из ее вариантов и смоделируем решение с помощью модуля random.

Представьте, что вам дали конверт A, вы открыли его и увидели сумму X. После этого вам предлагают заменить этот конверт на другой конверт B, сумма в котором с одинаковой вероятностью либо в два раза меньше, либо в два раза больше суммы X. Стоит ли вам заменить конверт с A на B?

задача о двух конвертах

Аналитическое решение

Решим задачу аналитически. Пусть суммы в конверте B могут быть равны X/2 и 2X, тогда в среднем результат замены конвертов составит

$$ \frac{1}{2} \left( \frac{X}{2} \right) + \frac{1}{2} \left( 2X \right) = \frac{5}{4}X $$

Это больше, чем изначальная сумма X в конверте A, а значит нам стоит взять конверт B.

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

$$ \frac{1}{2} \left( \frac{5000}{2} \right) + \frac{1}{2} \left( 2 \cdot 5000 \right) = \frac{5}{4} \cdot 5000 = 6250 $$

Моделирование с помощью модуля random

Проверим это решение с помощью Питона и Закона больших чисел.

Как мы видим, компьютер подтвердил, что в среднем (при достаточном количестве испытаний) будет выгоднее заменить конверт A на конверт B.

Парадокс двух конвертов

Добавлю, что если изменить условие задачи и не открывать конверт A, то задача превратится в парадокс двух конвертов (two envelopes paradox). Ведь если в среднем выгоднее заменить конверт A на B, то будет выгодно и обратное действие (заменить B на A).

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

Рассмотрение этого парадокса выходит за рамки сегодняшнего занятия.

Метод случайных блужданий

Генератор случайных чисел также применяется в так называемом методе случайных блужданий (random walk models, random walkers). С его помощью мы можем конструировать путь, состоящий из последовательности случайных шагов в определенном пространстве и, таким образом, моделировать различные процессы.

Случайное блуждание по числовой прямой

В самом простом варианте мы будем подбрасывать монету и смотреть, куда мы сдвинулись относительно предыдущего шага на числовой прямой. Посмотрим на алгоритм.

Выведем сделанные шаги.

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

Для наглядности давайте посмотрим на сделанные шаги на числовой прямой.

метод случайных блужданий

Начальное положение ($s_0$) находится в точке ноль. При первом подбрасывании монеты выпадает 0, и на первом шаге ($s_1$) мы остаемся на месте. На втором шаге ($s_2$) выпадает единица, и мы сдвигаемся в точку один. На третьем шаге ($s_3$) мы перемещаемся в точку два. После этого, на шагах четыре ($s_4$) и пять ($s_5$) мы остаемся на месте.

Доедем ли мы из Москвы до Санкт-Петербурга

метод случайных блужданий: путь от Москвы до Санкт-Петербурга

Теперь усложним задачу и посмотрим, сможем ли мы добраться из Москвы до Санкт-Петербурга на довольно необычном транспортном средстве, движение которого зависит от бросания игральной кости. Условия следующие:

  • Нам нужно преодолеть ровно 700 километров;
  • Если при бросании кости выпадет единица или двойка, мы сместимся на 5 километров вперед;
  • Если тройка или четверка, то на 5 километров назад;
  • В случае если выпадет пять или шесть, мы повторно бросим кость и сдвинемся на выпавшее число километров умноженное на пять;
  • Кроме того, есть вероятность внезапной поломки. Она составляет менее 0,001, то есть менее 0,1%. В этом случае на эвакуаторе нас отвезут обратно в Москву для ремонта.

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

Теперь давайте с помощью Питона построим такой механизм.

  • Как уже было сказано, для бросания игральной кости мы будем использовать функцию np.random.randint() с параметрами 1 и 7.
  • Для того чтобы не уйти в минус, используем функцию max(). Первым параметром передадим 0, вторым — километр, на котором мы окажемся, если сдвинемся назад (то есть если выпадет твойка или четверка). В случае если второе число окажется отрицательным, то функция max() выберет ноль.
  • Поломку мы будем имитировать с помощью функции np.random.rand(). Эта функция генерирует равномерное распределение (об этом ниже) чисел в интервале [0, 1) и, если выпавшее число меньше, чем 0,001, мы вернемся в начало пути.

Посмотрим на код.

Функция случайного блуждания

Создание одного случайного блуждания

Вызовем эту функцию для генерации одного случайного блуждания.

Выведем результат с помощью графика.

одно случайное блуждание

Посмотрим на расстояние, которое удалось преодолеть.

Очевидно в этот раз мы не смогли доехать до Санкт-Петербурга.

Моделирование нескольких случайных блужданий

Давайте посмотрим на то, какова вероятность в среднем преодолеть этот путь.

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

Для наглядности мы можем использовать DataFrame библиотеки Pandas.

Первый столбец соответствует первому случайному блужданию и на шаге 100 (с индексом 99) мы достигли 430-ого километра.

Выведем результат на графике.

несколько случайных блужданий (график)

Обратите внимание на «упавшие» кривые. Это как раз те случаи, когда наше воображаемое транспортное средство ломалось, и мы возвращались в Москву на ремонт.

Вероятность преодоления пути

Теперь давайте возьмем достигнутые конечные точки каждого из 1000 испытаний и посмотрим на их распределение с помощью гистограммы.

несколько случайных блужданий (гистограмма)

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

Кроме того, мы можем рассчитать вероятность преодоления 700 км. Для этого с помощью функции np.count_nonzero() рассчитаем количество блужданий, преодолевших расстояние между двумя столицами, и разделим получившийся результат на общее количество испытаний.

Транспортное средство не слишком надежно, лишь в 12,8% случаев мы доезжали до конечной точки.

Продолжим знакомиться с теорией вероятностей и параллельно изучать возможности модуля random.