Анализ данных и сексизм

Недавно в фейсбуке, Эдуард Бабушкин предложил решить кейс, выдал данные по опросу и поставил задачу:

«Жду от вас:
1) любое решение в области регрессии (например, что влияет на уровень ЗП)
2) любую задачу по классификации
3) любое решение в анализе дожития с помощью регрессии Кокса
Обратите внимание 1) в конце уже сделал за вас много работы: посчитал стаж и выставил event
2) ЗП дается в национальной валюте
Работаем по боевому
Лайкните, если задача понятна» (с)

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

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

Для начала загрузим данные:

# load our data
mydata <- read.csv("hrmcase.csv", header = T, sep = ",", dec = ".", 
                   na.strings = c("", "NA", "#Н/Д"))
# what about structure
str(mydata)

## 'data.frame':    2354 obs. of  80 variables:

Итого у нас 80 переменных и 2354 наблюдений. Данный на любой вкус: даты, категориальные, числовые, ещё и разделители разные, где-то десятые точка разделяет, где-то запятая.

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

В процессе просмотра переменных выяснилось, что есть категориальные переменные (факторы) с задвоенными уровнями, как например образование:

# fix levels of education
education <- education1

levels(education)

## [1] "высшее профессиональное"    "высшее профессиональное "  
## [3] "научная степень"            "начальное профессиональное"
## [5] "основное общее"             "основное общее "           
## [7] "среднее (полное) общее"     "среднее профессиональное"

levels(education) <- c("высшее профессиональное", "высшее профессиональное", 
                       "научная степень", "начальное профессиональное", "основное общее", 
                       "основное общее", "среднее (полное) общее", "среднее профессиональное")

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

После предобработки, соберём рабочий набор данных и посмотрим при помощи функции findLinearCombos{ caret} нет ли у нас проблемных переменных:

findLinearCombos(testdata)

## $linearCombos
## $linearCombos[[1]]
## [1] 25 23 24
## 
## 
## $remove
## [1] 25

25 – это год рождения, стоило ожидать, что они вместе с возрастом и датой трудоустройства не станут сидеть тихо, а устроят шум. Итого у нас 44 переменные.

Для начала построим обычную линейную регрессию:

## Residual standard error: 62320 on 1372 degrees of freedom
## Multiple R-squared:  0.2625, Adjusted R-squared:  0.2389 
## F-statistic:  11.1 on 44 and 1372 DF,  p-value: < 2.2e-16

Значимость устраивает, доля объяснённой дисперсии (R-squared) 0,2625 – не густо, с другой стороны мы ведь не с физическими процессами работаем, где запредельная стабильность и точность. Но всё равно маловато будет (с) значит, будем улучшать.

Посмотрим на графическую статистику модели:

Residuals vs Fitted показывает нам, что имеет место некоторая зависимость от Y.

Судя по Normal Q-Q распределение не нормально, кроме того есть выбросы.

Нужно что-то менять, начнём с логарифмирования переменной отклика, нашей Y – заработной платы.

Нормально.

Построим модель ещё раз, но уже с логарифмированной Y.

## Residual standard error: 0.7893 on 1372 degrees of freedom
## Multiple R-squared:  0.3878, Adjusted R-squared:  0.3682 
## F-statistic: 19.75 on 44 and 1372 DF,  p-value: < 2.2e-16

Намного лучше, R2: 0,3878.

Но мы только начали, теперь попробуем отобрать оптимальную модель по критерию Акаике (AIC) при помощи функции stepAIC {MASS}, и снова запустим линейную регрессию.

## Residual standard error: 0.7878 on 1385 degrees of freedom
## Multiple R-squared:  0.3843, Adjusted R-squared:  0.3705 
## F-statistic: 27.89 on 31 and 1385 DF,  p-value: < 2.2e-16

Унывать рано – впереди кроссс-валидация с помощью функции train {caret}, данные на тестовых выборках:
- для полной модели:

##   RMSE       Rsquared 
##   0.8031386  0.3547105

- для AIC модели:

##   RMSE       Rsquared 
##   0.7977403  0.3636119

Запускаем тяжёлую артиллерию - градиентный бустинг, результаты на тесте:

#RMSE 0.5988630 and R2 0.6406402 

Отбор оптимальной модели происходит по RMSE – квадратный корень из среднеквадратичной ошибки. Поэтому будем ориентироваться на неё.

На очереди случайный лес – Random Forest:

# RMSE 0.6413813 and R2 0.5889833 

Что ж, градиентный бустинг уверенно лидирует, а теперь экстримальный бустинг с XGBOOST, который даёт нам:

#  test-rmse:0.582180 

И XGBOOST вырывает победу в этом состязании. Но, нет предела совершенству и я вспоминаю, что в опросе принимали участия респонденты из разных стран, а значит они могли указывать зарплату в нац.валюте. Посмотрим:

levels(fixdata$region)

##  [1] "Азербайджан"                                             
##  [2] "Армения"                                                 
##  [3] "Беларусь"                                                
##  [4] "Грузия"                                                  
##  [5] "Дальневосточный федеральный округ"                       
##  [6] "другое"                                                  
##  [7] "Казахстан"                                               
##  [8] "Латвия"                                                  
##  [9] "Молдова"                                                 
## [10] "Москва"                                                  
## [11] "Приволжский федеральный округ"                           
## [12] "Санкт-Петербург"                                         
## [13] "Северо-Западный федеральный округ (без Санкт-Петербурга)"
## [14] "Северо-Кавказский федеральный округ"                     
## [15] "Сибирский федеральный округ"                             
## [16] "Узбекистан"                                              
## [17] "Украина"                                                 
## [18] "Уральский федеральный округ"                             
## [19] "Центральный федеральный округ (без Москвы)"              
## [20] "Эстония"                                                 
## [21] "Южный федеральный округ"

Переведём всё в рубли, по текущему курсу (грубовато, но пока так). Посмотрим, что произошло с зарплатой после перевода:

Либо в опросе принимал участие кто-то из списка форбс, либо иностранные респонденты указывали зарплату в рублях, по-крайней мере некоторые.

Заранее прошу прощенья у всех иностранных респондентов, но видимо стоит оставить только представителей  Российской Федерации.

Немного поработаем с параметрами XGBOOST и получаем:

# test-rmse:0.435905

Посмотрим на топ 10 предикторов:

И это:

1. Возраст.

2. Год трудоустройства.

3. Регион.

4. Позиция.

5. Соц.пакет, предоставляемый компанией.

6. Масштаб населённого пункта.

7. Соц.сеть, в которой респондент проявляет максимальную активность.

8. Пол.

9. Отрасль компании.

10. Самообучение.

Ради интереса посмотрим на детальный портрет потенциального высокооплачиваемого респондента:

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

На очереди вторая задача – классификация.

Проверим предположение о том, что мужчины лучше проходят тесты на интеллект (да да, сексизм, но только ради науки).

Респонденты, по желанию заполняли тесты из батареи Лаборатории «Гуманитарные технологии».

Возьмём следующие шкалы:

ВЕРБАЛЬНЫЙ.IQ, ЭРУДИЦИЯ, ЧИСЛОВОЙ.IQ, ОБРАБОТКА.ИНФОРМАЦИИ

Я не сексист и интерпретировать результаты могу следующим образом:

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

А мы идём дальше и третья задача – анализ дожития.

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

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