Устойчивая красота неприличных моделей

— Не могли бы вы построить нам статистическую модель?
— С удовольствием. Можно посмотреть на ваши исторические данные?
— Данных у нас ещё нет. Но модель всё равно нужна.

Знакомый диалог, не правда ли? Далее возможны два варианта развития событий:

A. «Тогда приходите, когда появятся данные.» Вариант рассматриваться не будет как тривиальный.
Б. «Расскажите, какие факторы по вашему мнению наиболее важны.» Остаток статьи про это.

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

Откуда такое странное название

Из статьи The Robust Beauty of Improper Linear Models in Decision Making товарища Robyn M. Dawes, который в свою очередь ссылается на работы предшественников.

Приличные линейные модели (proper linear models) — те, в которых веса назначаются таким образом, чтобы полученная линейная комбинация предсказывала интересующую величину оптимальным образом. Примером может служить обычная линейная регрессия, подогнанная методом наименьших квадратов. [...]

Неприличные линейные модели (improper linear models) — те, в которых веса определяются неоптимально. Например, назначаются на основании интуиции, предыдущего опыта или приравниваются единице.

Обычная линейная модель выглядит так:

ŷ = β₀ + β₁x₁ + β₂x₂ +… + β_mx_m

Здесь ŷ — оценка зависимой переменной, β₀ — свободный член, β₁...β_m — коэффициенты регрессии, x₁...x_m — независимые переменные.

По причине отсутствия данных для определения значений β₁...β_m мы поступим волюнтаристски: назначим им значения +1 (влияет положительно), -1 (влияет отрицательно) и 0 (не влияет).

Притянутый за уши пример

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

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

• Предсказание выживания пассажиров Титаника при помощи Azure Machine Learning
• Титаник на Kaggle: вы не дочитаете этот пост до конца

Участники состязаний

• ImproperModel — регрессия с единичными весами (тема этой статьи)
• ProperModel — логистическая регрессия с теми же независимыми переменными
• LogisticRegression — логистическая регрессия со всеми независимыми переменными
• RandomForest — один из широко распространённых современных алгоритмов

Как будем сравнивать

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

Чудес не бывает и совсем без каких-то знаний толковую модель не построить. В данном случае таким знанием будет галантность — женщин и детей при посадке в шлюпки пропускали вперёд. На основании этого априорного знания для регрессии с единичными весами назначаем переменной female вес +1, а переменной age вес -1. Кроме того, пассажирам с билетами второго и третьего класса, у которых каюты в трюме, дальше бежать до шлюпок, поэтому переменной pclass назначаем вес -1.

Чтобы два раза не вставать, тут же проведём нормализацию — поделим на стандартное отклонение в обучающей выборке. Это нужно для приведения независимых переменных к сравнимой шкале. Можно было поделить на диапазон и получить похожий результат.

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

Повторим эксперимент тысячу раз и посмотрим на результаты. Качество модели будем измерять с помощью коэффициента Gini (G1).

Код на R с пошаговыми комментариями

install.packages("pROC")
install.packages("randomForest")
library(pROC)
library(randomForest)

Читаем исходные данные, также доступны тут: titanic3.csv

data <- read.csv('~/habr/unit_weights/titanic3.csv')

Зависимая переменная survived

data$survived <- as.factor(data$survived)

Исключаем переменные, у которых слишком много разных значений для нормального анализа

data$name <- NULL
data$ticket <- NULL
data$cabin <- NULL
data$home.dest <- NULL
data$embarked <- NULL

Не используем переменные из будущего

data$boat <- NULL
data$body <- NULL

Заменяем пропущенные значения средним

data$age[is.na(data$age)] <- mean(data$age, na.rm=TRUE)
data$fare[is.na(data$fare)] <- mean(data$fare, na.rm=TRUE)

Преобразуем пол в индикаторную переменную

data$female <- 0
data$female[which(data$sex == 'female')] <- 1
data$sex <- NULL

Смотрим на то, что осталось:

 survived     pclass           age              sibsp       
 0:809    Min.   :1.000   Min.   : 0.1667   Min.   :0.0000  
 1:500    1st Qu.:2.000   1st Qu.:22.0000   1st Qu.:0.0000  
          Median :3.000   Median :29.8811   Median :0.0000  
          Mean   :2.295   Mean   :29.8811   Mean   :0.4989  
          3rd Qu.:3.000   3rd Qu.:35.0000   3rd Qu.:1.0000  
          Max.   :3.000   Max.   :80.0000   Max.   :8.0000  
     parch            fare             female     
 Min.   :0.000   Min.   :  0.000   Min.   :0.000  
 1st Qu.:0.000   1st Qu.:  7.896   1st Qu.:0.000  
 Median :0.000   Median : 14.454   Median :0.000  
 Mean   :0.385   Mean   : 33.295   Mean   :0.356  
 3rd Qu.:0.000   3rd Qu.: 31.275   3rd Qu.:1.000  
 Max.   :9.000   Max.   :512.329   Max.   :1.000

Повторяем эксперимент тысячу раз

im.gini = NULL
pm.gini = NULL
lr.gini = NULL
rf.gini = NULL
set.seed(42)
for (i in 1:1000) {

Разбиваем пассажиров на две выборки — 70% отправляем в обучающую выборку, на остальных 30% будем мерять качество модели.

  data$random_number <- runif(nrow(data),0,1)
  development <- data[ which(data$random_number > 0.3), ]
  holdout     <- data[ which(data$random_number <= 0.3), ]
  development$random_number <- NULL
  holdout$random_number <- NULL

Модель с единичными весами

  beta_pclass <- -1/sd(development$pclass)
  beta_age    <- -1/sd(development$age   )
  beta_female <-  1/sd(development$female)
  im.score <- beta_pclass*holdout$pclass + beta_age*holdout$age + beta_female*holdout$female 
  im.roc <- roc(holdout$survived, im.score)
  im.gini[i] <- 2*im.roc$auc-1

Обычная модель — логистическая регрессия с теми же независимыми переменными

  pm.model = glm(survived~pclass+age+female, family=binomial(logit), data=development)
  pm.score <- predict(pm.model, holdout, type="response")
  pm.roc <- roc(holdout$survived, pm.score)
  pm.gini[i] <- 2*pm.roc$auc-1

Логистическая регрессия со всеми переменными

  lr.model = glm(survived~., family=binomial(logit), data=development)
  lr.score <- predict(lr.model, holdout, type="response")
  lr.roc <- roc(holdout$survived, lr.score)
  lr.gini[i] <- 2*lr.roc$auc-1

Всеми любимый (и не без оснований) RandomForest

  rf.model <- randomForest(survived~., development)
  rf.score <- predict(rf.model, holdout, type = "prob")
  rf.roc <- roc(holdout$survived, rf.score[,1])
  rf.gini[i] <- 2*rf.roc$auc-1
}

Выводим результаты

bpd<-data.frame(ImproperModel=im.gini, ProperModel=pm.gini, LogisticRegression=lr.gini, RandomForest=rf.gini)
png('~/habr/unit_weights/auc_comparison.png', height=700, width=400, res=120, units='px')
boxplot(bpd, las=2, ylab="Gini", ylim=c(0,1), par(mar=c(9,5,1,1)+ 0.1), col=c("red","green","royalblue2","brown"))
dev.off()
mean(im.gini)
mean(pm.gini)
mean(lr.gini)
mean(rf.gini)
mean(im.gini)/mean(rf.gini)
mean(pm.gini)/mean(rf.gini)
mean(lr.gini)/mean(rf.gini)
mean(rf.gini)/mean(rf.gini)

Результаты

Модель	Gini	В процентах от лучшего
ImproperModel	0.639	90.4%
ProperModel	0.667	94.3%
LogisticRegression	0.679	96.0%
RandomForest	0.707	100%

Далекоидущие выводы

Пессимист скажет, что модель с единичными весами в этом сравнении заняла последнее место, и это верно. Но так же верно, что она набрала 90% от лучшего результата без использования исторических данных, и отстала от обычной логистической регрессии с теми же независимыми переменными всего на 4%.

Почему так происходит

В статье с необычным названием Estimating Coefficients in Linear Models: It Don't Make No Nevermind её автор (Ховард Вайнер) приводит такую теорему равных весов:

Если k линейно независимых переменных x_i (i = 1, ..., k) с нулевым средним значением и единичной дисперсией используются для предсказания переменной y, которая также отмасштабирована до нулевого среднего и единичной дисперсии, и стандартизованные значения коэффициентов регрессии по методу наименьших квадратов β_i (i = 1, ..., k) равномерно распределены на интервале [0.25, 0.75], то при переходе к равным весам (0.5) матожидание уменьшения доли дисперсии зависимой переменной, объясняемой моделью, будет меньше k/96. Потери ещё меньше, если x_i коррелированы между собой.

В приведённом выше примере регрессия логистическая, но эффект всё равно виден.

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