Як виконати множинну лінійну регресію в r
У цьому посібнику показано приклад виконання множинної лінійної регресії в R, зокрема:
- Вивчіть дані перед підгонкою моделі
- Коригування моделі
- Перевірка припущень моделі
- Інтерпретація вихідних даних моделі
- Оцінка відповідності моделі
- Використовуйте модель для прогнозування
Ходімо!
Об’єкт
Для цього прикладу ми використаємо вбудований набір даних R mtcars , який містить інформацію про різні атрибути 32 різних автомобілів:
#view first six lines of mtcars
head(mtcars)
# mpg cyl disp hp drat wt qsec vs am gear carb
#Mazda RX4 21.0 6 160 110 3.90 2.620 16.46 0 1 4 4
#Mazda RX4 Wag 21.0 6 160 110 3.90 2.875 17.02 0 1 4 4
#Datsun 710 22.8 4 108 93 3.85 2.320 18.61 1 1 4 1
#Hornet 4 Drive 21.4 6 258 110 3.08 3.215 19.44 1 0 3 1
#Hornet Sportabout 18.7 8 360 175 3.15 3.440 17.02 0 0 3 2
#Valiant 18.1 6 225 105 2.76 3.460 20.22 1 0 3 1
У цьому прикладі ми побудуємо модель множинної лінійної регресії, яка використовує mpg як змінну відповіді та disp , hp і drat як змінні прогнозу.
#create new data frame that contains only the variables we would like to use to data <- mtcars[, c("mpg", "disp", "hp", "drat")] #view first six rows of new data frame head(data) # mpg disp hp drat #Mazda RX4 21.0 160 110 3.90 #Mazda RX4 Wag 21.0 160 110 3.90 #Datsun 710 22.8 108 93 3.85 #Hornet 4 Drive 21.4 258 110 3.08 #Hornet Sportabout 18.7 360 175 3.15 #Valiant 18.1 225 105 2.76
Огляд даних
Перш ніж підбирати модель, ми можемо переглянути дані, щоб краще зрозуміти їх, а також візуально оцінити, чи може множинна лінійна регресія бути хорошою моделлю для підгонки цих даних.
Зокрема, нам потрібно перевірити, чи змінні предикторів мають лінійний зв’язок зі змінною відповіді, що вказує на те, що модель множинної лінійної регресії може бути придатною.
Для цього ми можемо використати функцію pairs() , щоб створити діаграму розсіювання кожної можливої пари змінних:
pairs(data, pch = 18, col = "steelblue")
З цього графіка пар ми можемо побачити наступне:
- mpg та доступність, здається, мають сильну негативну лінійну кореляцію
- mpg і hp, здається, мають сильну позитивну лінійну кореляцію
- MPG і drat, здається, мають помірну негативну лінійну кореляцію
Зауважте, що ми також можемо використати функцію ggpairs() із бібліотеки GGally , щоб створити подібний графік, що містить фактичні лінійні коефіцієнти кореляції для кожної пари змінних:
#install and load the GGally library install.packages("GGally") library(GGally) #generate the pairs plot ggpairs(data)
Здається, що кожна зі змінних предикторів має помітну лінійну кореляцію зі змінною відповіді mpg , тому ми продовжимо адаптувати модель лінійної регресії до даних.
Коригування моделі
Основний синтаксис для підгонки моделі множинної лінійної регресії в R:
lm(response_variable ~ predictor_variable1 + predictor_variable2 + ..., data = data)
Використовуючи наші дані, ми можемо адаптувати модель за допомогою наступного коду:
model <- lm(mpg ~ disp + hp + drat, data = data)
Перевірка припущень моделі
Перш ніж перейти до перевірки результатів моделі, ми повинні спочатку переконатися, що припущення моделі виконуються. А саме, нам потрібно перевірити наступне:
1. Розподіл модельних залишків має бути приблизно нормальним.
Ми можемо перевірити, чи виконується це припущення, створивши просту гістограму залишків:
hist(residuals(model), col = "steelblue")
Хоча розподіл дещо зміщений праворуч , він не є настільки ненормальним, щоб викликати серйозне занепокоєння.
2. Дисперсія залишків має бути узгодженою для всіх спостережень.
Цей бажаний стан відомий як гомоскедастичність. Порушення цього припущення відоме як гетероскедастичність .
Щоб перевірити, чи виконується це припущення, ми можемо створити графік скоригованої/залишкової вартості:
#create fitted value vs residual plot plot(fitted(model), residuals(model)) #add horizontal line at 0 abline(h = 0, lty = 2)
В ідеалі ми хотіли б, щоб залишки були рівномірно розподілені для кожного підігнаного значення. На графіку ми бачимо, що дисперсія має тенденцію до збільшення для більших підігнаних значень, але ця тенденція не є настільки екстремальною, щоб викликати зайве занепокоєння.
Інтерпретація вихідних даних моделі
Після того, як ми переконалися, що припущення моделі достатньо виконані, ми можемо перевірити вихід моделі за допомогою функції summary() :
summary(model) #Call: #lm(formula = mpg ~ disp + hp + drat, data = data) # #Residuals: # Min 1Q Median 3Q Max #-5.1225 -1.8454 -0.4456 1.1342 6.4958 # #Coefficients: #Estimate Std. Error t value Pr(>|t|) #(Intercept) 19.344293 6.370882 3.036 0.00513 ** #disp -0.019232 0.009371 -2.052 0.04960 * #hp -0.031229 0.013345 -2.340 0.02663 * #drat 2.714975 1.487366 1.825 0.07863 . #--- #Significant. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1 # #Residual standard error: 3.008 on 28 degrees of freedom #Multiple R-squared: 0.775, Adjusted R-squared: 0.7509 #F-statistic: 32.15 on 3 and 28 DF, p-value: 3.28e-09
З результату ми бачимо наступне:
- Загальна F-статистика моделі становить 32,15 , а відповідне значення p — 3,28e-09 . Це вказує на те, що загальна модель є статистично значущою. Іншими словами, модель регресії в цілому корисна.
- disp є статистично значущим на рівні значущості 0,10. Зокрема, коефіцієнт із результатів моделі вказує на те, що збільшення доступності на одну одиницю пов’язане зі зменшенням у середньому на -0,019 одиниці миль на галлон , припускаючи, що потужність і споживання палива залишаються постійними. .
- hp є статистично значущим на рівні значущості 0,10. Зокрема, коефіцієнт із результатів моделі вказує на те, що збільшення кінських сил на одну одиницю пов’язане зі зменшенням у середньому на -0,031 одиниць у милях на галлон , припускаючи, що disp і drat залишаються постійними.
- drat є статистично значущим на рівні значущості 0,10. Зокрема, коефіцієнт із результатів моделі вказує на те, що збільшення споживання бензину на одну одиницю пов’язане із середнім збільшенням на 2715 одиниць миль на галлон , припускаючи, що швидкість потоку та потужність залишаються постійними.
Оцінка відповідності моделі
Щоб оцінити, наскільки модель регресії відповідає даним, ми можемо розглянути кілька різних показників:
1. Кілька R-квадратів
Це вимірює силу лінійного зв’язку між змінними предиктора та змінною відповіді. R-квадрат, кратний 1, означає ідеальний лінійний зв’язок, тоді як R-квадрат, кратний 0, вказує на відсутність лінійного зв’язку.
Кратне R також є квадратним коренем із R у квадраті, що є часткою дисперсії у змінній відповіді, яку можна пояснити змінними предиктора. У цьому прикладі кратне R-квадрат дорівнює 0,775 . Отже, R у квадраті дорівнює 0,775 2 = 0,601 . Це вказує на те, що 60,1% дисперсії миль на галон можна пояснити предикторами моделі.
За темою: що таке хороше значення R-квадрат?
2. Залишкова стандартна помилка
Це вимірює середню відстань між спостережуваними значеннями та лінією регресії. У цьому прикладі спостережувані значення відхиляються в середньому на 3,008 одиниць від лінії регресії .
пов’язані: Розуміння стандартної помилки регресії
Використовуйте модель для прогнозування
З результатів моделі ми знаємо, що підігнане рівняння множинної лінійної регресії таке:
hat mpg = -19,343 – 0,019*disp – 0,031*hp + 2,715*drat
Ми можемо використовувати це рівняння, щоб робити прогнози про те, яким буде миль на галлон для нових спостережень . Наприклад, ми можемо знайти прогнозоване значення миль на галон для автомобіля, який має такі атрибути:
- дисплей = 220
- ch = 150
- драт = 3
#define the coefficients from the model output intercept <- coef(summary(model))["(Intercept)", "Estimate"] disp <- coef(summary(model))["disp", "Estimate"] hp <- coef(summary(model))["hp", "Estimate"] drat <- coef(summary(model))["drat", "Estimate"] #use the model coefficients to predict the value for mpg intercept + disp*220 + hp*150 + drat*3 #[1] 18.57373
Для автомобіля з disp = 220, hp = 150 і drat = 3 модель передбачає, що автомобіль отримає 18,57373 mpg .
Ви можете знайти повний код R, використаний у цьому посібнику, тут .
Додаткові ресурси
У наступних посібниках пояснюється, як підігнати інші типи моделей регресії в R:
Як виконати квадратичну регресію в R
Як виконати поліноміальну регресію в R
Як виконати експоненціальну регресію в R