Множественная линейная регрессия
Андрей Ловаков
17.03.2017
Установка пакетов и рабочей директории
Установим необходимые для занятия пакеты (это действие необходимо, только если пакеты ещё не установлены на ваш компьютер) (для установки выполните следующую строку без символа # вначале).
# install.packages(c('foreign', 'ggplot2', 'visreg'))Подгружаем загруженные пакеты для их использования в текущей сессии работы R.
library(foreign)
library(ggplot2)
library(visreg)Устанавливаем рабочую директорию: через меню или командой: Session -> Set working directory -> Choose directory
#setwd('...')
#(вместо '...' укажите путь к папке, в которой хранится файл)Открываем файл с данными child_data.sav. Он в формате sav (формат, используемый в SPSS), поэтому его можно открыть с помощью функции read.spss из пакета foreign.
df <- read.spss("child_data.sav", use.value.labels = T, to.data.frame = T, use.missings = T)Множественная регрессия (метод наименьших квадратов, МНК)
Умение читать может зависить не только от одного фактора. Для улучшение качества и предсказательной силы модели добавим в неё ещё один предиктор - возраст, т.к. можно предположить, что более взрослые дети также могут читать лучше.
model2 <- lm(READ_AB ~ MEM_SPAN + AGE, data=df)
summary(model2)##
## Call:
## lm(formula = READ_AB ~ MEM_SPAN + AGE, data = df)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.65555 -0.15167 0.01677 0.11466 0.73690
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 1.89714 0.57819 3.281 0.00441 **
## MEM_SPAN 0.52096 0.18353 2.839 0.01135 *
## AGE 0.33936 0.09881 3.435 0.00316 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.3192 on 17 degrees of freedom
## Multiple R-squared: 0.8071, Adjusted R-squared: 0.7844
## F-statistic: 35.57 on 2 and 17 DF, p-value: 8.414e-07p-value для F-статистики для модели равняется 8.414e-07, что меньше 0.05, значит модель не является бессмысленной и объясняет 81% дисперсии умения читать (Multiple R-squared: 0.8071).
Для переменной MEM_SPAN p-value меньше 0.05, следовательно мы должны отклонить нулевую гипотезу о равенстве нулю коэффициента при MEM_SPAN. Оценка этого коэффициента равна 0.52096.
Для переменной AGE p-value также меньше 0.05, следовательно мы должны отклонить нулевую гипотезу о равенстве этого коэффициента нулю. Оценка этого коэффициента равна 0.33936.
Запишем регрессионное уравнение полностью:
READ_AB = 1.89714 + 0.52096 * MEM_SPAN + 0.33936 * AGE
При увеличении MEM_SPAN на единицу ожидается увеличение умения читать на 0.52096 при постоянном значении AGE.
При увеличении AGE на единицу ожидается увеличение умения читать на 0.33936 при постоянном значении MEM_SPAN.
Но память также развивается с возрастом, поэтому взаимодействие памяти и возраста также может влиять на умение читать. Проверим эту гипотезу, добавив в модель ещё один предиктор.
model3 <- lm(READ_AB ~ MEM_SPAN + AGE + MEM_SPAN:AGE, data=df)
summary(model3)##
## Call:
## lm(formula = READ_AB ~ MEM_SPAN + AGE + MEM_SPAN:AGE, data = df)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.48254 -0.15311 0.00693 0.14414 0.62320
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -10.6151 4.8385 -2.194 0.04336 *
## MEM_SPAN 3.4686 1.1448 3.030 0.00797 **
## AGE 2.2162 0.7269 3.049 0.00766 **
## MEM_SPAN:AGE -0.4382 0.1685 -2.600 0.01935 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.2759 on 16 degrees of freedom
## Multiple R-squared: 0.8644, Adjusted R-squared: 0.839
## F-statistic: 34 on 3 and 16 DF, p-value: 3.578e-07p-value для фактора взаимодействия памяти и возраста также равняется 0.01935, что также меньше 0.05, значит мы должны отклонить нулевую гипотезу о равенстве нулю коэффициента при этом факторе. Оценка этого коэффициента равна -0.4382. Это означает, что на взаимосвязь между одной независимой переменной и зависимой влияют значения другой независимой переменной.
Запишем регрессионное уравнение полностью:
READ_AB = -10.6151 + 3.4686 * MEM_SPAN + 2.2162 * AGE + (-0.4382 * MEM_SPAN:AGE)
Вместе эти три предиктора объясняют уже 86% дисперсии умения читать (Multiple R-squared: 0.8644).
Визуализируем взаимодействие с помощью функции visreg() из пакета visreg.
visreg(model3, "MEM_SPAN", by = "AGE", overlay=TRUE, band=TRUE)
Попробуйте добавить в нашу модель ещё один предиктор - IQ.
Множественная регрессия с номинальным предиктором (метод наименьших квадратов, МНК)
В качестве предикторов можно использовать не только количественные переменные, но и качественные. Для это их надо трансформировать в дамми-переменные, принимающие только два значения: 1 или 0.
Во встроенном в R наборе данных mtcars есть данные о расходе топлива (переменная mpg, Miles/(US) gallon) автомобилей с двумя типами трансмиссии (переменная am, (0 = automatic, 1 = manual)) и разной мощностью (переменная hp, в лошадиных силах).
str(mtcars)## 'data.frame': 32 obs. of 11 variables:
## $ mpg : num 21 21 22.8 21.4 18.7 18.1 14.3 24.4 22.8 19.2 ...
## $ cyl : num 6 6 4 6 8 6 8 4 4 6 ...
## $ disp: num 160 160 108 258 360 ...
## $ hp : num 110 110 93 110 175 105 245 62 95 123 ...
## $ drat: num 3.9 3.9 3.85 3.08 3.15 2.76 3.21 3.69 3.92 3.92 ...
## $ wt : num 2.62 2.88 2.32 3.21 3.44 ...
## $ qsec: num 16.5 17 18.6 19.4 17 ...
## $ vs : num 0 0 1 1 0 1 0 1 1 1 ...
## $ am : num 1 1 1 0 0 0 0 0 0 0 ...
## $ gear: num 4 4 4 3 3 3 3 4 4 4 ...
## $ carb: num 4 4 1 1 2 1 4 2 2 4 ...Попробуем визуально оценить связь между этими переменными.
ggplot(mtcars, aes(mpg, hp)) + geom_point(size = 4) + theme_classic() + geom_smooth(method=lm)
ggplot(mtcars, aes(as.factor(am), mpg)) + geom_violin(aes(fill=as.factor(am))) + stat_summary(fun.data=mean_sdl, mult=1, geom="pointrange", color="black") + theme_classic()
С помощью регрессионного анализа проверим, влияют ли мощность и тип трансмиссии на расход топлива
model4 <- lm(mpg ~ hp + am, data=mtcars)
summary(model4)##
## Call:
## lm(formula = mpg ~ hp + am, data = mtcars)
##
## Residuals:
## Min 1Q Median 3Q Max
## -4.3843 -2.2642 0.1366 1.6968 5.8657
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 26.584914 1.425094 18.655 < 2e-16 ***
## hp -0.058888 0.007857 -7.495 2.92e-08 ***
## am 5.277085 1.079541 4.888 3.46e-05 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 2.909 on 29 degrees of freedom
## Multiple R-squared: 0.782, Adjusted R-squared: 0.767
## F-statistic: 52.02 on 2 and 29 DF, p-value: 2.55e-10p-value для F-статистики для модели равняется 2.55e-10, что меньше 0.05, значит модель не является бессмысленной и объясняет 78% дисперсии расхода топлива (Multiple R-squared: 0.782).
Для переменной hp p-value меньше 0.05, следовательно мы должны отклонить нулевую гипотезу о равенстве нулю коэффициента при hp. Оценка этого коэффициента равна -0.058888.
Для переменной am p-value также меньше 0.05, следовательно мы должны отклонить нулевую гипотезу о равенстве этого коэффициента нулю. Оценка этого коэффициента равна 5.277085.
Запишем регрессионное уравнение полностью. Поскольку один из предикторов является качественном переменной, у нас будет два уравнения.
- Уравнение для автомобилей с ручной трансмиссией (am=1):
mpg = 26.584914 + (-0.058888 * hp) + 5.277085 * am, поскольку am=1, то упрощаем уравнение ->
mpg = (26.584914 + 5.277085) - 0.058888 * hp
При увеличении мощности на одну лошадиную силу ожидается уменьшение растояние, которое автомобиль с ручной трансмиссией может проехать, на 0.058888 миль на галон (т.е. расход топлива растёт).
- Уравнение для автомобилей с автоматической трансмиссией (am=0):
mpg = 26.584914 + (-0.058888 * hp), поскольку am=0, коэфициент при этой переменной обнуляется