Проверка гипотез: откуда берётся p-value

Фрея · риск-менеджер необанка «Монетка»

Вы пишете то p = 0,053, то p = 0,00004, и в одном случае «не доказано», в другом «эффект есть». Я не подписываю того, чего не понимаю. Что это за число? И как часто вы с ним ошибаетесь?

В уроке про статистику A/B на Python мы получили p-value рецептом: 2 * stats.norm.sf(abs(z)). В прошлом уроке разобрались, почему средние и доли гуляют по колоколу и что такое стандартная ошибка. Сегодня соберём p-value руками, без формул, и Фрея увидит, откуда оно берётся.

Данные те же, что в уроке 68: в cafe_clients.csv по строке на клиента теста кешбэка, клиенты из обеих групп сразу уже исключены. Главная метрика урока 54: доля клиентов, которые с 1 апреля хоть раз заплатили в кафе. В тесте 90,27%, в контроле 87,53%, разница 2,74 пункта.

Суд над кешбэком

В суде подсудимый невиновен, пока вина не доказана. Прокурор приносит улики, и суд решает: улик хватает, чтобы осудить, или нет. Если не хватает, приговор звучит «оправдан», а не «доказано, что невиновен».

Проверка гипотез устроена так же.

• Нулевая гипотеза H0 это презумпция невиновности: кешбэк ничего не меняет, а 2,74 пункта нагнал случай.
• Альтернативная гипотеза H1: кешбэк меняет долю клиентов в кафе, в плюс или в минус.
• Улика это статистика, одно число по данным: у нас разница долей, 2,74 пункта.

Вопрос суда: насколько странна такая улика, если подсудимый невиновен?

Мир, где кешбэк не работает

Если H0 верна, метки test и control ничего не значат. Клиент сходил бы в кафе или нет независимо от того, какая буква стоит у него в строке. Значит, метки можно перемешать, и мир не изменится. Перемешаем 10 000 раз и каждый раз посчитаем разницу:

import numpy as np
import pandas as pd

clients = pd.read_csv("cafe_clients.csv")
is_test = (clients["test_group"] == "test").to_numpy()
converted = (clients["cafe_visits"] > 0).to_numpy()
observed = converted[is_test].mean() - converted[~is_test].mean()   # 0.0274

rng = np.random.default_rng(68)
null = np.empty(10_000)
for i in range(10_000):
    shuffled = rng.permutation(is_test)
    null[i] = converted[shuffled].mean() - converted[~shuffled].mean()

rng.permutation(is_test) переставляет метки: в группе test по-прежнему 976 клиентов, но уже других. Получилось распределение статистики при H0: какие разницы рождает случай. Вот оно по полпункта:

разница, п.п.   перестановок
  от -5,0        6
  от -4,5       14
  от -4,0       60  █
  от -3,5       91  ██
  от -3,0      297  ██████
  от -2,5      342  ███████
  от -2,0      788  ████████████████
  от -1,5      818  ████████████████
  от -1,0     1509  ██████████████████████████████
  от -0,5     1098  ██████████████████████
  от  0,0     1121  ██████████████████████
  от +0,5     1510  ██████████████████████████████
  от +1,0      764  ███████████████
  от +1,5      838  █████████████████
  от +2,0      315  ██████
  от +2,5      266  █████     ← наши 2,74
  от +3,0       81  ██
  от +3,5       52  █
  от +4,0       18
  от +4,5        8

Колокол вокруг нуля, как и обещала ЦПТ. Его разброс null.std() равен 0,0143, а стандартная ошибка по формуле из урока 54 была 0,0141: перестановки и формула описывают один и тот же мир. Зубчики на колоколе не ошибка. Доля в группе из 976 человек меняется ступеньками по одному клиенту, и некоторые ступеньки попадают в столбец чаще соседних.

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

p-value

Теперь вопрос суда в числах: в какой доле перестановок разница вышла не меньше нашей по модулю?

np.mean(np.abs(null) >= abs(observed))   # 0.0654

В 654 перестановках из 10 000 случай нагнал разницу в 2,74 пункта или больше, в любую сторону. Это и есть p-value: 0,065. В мире, где кешбэк не работает, такие улики встречаются в 6,5% случаев. Не каждый день, но и не чудо.

Формула урока 68 дала 0,053. Разница в ступеньках: формула заменяет их гладкой кривой, а перестановки считают ступеньку ровно на 2,74 пункта, и на неё одну приходится 0,9% перестановок. Готовая scipy.stats.permutation_test на тех же данных даёт 0,059: двусторонний p-value там считают как удвоенный меньший хвост. Три рецепта, три числа в третьем знаке, и все больше 0,05. Отсюда урок: p = 0,049 и p = 0,051 это одна и та же улика, а не два разных мира.

Если ни одна перестановка не дошла до наблюдаемой разницы, p-value не ноль. Пиши «p < 1/N», где N число перестановок: мы знаем только, что он меньше шага нашей симуляции. На 10 000 перестановках это «p < 1/10 000». Так выходит для частоты походов в кафе: из 5000 перестановок ни одна не дотянула до 0,63 похода на клиента, значит, p < 1/5000.