Привет, Хабр Недавно посмотрел выступление Валерия Бабушкина, которое было опубликовано в далеком ковидном 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 = 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])

Ну тут уже в любом случае проигрыш: необходимо 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()

Что-то пошло не так...Доля наблюдаемых ошибок первого рода составляет порядка 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()

А теперь всё встало на свои места, ожидания и реальность совпали. При таком разбиении можно спокойно запускать 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()

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

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

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

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

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

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

В заключение

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

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