Сколько ждать успеха: геометрическое распределение и задача о купонах

Если ты выведешь 1/p одной строкой, посчитаешь вероятность ждать дольше k попыток и объяснишь, откуда в задаче о купонах берётся логарифм, - листай до «Теперь сам».

Как думаешь, сколько запросов нужно, чтобы задеть все серверы

Балансировщик отправляет каждый запрос на случайный из 50 серверов, все серверы равновероятны. Ты выкатываешь проверку и хочешь, чтобы через каждый сервер прошёл хотя бы один запрос.

Попробуй угадать до чтения, сколько запросов в среднем для этого нужно: 50, около 100 или больше 200? Ответ - 225. Первые серверы закрываются быстро, а последний приходится ждать в среднем 50 запросов: каждый запрос попадает в него с вероятностью 1/50.

На что это похоже

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

Сбор коллекции - это цепочка ожиданий, и каждое следующее длиннее предыдущего. Сначала ждут события с вероятностью почти единица, в конце - с вероятностью 1/n.

Механизм: одно ожидание и цепочка ожиданий

Испытание Бернулли - попытка с двумя результатами: успех с вероятностью p, неудача с вероятностью 1 − p. Попытки независимы, p не меняется.

Пусть X - номер попытки, на которой случился первый успех. Чтобы первый успех пришёлся на попытку k, нужны k − 1 неудач подряд и затем успех: P(X = k) = (1 − p)ᵏ⁻¹ · p. Это геометрическое распределение: вероятности идут геометрической прогрессией, каждая следующая в 1 − p раз меньше.

Ожидание проще всего найти, посмотрев на первую попытку. Она потрачена в любом случае. С вероятностью p на ней всё закончилось. С вероятностью 1 − p случилась неудача, и дальше ситуация ровно та же, что была в начале: прошлые неудачи ничего не изменили. Поэтому

E = 1 + (1 − p) · E
E · p = 1
E = 1 / p

Почему деление: p - доля удачных попыток в длинном ряду. Если удаётся одна из шести, на один успех приходится шесть попыток.

Вероятность ждать дольше k попыток - это k неудач подряд: P(X > k) = (1 − p)ᵏ.

Теперь купоны. Вариантов n, каждая попытка даёт случайный из них. Пусть уже собрано i разных. Новый вариант выпадает с вероятностью (n − i) / n, и ждать его в среднем n / (n − i) попыток. Полное время - сумма ожиданий по этапам, по линейности:

этап:        0 собрано   1 собран   2 собрано   …   n − 1 собрано
p нового:    n/n         (n−1)/n    (n−2)/n     …   1/n
ждать:       n/n         n/(n−1)    n/(n−2)     …   n/1

E = n · (1/n + 1/(n−1) + … + 1/2 + 1/1) = n · Hₙ

Hₙ - гармоническое число, примерно ln n + 0,577. Логарифм здесь - сумма долей 1/k: каждый следующий этап дороже, и цены этапов образуют гармонический ряд, записанный с конца.

Правило

Число независимых попыток до первого успеха при вероятности успеха p имеет ожидание 1/p, а вероятность прождать больше k попыток равна (1 − p)ᵏ. Ожидаемое число попыток, за которое выпадут все n равновероятных вариантов, равно n · Hₙ ≈ n · ln n + 0,577 · n. Число успехов в n попытках имеет ожидание n · p.

Разбор: ожидание шестёрки

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

import random

random.seed(4)
opytov = 5000
vsego = 0
dolgo = 0
for opyt in range(opytov):
    broskov = 1
    while random.randint(1, 6) != 6:
        broskov = broskov + 1
    vsego = vsego + broskov
    if broskov > 12:
        dolgo = dolgo + 1
print(round(vsego / opytov, 3), round(dolgo / opytov, 4), round((5 / 6) ** 12, 4))

5.961 0.1158 0.1122

Среднее 5,96 при ожидании 6. Дольше двенадцати бросков - вдвое дольше среднего - ждать пришлось в 11,6% опытов; формула (5/6)¹² даёт 11,2%. Ждать втрое дольше среднего придётся с вероятностью (5/6)¹⁸ ≈ 3,8%, вчетверо - 1,3%. Каждые шесть лишних бросков уменьшают вероятность примерно в три раза: хвост геометрического распределения тает быстро, и это свойство пригодится, когда повторы станут способом снижать ошибку алгоритма.

Разбор: полная коллекция из 50 вариантов

import random

def sobrat(n):
    est = [False] * n
    raznyh = 0
    pokupok = 0
    while raznyh < n:
        k = random.randrange(n)
        pokupok = pokupok + 1
        if not est[k]:
            est[k] = True
            raznyh = raznyh + 1
    return pokupok

random.seed(21)
n = 50
opytov = 200
vsego = 0
hudshij = 0
for opyt in range(opytov):
    p = sobrat(n)
    vsego = vsego + p
    hudshij = max(hudshij, p)

garm = 0
for i in range(1, n + 1):
    garm = garm + 1 / i
print(round(vsego / opytov, 1), round(n * garm, 1), hudshij)

227.7 225.0 518

random.randrange(n) возвращает случайное целое от 0 до n − 1. Среднее по двумстам опытам - 227,7, по формуле 50 · H₅₀ = 225,0. Худший из двухсот опытов занял 518 попыток, больше двух средних: время сбора определяют последние варианты, а у них самое длинное и самое изменчивое ожидание.

Как распределено время по этапам, видно из таблицы ожиданий n / (n − i):

первые 25 вариантов из 50:    50/50 + 50/49 + … + 50/26  ≈  34 попытки
следующие 20:                 50/25 + … + 50/6           ≈  77 попыток
последние 5:                  50/5 + 50/4 + … + 50/1     ≈ 114 попыток

Половина коллекции собирается за 34 попытки, а последние пять вариантов стоят половину всего времени.

Граница формулы - равные вероятности. Если один вариант выпадает с вероятностью 1/1000, а остальные 49 делят между собой остальное, ждать редкий придётся в среднем тысячу попыток, и к этому моменту все прочие давно собраны. При неравных вероятностях время сбора определяют редкие варианты: оно не меньше 1 / p_min, а если одинаково редких вариантов несколько, ждать приходится их все. К n · ln n это время отношения не имеет.

Что ломается без этого

Память у независимых попыток. После десяти неудач до успеха в среднем остаётся те же 1/p попыток, что и в начале. Счётчик «уже десять раз не повезло, сейчас повезёт» ничего не предсказывает.

Формула 1/p для зависимых попыток. Клиент повторяет запрос к перегруженному серверу. Если повторяют все клиенты сразу, нагрузка растёт, вероятность успеха падает, и ожидание уже не 1/p: каждая неудача ухудшает следующую попытку. Поэтому повторы делают с растущей паузой и случайной добавкой к ней - чтобы попытки разных клиентов не совпадали во времени.

Оценка «n попыток на n вариантов». Пятьдесят случайных запросов заденут в среднем только 32 сервера из пятидесяти: каждый остаётся незадетым с вероятностью (49/50)⁵⁰ ≈ 0,364.

Тест «до первой удачи» без предела. Цикл «генерируй случайный вариант, пока не подойдёт» при p = 0 не кончается никогда. Ожидание 1/p существует только при p > 0, и перед таким циклом нужно убедиться, что подходящий вариант вообще есть.

Зачем это в работе

Открытая адресация в хеш-таблице. Если каждая проба попадает в случайную ячейку, то при заполнении α она оказывается свободной с вероятностью 1 − α, и неудачный поиск делает в среднем 1 / (1 − α) проб - это и есть ожидание геометрической величины. При линейном пробировании занятые ячейки слипаются в серии, пробы перестают быть независимыми, и проб в среднем больше. При α = 0,5 две пробы, при α = 0,9 - десять; отсюда пороги расширения таблиц.

Генерация случайного объекта отбором: «возьми случайное число нужной длины, проверь, простое ли». Простых среди 1024-битных чисел примерно одно на 710, значит до успеха в среднем 710 кандидатов, или 355, если брать только нечётные. Так выбирают числа для ключа RSA.

Случайное тестирование: если у функции n веток и каждая случайная проверка попадает в одну из них, задача о купонах говорит, сколько проверок нужно, чтобы пройти все, и предупреждает, что редкая ветка потребует 1 / p проверок сама по себе.

Теперь сам

1. Пакет теряется с вероятностью 0,2, передачи независимы. Сколько передач в среднем нужно для доставки?
2. Какова вероятность, что тот же пакет не доставлен после четырёх передач?
3. Кубик бросают, пока не выпадут все шесть граней. Сколько бросков в среднем?
4. Собрано 45 вариантов из 50. Сколько попыток в среднем нужно до следующего нового?
5. Хеш-таблица с открытой адресацией заполнена на 75%, каждая проба попадает в случайную ячейку. Сколько проб в среднем делает поиск ключа, которого в таблице нет?