Raport 5
Maja i Zuzia
2024-12-03
Raport 5 - Regresja grzbietowa, LASSO i elastic net
library(ISLR)
library(dplyr)
library(glmnet)
library(knitr)Dane
Dane pochodzą z pakietu ISLR - Wage. Baza posiada 3000 obserwacji, opisujących pracowników: ich wiek, stan cywilny, rasę, poziom edukacji, miejsce zamieszkania, rodzaj wykonywanego zawodu, poziom zdrowia oraz wielkość wynagrodzenia.
dane <- Wage
dane = na.omit(dane)
kable(head(dane, 5),format="html",caption = "Dane z pakietu Wage")| year | age | maritl | race | education | region | jobclass | health | health_ins | logwage | wage | |
|---|---|---|---|---|---|---|---|---|---|---|---|
| 231655 | 2006 | 18 |
|
|
|
|
|
|
|
4.318063 | 75.04315 |
| 86582 | 2004 | 24 |
|
|
|
|
|
|
|
4.255273 | 70.47602 |
| 161300 | 2003 | 45 |
|
|
|
|
|
|
|
4.875061 | 130.98218 |
| 155159 | 2003 | 43 |
|
|
|
|
|
|
|
5.041393 | 154.68529 |
| 11443 | 2005 | 50 |
|
|
|
|
|
|
|
4.318063 | 75.04315 |
Zbiór treningowy i testowy
Najpierw przekształciliśmy dane, tak, aby pasowały do modelowania
glmnet(). Naszą zmienną objaśnianą będzie wielkość
wynagrodzenia (wage), a macierz zmiennych objasniajacych utworzy zbiór
pozostałych dostepnych zmiennych.
x = model.matrix(wage~age+maritl+race+education+
region+jobclass+health+health_ins, dane)
y = dane %>%
select(wage) %>%
unlist() %>%
as.numeric()W raporcie porównujemy regresję grzbietową i LASSO w odniesieniu do zwykłej regresji OLS, co jest możliwe dzięki porównaniu błędu średniokwadratowego MSE. Podzielenie próbki na zbiór treningowy i testowy umozliwi oszacowanie błędu testu.
set.seed(1)
train = dane %>%
sample_frac(0.5)
test = dane %>%
setdiff(train)
x_train = model.matrix(wage~age+maritl+race+education+
region+jobclass+health+health_ins, train)
x_test = model.matrix(wage~age+maritl+race+education+
region+jobclass+health+health_ins, test)
y_train = train %>%
select(wage) %>%
unlist() %>%
as.numeric()
y_test = test %>%
select(wage) %>%
unlist() %>%
as.numeric()Regresja Grzbietowa
Skalowanie:
Wybrany zakres wartości to od lambda = 10^10 do lambda = 10^-2. Zmienne poddane zostały standaryzacji.
grid = 10^seq(10, -2, length = 100)
ridge_mod = glmnet(x, y, alpha = 0, lambda = grid) # alpha=0 - regresja grzbietowaplot(ridge_mod, xlab="")
title(sub="Wykres współczynników L2")
Oś X reprezentuje wartości parametru regularizacji λ, wyższe wartości odpowiadają silniejszej regularyzacji. Na tej podstawie, możemy stwierdzić, że zmienne objaśniające w modelu są zróżnicowane: niektóre są bardziej odporne na regularyzację (wraz ze wzrostem wartości λ, wartość współczynnika oscyluje koło 0), a inne są szybko redukowane.
Optymalna wartość parametru λ
W celu znalezienia optymalnej wartości parametru λ, która minimalizuje wartość błedu MSE, użyjemy walidacji krzyżowej:
set.seed(123) # dla powtarzalności wyników
cv_ridge <- cv.glmnet(x_train, y_train, alpha = 0)plot(cv_ridge)
title(sub="Walidacja krzyżowa dla regresji grzbietowej")
# optymalna wartość lambdy
best_lambda <- cv_ridge$lambda.min
# przewidywanie danych testowych
pred_ridge = predict(ridge_mod, s = best_lambda, newx = x_test)
# błąd średniokwadratowy MSE
mse_pred_rigde = mean((pred_ridge - y_test)^2)
# błąd RMSE i RSMPE
rmse_pred_rigde = sqrt(mse_pred_rigde)
rmspe_pred_rigde = rmse_pred_rigde / mean(y_test)| Najlepszaλ | MSE | RMSE | RSMPE |
|---|---|---|---|
| 1.699209 | 1163.668 | 34.11258 | 0.3046133 |
Wartość λ na poziomie 0.1233 daje najbardziej precyzyjny model. Błąd średniokwadratowy tego modelu wyniósł 1159. Zatem możemy powiedzieć, że prognoza dla zbioru testowego wiekości wynagrodzenia pracowników, odchyla się od wartości rzeczywistych średnio o 34 jednostki, co w ujęciu procentowym wynosi 30.4%.
Oszacowanie modelu
model_ridge = glmnet(x, y, alpha = 0)
predict(model_ridge, type = "coefficients", s = best_lambda)[1:25,]## (Intercept) (Intercept)
## 71.71642208 0.00000000
## age maritl2. Married
## 0.29155099 16.04786943
## maritl3. Widowed maritl4. Divorced
## 0.08833332 2.81324269
## maritl5. Separated race2. Black
## 9.18857839 -5.11495712
## race3. Asian race4. Other
## -2.01741271 -7.21275543
## education2. HS Grad education3. Some College
## 1.94636815 11.97136910
## education4. College Grad education5. Advanced Degree
## 24.46796238 46.35258573
## region2. Middle Atlantic region3. East North Central
## 0.00000000 0.00000000
## region4. West North Central region5. South Atlantic
## 0.00000000 0.00000000
## region6. East South Central region7. West South Central
## 0.00000000 0.00000000
## region8. Mountain region9. Pacific
## 0.00000000 0.00000000
## jobclass2. Information health2. >=Very Good
## 4.11561881 6.77240059
## health_ins2. No
## -17.47962728Najważniejsze predyktory:
age - wiek pracownika
maritl - stan cywilny
race - rasa
education - poziom wykształcenia
jobclass - praca fizyczna/umysłowa
health - stan zdrowia
Współczynniki zerowe dla wszystkich zmiennych regionu oraz sugerują, że lokalizacja geograficzna nie wpływa znacząco na wynagrodzenie w modelu tegresji grzbietowej.
Regresja LASSO
W celu stworzenia modelu lasso użyta została ponownie funkcja glmnet(), tak jak w przypadku regresji grzbietowej, jednak w tym przypadku alpha wynosi 1.
lasso_mod = glmnet(x_train, y_train, alpha = 1, lambda = grid)
plot(lasso_mod)## Warning in regularize.values(x, y, ties, missing(ties), na.rm = na.rm):
## zwijanie do unikalnych wartości 'x'title(sub="Wykres współczynników L1")
Na wykresie współczynników można zauważyć, że obserwowane wartości się bardziej zbliżone do 0, niż w przypadku regresji grzbietowej. Dla małych wartości parametru λ większość współczynników osiąga wartość 0 lub bardzo zbliżone.
Optymalna wartość parametru λ
set.seed(1)
cv_lasso = cv.glmnet(x_train, y_train, alpha = 1)
plot(cv_lasso)
title(sub="Walidacja krzyżowa dla regresji lasso")
Porównując wyniki walidacji krzyżowej dla obu regresji można zauważyć, że dla regresji lasso, dla małych λ błąd MSE jest niski, ale gwałtownie rośnie przy większych 𝜆, natomiast w przypadku regresji grzbietowej rośnie stopniowo. Wynika to z tego, że LASSO zeruje mniej istotne współczynniki, co upraszcza model, ale przy dużych𝜆eliminuje zbyt wiele informacji, powodując wzrost błędu.
# optymalna wartość lambdy
best_lambda1 <- cv_lasso$lambda.min
# przewidywanie danych testowych
pred_lasso = predict(lasso_mod, s = best_lambda1, newx = x_test)
# błąd średniokwadratowy MSE
mse_pred_lasso = mean((pred_lasso - y_test)^2)
# błąd RMSE i RSMPE
rmse_pred_lasso = sqrt(mse_pred_lasso)
rmspe_pred_lasso = rmse_pred_lasso / mean(y_test)| Najlepszaλ | MSE | RMSE | RSMPE |
|---|---|---|---|
| 0.0229911 | 1187.061 | 34.45376 | 0.3076599 |
Błędy uzyskane za pomocą regresji lasso są nieznacznie, lecz jednak wyższe niż w przypadku regresji grzbietowej. Optymalna wartość λ wynosi ok. 0,023 i jest znacząco niższa niż w przypadku regresji grzbietowej, co oznacza, że lasso już przy niskich wartościach parametru λ skutecznie eliminuje mniej istotne cechy. Dlatego jej mniejsze wartości wystarczają, by uzyskać model, który dobrze równoważy dokładność i prostotę. W regresji grzbietowej λ musi być większe, aby znacząco zmniejszyć współczynniki i zapobiec nadmiernemu dopasowaniu.
Oszacowanie modelu
out = glmnet(x, y, alpha = 1, lambda = grid)
lasso_coef = predict(out, type = "coefficients", s = best_lambda1)[1:25,]
lasso_coef## (Intercept) (Intercept)
## 66.4319729 0.0000000
## age maritl2. Married
## 0.2839197 16.8198441
## maritl3. Widowed maritl4. Divorced
## 0.4868734 3.4601739
## maritl5. Separated race2. Black
## 11.1831388 -4.8267626
## race3. Asian race4. Other
## -2.3911882 -5.8587670
## education2. HS Grad education3. Some College
## 7.3046667 17.7500675
## education4. College Grad education5. Advanced Degree
## 30.7800546 53.6109057
## region2. Middle Atlantic region3. East North Central
## 0.0000000 0.0000000
## region4. West North Central region5. South Atlantic
## 0.0000000 0.0000000
## region6. East South Central region7. West South Central
## 0.0000000 0.0000000
## region8. Mountain region9. Pacific
## 0.0000000 0.0000000
## jobclass2. Information health2. >=Very Good
## 3.4716311 6.5214179
## health_ins2. No
## -17.4677452lasso_coef[lasso_coef != 0]## (Intercept) age
## 66.4319729 0.2839197
## maritl2. Married maritl3. Widowed
## 16.8198441 0.4868734
## maritl4. Divorced maritl5. Separated
## 3.4601739 11.1831388
## race2. Black race3. Asian
## -4.8267626 -2.3911882
## race4. Other education2. HS Grad
## -5.8587670 7.3046667
## education3. Some College education4. College Grad
## 17.7500675 30.7800546
## education5. Advanced Degree jobclass2. Information
## 53.6109057 3.4716311
## health2. >=Very Good health_ins2. No
## 6.5214179 -17.4677452Najważniejsze predyktory:
age - wiek pracownika
maritl - stan cywilny
race - rasa
education - poziom wykształcenia
jobclass - praca fizyczna/umysłowa
health - stan zdrowia
W przypadku obu modeli regresji te same predyktory okazały się być najważniejsze
Regresja grzbietowa vs Regresja Lasso
Wybór między lasso a regresją grzbietową zależy od celu analizy: selekcji cech (LASSO) czy kontroli nadmiaru dopasowania (grzbietowa). Regresja grzbietowa jest bardziej konserwatywna, a jej błąd rośnie łagodniej, co czyni ją stabilniejszą przy większych𝜆, natomiast lasso jest bardziej agresywne w eliminacji cech, co prowadzi do większego wzrostu błędu przy dużych 𝜆. Metoda lasso wydaje się być lepsza do wybierania najważniejszych cech, a regresja grzbietowa sprawdza się, gdy wszystkie cechy są istotne.
Regresja Elastic Net
Elastic Net wykorzystuje parametr α, który kontroluje równowagę między regularyzacją L1 i L2. Gdy α=1, model jest równoważny LASSO, a gdy α=0, odpowiada regresji grzbietowej(Ridge Regression). Wybrana wartość α=0.5 równoważy oba typy regularyzacji.
elastic_mod = glmnet(x_train, y_train, alpha = 0.5, lambda = grid)
plot(elastic_mod, xlab='')## Warning in regularize.values(x, y, ties, missing(ties), na.rm = na.rm):
## zwijanie do unikalnych wartości 'x'title(sub="Wykres współczynników - elastic net")
Optymalna wartość parametru λ:
cv_elastic = cv.glmnet(x_train, y_train, alpha = 0.5)
plot(cv_elastic)
title(sub="Walidacja krzyżowa dla regresji Elastic Net")
# optymalna wartość lambdy
best_lambda2 <- cv_elastic$lambda.min
# przewidywanie danych testowych
pred_elastic = predict(elastic_mod, s = best_lambda2, newx = x_test)
# błąd średniokwadratowy MSE
mse_pred_elastic = mean((pred_elastic - y_test)^2)
# błąd RMSE i RSMPE
rmse_pred_elastic = sqrt(mse_pred_elastic)
rmspe_pred_elastic = rmse_pred_elastic / mean(y_test)| Najlepszaλ | MSE | RMSE | RSMPE |
|---|---|---|---|
| 0.0459822 | 1187.093 | 34.45421 | 0.3076639 |
Elastic Net, przy optymalnej wartości λ, osiąga błąd średniokwadratowy na poziomie MSE=1187, co odpowiada prognozom obarczonym przeciętnym błędem na poziomie RSMPE=30.77%. Wyniki tej regresji są bardzo zbliżone do regresji z wykorzystamiem regularyzacji LASSO.
model_elastic = glmnet(x, y, alpha = 0.5, lambda=grid)
predict(model_elastic, type = "coefficients", s = best_lambda2)[1:25,]## (Intercept) (Intercept)
## 66.5361212 0.0000000
## age maritl2. Married
## 0.2840374 16.8073401
## maritl3. Widowed maritl4. Divorced
## 0.4746751 3.4483366
## maritl5. Separated race2. Black
## 11.1444592 -4.8307941
## race3. Asian race4. Other
## -2.3837220 -5.8822084
## education2. HS Grad education3. Some College
## 7.1977962 17.6361652
## education4. College Grad education5. Advanced Degree
## 30.6580037 53.4744436
## region2. Middle Atlantic region3. East North Central
## 0.0000000 0.0000000
## region4. West North Central region5. South Atlantic
## 0.0000000 0.0000000
## region6. East South Central region7. West South Central
## 0.0000000 0.0000000
## region8. Mountain region9. Pacific
## 0.0000000 0.0000000
## jobclass2. Information health2. >=Very Good
## 3.4828183 6.5260460
## health_ins2. No
## -17.4707279Podobnie jak w przypadku LASSO i regredji grzbietowej, najważniejszymi predyktorami okazały się:
age - wiek pracownika
maritl - stan cywilny
race - rasa
education - poziom wykształcenia
jobclass - praca fizyczna/umysłowa
health - stan zdrowia
Regresja liniowa
OLS <- lm(y~x)
summary(OLS)##
## Call:
## lm(formula = y ~ x)
##
## Residuals:
## Min 1Q Median 3Q Max
## -103.663 -18.706 -3.473 13.853 211.966
##
## Coefficients: (9 not defined because of singularities)
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 65.8356 3.5407 18.594 < 2e-16 ***
## x(Intercept) NA NA NA NA
## xage 0.2837 0.0623 4.554 5.47e-06 ***
## xmaritl2. Married 16.9253 1.7234 9.821 < 2e-16 ***
## xmaritl3. Widowed 0.9009 8.0206 0.112 0.910578
## xmaritl4. Divorced 3.6329 2.8929 1.256 0.209287
## xmaritl5. Separated 11.5439 4.8563 2.377 0.017512 *
## xrace2. Black -4.8977 2.1505 -2.277 0.022830 *
## xrace3. Asian -2.5041 2.6087 -0.960 0.337193
## xrace4. Other -5.9525 5.6809 -1.048 0.294809
## xeducation2. HS Grad 7.8432 2.3754 3.302 0.000972 ***
## xeducation3. Some College 18.3040 2.5265 7.245 5.49e-13 ***
## xeducation4. College Grad 31.3257 2.5547 12.262 < 2e-16 ***
## xeducation5. Advanced Degree 54.1677 2.8180 19.222 < 2e-16 ***
## xregion2. Middle Atlantic NA NA NA NA
## xregion3. East North Central NA NA NA NA
## xregion4. West North Central NA NA NA NA
## xregion5. South Atlantic NA NA NA NA
## xregion6. East South Central NA NA NA NA
## xregion7. West South Central NA NA NA NA
## xregion8. Mountain NA NA NA NA
## xregion9. Pacific NA NA NA NA
## xjobclass2. Information 3.4806 1.3273 2.622 0.008775 **
## xhealth2. >=Very Good 6.5454 1.4244 4.595 4.51e-06 ***
## xhealth_ins2. No -17.4482 1.4069 -12.402 < 2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 34.09 on 2984 degrees of freedom
## Multiple R-squared: 0.336, Adjusted R-squared: 0.3327
## F-statistic: 100.7 on 15 and 2984 DF, p-value: < 2.2e-16Najważniejsze predyktory:
age - wiek pracownika
maritl - stan cywilny
race - rasa
education - poziom wykształcenia
jobclass - praca fizyczna/umysłowa
health - stan zdrowia
Najważniejsze predykotry są niezmienne przy każdej użytej w raporcie metodzie
Porównanie regresji liniowej z regresją z wykorzystaniem regularyzacji
bledy <- data.frame(
"Błąd"=c("MSE", "RMSE", "RMSPE"),
"Grzbietowa"=c(mse_pred_rigde, rmse_pred_rigde, rmspe_pred_rigde),
"Lasso"=c(mse_pred_lasso, rmse_pred_lasso, rmspe_pred_lasso),
"Elastic Net"=c(mse_pred_elastic, rmse_pred_elastic, rmspe_pred_elastic),
"OLS"=c(mse_pred_ols, rmse_pred_ols, rmspe_pred_ols)
)
kable(bledy)| Błąd | Grzbietowa | Lasso | Elastic.Net | OLS |
|---|---|---|---|---|
| MSE | 1163.6682155 | 1187.0612857 | 1187.0926676 | 2373.262817 |
| RMSE | 34.1125815 | 34.4537558 | 34.4542112 | 48.716145 |
| RMSPE | 0.3046133 | 0.3076599 | 0.3076639 | 0.435018 |
Porównując wyniki modeli, regresja OLS wypada znacząco gorzej pod względem wszystkich miar błędu. Średni błąd kwadratowy dla OLS jest ponad dwukrotnie większy niż regresji wykorzystujących regularyzację. Błędy RMSE RMSPE dla OLS są wyraźnie wyższe niż dla regresji grzbietowej lasso i elastic net, świadcząc o słabszym dopasowaniu.
Powyższe wyniki sugerują, że OLS jest mniej skuteczny w modelowaniu danych, zwłaszcza w kontekście, gdzie regularyzacja pomaga lepiej radzić sobie z problemami, takimi jak wielowymiarowość czy nadmierne dopasowanie.