Когнитивные технологии

Алла Тамбовцева, НИУ ВШЭ

Фиктивные (дамми) переменные и эффекты взаимодействия в линейных моделях

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

In [1]:
from statsmodels.formula.api import ols 
import pandas as pd
In [2]:
df = pd.read_csv("https://vincentarelbundock.github.io/Rdatasets/csv/carData/Salaries.csv")

# опять переименуем столбцы с точкой
cols = list(df.columns)
cols[3], cols[4] = 'phd', 'service'
df.columns = cols

# посмотрим
df.head()
Out[2]:
Unnamed: 0 rank discipline phd service sex salary
0 1 Prof B 19 18 Male 139750
1 2 Prof B 20 16 Male 173200
2 3 AsstProf B 4 3 Male 79750
3 4 Prof B 45 39 Male 115000
4 5 Prof B 40 41 Male 141500

Попробуем учесть в модели пол преподавателя. В нашей таблице переменная sex является текстовой (два значения female и male), поэтому не очень ясно, как включать её в модель. Обычно в таких случаях переменную превращают в бинарную. Но в данном случае этого делать не нужно, Python выполнит преобразования автоматически: упорядочит два значения по алфавиту, первому присвоит значение $0$, второму – значение $1$.

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

In [3]:
df['salary_th'] = df['salary'] / 1000

Построим модель, где предикторами являются переменные service, phd, и sex:

In [4]:
model1 = ols('salary_th ~ service + phd + sex', df).fit()
print(model1.summary())
                            OLS Regression Results                            
==============================================================================
Dep. Variable:              salary_th   R-squared:                       0.195
Model:                            OLS   Adj. R-squared:                  0.189
Method:                 Least Squares   F-statistic:                     31.75
Date:                Fri, 04 Jan 2019   Prob (F-statistic):           2.13e-18
Time:                        23:42:59   Log-Likelihood:                -1873.8
No. Observations:                 397   AIC:                             3756.
Df Residuals:                     393   BIC:                             3772.
Df Model:                           3                                         
Covariance Type:            nonrobust                                         
===============================================================================
                  coef    std err          t      P>|t|      [0.025      0.975]
-------------------------------------------------------------------------------
Intercept      82.8759      4.801     17.264      0.000      73.438      92.314
sex[T.Male]     8.4571      4.656      1.816      0.070      -0.697      17.611
service        -0.6498      0.254     -2.558      0.011      -1.149      -0.150
phd             1.5528      0.256      6.062      0.000       1.049       2.056
==============================================================================
Omnibus:                       14.548   Durbin-Watson:                   1.888
Prob(Omnibus):                  0.001   Jarque-Bera (JB):               15.448
Skew:                           0.425   Prob(JB):                     0.000442
Kurtosis:                       3.460   Cond. No.                         156.
==============================================================================

Warnings:
[1] Standard Errors assume that the covariance matrix of the errors is correctly specified.

В модели добавился предиктор sex[T.Male], коэффициент при котором показывает разницу в заработной плате мужчин и женщин. Так как в качестве базовой категории (значения $0$, категории, с которой мы сравниваем остальные) используется значение female, коэффициент при sex показывает, на сколько больше/меньше заработная плата у мужчин по сравнению с женщинами.

Запишем уравнение с учётом коэффициентов модели:

$$ salary = 82.9 + 8.46 \cdot sex - 0.65 \cdot service + 1.56 \cdot phd $$

Известно, что для преподавателей женского пола $sex=0$, для преподавателей мужского пола $sex=1$. Учтём этот факт и запишем уравнения отдельно для каждого пола:

$$ salary_{female} = 82.9 - 0.65 \cdot service + 1.56 \cdot phd $$$$ salary_{male} = 82.9 + 8.46 - 0.65 \cdot service + 1.56 \cdot phd $$$$ salary_{male} = 91.36 - 0.65 \cdot service + 1.56 \cdot phd $$

Видно, что в среднем, при прочих равных, заработная плата у мужчин на 8.46 тысяч выше, чем заработная плата у женщин. Это и есть коэффициент при sex[T.Male] в таблице. Теперь пойдём дальше и попробуем учесть не только влияние пола преподавателя, но и тот факт, что влияние опыта работы у мужчин и женщин неодинаково. Добавим в модель эффект взаимодействия – предиктор, представляющий собой произведение двух переменных, в данном случае произведение sex и service. В общем виде уравнение модели выглядит так:

$$ salary = \beta_0 + \beta_1 \cdot sex + \beta_2 \cdot service + \beta_3 \cdot sex \cdot service + \beta_4 \cdot phd. $$

Важно: чтобы избежать смещения оценок коэффициентов и проблемы пропущенных переменных, переменные которые входят в эффект взаимодействия, должны быть включены в модель по отдельности. Так, здесь мы не можем убрать слагаемые $\beta_1 \cdot sex$ и $\beta_2 \cdot service$.

В statsmodels переменная взаимодействия записывается через ::

In [5]:
model2 = ols('salary_th ~ service + phd + sex +  sex:service', df).fit()
print(model2.summary())
                            OLS Regression Results                            
==============================================================================
Dep. Variable:              salary_th   R-squared:                       0.201
Model:                            OLS   Adj. R-squared:                  0.193
Method:                 Least Squares   F-statistic:                     24.60
Date:                Fri, 04 Jan 2019   Prob (F-statistic):           3.47e-18
Time:                        23:42:59   Log-Likelihood:                -1872.4
No. Observations:                 397   AIC:                             3755.
Df Residuals:                     392   BIC:                             3775.
Df Model:                           4                                         
Covariance Type:            nonrobust                                         
=======================================================================================
                          coef    std err          t      P>|t|      [0.025      0.975]
---------------------------------------------------------------------------------------
Intercept              73.5908      7.385      9.965      0.000      59.072      88.110
sex[T.Male]            18.5159      7.659      2.418      0.016       3.458      33.574
service                 0.1696      0.557      0.305      0.761      -0.925       1.265
sex[T.Male]:service    -0.8473      0.513     -1.652      0.099      -1.856       0.161
phd                     1.5412      0.256      6.028      0.000       1.039       2.044
==============================================================================
Omnibus:                       15.188   Durbin-Watson:                   1.888
Prob(Omnibus):                  0.001   Jarque-Bera (JB):               16.321
Skew:                           0.430   Prob(JB):                     0.000286
Kurtosis:                       3.496   Cond. No.                         302.
==============================================================================

Warnings:
[1] Standard Errors assume that the covariance matrix of the errors is correctly specified.

Уравнение модели:

$$ salary = 73.6 + 18.5 \cdot sex + 0.17 \cdot service - 0.85 \cdot sex \cdot service + 1.54 \cdot phd $$

Запишем уравнения модели отдельно для женщин и мужчин.

sex = 0

$$ salary_{female} = 73.6 + 0.17 \cdot service + 1.54 \cdot phd $$

sex = 1

$$ salary_{male} = 73.6 + 18.5 + 0.17 \cdot service - 0.85 \cdot service + 1.54 \cdot phd $$$$ salary_{male} = 92.1 - 0.68 \cdot service + 1.54 \cdot phd $$

Что мы здесь видим? Во-первых, заработная плата мужчин выше, чем заработная плата женщин. Во-вторых, опыт работы у женщин оказывает положительное влияние на заработную плату, а у мужчин – отрицательное.

Можем записать уравнение, описывающее влияние переменной service в зависимости от пола, предельный эффект опыта работа на заработную плату. Предельный эффект какой-либо переменной определяется как частная производная уравнения регрессии по этой переменной:

$$ \frac{\partial salary}{\partial service} = 0.17 \cdot service - 0.85 \cdot sex $$

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

Посмотрим на значимость коэффициентов. Интересно, что при добавлении эффекта взаимодействия, соответствующий ему предиктор «оттянул» на себя значимость переменной service. Теперь сам по себе опыт работы не оказывает статистически значимого влияния на заработную плату, а вот с учётом пола преподавателя – оказывает (значим на 10% уровня значимости). Кроме того, значимое влияние на заработную плату оказывает пол преподавателя и число лет после защиты диссертации.

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

In [6]:
df['rank'].unique()
Out[6]:
array(['Prof', 'AsstProf', 'AssocProf'], dtype=object)

В отличие от переменной «пол» здесь уже не два значения, а три. Как быть? Включить в модель не саму переменную rank, а набор фиктивных переменных (дамми-переменных), которые являются бинарными.

Так, в нашем случае вместо rank будут созданы три дамми-переменных: rank[AssocProf], rank[T.AsstProf] и rank[T.Prof]. Проиллюстрируем их смысл на небольшом фрагменте таблицы:

            rank    rank[AssocProf]    rank[T.AsstProf]    rank[T.Prof]
            Prof                  0                   0               1  
            AssocProf             1                   0               0
            AsstProf              0                   1               0


На самом деле, нет необходимости создавать дамми-переменные самим, Python опять это сделает самостоятельно, но важно понимать, что именно происходит. Посмотрим на модель:

In [7]:
model3 = ols('salary_th ~ service + phd + sex + rank', df).fit()
print(model3.summary())
                            OLS Regression Results                            
==============================================================================
Dep. Variable:              salary_th   R-squared:                       0.402
Model:                            OLS   Adj. R-squared:                  0.394
Method:                 Least Squares   F-statistic:                     52.51
Date:                Fri, 04 Jan 2019   Prob (F-statistic):           1.29e-41
Time:                        23:42:59   Log-Likelihood:                -1814.9
No. Observations:                 397   AIC:                             3642.
Df Residuals:                     391   BIC:                             3666.
Df Model:                           5                                         
Covariance Type:            nonrobust                                         
====================================================================================
                       coef    std err          t      P>|t|      [0.025      0.975]
------------------------------------------------------------------------------------
Intercept           89.5964      4.891     18.318      0.000      79.980      99.212
sex[T.Male]          5.4991      4.035      1.363      0.174      -2.433      13.431
rank[T.AsstProf]   -13.8866      4.333     -3.205      0.001     -22.406      -5.367
rank[T.Prof]        33.0534      3.701      8.932      0.000      25.778      40.329
service             -0.3738      0.221     -1.693      0.091      -0.808       0.060
phd                  0.2658      0.248      1.072      0.284      -0.222       0.753
==============================================================================
Omnibus:                       39.686   Durbin-Watson:                   1.769
Prob(Omnibus):                  0.000   Jarque-Bera (JB):               60.293
Skew:                           0.668   Prob(JB):                     8.08e-14
Kurtosis:                       4.363   Cond. No.                         178.
==============================================================================

Warnings:
[1] Standard Errors assume that the covariance matrix of the errors is correctly specified.

Видно, что в модель были включены не три дамми-переменные для должности, а только две. Так будет всегда: чтобы избежать строгой мультиколлинеарности (два предиктора абсолютно скоррелированы, коэффициент корреляции равен 1), в модель будет включаться на одну дамми-переменную меньше. В модель не включена переменная rank[AssocProf], следовательно, в качестве базового уровня выбрана должность доцента (AssocProf). С ней и будем сравнивать заработную плату других категорий преподавателей. Проинтерпретируем полученные результаты.

  1. Заработная плата профессора статистически значимо отличается от заработной платы доцента: при прочих равных зарплата профессора в среднем на 33 тысячи выше зарплаты доцента.

  2. Заработная плата обычного преподавателя значимо отличается от заработной платы доцента: при прочих равных зарплата преподавателя в среднем на 14 тысяч ниже зарплаты доцента.

  3. Заработная плата преподавателей-мужчин в среднем на 5.5 тысяч выше, но если принимать во внимание должность преподавателя, эта разница не является статистически значимой.

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

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

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