Продолжим работать с базой по заработной плате: загрузим необходимые библиотеки и сам файл, переименуем столбцы с точкой в названии:
from statsmodels.formula.api import ols
import pandas as pd
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()
Попробуем учесть в модели пол преподавателя. В нашей таблице переменная sex является текстовой (два значения female и male), поэтому не очень ясно, как включать её в модель. Обычно в таких случаях переменную превращают в бинарную. Но в данном случае этого делать не нужно, Python выполнит преобразования автоматически: упорядочит два значения по алфавиту, первому присвоит значение $0$, второму – значение $1$.
В предыдущих моделях заработная плата у нас измерялась как есть, в условных единицах, но для качества модели и для интерпретации будет лучше, если разброс значений зависимой переменной будет меньше. Давайте добавим переменную salary_th – заработную плату в тысячах, а затем будем её использовать в качестве зависимой переменной.
df['salary_th'] = df['salary'] / 1000
Построим модель, где предикторами являются переменные service, phd, и sex:
model1 = ols('salary_th ~ service + phd + sex', df).fit()
print(model1.summary())
В модели добавился предиктор sex[T.Male], коэффициент при котором показывает разницу в заработной плате мужчин и женщин. Так как в качестве базовой категории (значения $0$, категории, с которой мы сравниваем остальные) используется значение female, коэффициент при sex показывает, на сколько больше/меньше заработная плата у мужчин по сравнению с женщинами.
Запишем уравнение с учётом коэффициентов модели:
Известно, что для преподавателей женского пола $sex=0$, для преподавателей мужского пола $sex=1$. Учтём этот факт и запишем уравнения отдельно для каждого пола:
Видно, что в среднем, при прочих равных, заработная плата у мужчин на 8.46 тысяч выше, чем заработная плата у женщин. Это и есть коэффициент при sex[T.Male] в таблице. Теперь пойдём дальше и попробуем учесть не только влияние пола преподавателя, но и тот факт, что влияние опыта работы у мужчин и женщин неодинаково. Добавим в модель эффект взаимодействия – предиктор, представляющий собой произведение двух переменных, в данном случае произведение sex и service. В общем виде уравнение модели выглядит так:
Важно: чтобы избежать смещения оценок коэффициентов и проблемы пропущенных переменных, переменные которые входят в эффект взаимодействия, должны быть включены в модель по отдельности. Так, здесь мы не можем убрать слагаемые $\beta_1 \cdot sex$ и $\beta_2 \cdot service$.
В statsmodels переменная взаимодействия записывается через ::
model2 = ols('salary_th ~ service + phd + sex + sex:service', df).fit()
print(model2.summary())
Уравнение модели:
$$ salary = 73.6 + 18.5 \cdot sex + 0.17 \cdot service - 0.85 \cdot sex \cdot service + 1.54 \cdot phd $$Запишем уравнения модели отдельно для женщин и мужчин.
sex = 0
sex = 1
Что мы здесь видим? Во-первых, заработная плата мужчин выше, чем заработная плата женщин. Во-вторых, опыт работы у женщин оказывает положительное влияние на заработную плату, а у мужчин – отрицательное.
Можем записать уравнение, описывающее влияние переменной service в зависимости от пола, предельный эффект опыта работа на заработную плату. Предельный эффект какой-либо переменной определяется как частная производная уравнения регрессии по этой переменной:
И здесь мы опять видим, что влияние опыта работы у женщин более значимо отражается на заработной плате.
Посмотрим на значимость коэффициентов. Интересно, что при добавлении эффекта взаимодействия, соответствующий ему предиктор «оттянул» на себя значимость переменной service. Теперь сам по себе опыт работы не оказывает статистически значимого влияния на заработную плату, а вот с учётом пола преподавателя – оказывает (значим на 10% уровня значимости). Кроме того, значимое влияние на заработную плату оказывает пол преподавателя и число лет после защиты диссертации.
Теперь давайте учтём в модели статус преподавателя, его должность. Посмотрим на уникальные значения переменной rank:
df['rank'].unique()
В отличие от переменной «пол» здесь уже не два значения, а три. Как быть? Включить в модель не саму переменную 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 опять это сделает самостоятельно, но важно понимать, что именно происходит. Посмотрим на модель:
model3 = ols('salary_th ~ service + phd + sex + rank', df).fit()
print(model3.summary())
Видно, что в модель были включены не три дамми-переменные для должности, а только две. Так будет всегда: чтобы избежать строгой мультиколлинеарности (два предиктора абсолютно скоррелированы, коэффициент корреляции равен 1), в модель будет включаться на одну дамми-переменную меньше. В модель не включена переменная rank[AssocProf], следовательно, в качестве базового уровня выбрана должность доцента (AssocProf). С ней и будем сравнивать заработную плату других категорий преподавателей. Проинтерпретируем полученные результаты.
Заработная плата профессора статистически значимо отличается от заработной платы доцента: при прочих равных зарплата профессора в среднем на 33 тысячи выше зарплаты доцента.
Заработная плата обычного преподавателя значимо отличается от заработной платы доцента: при прочих равных зарплата преподавателя в среднем на 14 тысяч ниже зарплаты доцента.
Заработная плата преподавателей-мужчин в среднем на 5.5 тысяч выше, но если принимать во внимание должность преподавателя, эта разница не является статистически значимой.
Если принимать во внимание должность преподавателя, то число лет после получения степени и опыт работы не оказывают значимого влияния на заработную плату. Должность «оттягивает» на себя значимость других предикторов и оказывает решающее влияние на размер заработной платы.
Значение константы можно проинтерпретировать так: средняя заработная плата преподавателей равна 88 тысяч (без учёта пола преподавателя, числа лет опыта работы и других признаков, считаем все предикторы равными нулю).
Качество модели не очень высокое, но гораздо лучше, чем в предыдущих моделях (в двух предыдущих ноутбуках), наша последняя модель описывает 40% вариации заработной платы. Самостоятельно можете проверить модель на наличие мультиколлинеарности, гетероскедастичности и влиятельные наблюдения.