Дискретные распределения

Большинство вопросов о случайности — вопросы о числах: сколько орлов выпадет в 10 бросках, сколько кликов придётся на 1000 показов, сколько бросков понадобится до первой шестёрки, сколько запросов поступит на сервер за секунду. Случайная величина — это такое число, а её распределение говорит, как вероятность разнесена по возможным значениям. Несколько распределений, у каждого из которых есть короткая порождающая его история, описывают большинство встречающихся на практике счётных величин. В этой заметке мы вводим случайные величины и пять таких историй; следующая заметка делает то же самое для измерений, которые меняются непрерывно.

Обозначения. $X, Y$ — случайные величины. $p_X(k) = P(X = k)$ — функция вероятности (probability mass function, PMF; закон распределения) дискретной величины $X$, а $F_X(x) = P(X \le x)$ — её функция распределения (cumulative distribution function, CDF). Запись $X \sim \mathrm{Bin}(n, p)$ читается как «$X$ имеет биномиальное распределение с параметрами $n$ и $p$». $\mathbb 1_A$ — индикатор события $A$: 1, если $A$ происходит, и 0 иначе.

Случайные величины#

Формально случайная величина — это функция $X : \Omega \to \mathbb R$, которая сопоставляет каждому исходу число. Бросим две кости; сумма $S$ переводит исход $(3, 5)$ в 8. События, связанные с $X$, — это события, связанные с исходами: $\{S = 8\}$ — множество из пяти пар с суммой 8, поэтому $p_S(8) = 5/36$. Вся функция вероятности $S$:

$$p_S(2), \dots, p_S(12) = \tfrac{1}{36}, \tfrac{2}{36}, \tfrac{3}{36}, \tfrac{4}{36}, \tfrac{5}{36}, \tfrac{6}{36}, \tfrac{5}{36}, \tfrac{4}{36}, \tfrac{3}{36}, \tfrac{2}{36}, \tfrac{1}{36} .$$

Функция вероятности неотрицательна и в сумме по всем значениям даёт 1, а любая вероятность, связанная с $X$, — сумма её значений: $P(S \ge 10) = \frac{3 + 2 + 1}{36} = \frac16$. Когда распределение известно, о пространстве исходов можно забыть. Два разных эксперимента — бросок кости и вращение колеса с шестью равными секторами — дают случайные величины с одним и тем же распределением, и все вероятностные утверждения о них совпадают.

Функция распределения $F_X(x) = P(X \le x)$ несёт ту же информацию в накопленном виде. Для дискретной величины это ступенчатая функция, которая подскакивает на $p_X(k)$ в каждом значении $k$, и $P(a < X \le b) = F_X(b) - F_X(a)$.

Распределения Бернулли и категориальное#

Простейшая случайная величина принимает значение 1 с вероятностью $p$ и 0 с вероятностью $1 - p$: $X \sim \mathrm{Bernoulli}(p)$ (распределение Бернулли). Каждое событие типа «да/нет» даёт такую величину — свой индикатор $\mathbb 1_A \sim \mathrm{Bernoulli}(P(A))$. Спам-фильтр, модель кликов и бинарный классификатор — все они выдают параметр $p$ распределения Бернулли для метки.

Если возможных значений не два, а $K$, категориальное распределение (categorical) приписывает значению $k$ вероятность $p_k$, где $p_1 + \dots + p_K = 1$. Выход softmax классификатора — категориальное распределение по классам, а выход языковой модели — категориальное распределение по словарю. Обучать их по кросс-энтропии — значит подбирать эти параметры методом максимального правдоподобия.

Биномиальное распределение#

Повторим испытание Бернулли $n$ раз независимо, с одним и тем же $p$, и подсчитаем успехи: $X \sim \mathrm{Bin}(n, p)$. Конкретная последовательность с $k$ успехами и $n - k$ неудачами в силу независимости имеет вероятность $p^k (1-p)^{n-k}$, а таких последовательностей $\binom nk$, поэтому

$$P(X = k) = \binom nk p^k (1-p)^{n-k}, \qquad k = 0, 1, \dots, n .$$

Десять бросков честной монеты дают ровно 5 орлов с вероятностью $\binom{10}{5} / 2^{10} = 252/1024 = 0.246$, а не меньше 8 орлов — с вероятностью $56/1024 = 0.055$. Объявление, показанное 1000 раз с вероятностью клика 2%, получает $\mathrm{Bin}(1000, 0.02)$ кликов: в среднем 20, а 30 и больше — с вероятностью всего $0.021$. Говорят ли 30 кликов о том, что объявление лучше, или это везение, — вопрос для проверки гипотез; вычисления в ней опираются на биномиальное распределение. Число правильных ответов классификатора на тестовой выборке тоже имеет биномиальное распределение, если тестовые примеры независимы.

Геометрическое распределение#

Будем повторять испытание до первого успеха и подсчитаем число испытаний: $X \sim \mathrm{Geom}(p)$. Первый успех наступает в испытании $k$, если первые $k - 1$ испытаний неудачны, а $k$-е успешно:

$$P(X = k) = (1-p)^{k-1} p, \qquad P(X > k) = (1-p)^k, \qquad k = 1, 2, \dots$$

Ожидание шестёрки занимает больше 10 бросков с вероятностью $(5/6)^{10} = 0.16$. (В некоторых книгах вместо этого считают неудачи до первого успеха, что сдвигает значения на единицу.)

Геометрическое распределение не имеет памяти (memoryless; свойство отсутствия последействия): то, что уже было $m$ неудач, ничего не меняет в предстоящем ожидании:

$$P(X > m + k \mid X > m) = \frac{(1-p)^{m+k}}{(1-p)^m} = (1-p)^k = P(X > k) .$$

После 10 бросков без шестёрки вероятность ещё 10 бросков без неё по-прежнему $0.16$. У кости нет памяти; убеждение, что шестёрка «должна вот-вот выпасть», — это ошибка игрока (gambler’s fallacy).

Распределение Пуассона#

Многие счётные величины — это числа редких событий среди очень большого числа возможностей: опечатки на странице, запросы к серверу за секунду, мутации в геноме, голы в футбольном матче. Если возможностей $n$ и у каждой малая вероятность $p = \lambda / n$, то число событий распределено как $\mathrm{Bin}(n, \lambda/n)$, и при росте $n$ с фиксированным средним $\lambda$

$$\binom nk \Big(\frac{\lambda}{n}\Big)^k \Big(1 - \frac{\lambda}{n}\Big)^{n-k} = \underbrace{\frac{n (n-1) \cdots (n-k+1)}{n^k}}_{\to 1} \frac{\lambda^k}{k!} \underbrace{\Big(1 - \frac{\lambda}{n}\Big)^{n}}_{\to e^{-\lambda}} \underbrace{\Big(1 - \frac{\lambda}{n}\Big)^{-k}}_{\to 1} \to e^{-\lambda} \frac{\lambda^k}{k!} .$$

Словами: первый множитель сравнивает число способов расположить $k$ событий, $\binom nk$, с $n^k / k!$ и стремится к 1, поскольку $k$ фиксировано, а $n$ растёт; так же ведёт себя последний множитель; а средний — вероятность того, что ни одна из $n$ возможностей не сработала, — стремится к $e^{-\lambda}$ по известному из анализа пределу $(1 + x/n)^n \to e^x$, тому самому, что стоит за непрерывным начислением процентов.

Предел — это распределение Пуассона, $X \sim \mathrm{Poisson}(\lambda)$, с $P(X = k) = e^{-\lambda} \lambda^k / k!$ при $k = 0, 1, 2, \dots$ У него один параметр — среднее число событий, и ему не нужны $n$ и $p$ по отдельности, поэтому оно подходит для счётных величин, у которых число возможностей никто не знает. Приближение хорошее уже при умеренных $n$:

Рисунок

Биномиальное распределение с малым p приближается к пуассоновскому

Классический набор данных — подсчёт Ладислава Борткевича: число солдат прусской кавалерии, погибших от удара копытом лошади. За 200 корпусо-лет погибли 122 человека, в среднем $0.61$ на корпус в год. Распределение Пуассона с $\lambda = 0.61$, подобранное по одному этому числу, предсказывает, сколько корпусо-лет должно было прийтись на 0, 1, 2, 3 и 4 смерти:

Рисунок

Гибель от удара копытом в прусской кавалерии и распределение Пуассона

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

Гипергеометрическое распределение#

Вынем $n$ элементов без возвращения из совокупности из $N$ элементов, $K$ из которых особые, и подсчитаем особые. Заметка о подсчёте уже вычислила это распределение для двух тузов в раздаче из пяти карт:

$$P(X = k) = \frac{\binom Kk \binom{N - K}{n - k}}{\binom Nn} ,$$

и $P(X = 2) = 0.0399$ при $N = 52$, $K = 4$, $n = 5$. Число особых элементов меняется от $\max(0, n - (N - K))$ — когда выборка так велика, что избежать всех особых элементов нельзя, — до $\min(n, K)$; вне этого диапазона один из биномиальных коэффициентов в формуле равен нулю. Извлечения зависимы: каждый вынутый туз делает следующий менее вероятным. Если считать их независимыми испытаниями с $p = 4/52$, получится биномиальный ответ $0.0465$ — на 17% больше. Когда совокупность гораздо больше выборки, удаление нескольких элементов почти не меняет пропорций, и гипергеометрическое распределение близко к биномиальному; поэтому в опросах 1000 человек в многомиллионной стране пользуются биномиальными формулами.

Какая история даёт какое распределение#

История	Распределение	Значения	Функция вероятности
Одно испытание «да/нет»	$\mathrm{Bernoulli}(p)$	0, 1	$p^k (1-p)^{1-k}$
Одна из $K$ категорий	$\mathrm{Categorical}(p_1, \dots, p_K)$	$1, \dots, K$	$p_k$
Успехи в $n$ независимых испытаниях	$\mathrm{Bin}(n, p)$	$0, \dots, n$	$\binom nk p^k (1-p)^{n-k}$
Испытания до первого успеха	$\mathrm{Geom}(p)$	$1, 2, \dots$	$(1-p)^{k-1} p$
Редкие события среди многих возможностей	$\mathrm{Poisson}(\lambda)$	$0, 1, \dots$	$e^{-\lambda} \lambda^k / k!$
Особые элементы в выборке без возвращения	$\mathrm{HGeom}(N, K, n)$	$\max(0, n - N + K), \dots, \min(n, K)$	$\binom Kk \binom{N-K}{n-k} / \binom Nn$

Что изменилось#

	События	Случайные величины
Что описывается	одно утверждение «да/нет» за раз	число и сразу все утверждения о нём
Как	$P(A)$	функция вероятности $p_X(k)$ или функция распределения $F_X(x)$
Повторное использование	новое вычисление для каждой задачи	именованное распределение для каждой истории: биномиальное, геометрическое, пуассоновское, …

Ограничения#

• Распределение хорошо ровно настолько, насколько хороша его история. Биномиальному распределению нужны независимые испытания с одним и тем же $p$; клики одного и того же пользователя или тестовая выборка с почти-дубликатами нарушают оба условия, и реальный разброс оказывается шире биномиального.
• Реальные счётные данные часто сверхдисперсны (overdispersed). Числа событий, которые приходят всплесками (запросы от ботов, слова в документах), меняются сильнее, чем допускает распределение Пуассона с тем же средним; тогда используют отрицательное биномиальное распределение или модели смесей.
• Значений должно быть не более чем счётное число. Измерениям вроде времени ожидания или роста нужны плотности — тема следующей заметки.

Что нужно запомнить#

• Случайная величина — число, определяемое исходом; её функция вероятности $p_X(k) = P(X = k)$ или функция распределения $F_X(x) = P(X \le x)$ отвечает на любой вероятностный вопрос о ней.
• Распределения Бернулли и категориальное: одно испытание с двумя или $K$ исходами, выход классификатора.
• Биномиальное: число успехов в $n$ независимых испытаниях, $\binom nk p^k (1-p)^{n-k}$.
• Геометрическое: число испытаний до первого успеха, $P(X > k) = (1-p)^k$, нет памяти.
• Пуассоновское: предел $\mathrm{Bin}(n, \lambda/n)$, $e^{-\lambda} \lambda^k / k!$, для редких событий среди многих возможностей.
• Гипергеометрическое: выборка без возвращения; близко к биномиальному, когда совокупность гораздо больше выборки.

Код#

Все числа заметки с помощью scipy.stats:

import numpy as np
from scipy import stats

# the sum of two dice: a random variable and its PMF (counts out of 36)
s = np.add.outer(np.arange(1, 7), np.arange(1, 7)).ravel()
print(np.bincount(s)[2:])                                         # [1 2 3 4 5 6 5 4 3 2 1]

print(round(stats.binom(10, 0.5).pmf(5), 4), round(stats.binom(10, 0.5).sf(7), 4))   # 0.2461 0.0547
clicks = stats.binom(1000, 0.02)
print(clicks.mean(), round(clicks.sf(29), 4))                     # 20.0 0.0207

print(round((5 / 6) ** 10, 4))                                    # 0.1615: more than 10 rolls to the first six

# binomial(n, 3/n) against Poisson(3): total variation distance
for n in [5, 10, 100]:
    k = np.arange(n + 1)
    tv = 0.5 * (np.abs(stats.binom(n, 3 / n).pmf(k) - stats.poisson(3).pmf(k)).sum() + stats.poisson(3).sf(n))
    print(n, round(tv, 3))                                        # 5 0.219, 10 0.086, 100 0.008

# Bortkiewicz's horse kicks: corps-years with 0..4 deaths, and the Poisson fit
observed = np.array([109, 65, 22, 3, 1])
lam = (observed * np.arange(5)).sum() / observed.sum()
print(lam, (observed.sum() * stats.poisson(lam).pmf(np.arange(5))).round(1))   # 0.61 [108.7  66.3  20.2   4.1   0.6]

# two aces in five cards: hypergeometric, and the binomial that ignores the missing cards
print(round(stats.hypergeom(52, 4, 5).pmf(2), 4), round(stats.binom(5, 4 / 52).pmf(2), 4))   # 0.0399 0.0465

Дальше#

Счётные величины принимают отдельные значения, и их вероятности — массы, сосредоточенные в этих значениях. Время ожидания, рост или положение случайной точки могут принимать любое значение из интервала, и каждое отдельное значение имеет вероятность 0, как заметка об аксиомах показала для точки на $[0, 1]$. Следующая заметка заменяет функцию вероятности плотностью, вводит равномерное, показательное и нормальное распределения и показывает, что показательное распределение — непрерывный аналог геометрического: единственное время ожидания без памяти.

Литература#

• J. K. Blitzstein, J. Hwang. Introduction to Probability, 2nd ed. CRC Press, 2019. Глава 3 (случайные величины, распределения Бернулли, биномиальное и гипергеометрическое), глава 4 (геометрическое распределение, распределение Пуассона, пуассоновское приближение).
• W. Feller. An Introduction to Probability Theory and Its Applications, Vol. 1, 3rd ed. Wiley, 1968. Глава VI (биномиальное распределение и распределение Пуассона).
• L. von Bortkiewicz. Das Gesetz der kleinen Zahlen. Teubner, 1898. Данные о гибели от ударов копытом.