Увеличиваем размер выборки и прокрашиваем серые метрики: неочевидная ошибка при проведении А/B — тестов

от автора

Привет, Хабр Недавно посмотрел выступление Валерия Бабушкина, которое было опубликовано в далеком ковидном 2020 году, но тем не менее основная часть информации из этого видео остается актуальной и по сей день (всем новичкам рекомендуется к просмотру). В его докладе был один единственный момент, который меня немного смутил: в А/B‑тестах компании X5 Tech на тот момент в качестве единицы анализа рассматривается «Пятёрочка‑день». Если так подумать, то действительно, такой подход звучит как легальный способ увеличить количество наблюдений в группах, что, казалось бы, влечет за собой массу положительных моментов. Увеличение размеров групп положительно влияет на мощность тестов и MDE.

В этой статье я постараюсь кратко подчеркнуть негативные последствия такого подхода, покажу всё с помощью кода и графиков.

Генерация данных

К сожалению, нырнуть в базы данных X5 мы не имеем возможности, но провести синтезацию похожих по распределению данных мы можем легко с помощью инструментов библиотек для анализа данных.

Возьмем 1000 магазинов (shop_id), для каждого магазина создадим равное количеству дней эксперимента число строк (пусть наш синтетический тест будет длиться 2 недели), каждой записи shop_id-date присвоим значение метрики. Метрика может быть любой: количество доставок, среднее время доставки заказа, количество отмен, среднее количество заказов на курьера в этот день и др. Для наглядности предлагаю взять курьерскую доставку и обозначить в качестве целевой метрики «Среднее количество доставленных заказов за день».

import pandas as pd import numpy as np from datetime import datetime, timedelta  # Параметры num_shops = 1000  # Количество магазинов num_days = 14    # Количество дней  # Генерация случайных параметров для магазинов shop_params = {     shop_id: {         'mean': np.random.normal(loc=200, scale=25),  # Случайное выборочное среднее из ГС         'std_dev': np.random.normal(loc=25, scale=2)  # Случайное стандартное отклонение для выборки     }     for shop_id in range(1, num_shops + 1) }  # Генерация данных data = []  for shop_id in range(1, num_shops + 1):     for day in range(num_days):         # Вычисляем дату (только день)         timestamp = (datetime.now() - timedelta(days=num_days - day)).date()  # Получаем дату         # Получаем параметры метрики для текущего магазина         mean = shop_params[shop_id]['mean']         std_dev = shop_params[shop_id]['std_dev']         # Генерируем значение метрики из нормального распределения         metric_value = np.random.normal(loc=mean, scale=std_dev)  # Используем случайные mean и std_dev         # Добавляем запись в список         data.append({             'shop_id': shop_id,             'date': timestamp,             'metric': metric_value         })  df = pd.DataFrame(data)  df.sort_values(by=['shop_id', 'date'], inplace=True)  df.head(20)

Данные выглядят следующим образом:

Фрагмент данных

Фрагмент данных

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

Гипотеза

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

  • новый алгоритм построения маршрута для курьеров

  • оптимизация распределения заказов между вело/пешими/авто-курьерами

  • новые мотивирующие факторы, которые уменьшают время простоя курьера на пункте сбора

Как будем проверять

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

Распределение нашей метрики

Распределение нашей метрики

Супер, на картинке знакомый колокол, можем применять t-тест для наших данных!

А правда ли выиграем в размере выборки?

С помощью инструментов Python посчитаем сколько наблюдений нам потребуется:

from typing import Optional, Union from pingouin import power_ttest import numpy as np import pandas as pd  def get_stat(     array: Union[pd.Series, np.array],     effect_percent: Optional[float] = None,     n: Optional[int] = None,     alpha: Optional[float] = None,     power: Optional[float] = None,     alternative: str = 'two-sided',     paired: bool = False, ) -> float:     mean = array.mean()     std = array.std()     if effect_percent:         cohen_effect = effect_percent * (mean / std)     else:         cohen_effect = None     contrast = 'paired' if paired else 'two-samples'     result = power_ttest(         d=cohen_effect,         n=n,         power=power,         alpha=alpha,         contrast=contrast,         alternative=alternative,     )     if not effect_percent:         result = round(result * (std / mean) * 100, 2)     return result  def get_stat_table(     array: Union[pd.Series, np.array],     effect_percent: Union[np.array, list, float, None] = None,     n: Union[np.array, list, float, None] = None,     alpha: Union[np.array, list, float, None] = None,     power: Union[np.array, list, float, None] = None,     alternative: str = 'two-sided',     paired: bool = False, ) -> pd.DataFrame:     stats = n if n else effect_percent     result = pd.DataFrame(index=stats)     for p in power:         for a in alpha:             outputs = []             for stat in stats:                 if n:                     output = get_stat(                         array=array,                         effect_percent=None,                         n=stat,                         alpha=a,                         power=p,                         alternative=alternative,                         paired=paired,                     )                 else:                     output = get_stat(                         array=array,                         effect_percent=stat,                         n=None,                         alpha=a,                         power=p,                         alternative=alternative,                         paired=paired,                     )                 outputs.append(round(output, 2))             # Создаем многоуровневый индекс для заголовков             result[(f'α={round(a, 3)}', f'β={round(1-p, 3)}')] = outputs             result.index.name = 'n' if n else 'MDE'          # Переименовываем столбцы, чтобы они отображались в виде двух строк     result.columns = pd.MultiIndex.from_tuples(result.columns)          return result
get_stat_table(np.array(df['metric']),                 effect_percent=[0.01, 0.05, 0.06, 0.07, 0.08, 0.09, 0.1], #эффект в долях                n=None,                 alpha=[0.01, 0.05, 0.1],                 power=[0.8, 0.9])
Таблица зависимости кол-ва наблюдений в одной группе от MDE, α и β для случая, когда единицей анализа является "Пятёрочка-день"

Таблица зависимости кол-ва наблюдений в одной группе от MDE, α и β для случая,
когда единицей анализа является «Пятёрочка-день»

Предположим, что мы уже имеем успешный продукт, сорвали все «низковисящие фрукты», и теперь хотим улавливать очень маленькие эффекты: установим MDE = 5%. Уровень значимости и ошибку второго рода установим стандартными α = 0.05 и β = 0.2. В таблице видим — необходимо всего лишь 202 наблюдения в группе для выбранных параметров. Это что же получается, мы можем держать тест 2 недели на всего лишь 15 магазинах? А если еще поверх всего этого применять всякие стратификации, CUPED’ы и остальные хайповые аналитические инструменты, то профита от тестов будет еще больше.
Звучит заманчиво!

Так, хорошо…А если бы мы анализировали тест стандартным подходом, где каждый магазин будет единственно представлен в свой группе? Для этого усредним метрики для каждого магазина:

get_stat_table(np.array(df.groupby('shop_id')['metric'].mean().reset_index()['metric']),                 effect_percent=[0.01, 0.05, 0.06, 0.07, 0.08, 0.09, 0.1],                 n=None,                 alpha=[0.01, 0.05, 0.1],                 power=[0.8, 0.9])
Таблица зависимости кол-ва наблюдений в одной группе от MDE, α и β для случая, когда единицей анализа является "Пятёрочка"

Таблица зависимости кол-ва наблюдений в одной группе от MDE, α и β для случая,
когда единицей анализа является «Пятёрочка»

Ну тут уже в любом случае проигрыш: необходимо 109 магазинов в каждой группе, проигрыш в 7 раз!

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

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

Ну что, давайте тестировать?

Всё начинается с проверки критерия

Прежде чем приступать к запуску А/В-тестов необходимо всегда проверять работу статкритериев на А/А-тестах. От А/А-теста ожидаем увидеть отсутствие различий между группами. Как это можно проверить? Очень просто: 1000 раз проводим сравнение статзначимости между группами, в которых нет изменений, попутно каждый раз записывая значение p-value для текущих выборок.

Что мы ожидаем увидеть? Мы ожидаем, что распределение значений p-value будет равномерным, поскольку в A/A-тесте группы не должны отличаться друг от друга. Это означает, что каждое значение p-value должно встречаться с равной вероятностью.

Давайте теперь это проверим: в группы А и В будем брать по 100 случайных shop_id.

from scipy import stats from tqdm import tqdm  pvalues_aa = []  for _ in tqdm(range(1000), desc="Conducting A/A Tests"):     # Случайное деление магазинов на две группы     shops = df['shop_id'].unique()          # Случайно выбираем 100 уникальных магазинов     selected_shops = np.random.choice(shops, size=200, replace=False)  # 200 магазинов всего     group_a = selected_shops[:100]  # Первая группа     group_b = selected_shops[100:]  # Вторая группа          # Получаем метрики для каждой группы     metrics_a = df[df['shop_id'].isin(group_a)]['metric']     metrics_b = df[df['shop_id'].isin(group_b)]['metric']          # Проводим t-тест     pvalue = stats.ttest_ind(metrics_a, metrics_b).pvalue     pvalues_aa.append(pvalue)  plt.figure(figsize=(10, 6)) plt.hist(pvalues_aa, bins=50, color='skyblue', edgecolor='black', density=False) plt.title('Распределение p-значений A/A-тестов') plt.xlabel('p-value') plt.ylabel('Кол-во') plt.grid(axis='y', alpha=0.75) plt.axvline(x=0.05, color='red', linestyle='--', label='Уровень значимости α = 0.05') plt.legend() plt.show()
Распределение p-value для изначальной метрики

Распределение p-value для изначальной метрики

Что-то пошло не так…Доля наблюдаемых ошибок первого рода составляет порядка 40%! Иными словами, при запуске A/B-теста мы бы наблюдали статзначимые различия там, где их на самом деле нет. Метрики бы красились практически в каждом втором эксперименте. Думаю, что многие читатели еще до запуска эксперимента предвидели такой результат. Проблема заключается в следующем: t-тест для независимых наблюдений был проведен на зависимых данных, от каждого магазина в выборке оказывается по 14 наблюдений. Единица разбиения (shop_id) не совпадает с единицей анализа (shop_id-datetime, т.е. «Пятёрочка-день»).

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

pvalues_aa_fixed = []  # Общее количество записей total_records = len(df)  for _ in tqdm(range(1000), desc="Conducting A/A Tests"):     # Случайное деление всех записей на две группы     # Перемешиваем индексы записей     indices = np.arange(total_records)     np.random.shuffle(indices)  # Перемешиваем индексы      # Делим индексы на две группы     split_index = total_records // 2     group_a_indices = indices[:split_index]     group_b_indices = indices[split_index:]      # Получаем метрики для каждой группы     metrics_a = df.iloc[group_a_indices]['metric']     metrics_b = df.iloc[group_b_indices]['metric']      # Проводим t-тест     pvalue = stats.ttest_ind(metrics_a, metrics_b).pvalue     pvalues_aa_fixed.append(pvalue)  plt.figure(figsize=(10, 6)) plt.hist(pvalues_aa_fixed, bins=50, color='skyblue', edgecolor='black', density=False) plt.title('Распределение p-value A/A-тестов') plt.xlabel('p-value') plt.ylabel('Кол-во') plt.grid(axis='y', alpha=0.75) plt.axvline(x=0.05, color='red', linestyle='--', label='Уровень значимости α = 0.05') plt.legend() plt.show()
Распределение p-value для усредненной метрики

Распределение p-value для усредненной метрики

А теперь всё встало на свои места, ожидания и реальность совпали. При таком разбиении можно спокойно запускать A/B.

Switchback-подход

На самом деле использовать в качестве единицы разбиения «Пятерочку-день» можно, в таком случае каждый магазин имеет шанс оказаться одновременно и в тестовой, и в контрольной группе. Что я под этим имею в виду: на 1-ый день эксперимента магазин попадает в контроль, на 2-й — в тест, на 3-й — снова в тест, 4-й — контроль, и так далее. Да, это алгоритм switchback-тестирования (подробно об этом у ребят из Ситимобил тут), разве что с более широким окном перемешивания групп.

Продемонстрируем на А/А-тесте, что такой подход не противоречит условиям применения
t-теста, и снова построим распределение p-value:

import numpy as np import pandas as pd from scipy import stats import matplotlib.pyplot as plt from tqdm import tqdm  pvalues_aa_switchback = []  for _ in tqdm(range(1000), desc="Conducting A/A Tests"):     # Случайный выбор 100 уникальных магазинов     shops = df['shop_id'].unique()     selected_shops = np.random.choice(shops, size=100, replace=False)  # 100 магазинов всего      # Получаем все записи для выбранных магазинов     selected_records = df[df['shop_id'].isin(selected_shops)]      # Перемешиваем записи     shuffled_records = selected_records.sample(frac=1).reset_index(drop=True)      # Делим записи на две группы     split_index = len(shuffled_records) // 2     group_a_records = shuffled_records.iloc[:split_index]     group_b_records = shuffled_records.iloc[split_index:]      # Получаем метрики для каждой группы     metrics_a = group_a_records['metric']     metrics_b = group_b_records['metric']      # Проводим t-тест     pvalue = stats.ttest_ind(metrics_a, metrics_b).pvalue     pvalues_aa_switchback.append(pvalue)  plt.figure(figsize=(10, 6)) plt.hist(pvalues_aa_switchback, bins=50, color='skyblue', edgecolor='black', density=False) plt.title('Распределение p-value (A/A, switchback)') plt.xlabel('p-value') plt.ylabel('Кол-во') plt.grid(axis='y', alpha=0.75) plt.axvline(x=0.05, color='red', linestyle='--', label='Уровень значимости α = 0.05') plt.legend() plt.show()
Распределение p-value с использованием switchback'a

Распределение p-value с использованием switchback’a

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

Несомненно, можно долго обсуждать главное преимущества switchback’a, такое как снижение сетевого эффекта (курьеры с разных магазинов могут взаимодействовать между собой и влиять таким образом на результаты теста), но всё-таки такой подход может потребовать времени на разработку систем автоматизации для переключения магазина между группами (или нам вообще придется заставить человека вручную шаффлить группы).

А что делать, если мы хотим использовать такой подход, но при этом не можем себе позволить постоянно перемешивать группы? Например, включение и выключение изо дня в день какой-нибудь функции в приложении курьера может смутить его, что в теории может привести к изменению его поведения (или он вообще подумает, что приложение работает с багами, и уйдет к конкуренту…но это я, конечно, утрирую).

Здесь, на мой взгляд, может помочь применение моделей машинного обучения. Для запуска теста нам всего понадобится одна группа, на которую мы раскатываем изменение на период эксперимента. Предварительно перед экспериментом отбираем признаки, на основе которых будем обучать модель, проводим тестирование, валидацию, тюним гиперпараметры при необходимости (в общем, стараемся добиться улучшения метрик качества модели).
В течение эксперимента мы параллельно предсказываем значение метрики с помощью построенной регрессионной модели машинного обучения, и эти значения будут составлять нашу контрольную группу (уникальная связка shop_id-date попадет только в одну из групп).

Схема формирования групп

Схема формирования групп

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

В общем, здесь создается широкое поле для размышлений и анализа.

В заключение

Целью этой статьи было очередное напоминание о важности валидирования ваших статкритериев перед запуском А/B‑тестов. Всегда нужно учитывать, что использование t‑теста для проверки значимости разницы средних не противоречит имеющимся данным: зависимость в них сильно искажает результаты.

Это моя первая статья на Хабре, надеюсь, что для кого-то она была полезной. Буду рад обсудить все вопросы в комментариях!) Тг на всякий @dezluvv


ссылка на оригинал статьи https://habr.com/ru/articles/859088/


Комментарии

Добавить комментарий

Ваш адрес email не будет опубликован. Обязательные поля помечены *