Jak wykonać regresję ols w r (z przykładem)

Regresja metodą najmniejszych kwadratów (OLS) to metoda pozwalająca znaleźć linię najlepiej opisującą związek między jedną lub większą liczbą zmiennych predykcyjnych azmienną odpowiedzi .

Metoda ta pozwala nam znaleźć następujące równanie:

ŷ = b ₀ + b ₁ x

Złoto:

• ŷ : Szacowana wartość odpowiedzi
• b ₀ : Początek linii regresji
• b ₁ : Nachylenie linii regresji

Równanie to może pomóc nam zrozumieć związek między predyktorem a zmienną odpowiedzi, a także można je wykorzystać do przewidywania wartości zmiennej odpowiedzi, biorąc pod uwagę wartość zmiennej predykcyjnej.

Poniższy przykład pokazuje krok po kroku, jak przeprowadzić regresję OLS w R.

Krok 1: Utwórz dane

Na potrzeby tego przykładu utworzymy zbiór danych zawierający dwie zmienne dla 15 uczniów:

• Łączna liczba przepracowanych godzin
• Wynik egazminu

Przeprowadzimy regresję OLS, wykorzystując godziny jako zmienną predykcyjną i wynik egzaminu jako zmienną odpowiedzi.

Poniższy kod pokazuje, jak utworzyć ten fałszywy zbiór danych w R:

 #create dataset
df <- data. frame (hours=c(1, 2, 4, 5, 5, 6, 6, 7, 8, 10, 11, 11, 12, 12, 14),
                 score=c(64, 66, 76, 73, 74, 81, 83, 82, 80, 88, 84, 82, 91, 93, 89))

#view first six rows of dataset
head(df)

  hours score
1 1 64
2 2 66
3 4 76
4 5 73
5 5 74
6 6 81

Krok 2: Wizualizuj dane

Przed wykonaniem regresji OLS utwórzmy wykres rozrzutu, aby zwizualizować związek między liczbą godzin a wynikiem egzaminu:

 library (ggplot2)

#create scatterplot
ggplot(df, aes(x=hours, y=score)) +
  geom_point(size= 2 )

Jedno z czterech założeń regresji liniowej zakłada, że istnieje liniowa zależność między predyktorem a zmienną odpowiedzi.

Z wykresu widać, że zależność wydaje się być liniowa. Wraz ze wzrostem liczby godzin wynik również ma tendencję do zwiększania się liniowo.

Następnie możemy utworzyć wykres pudełkowy, aby zwizualizować rozkład wyników egzaminu i sprawdzić wartości odstające.

Uwaga : R definiuje obserwację jako wartość odstającą, jeśli jest 1,5-krotnością rozstępu międzykwartylowego powyżej trzeciego kwartyla lub 1,5-krotności rozstępu międzykwartylowego poniżej pierwszego kwartyla.

Jeśli obserwacja jest odstająca, na wykresie pudełkowym pojawi się małe kółko:

 library (ggplot2)

#create scatterplot
ggplot(df, aes(y=score)) +
  geom_boxplot() 

Na wykresie pudełkowym nie ma małych okręgów, co oznacza, że w naszym zbiorze danych nie ma wartości odstających.

Krok 3: Wykonaj regresję OLS

Następnie możemy użyć funkcji lm() w R, aby przeprowadzić regresję OLS, używając godzin jako zmiennej predykcyjnej i wyniku jako zmiennej odpowiedzi:

 #fit simple linear regression model
model <- lm(score~hours, data=df)

#view model summary
summary(model)

Call:
lm(formula = score ~ hours)

Residuals:
   Min 1Q Median 3Q Max 
-5,140 -3,219 -1,193 2,816 5,772 

Coefficients:
            Estimate Std. Error t value Pr(>|t|)    
(Intercept) 65,334 2,106 31,023 1.41e-13 ***
hours 1.982 0.248 7.995 2.25e-06 ***
---
Significant. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 3.641 on 13 degrees of freedom
Multiple R-squared: 0.831, Adjusted R-squared: 0.818 
F-statistic: 63.91 on 1 and 13 DF, p-value: 2.253e-06

Z podsumowania modelu widzimy, że dopasowane równanie regresji ma postać:

Wynik = 65,334 + 1,982*(godziny)

Oznacza to, że każda dodatkowa godzina nauki wiąże się ze średnim wzrostem wyniku z egzaminu o 1982 punkty.

Oryginalna wartość 65 334 mówi nam o średnim oczekiwanym wyniku egzaminu dla studenta studiującego przez zero godzin.

Możemy również użyć tego równania, aby znaleźć oczekiwany wynik egzaminu na podstawie liczby godzin spędzonych przez studenta.

Przykładowo student studiujący 10 godzin powinien uzyskać z egzaminu wynik 85,15 :

Wynik = 65,334 + 1,982*(10) = 85,15

Oto jak zinterpretować pozostałą część podsumowania modelu:

• Pr(>|t|): Jest to wartość p powiązana ze współczynnikami modelu. Ponieważ wartość p dla godzin (2,25e-06) jest znacznie mniejsza niż 0,05, możemy powiedzieć, że istnieje statystycznie istotny związek między godzinami a wynikiem .
• Wielokrotne R-kwadrat: Liczba ta mówi nam, że procent zróżnicowania wyników egzaminu można wyjaśnić liczbą godzin nauki. Ogólnie rzecz biorąc, im większa wartość R-kwadrat modelu regresji, tym lepiej zmienne predykcyjne przewidują wartość zmiennej odpowiedzi. W tym przypadku 83,1% różnic w wynikach można wyjaśnić liczbą przepracowanych godzin.
• Resztkowy błąd standardowy: jest to średnia odległość pomiędzy obserwowanymi wartościami a linią regresji. Im niższa jest ta wartość, tym bardziej linia regresji może odpowiadać obserwowanym danym. W tym przypadku średni wynik zaobserwowany na egzaminie odbiega o 3641 punktów od wyniku przewidywanego przez linię regresji.
• Statystyka F i wartość p: Statystyka F ( 63,91 ) i odpowiadająca jej wartość p ( 2.253e-06 ) mówią nam o ogólnym znaczeniu modelu regresji, tj. czy zmienne predykcyjne w modelu są przydatne do wyjaśnienia zmienności . w zmiennej odpowiedzi. Ponieważ wartość p w tym przykładzie jest mniejsza niż 0,05, nasz model jest istotny statystycznie i godziny uważa się za przydatne do wyjaśnienia zmienności wyniku .

Krok 4: Utwórz działki resztkowe

Na koniec musimy utworzyć wykresy resztowe, aby sprawdzić założenia homoskedastyczności i normalności .

Założeniem homoskedastyczności jest to, że reszty modelu regresji mają w przybliżeniu równą wariancję na każdym poziomie zmiennej predykcyjnej.

Aby sprawdzić, czy to założenie jest spełnione, możemy utworzyć wykres reszt w funkcji dopasowań .

Oś x wyświetla dopasowane wartości, a oś y wyświetla reszty. Dopóki reszty wydają się być losowo i równomiernie rozmieszczone na wykresie wokół wartości zerowej, możemy założyć, że homoskedastyczność nie jest naruszona:

 #define residuals
res <- resid(model)

#produce residual vs. fitted plot
plot(fitted(model), res)

#add a horizontal line at 0 
abline(0,0)

Reszty wydają się być losowo rozproszone wokół zera i nie wykazują zauważalnego wzoru, więc to założenie jest spełnione.

Założenie normalności stwierdza, że reszty modelu regresji mają w przybliżeniu rozkład normalny.

Aby sprawdzić, czy to założenie jest spełnione, możemy utworzyć wykres QQ . Jeśli punkty wykresu leżą w przybliżeniu na linii prostej tworzącej kąt 45 stopni, wówczas dane mają rozkład normalny:

 #create QQ plot for residuals
qqnorm(res)

#add a straight diagonal line to the plot
qqline(res) 

Wartości resztkowe odbiegają nieco od linii 45 stopni, ale nie na tyle, aby powodować poważne obawy. Można założyć, że założenie normalności jest spełnione.

Ponieważ reszty mają rozkład normalny i homoskedastyczny, sprawdziliśmy, czy spełnione są założenia modelu regresji OLS.

Zatem wynik naszego modelu jest niezawodny.

Uwaga : jeśli jedno lub więcej założeń nie zostało spełnionych, możemy podjąć próbę przekształcenia naszych danych.

Dodatkowe zasoby

Poniższe samouczki wyjaśniają, jak wykonywać inne typowe zadania w języku R:

Jak wykonać wielokrotną regresję liniową w R
Jak przeprowadzić regresję wykładniczą w R
Jak wykonać ważoną regresję metodą najmniejszych kwadratów w R