Привет! Меня зовут Николай Олигеров, я работаю продуктовым аналитиком в Яндекс Путешествиях. В этой статье я расскажу, как мы применяли PSM (Propensity Score Matching) — статистический метод, который позволяет корректно сравнивать группы, уменьшая систематические различия между ними. Подробно разберу, как выровнять группы теста и контроля с помощью PSM, расскажу о типичных ошибках (например, утечке признаков), дам практические рекомендации по сбору и выбору фич для мэтчинга, а также покажу, как валидировать полученные результаты и оценить их достоверность.

Почему A/B-тесты не всегда применимы

Если мы наблюдаем корреляцию между двумя явлениями, нам нужно знать направление корреляции. Самый лучший способ узнать его — провести A/B‑тест. Мы ставим гипотезу, изолируем тестовую группу от контрольной, воздействуем эффектом (далее по тексту — тритмент, от treatment) на тестовую группу и сравниваем результаты. Тут особо останавливаться не буду, потому что про A/B‑тесты уже много всего написано.

Но не всегда есть возможность провести честный A/B‑тест. Иногда его нельзя провести из‑за его высокой стоимости или высоких рисков. 

Одна из таких гипотез, которую мы не могли проверить в эксперименте, — какой инкрементальный GBV (gross booking value, общая стоимость бронирования) приносят новые подключённые отели. Предположим, у нас в городе N в радиусе 5 км доступны пять отелей для бронирования на Яндекс Путешествиях. При этом есть еще два отеля, которые не подключены к нашему сервису. 

Вопрос: увеличим ли мы общее количество бронирований, если подключим эти отели? Да, через них пойдёт какое‑то количество заказов. Но какое количество заказов пользователи бы всё равно оформили при отсутствии этих двух отелей, а какое количество заказов мы бы потеряли? Те заказы, которые пользователи бы не оформили при отсутствии этих объектов, — это как раз инкрементальные заказы.

Этот вопрос критически важен, потому что если новые объекты не дают инкремента, нет смысла выделять ресурсы на их подключение и нужно фокусироваться на удержании существующих, и наоборот.

Почему эту задачу нельзя было решить в классическом A/B?

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

Во‑вторых, мы не хотим, чтобы даже в рамках эксперимента подключённые к нам отели теряли часть заказов.

В‑третьих, если отели‑новички дают высокий инкремент, то проведение такого эксперимента стоило бы очень дорого. 

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

Поэтому для решения этой задачи мы и воспользовались альтернативой A/B‑тесту — PSM. 

Как работает PSM

Сразу дам дисклеймер: некоторые вещи, связанные с A/B‑тестами, буду намеренно упрощать, чтобы текст был доступен широкому кругу читателей. Например, когда говорю, что в A/B‑тестах есть соответствие между юзерами тестовой группы и контрольной (в реальности всё сложнее и всегда есть погрешность в разбиении, из‑за которой приходится смотреть на p‑value). В этой статье я фокусируюсь именно на технической реализации PSM и подводных камнях. 

Итак, в классическом A/B‑тесте мы случайным образом распределяем пользователей на тест и контроль. И A/B‑тесты работают, потому что наиболее вероятным является равномерное перемешивание пользователей. Предположим, они отличаются каким‑то одним параметром. Возьмём восемь пользователей (четырёх представим в форме треугольника и других четырёх — в форме кружка). Тогда число способов их случайного разбиения будет выглядеть таким образом: 

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

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

Изобразим это графически применимо к Яндекс Путешествиям. Предположим, что жёлтые треугольники — это юзеры с более высоким средним чеком, которым более важен комфорт. Зелёные кружки — это юзеры, у которых средний диапазон среднего чека. И красные квадраты — это те, кому важна цена в первую очередь. 

В классическом A/B‑тесте наиболее вероятное состояние после разбивки на тест/контроль было бы, как на схеме ниже. Юзерам из тестовой группы найдётся примерно такое же количество аналогичных юзеров из контрольной группы. 

Теперь посмотрим, что будет, если сравнить юзеров, которые заходили на страницы отелей‑новичков без A/B‑теста. Для этого оставим в тесте всех юзеров, а в контроле только тех, кто заходил на страницы этих отелей. И предположим, что мы подключили отели с низким средним чеком, поэтому ими будут интересоваться соответствующие юзеры.

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

В этом случае группы будут отличаться между собой наличием тритмента (в одной группе юзеры видели отели‑новички, в другой — нет), и для каждого юзера из теста найдётся эквивалентный в контроле. Поэтому эти группы, как и в A/B‑тесте, можно сравнивать между собой. А поскольку юзеры эквивалентны в группах и отличаются только наличием или отсутствием тритмента (видели отели‑новички / не видели отели‑новички), в случае получения статистически значимой разницы в метриках между тестом и контролем можно утверждать, что эта разница вызвана тритментом (в рамках p‑value).

Как это выглядит в коде

Сгенерируем тестовый семпл из 100 000 юзеров.

# Импортируем нужные библиотеки
import pandas as pd
import numpy as np
from sklearn.linear_model import LogisticRegression
from sklearn.neighbors import NearestNeighbors
import seaborn as sns
import scipy.stats as sps
from scipy import stats
import matplotlib.pyplot as plt

# Генерируем данные
n_users = 100000

# Создаём DataFrame с уникальными ID
df = pd.DataFrame({
    'user_id': range(1, n_users + 1)
})

Пусть 90% юзеров — пользователи устройств с iOS.

# Фиксируем seed
np.random.seed(15)

# 90% пользователей –— iOS
n_ios = int(0.9 * n_users)
n_android = n_users - n_ios

# Создаём массив: [1, 1, ..., 0, 0] (90% единиц, 10% нулей)
ios_flg = np.array([1] * n_ios + [0] * n_android)

# Перемешиваем, чтобы распределение было случайным
np.random.shuffle(ios_flg)

# Android –— обратное значение
android_flg = 1 - ios_flg

df['ios_flg'] = ios_flg
df['android_flg'] = android_flg

Теперь сгенерируем фичи. Часть будет случайно распределена с заданной вероятностью, другая часть будет зависеть от платформы. 

# Функция для генерации флагов фичи, если фича имеет разное распределение между iOS и Android
def set_dependent_value(df, column_name, ios_prob, android_prob):
    df[column_name] = 0
    
    # Выставляем вероятность фичи для юзеров устройств c Android в соответствии с android_prob
    mask_android = df['android_flg'] == 1
    n_android = mask_android.sum()
    df.loc[mask_android, column_name] = np.random.binomial(1, android_prob, size=n_android)
    
    # Аналогично для юзеров устройств с iOS
    mask_ios = df['ios_flg'] == 1
    n_ios = mask_ios.sum()
    df.loc[mask_ios, column_name] = np.random.binomial(1, ios_prob, size=n_ios)
    return df

# Функция для генерации рандомных флагов фичи
def set_independent_value(df, column_name, prob):
    
    n = int(len(df) * prob)
    value_flg = np.array([1] * n + [0] * (len(df) - n))
    np.random.shuffle(value_flg)

    df[column_name] = value_flg
    return df

independent_values_list = ['independent_value_1_flg', 
                         'independent_value_2_flg', 
                         'independent_value_3_flg', 
                         'independent_value_4_flg', 
                         'independent_value_5_flg', 
                         'independent_value_6_flg', 
                         'independent_value_7_flg', 
                         'independent_value_8_flg', 
                         'independent_value_9_flg', 
                         'independent_value_10_flg'
                         ]
# Генерируем флаги фич с вероятностью 10%
for value in independent_values_list:
    df = set_independent_value(df, value, 0.1)

# Генерируем флаги фич, которые зависят от iOS/Android
df = set_dependent_value(df, 'dependent_value_1_flg', 0.5, 0.7)
df = set_dependent_value(df, 'dependent_value_2_flg', 0.9, 0.1)
df = set_dependent_value(df, 'dependent_value_3_flg', 0.3, 0.6)
df = set_dependent_value(df, 'dependent_value_4_flg', 0.1, 0.2)
df = set_dependent_value(df, 'dependent_value_5_flg', 0.3, 0.05)
df = set_dependent_value(df, 'dependent_value_6_flg', 0.2, 0.4)
df = set_dependent_value(df, 'dependent_value_7_flg', 0.4, 0.1)
df = set_dependent_value(df, 'dependent_value_8_flg', 0.7, 0.2)
df = set_dependent_value(df, 'dependent_value_9_flg', 0.05, 0.1)
df = set_dependent_value(df, 'dependent_value_10_flg', 0.8, 0.9)
df = set_dependent_value(df, 'hotel_newbie_flg', 0.1, 0.2)

# Генерируем флаг того, что юзер заказал отели-новички
# Среди юзеров устройств с iOS таких 10%, среди юзеров устройств с Android — 20%
df = set_dependent_value(df, 'hotel_newbie_flg', 0.1, 0.2)
df['test_group_flg'] = df['hotel_newbie_flg']

user_features = ['dependent_value_1_flg', 
                         'dependent_value_2_flg', 
                         'dependent_value_3_flg', 
                         'dependent_value_4_flg', 
                         'dependent_value_5_flg', 
                         'dependent_value_6_flg', 
                         'dependent_value_7_flg', 
                         'dependent_value_8_flg', 
                         'dependent_value_9_flg', 
                         'dependent_value_10_flg',
                'independent_value_1_flg', 
                         'independent_value_2_flg', 
                         'independent_value_3_flg', 
                         'independent_value_4_flg', 
                         'independent_value_5_flg', 
                         'independent_value_6_flg', 
                         'independent_value_7_flg', 
                         'independent_value_8_flg', 
                         'independent_value_9_flg', 
                         'independent_value_10_flg'] # 'ios_flg', 'android_flg', 
# Поскольку фичи и флаг тестовой группы напрямую зависят от флага iOS/Android, они в датасете при обучении не будут использоваться 
# Эти флаги дают явную подсказку, куда определить юзеров — в тест или контроль, а это плохо. Ниже будет объяснение, почему

# Теперь посмотрим, сколько у нас юзеров в тестовой группе, а сколько — в контрольной:

print(f"""
Юзеров в контрольной группе: {df[df['test_group_flg'] == 0]['user_id'].count()}
Юзеров в тестовой группе: {df[df['test_group_flg'] == 1]['user_id'].count()}""")

# Результаты:
# Юзеров в контрольной группе: 89 176
# Юзеров в тестовой группе: 10 824

Здесь получилось соотношение примерно 1 к 10 (юзеры в тестовой группе / юзеры в контрольной группе). Это соотношение должно быть не меньше 1 к 5. Если оно будет меньше, то алгоритм будет хуже мэтчить пары тест/контроль. Если соотношение меньше 1 к 5, нужно сделать случайную выборку из тестовой группы. Например, если мы хотим в два раза уменьшить размер тестовой группы, это будет выглядеть вот так:

seed = 15 
subset = df[df['test_group_flg'] == 1]
sampled = subset.sample(n=len(subset) // 2, random_state=seed)
df_reduced = pd.concat([df[df['test_group_flg'] == 0], sampled])

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

Важно отслеживать именно нижнюю границу. Чем выше соотношение теста и контроля, тем лучше. То есть если на одного юзера из тестовой группы приходится 100 юзеров в контроле, тестовую группу не нужно уменьшать.

Теперь создадим датасеты для обучения:

x = df[user_features]
y = df['test_group_flg']

В y — флаг тестовой группы (тест: заходил ли юзер на страницу отеля‑новичка → y = 1; контроль: не заходил → y = 0). Выборки тех, кто заходил и не заходил на страницы отелей‑новичков, смещены друг относительно друга. В процессе обработки нужно будет получить две группы эквивалентных между собой юзеров, в одной из которых будут только те, которые заходили на страницы новичков, в другой — которые не заходили на страницы новичков.

Фичи — это то, чем пользователи могут отличаться друг от друга. В случае Яндекс Путешествий это регионы пользователя и отеля, платформа устройства (iOS / Android / web / mobile web), канал привлечения юзера (через рекламу/органику) и тому подобное В финальных версиях у меня было больше 140 фич, при этом всё хорошо работало. До того чтобы от количества фич ломался PSM, я не доходил. 

Но важно качество фич. Они не должны явно подсказывать алгоритму, кого отправить в тест, а кого — в контроль. Например, флаг того, оставлял ли юзер отзыв на отель‑новичок, — плохая фича, потому что она явно подсказывает алгоритму, какие юзеры точно будут в тестовой группе. Более того, они не должны включать в себя те метрики, которые у тестовой группы могут отличаться от контрольной (например, конверсия в заказ, средний чек или доля отменённых заказов). Если попробовать подобрать юзерам из теста пару из контроля по среднему чеку, то мы подберём в тесте тех юзеров, у которых изначально смещены метрики в ту же сторону, что и в эксперименте.

В моём применении почти все фичи были бинарными — принимали значение 1 или 0. Если фича может принимать больше двух значений, то я выставлял флаги 1 или 0 для каждого значения. Например, если нужно задать флаг интереса отелей в Москве, Санкт‑Петербурге, Нижегородской области и Татарстане, то я создавал четыре колонки: is_hotel_mow, is_hotel_spb, is_hotel_nn, is_hotel_tat. Причина этого в том, что один юзер может интересоваться отелями сразу в нескольких регионах.

Дальше нужно применить алгоритм классификации и получить propensity score — вероятность оказаться в тестовой или контрольной группе. Какого‑то строгого требования нет, можно использовать любой алгоритм классификации: catboost, random tree forest, логистическую регрессию. По моему опыту, существенной разницы от выбранного алгоритма в результате нет, гораздо важнее правильно подобрать фичи. Я использовал логистическую регрессию.

lr = LogisticRegression()
lr.fit(X, y)

pred_prob = lr.predict_proba(X)  
df['ps'] = pred_prob[:, 1] # вероятность оказаться в тесте или контроле

Качество выбранных метрик можно оценить, построив гистограмму распределения ps для тестовой группы и контрольной группы:

sns.histplot(data=df, x='ps', hue='test_group_flg')

По факту всё то множество фич для каждого юзера, которое выбрали выше, мы сводим к одному числу — ps. И чем ближе друг к другу значения ps, тем более похожи юзеры между собой. 

График выше нужен, чтобы оценить, насколько похожи между собой юзеры из теста и контроля, есть ли сильное смещение между ними. Если юзеры похожи, мы будем наблюдать картину, как на графике: ps для тестовой группы пересекается с ps контрольной группы. Это значит, что мы сможем хорошо подобрать для тестовой группы юзеров в контроле. 

А вот примеры плохого пересечения (я отдельно искусственно сгенерировал данные):

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

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

Вот ещё один пример плохого распределения: 

Тут утечки нет, но есть сильное смещение юзеров друг относительно друга (несовпадающие пики у тестовой и контрольной группы между ps 0,4 и 0,45; для юзеров из тестовой группы с ps около 0,55 нет соответствующих юзеров из контрольной группы). Если с утечкой мы можем работать, добавляя или убирая какие‑то фичи, то такая форма распределения, по моему опыту, вызвана именно сильным смещением юзеров в тестовой и контрольной группах. 

Например, в тестовой группе оказалось 50% юзеров, которые едут в Санкт‑Петербург, и 25% юзеров, которые едут в Москву, а в контрольной группе оказалось 50% юзеров, которые едут в Москву, и 25% юзеров, которые едут в Санкт‑Петербург. Причина такого распределения может быть такая: в Москве очень мало отелей‑новичков, а в Санкт‑Петербурге — много. Тут добавлением и удалением фич мы вряд ли сможем как‑то улучшить качество выборки. Вариантов в этом случае два: 

• Считать as is, приняв тот факт, что результат у нас будет очень сильно смещённым. 

• Взять на проверку другую гипотезу. Может, стоит посмотреть на юзеров только по тем гео, где нет такого смещения по просмотру отелей‑новичков, например Краснодарский край.

Дальше нам нужно каждому юзеру из тестовой группы подобрать наиболее похожего из контрольной группы. Степень похожести определяем по ps: чем ближе у юзеров ps, тем больше выбранных фич у них совпадает. Я использовал алгоритм K‑ближайших соседей (в примере подбираем любых, ещё не мэтчим пары тест/контроль):

caliper = np.std(df.ps) * 0.25 # максимальный радиус, внутри которого ищем ближайшего соседа — 0,25 стандартного отклонения ps
n_neighbors = 10 # ищем 10 ближайших соседей внутри радиуса

knn = NearestNeighbors(n_neighbors=n_neighbors, radius=caliper)

ps = df[['ps']]  
knn.fit(ps)

distances, neighbor_indexes = knn.kneighbors(ps)

Размер радиуса выбирается исходя из того, насколько хорошо пересекаются ps для тестовой и контрольной группы. Если пересечение хорошее, достаточно 0,25 стандартного отклонения. Если пересечение плохое, можно взять одно стандартное отклонение (но от этого может упасть качество мэтчинга).

По поводу того, сколько брать ближайших соседей (n_neighbors): чем больше, тем лучше, но с увеличением количества соседей увеличивается время расчёта. От 10 до 20 соседей — оптимальное значение. Если берём меньше 10, то это сильно снижает количество пар теста/контроля на выходе. Если взять больше 20, количество пар теста/контроля при увеличении n особо не будет увеличиваться.

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

Выше я подобрал K‑ближайших соседей для каждого юзера по ps. Теперь подберём для каждого юзера из тестовой группы юзера в контроле:

matched_index = set()

for current_index, row in df.iterrows():  # проходимся по всем строкам датафрейма
    if row.test_group_flg == 0:  # этот юзер в контрольной группе
        df.loc[current_index, 'matched'] = np.nan # так как ищем пару для юзеров из теста, просто указываем, что пары для юзеров из контроля нет
    else: # если юзер в контрольной группе
        for idx in neighbor_indexes[current_index, :]: # итеративно проходимся по всем подобранным соседям
            if (current_index != idx) and (df.loc[idx].test_group_flg == 0) and (idx not in matched_index):
            # Первое условие — проверяем, что это не тот же юзер, которому мы ищем пару
            # Второе условие — это юзер из контроля
            # Третье условие — что для данного юзера из контроля мы ещё не нашли пару из теста
                df.loc[current_index, 'matched'] = idx  # фиксируем индекс смэтченного юзера
                matched_index.add(idx)  # добавляем в множество смэтченного юзера, чтобы не допустить выборки с повторениями
                break # останавливаем поиск пары для текущего юзера из теста, т. к. пару мы ему нашли

Посмотрим, сколько юзеров получилось смэтчить:

print('Кол-во юзеров в тестовой группе:', len(df[df.test_group_flg==1]))
print('Кол-во смэтченных юзеров:', len(matched_index))

Кол-во юзеров в тестовой группе: 10 824
Кол-во смэтченных юзеров: 8502

То есть 2322 юзерам алгоритм не нашёл пары в контроле, и мы их удаляем из рассмотрения. 

Теперь объединим всё в один датафрейм: 

treatment_matched = df.dropna(subset=['matched'])  # берём только смэтченных юзеров из теста

control_matched_idx = treatment_matched.matched.astype(int) # получаем индексы пар контроля для тестовой группы
control_matched = df.loc[control_matched_idx, :]  # получаем датафрейм со смэтченными юзерами из контроля

df_matched = pd.concat([treatment_matched, control_matched]).reset_index(drop=True) # объединяем датафреймы с юзерами из контроля и теста

Дальше по этим юзерам можно считать метрики, как в обычном A/B‑тесте. 

Важный момент: алгоритм выше — детерминированный. То есть если мы запустим его 100 раз, получим 100 одинаковых результатов. Это отличает его от обычных A/B‑тестов, где в случае многократных случайных разбиений на тест/контроль мы периодически будем получать ложноположительные результаты. 

Валидация результатов

После мэтчинга надо обязательно убедиться, что пары теста/контроля подобраны верно (между юзерами из подобранного теста и подобранного контроля нет смещения). Ниже приведены способы валидации, и я рекомендую при применении метода использовать их все. Иначе есть риск получить неверные результаты. 

Баланс ковариантов

Если мы случайным образом поделим юзеров на две группы и посчитаем, как эти две группы отличаются по сотне метрик, то с p‑value = 5% примерно в 5% случаев мы получим ложноположительные прокрасы. Это происходит потому, что при правильном случайном разбиении на две группы p‑value метрик должен быть распределён равномерно. 

Посмотрим, как были распределены p‑value метрик (коварианты) до балансировки групп через PSM и после балансировки:

# Функция для расчёта p-value групп теста/контроля
def calculate_t_test(test, control):
    return {"p value t-test": sps.ttest_ind(np.asarray(list(control)), np.asarray(list(test))).pvalue,
            "test value": np.mean(np.asarray(list(test))),
            "ctrl value": np.mean(np.asarray(list(control)))}


# Функция для проведения множественного тестирования
def multi_calculate_t_test(df, field_names, test_group_name):
    output_list = []
    for _ in range(len(field_names)):
        stat_values = calculate_t_test(df[(df[test_group_name] == 1)][field_names[_]],
                                       df[(df[test_group_name] == 0)][field_names[_]])
        output_list.append([field_names[_], 
                            stat_values['p value t-test'],
                            stat_values['test value'],
                            stat_values['ctrl value'],])

    return pd.DataFrame(output_list, columns=['Metric names', 
                                              'p value t-test',
                                              'test value', 
                                              'ctrl value'])

# Распределение p-value метрик до выравнивания групп через PSM
output_df_before_psm = multi_calculate_t_test(df, 
                                              user_features, 
                                              'test_group_flg').sort_values(by='p value t-test')
fig, ax = plt.subplots(figsize=(9, 6))
output_df_before_psm['p value t-test'].hist(bins=10, alpha=0.7, label='p value t-test', ax=ax)

ax.set_xlabel('Значение p-value', fontsize=12)
ax.set_ylabel('Частота', fontsize=12)
ax.set_title('Распределение p-value до PSM', fontsize=14)
ax.legend()

plt.tight_layout()
plt.show()

А это распределение p‑value метрик после выравнивания групп через PSM:

output_df_after_psm = multi_calculate_t_test(df_matched, 
                                             user_features, 
                                             'hotel_newbie_flg').sort_values(by='p value t-test')
fig, ax = plt.subplots(figsize=(9, 6))
output_df_after_psm['p value t-test'].hist(bins=10, alpha=0.7, label='p value t-test', ax=ax)

ax.set_xlabel('Значение p-value', fontsize=12)
ax.set_ylabel('Частота', fontsize=12)
ax.set_title('Распределение p-value после PSM', fontsize=14)
ax.legend()

plt.tight_layout()
plt.show()

То, что распределение было неравномерным до PSM, говорит о том, что было смещение между выбранными группами. А то, что распределение p‑value после PSM стало равномерным, — показатель того, что мы сбалансировали группы, оставив похожих юзеров в тесте и контроле. 

Чтобы убедиться, что распределение p‑value действительно равномерно, можно выполнить тест на соответствие эмпирического распределения теоретическому равномерному. Проверять будем через хи‑квадрат. Ставим нулевую гипотезу: распределение p‑value равномерно. Альтернативная гипотеза — распределение p‑value ковариантов неравномерно. P‑value для статистически значимого результата — 5%.

Примечание: тут фигурирует два разных p‑value. Одно — p‑value ковариантов, которое мы посчитали для метрик, второе — которое получается при проверке распределения p‑value ковариантов на равномерность.

Сначала проверим для распределения p‑value метрик до PSM:

# Разбиваем данные на бины
num_bins = 10
counts, bin_edges = np.histogram(output_df_before_psm['p value t-test'], bins=num_bins)

# Ожидаемые частоты при равномерном распределении
expected_counts = sum(counts) / num_bins

# Тест хи‑квадрат
chi2_stat, p_value = stats.chisquare(counts, f_exp=expected_counts)

print(f"Хи‑квадрат статистика: {chi2_stat:.4f}")
print(f"p_value: {p_value:.4f}")

Хи‑квадрат статистика: 50.0000
p‑value: 0.0000

Получили, что до выравнивания групп через PSM p‑value метрик распределение неравномерно с p‑value = 0. Делаем то же самое для групп после выравнивания через PSM.

num_bins = 10
counts, bin_edges = np.histogram(output_df_after_psm['p value t-test'], bins=num_bins)

# Ожидаемые частоты при равномерном распределении
expected_counts = sum(counts) / num_bins

# Тест хи‑квадрат
chi2_stat, p_value = stats.chisquare(counts, f_exp=expected_counts)

print(f"Хи‑квадрат статистика: {chi2_stat:.4f}")
print(f"p_value: {p_value:.4f}")

Хи‑квадрат статистика: 4.0000
p‑value: 0.9114

Здесь мы получили p‑value 0,91. Мы не можем отвергнуть нулевую гипотезу, что распределение p‑value ковариантов равномерно. Значит, мы убрали смещение в группах теста/контроля по выбранным метрикам с помощью PSM. 

A/A-тесты

Второй способ проверить, что PSM подбирает пары теста/контроля верно, — провести A/A‑тест. Для этого нужно взять случайную выборку юзеров (5–10% от общего количества) и подобрать им пару через PSM. Затем для этих групп теста/контроля провести A/A‑тест: через бутстреп делать многократные случайные выборки с повторениями, для каждой выборки считать p‑value по какой‑то из метрик. Потом смотреть на распределение полученных p‑value: распределение должно получиться равномерным. В Яндекс Путешествиях мы проверяли распределение p‑value для таких метрик, как GBV неотменённых заказов / средняя сумма неотменённого заказа / конверсия в заказ / конверсия в отмену заказа. 

Во время проведения A/A‑тестов мы заметили особенность: из изначального датасета нужно убирать юзеров, у которых много незаполненных фич. Если этого не делать, то итоговый датасет со смэтченными юзерами не проходил валидацию по балансу ковариантов. В задаче по инкременту отелей‑новичков, например, оказалось, что нужно иметь минимум 10 ненулевых фич (значение фичи > 0) для каждого юзера.

Также можно взять изначально смещённых юзеров. Например, берём юзеров с высоким средним чеком относительно генеральной совокупности, подбираем для них пары с помощью PSM и проводим A/A‑тест. Это покажет, насколько алгоритм с выбранными фичами может корректно подбирать смещённых юзеров из генеральной совокупности. 

Сравнение с результатами уже проведённых A/B-тестов

Также в Яндекс Путешествиях мы сравнили результаты уже проведённых A/B‑тестов с PSM. Алгоритм проверки следующий:

1. Берём проведённый A/B‑тест с сильным эффектом (ключевые метрики должны статзначимо прокраситься более чем на 1%).

2. Берём случайную выборку из тестовой группы.

3. Берём контрольную группу.

4. Ищем пары через PSM.

5. Считаем метрики для подобранных пар теста/контроля.

6. Сравниваем с результатами A/B‑теста.

В Яндекс Путешествиях мы сравнили результаты трёх экспериментов с результатами PSM, и везде была одинаковая картина. Если метрики прокрашиваются более чем на 1%, то в PSM метрики тоже статистически значимо прокрашиваются. При этом в PSM дельта изменения метрик не совпадает (потому что мы берём часть юзеров, которые могут быть смещены относительно генеральной совокупности), но направление / факт того, прокрасится или не прокрасится метрика, совпадает.

Каким бы хорошим ни был PSM, он не заменяет A/B‑тесты. У него есть несколько существенных ограничений:

1. Поскольку PSM считает изменение метрик не по всей генеральной совокупности, а по смещённой выборке, то результаты для всей генеральной совокупности могут быть иными. Поэтому PSM стоит применять только для тех экспериментов, где эффект сильный. Если мы получим сильный эффект на смещённой выборке, тогда, скорее всего, направление эффекта на генеральной совокупности будет тем же (этот результат получили на сравнении результатов реальных экспериментов и PSM). 

2. Группы должны быть достаточно большими, и эффект должен быть сильным. Группа юзеров, которых мы оставляем в тесте/контроле, очень сильно уменьшается из‑за того, что:

   • мы выкидываем юзеров, про которых мало знаем / которые мало активны;

   • для балансировки тех юзеров, которых знаем (чтобы корректно подобрать тест/контроль) изначальное соотношение теста/контроля должно быть 1 к 10, то есть даже количество тех юзеров, которых хорошо знаем, уменьшается в среднем в 10 раз. 

В реальности у меня от изначальных групп оставалось в 50–80 раз меньше юзеров, чем в начале. Поэтому, как указано в первом пункте, эффект получается сильно смещённым на тех юзеров, которых мы хорошо знаем и которым мы смогли подобрать пару. Также из‑за уменьшения количества юзеров в итоговой выборке хуже красятся метрики, потому что в знаменателе стандартной ошибки стоит корень из числа юзеров. Меньше юзеров — выше стандартная ошибка.

Из положительных моментов отмечу то, что в PSM невозможен p‑хакинг (если использовать детерминистический алгоритм классификации, как логистическая регрессия). Если в A/A‑ или A/B‑тестах можно много раз перемешивать группы, и рано или поздно по ряду метрик мы получим статистически значимую разницу просто из‑за ошибки первого рода, то в PSM такое не произойдёт. В нём мы не делим пользователей случайно, а подбираем им пары по детерминистическому алгоритму.

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

Ввиду всех этих ограничений PSM не может служить заменой честным A/B‑тестам. Однако в ситуациях, где невозможно провести честный A/B‑тест, но нужно хотя бы приблизительно узнать направление эффекта, PSM может сильно выручить.