Nieklasyczne metody statystyki

Regularyzacja

To laboratorium na temat Regresji grzbietowej (Ridge Regression - RR) i Lasso w R pochodzi ze stron 251-255 książki “Introduction to Statistical Learning with Applications in R” autorstwa Garetha Jamesa, Danieli Witten, Trevora Hastie i Roberta Tibshirani. Zostało ono ponownie zaimplementowane jesienią 2016 roku w formacie tidyverse przez Amelię McNamarę i R. Jordana Crousera w Smith College.

W tym tygodniu omówimy dwie alternatywne formy regresji liniowej zwane regresją grzbietową i regresją LASSO. Te dwie metody są przykładami metod regularyzacji lub zmniejszania, w których zachęca się do tego, aby parametry modelu były małe.

Regresja Grzbietowa i Lasso

Wykorzystamy pakiet glmnet w celu przeprowadzenia regresji ridge i lasso. Główną funkcją w tym pakiecie jest glmnet(), która może być użyta do dopasowania modeli regresji grzbietowej, modeli lasso i innych.

Funkcja ta ma nieco inną składnię niż inne funkcje dopasowujące modele, z którymi zetknęliśmy się do tej pory. W szczególności, musimy przekazać macierz \(x\) jak również wektor \(y\) i nie używamy składni \(y \sim x\).

Zanim przejdziemy dalej, upewnijmy się najpierw, że brakujące wartości zostały zostały usunięte z danych, jak opisano w poprzednim laboratorium.

Hitters = na.omit(Hitters)

W raporcie tym przeprowadzimy regresję grzbietową i lasso, aby przewidzieć Salary na danych Hitters.

Skonfigurujmy nasze dane:

x = model.matrix(Salary~., Hitters)[,-1] # przycinam pierwszą kolumnę
                                         # zostawiam predyktory
y = Hitters %>%
  select(Salary) %>%
  unlist() %>%
  as.numeric()

Funkcja model.matrix() jest szczególnie przydatna do tworzenia \(x\); nie tylko nie tylko tworzy macierz odpowiadającą 19 predyktorom, ale również automatycznie przekształca wszelkie zmienne jakościowe w zmienne dummy.

Ta ostatnia właściwość jest ważna, ponieważ glmnet() może przyjmować tylko numeryczne, ilościowe dane wejściowe.

Bias vs Variance

Wybór modelu w problemach uczenia nadzorowanego wiąże się z realizacją dwóch sprzecznych celów:

1.) Model powienien być dobrze dopasowany do danych uczących, aby uchwycić zależność pomiędzy danymi.

2.) Model powinien dobrze przybliżać nieznane dane (zapewniać mały błąd generalizacji).

Modele złożone dobdrze dopasowują się do danych wyjściowych, ale charakteryzują się dużą zmiennością wartości wyjściowych. Ryzykiem jest nadmierne dopasowanie = overfitting!

Modele prostsze są obciążone dużym błędem systematyczny (bias) i ich zastosowanie niesie ryzyko niewystarczającego dopasowania (underfitting)!

Składnikiem błędów generalizacji jest nieredukowalny błąd związany ze zmiennością danych.

Regularyzacja

Duża liczna zmiennych objaśniających (predyktorów): Metoda OLS nie daje jednoznacznego rozwiązania, gdy macierz XTX nie jest odwracalna (tzn. gdy zmienne objaśniające są liniowo zależne).

Taka sytuacja może mieć miejsce gdy zmiennych objaśniających jest tyle samo lub więcej niż obserwacji.

Duża wartość θi oznacza dużą wrażliwość funkcji regresji na drobne fluktuacje cechy!

Lepszym rozwiązaniem jest gorsze dopasowanie do danych uczących przy równoczesnym ograniczeniu parametrów świadczących o potencjalnie dużym błędzie generalizacji.

Regresja Grzbietowa

Wprowadzenie

Regresja grzbietowa (ang. Ridge regression) to technika regresji liniowej, która wprowadza regularyzację \(L_2\) do estymacji współczynników modelu. Regularyzacja \(L_2\) polega na dodaniu do funkcji celu kary proporcjonalnej do kwadratu wartości współczynników regresji.

Podstawową ideą regresji grzbietowej jest minimalizacja funkcji celu, która składa się z dwóch składników: błędu dopasowania (sumy kwadratów różnic pomiędzy rzeczywistymi wartościami odpowiedzi a przewidywanymi wartościami modelu) i kary regularyzacyjnej \(L_2\).

Wzór funkcji celu dla regresji grzbietowej można przedstawić jako: Minimize: RSS + \(\lambda \|\beta\|_2^2\), gdzie:

• RSS to suma kwadratów różnic pomiędzy rzeczywistymi wartościami odpowiedzi a przewidywanymi wartościami modelu (błąd dopasowania),

• \(\lambda\) (lambda) to parametr regularyzacji, który kontroluje siłę regularyzacji,

• \(\|\beta\|_2^2\) to norma \(L_2\) współczynników regresji podniesiona do kwadratu.

Dodanie kary regularyzacyjnej \(L_2\) powoduje, że współczynniki regresji są skupione wokół zera, ale nie dokładnie równe zeru (chyba że \(\lambda\)=0).

Regresja grzbietowa zmniejsza wartości współczynników, ale nie powoduje, że stają się one równe zero. Im większa wartość \(\lambda\), tym bardziej są “sciskane” współczynniki regresji.

Regresja grzbietowa jest szczególnie przydatna, gdy mamy do czynienia z modelem, w którym występuje nadmierna wielowymiarowość lub wysokie korelacje między zmiennymi niezależnymi.

Poprzez zmniejszanie wartości współczynników, regresja grzbietowa może pomóc w redukcji wpływu mało istotnych cech, poprawić stabilność modelu i zmniejszyć ryzyko przeuczenia (overfitting).

Jednym ze sposobów kontroli złożoności modelu jest penalizacja jego wielkości. Na przykład, w problemie regresji liniowej:

\[ \min_{\beta \in \mathbb{R}^p} \sum_{i=1}^n (y_i - x_i^\top \beta)^2, \]

możemy kontrolować wielkość współczynników \(\beta\). Oczywiście wielkość \(\beta\) można zdefiniować na różne sposoby, np. norma-2: \(\|\beta\|_2\), norma-1: \(\|\beta\|_1\) czy norma-nieskończoność: \(\|\beta\|_{\infty}\). Regresja grzbietowa wiąże się z karą dwóch norm:

\[ \min_{\beta \in \mathbb{R}^p} \sum_{i=1}^n (y_i - x_i^\top \beta)^2 + \lambda \|\beta\|_2^2 \]

gdzie \(\lambda\) jest parametrem kontrolującym poziom regularyzacji. Zauważ, że \(X\) to macierz \(n\) na \(p\) wymiarów z wierszami: \(x_i^\top\), oraz \(Y\) to \(n\) na 1 wektor \(y_i\). Załóżmy, że \(X^\top X + \lambda I\) jest odwracalna, mamy dokładne rozwiązanie problemu regresji grzbietowej:

\[ \hat \beta_{`ridge`} = (X^\top X + \lambda I)^{-1}X^\top Y. \]

Przypomnijmy, że rozwiązaniem zwykłej regresji najmniejszych kwadratów jest (zakładając odwracalność macierzy \(X^\top X\)):

\[ \hat \beta_{ols} = (X^\top X)^{-1}X^\top Y. \]

Dwa fakty: kiedy \(\lambda \to 0\), \(\hat \beta_{`ridge`} \to \hat \beta_{ols}\); kiedy \(\lambda \to \infty\), \(\hat \beta_{`ridge`} \to 0\).

W szczególnych przypadkach \(X\) jest ortogonalna (tzn. kolumny \(X\) są ortogonalne), mamy:

\[ \hat \beta_{`ridge`} = \frac{\hat \beta_{ols}}{1 + \lambda}. \]

Widzimy więc, że estymator grzbietowy ma dodatkowo \(1/(1 + \lambda)\) tzw. “shrinkage factor”. W związku z tym na estymatorze grzbietowym występuje obciążliwość (bias).

Przykład

Funkcja glmnet() posiada argument alfa, który określa, jaki typ modelu jest dopasowywany.

Jeśli alfa = 0 to dopasowywany jest model regresji grzbietowej, a jeśli alfa = 1 to dopasowywany jest model lasso.

Najpierw dopasowujemy model regresji grzbietowej:

grid <- 10^seq(10, -2, length = 100)  # Tworzenie siatki lambda
ridge_mod <- glmnet(x, y, alpha = 0, lambda = grid)  # Dopasowanie modelu ridge regression

Domyślnie funkcja glmnet() wykonuje regresję grzbietową dla automatycznie wybranego wybranego zakresu wartości \(\lambda\). Jednakże, tutaj wybraliśmy implementację funkcję w zakresie wartości od \(\lambda = 10^{10}\) do \(\lambda = 10^{-2}\), zasadniczo pokrywając pełen zakres scenariuszy od modelu zerowego zawierającego tylko przechwyt, do dopasowania najmniejszego kwadratu.

Jak widać, możemy również obliczyć dopasowanie modelu dla konkretnej wartości \(\lambda\), która nie jest jedną z oryginalnych wartości siatki.

Zauważ, że domyślnie funkcja glmnet() standaryzuje zmienne tak, by były w tej samej skali. Aby wyłączyć to domyślne ustawienie, użyj argumentu standardize = FALSE.

Z każdą wartością \(\lambda\) związany jest wektor współczynników regresji grzbietowej, przechowywany w macierzy, do której można uzyskać dostęp przez coef(). W tym przypadku jest to macierz \(20 \times 100\), z 20 wierszami (po jednym dla każdego predyktora, plus intercept) i 100 kolumnami (po jednej dla każdej wartości \(\lambda\)).

library(glmnet)
dim(coef(ridge_mod))

## [1]  20 100

plot(ridge_mod)    # wykres współczynników

Spodziewamy się, że oszacowania współczynników będą znacznie mniejsze, w sensie normy \(l_2\), gdy używana jest duża wartość \(\lambda\), w porównaniu z małą wartością \(\lambda\).

Oto współczynniki, gdy \(\lambda = 11498\), wraz z ich normą \(l_2\):

ridge_mod$lambda[50] # Wyświetl 50-tą wartość lambdy

## [1] 11497.57

coef(ridge_mod)[,50] # Wyświetl współczynniki związane z 50-tą wartością lambdy

##   (Intercept)         AtBat          Hits         HmRun          Runs 
## 407.356050200   0.036957182   0.138180344   0.524629976   0.230701523 
##           RBI         Walks         Years        CAtBat         CHits 
##   0.239841459   0.289618741   1.107702929   0.003131815   0.011653637 
##        CHmRun         CRuns          CRBI        CWalks       LeagueN 
##   0.087545670   0.023379882   0.024138320   0.025015421   0.085028114 
##     DivisionW       PutOuts       Assists        Errors    NewLeagueN 
##  -6.215440973   0.016482577   0.002612988  -0.020502690   0.301433531

sqrt(sum(coef(ridge_mod)[-1,50]^2)) # Oblicz normę l2

## [1] 6.360612

Dla kontrastu, oto współczynniki, gdy \(\lambda = 705\), wraz z ich \(l_2\) normą. Zwróć uwagę na znacznie większą normę \(l_2\) współczynników związanych z tą mniejszą wartością \(\lambda\).

ridge_mod$lambda[60] # Wyświetl 60-tą wartość lambdy

## [1] 705.4802

coef(ridge_mod)[,60] # Wyświetl współczynniki powiązane z 60-tą wartość lambdy

##  (Intercept)        AtBat         Hits        HmRun         Runs          RBI 
##  54.32519950   0.11211115   0.65622409   1.17980910   0.93769713   0.84718546 
##        Walks        Years       CAtBat        CHits       CHmRun        CRuns 
##   1.31987948   2.59640425   0.01083413   0.04674557   0.33777318   0.09355528 
##         CRBI       CWalks      LeagueN    DivisionW      PutOuts      Assists 
##   0.09780402   0.07189612  13.68370191 -54.65877750   0.11852289   0.01606037 
##       Errors   NewLeagueN 
##  -0.70358655   8.61181213

sqrt(sum(coef(ridge_mod)[-1,60]^2)) # Oblicz normę l2

## [1] 57.11001

Funkcję predict() możemy wykorzystać do wielu celów. Na przykład, możemy uzyskać współczynniki regresji grzbietowej dla nowej wartości \(\lambda\), powiedzmy 50:

predict(ridge_mod, s = 50, type = "coefficients")[1:20,]

##   (Intercept)         AtBat          Hits         HmRun          Runs 
##  4.876610e+01 -3.580999e-01  1.969359e+00 -1.278248e+00  1.145892e+00 
##           RBI         Walks         Years        CAtBat         CHits 
##  8.038292e-01  2.716186e+00 -6.218319e+00  5.447837e-03  1.064895e-01 
##        CHmRun         CRuns          CRBI        CWalks       LeagueN 
##  6.244860e-01  2.214985e-01  2.186914e-01 -1.500245e-01  4.592589e+01 
##     DivisionW       PutOuts       Assists        Errors    NewLeagueN 
## -1.182011e+02  2.502322e-01  1.215665e-01 -3.278600e+00 -9.496680e+00

Podzielimy teraz próbki na zbiór treningowy i testowy w celu oszacować błąd testu regresji grzbietowej i lasso.

set.seed(1)

train = Hitters %>%
  sample_frac(0.5)

test = Hitters %>%
  setdiff(train)

x_train = model.matrix(Salary~., train)[,-1]
x_test = model.matrix(Salary~., test)[,-1]

y_train = train %>%
  select(Salary) %>%
  unlist() %>%
  as.numeric()

y_test = test %>%
  select(Salary) %>%
  unlist() %>%
  as.numeric()

Następnie dopasowujemy model regresji grzbietowej na zbiorze treningowym i oceniamy jego MSE na zbiorze testowym, używając \(\lambda = 4\). Zwróć uwagę na użycie funkcji predict(). Ponownie: tym razem otrzymujemy przewidywania dla zbioru testowego, zastępując type="coefficients" argumentem newx.

ridge_mod = glmnet(x_train, y_train, alpha=0, lambda = grid, thresh = 1e-12)
ridge_pred = predict(ridge_mod, s = 4, newx = x_test)
mean((ridge_pred - y_test)^2)

## [1] 139858.6

Testowe MSE wynosi 139858. Zauważ, że gdybyśmy zamiast tego dopasowali po prostu model tylko z wyrazem wolnym, przewidywalibyśmy każdą obserwację testową używając średniej z obserwacji zbioru treningowego. W takim przypadku moglibyśmy obliczyć MSE zestawu testowego w ten sposób:

mean((mean(y_train) - y_test)^2)

## [1] 224692.1

Moglibyśmy również uzyskać ten sam wynik, dopasowując model regresji grzbietowej z bardzo dużą wartością \(\lambda\). Zauważ, że 1e10 oznacza \(10^{10}\).

ridge_pred = predict(ridge_mod, s = 1e10, newx = x_test)
mean((ridge_pred - y_test)^2)

## [1] 224692.1

Tak więc dopasowanie modelu regresji grzbietowej z \(\lambda = 4\) prowadzi do znacznie niższego testu MSE niż dopasowanie modelu z samym przechwytem.

Sprawdzimy teraz, czy jest jakaś korzyść z wykonania regresji grzbietowej z \(\lambda = 4\) zamiast po prostu wykonać regresję najmniejszych kwadratów.

Przypomnijmy, że najmniejsza kwadratura to po prostu regresja grzbietowa z \(\lambda = 0\).

* Uwaga: Aby glmnet() dawał dokładne (exact) współczynniki najmniejszego kwadratu, gdy \(\lambda = 0\), używamy argumentu exact=T przy wywołaniu funkcji predict(). W przeciwnym razie, funkcja predict() będzie interpolować nad siatką wartości \(\lambda\) użytą w dopasowaniu modelu glmnet(), dając przybliżone wyniki. Nawet gdy użyjemy exact = T, pozostaje niewielka rozbieżność na trzecim miejscu po przecinku między wynikami glmnet(), gdy \(\lambda = 0\) i wyjściem z lm(); jest to spowodowane numerycznym przybliżeniem ze strony glmnet().

ridge_pred = predict(ridge_mod, s = 0, newx = x_test)
mean((ridge_pred - y_test)^2)

## [1] 174060

lm(Salary~., data = train)

## 
## Call:
## lm(formula = Salary ~ ., data = train)
## 
## Coefficients:
## (Intercept)        AtBat         Hits        HmRun         Runs          RBI  
##   2.398e+02   -1.639e-03   -2.179e+00    6.337e+00    7.139e-01    8.735e-01  
##       Walks        Years       CAtBat        CHits       CHmRun        CRuns  
##   3.594e+00   -1.309e+01   -7.136e-01    3.316e+00    3.407e+00   -5.671e-01  
##        CRBI       CWalks      LeagueN    DivisionW      PutOuts      Assists  
##  -7.525e-01    2.347e-01    1.322e+02   -1.346e+02    2.099e-01    6.229e-01  
##      Errors   NewLeagueN  
##  -4.616e+00   -8.330e+01

predict(ridge_mod, s = 0, type="coefficients")[1:20,]

##   (Intercept)         AtBat          Hits         HmRun          Runs 
##  239.89368111   -0.01946204   -2.07305757    6.44254692    0.64610179 
##           RBI         Walks         Years        CAtBat         CHits 
##    0.82179888    3.62448842  -13.28142313   -0.70314292    3.26064805 
##        CHmRun         CRuns          CRBI        CWalks       LeagueN 
##    3.33170237   -0.54000590   -0.72015101    0.22582579  131.41324242 
##     DivisionW       PutOuts       Assists        Errors    NewLeagueN 
## -134.76073238    0.20949301    0.61942855   -4.58545824  -82.35090554

Wygląda na to, że rzeczywiście poprawiamy się w stosunku do zwykłego najmniejszego kwadratu!

Uwaga: ogólnie, jeśli chcemy dopasować (niespenalizowany) model najmniejszych kwadratów, to powinniśmy użyć funkcji lm(), ponieważ ta funkcja dostarcza bardziej użytecznych wyjścia, takie jak błędy standardowe i wartości \(p\) dla współczynników.

Zamiast arbitralnie wybierać \(\lambda = 4\), lepiej byłoby użyć walidacji krzyżowej do wyboru parametru dostrojenia \(\lambda\). Możemy to zrobić używając wbudowanej funkcji walidacji krzyżowej, cv.glmnet(). Domyślnie funkcja ta wykonuje 10-krotną walidację krzyżową, choć można to zmienić używając argumentu argumentu folds. Zauważ, że najpierw ustawiamy losowe ziarno, aby nasze wyniki były powtarzalne, ponieważ wybór krotności walidacji krzyżowej jest losowy.

set.seed(1)
cv.out = cv.glmnet(x_train, y_train, alpha = 0) # Dopasuj model regresji grzbietowej na danych treningowych
bestlam_ridge = cv.out$lambda.min  # Wybierz lamdę, która minimalizuje treningowy MSE 
bestlam_ridge

## [1] 326.1406

Widzimy zatem, że wartość \(\lambda\), która powoduje najmniejszy błąd walidacji krzyżowej to 326. Możemy również wykreślić MSE jako funkcję \(\lambda\):

plot(cv.out) # Narysuj wykres treningowego MSE jako funkcję lambda

Jaki jest testowy MSE związany z tą wartością \(\lambda\)?

ridge_pred = predict(ridge_mod, s = bestlam_ridge, newx = x_test) # Użyj najlepszej lambdy do przewidywania danych testowych
mean((ridge_pred - y_test)^2) # Oblicz testowe MSE

## [1] 140056.2

Stanowi to dalszą poprawę w stosunku do testowego MSE, które uzyskaliśmy używając \(\lambda = 4\). Ostatecznie, ponownie wyznaczamy nasz model regresji grzbietowej na pełnym zestawie danych, używając wartości \(\lambda\) wybranej w walidacji krzyżowej, i sprawdzamy oszacowania współczynników.

out = glmnet(x, y, alpha = 0) # Dopasuj model regresji grzbietowej do pełnego zbioru danych
predict(out, type = "coefficients", s = bestlam_ridge)[1:20,] # Wyświetlanie współczynników przy użyciu lambda wybranego przez CV

##  (Intercept)        AtBat         Hits        HmRun         Runs          RBI 
##  15.44834992   0.07716945   0.85906253   0.60120338   1.06366687   0.87936073 
##        Walks        Years       CAtBat        CHits       CHmRun        CRuns 
##   1.62437580   1.35296285   0.01134998   0.05746377   0.40678422   0.11455696 
##         CRBI       CWalks      LeagueN    DivisionW      PutOuts      Assists 
##   0.12115916   0.05299953  22.08942756 -79.03490992   0.16618830   0.02941513 
##       Errors   NewLeagueN 
##  -1.36075645   9.12528397

Zgodnie z oczekiwaniami, żaden ze współczynników nie jest dokładnie zerowy - regresja grzbietowa nie dokonuje selekcji zmiennych!

Regresja Lasso

Wprowadzenie

Zamiast regularyzacji \(L_2\), LASSO używa penalizacji \(L_1\), to znaczy:

\[ \min_{\beta \in \mathbb{R}^p} \sum_{i=1}^n (y_i - x_i^\top \beta)^2 + \lambda \|\beta\|_1. \]

Ze względu na charakter normy \(L_1\), LASSO ma tendencję do dawania bardziej rzadkich rozwiązań niż regresja grzbietowa. Jest to typowo użyteczne w ustawieniach wielowymiarowych, gdy prawdziwy model jest w rzeczywistości niskowymiarowym osadzeniem.

Model regresji lasso został pierwotnie opracowany w 1989 roku. Jest to alternatywa dla klasycznego oszacowania metodą najmniejszych kwadratów, która unika wielu problemów z nadmiernym dopasowaniem (overfittingiem), gdy mamy dużą liczbę niezależnych zmiennych.

Regresja Lasso (Least Absolute Shrinkage and Selection Operator) to technika regresji liniowej stosowana do oszacowania współczynników modelu, która wprowadza regularyzację \(L_1\). Regularyzacja L1 polega na dodaniu do funkcji celu kary proporcjonalnej do wartości bezwzględnej współczynników regresji.

Regresja Lasso ma zdolność do jednoczesnego wykonania selekcji cech i regularyzacji, co oznacza, że może pomóc w identyfikacji najbardziej istotnych cech modelu, a także zmniejszyć wpływ mniej istotnych cech.

Podstawowym celem regresji Lasso jest minimalizacja funkcji celu, która składa się z dwóch składników: błędu dopasowania (sumy kwadratów różnic pomiędzy rzeczywistymi wartościami odpowiedzi a przewidywanymi wartościami modelu) i kary regularyzacyjnej \(L_1\).

Wzór funkcji celu dla regresji Lasso może być przedstawiony jako: Minimize: RSS + \(\lambda \|\beta\|_1\), gdzie:

• RSS to suma kwadratów różnic pomiędzy rzeczywistymi wartościami odpowiedzi a przewidywanymi wartościami modelu (błąd dopasowania),

• \(\lambda\) (lambda) to parametr regularyzacji, który kontroluje siłę regularyzacji, a \(\|\beta\|_1\) to norma \(L_1\) współczynników regresji.

Dodanie kary regularyzacyjnej \(L_1\) powoduje, że niektóre współczynniki regresji stają się równe zero, co prowadzi do selekcji cech. Im większa wartość \(\lambda\), tym większa jest tendencja do redukcji współczynników do zera, prowadząc do bardziej rzadkiego modelu z mniejszą liczbą cech.

Regresja Lasso jest przydatna w przypadkach, gdy mamy do czynienia z wieloma cechami, z których niektóre mogą być nieistotne. Może pomóc w identyfikacji istotnych cech, redukcji nadmiaru danych i zwiększeniu interpretowalności modelu.

Przykład

Zobaczyliśmy, że regresja grzbietowa z mądrym wyborem \(\lambda\) może przewyższać metodę najmniejszych kwadratów, jak również model zerowy na zbiorze danych Hitters.

Teraz zobaczmy, czy lasso może dać albo dokładniejszy, albo bardziej interpretowalny model niż regresja grzbietowa.

W celu dopasowania modelu lasso, po raz kolejny używamy funkcji glmnet(), jednak tym razem używamy argumentu alpha=1. Poza tą zmianą postępujemy tak samo jak w przypadku dopasowywania modelu regresji grzbietowej:

lasso_mod = glmnet(x_train, 
                   y_train, 
                   alpha = 1, 
                   lambda = grid) # Dopasuj model lasso do danych treningowych

plot(lasso_mod)    # Wykreśl współczynniki

Zauważmy, że na wykresie współczynników, w zależności od wyboru dostrojenia parametru, niektóre ze współczynników są dokładnie równe zeru. Teraz przeprowadzimy walidację krzyżową i obliczymy związany z nią błąd testu:

set.seed(1)
cv.out = cv.glmnet(x_train, y_train, alpha = 1) # Dopasuj model lasso do danych treningowych
plot(cv.out) # Narysuj wykres MSE dla próby uczącej jako funkcję lambda

bestlam_lasso = cv.out$lambda.min # Wybierz lamdę, która minimalizuje MSE w próbie uczącej
lasso_pred = predict(lasso_mod, s = bestlam_lasso, newx = x_test) # Użyj najlepszej lambdy do przewidywania danych testowych
mean((lasso_pred - y_test)^2) # Oblicz MSE w próbie testowej

## [1] 143273

Jest to znacznie niższe MSE zbioru testowego niż modelu zerowego i modelu najmniejszych kwadratów, i bardzo podobny do MSE testu regresji grzbietowej z \(\lambda\) wybranej przez walidację krzyżową.

Jednakże lasso ma istotną przewagę nad regresją grzbietową w tym, że wynikowe oszacowania współczynników są rzadkie. Tutaj widzimy, że 12 z 19 oszacowań współczynników jest dokładnie zerowych:

out = glmnet(x, y, alpha = 1, lambda = grid) # Dopasuj model lasso do pełnego zbioru danych
lasso_coef = predict(out, type = "coefficients", s = bestlam_lasso)[1:20,] # Wyświetlanie współczynników przy użyciu lambda wybranego przez CV
lasso_coef

##   (Intercept)         AtBat          Hits         HmRun          Runs 
##    1.27429897   -0.05490834    2.18012455    0.00000000    0.00000000 
##           RBI         Walks         Years        CAtBat         CHits 
##    0.00000000    2.29189433   -0.33767315    0.00000000    0.00000000 
##        CHmRun         CRuns          CRBI        CWalks       LeagueN 
##    0.02822467    0.21627609    0.41713051    0.00000000   20.28190194 
##     DivisionW       PutOuts       Assists        Errors    NewLeagueN 
## -116.16524424    0.23751978    0.00000000   -0.85604181    0.00000000

Wybierając tylko predyktory o niezerowych współczynnikach widzimy, że model lasso z \(\lambda\) wybranym przez walidację krzyżową zawiera tylko siedem zmiennych:

lasso_coef[lasso_coef != 0] # Wyświetlanie tylko niezerowych współczynników

##   (Intercept)         AtBat          Hits         Walks         Years 
##    1.27429897   -0.05490834    2.18012455    2.29189433   -0.33767315 
##        CHmRun         CRuns          CRBI       LeagueN     DivisionW 
##    0.02822467    0.21627609    0.41713051   20.28190194 -116.16524424 
##       PutOuts        Errors 
##    0.23751978   -0.85604181

Twoja kolej!

Teraz nadszedł czas na przetestowanie tych metod (regresja grzbietowa i lasso) oraz metod oceny (zestaw walidacyjny, walidacja krzyżowa) na innych zbiorach danych. Możesz pracować z zespołem nad tą częścią laboratorium.

Możesz użyć dowolnego zbioru danych zawartego w ISLR lub wybrać jeden z pakietów danych na Kaggle/Data World itp. (zmienna zależna musi być ciągła).

Pobierz zbiór danych i spróbuj określić optymalny zestaw parametrów, które należy użyć do jego modelowania!

#Przygotowanie zbioru danych
College = na.omit(College)

• Który zbiór danych wybrałeś?

Wybrałam zbiór danych -> College (ISLR).

• Jaka była Twoja zmienna zależna (tzn. co próbowałeś modelować)?

Zmienną zależną w analizie jest Grad.Rate, czyli odsetek studentów, którzy ukończyli studia, będący wskaźnikiem sukcesu zarówno uczelni, jak i studentów. Większość zmiennych niezależnych w zbiorze danych może być uznana za istotne predyktory tej zmiennej. Na Grad.Rate mogą wpływać takie czynniki, jak liczba aplikacji (Apps), liczba przyjętych studentów (Accept), liczba studentów rozpoczynających naukę (Enroll), wysokość czesnego (Outstate), poziom kwalifikacji kadry akademickiej (PhD, Terminal) oraz dostępność zasobów uczelni, takich jak wydatki na studenta (Expend) i koszty zakwaterowania oraz wyżywienia (Room.Board).

#Tworzenie macierzy projektowej
x1 = model.matrix(Grad.Rate~., College)[,-1] # przycinam pierwszą kolumnę
                                              # zostawiam predyktory

#Wybór zmiennej objaśnianej
y1 = College %>%
  select(Grad.Rate) %>%
  unlist() %>%
  as.numeric()

Utworzenie macierzy projektowej x1, zawierającej wszystkie zmienne predykcyjne ze zbioru danych College, z pominięciem pierwszej kolumny (intercept). Następnie wyodrębniono kolumnę Grad.Rate jako zmienną objaśnianą y1, przedstawioną w formie wektora liczbowego.

Regresja grzbietowa

#Ustawienie siatki parametrów lambda
grid = 10^seq(10, -2, length = 100)
#Dopasowanie modelu regresji grzbietowej
ridge_mod2 = glmnet(x1, y1, alpha = 0, lambda = grid)

Utworzenie siatki potencjalnych wartości parametru regularizacji (lambda) oraz dopasowanie modelu ridge regression przy użyciu funkcji glmnet, gdzie ustawienie parametru alpha = 0 oznacza zastosowanie regresji grzbietowej.

#Wymiary współczynników
dim(coef(ridge_mod2))

## [1]  18 100

plot(ridge_mod2)    # wykres współczynników funkcji normy L1

Wyświetlenie wymiarów macierzy współczynników, gdzie pierwsza wartość (18) reprezentuje liczbę współczynników dla każdego poziomu parametru lambda. Dodatkowo, wizualizacja zmian współczynników w funkcji normy 𝐿1 ilustruje wpływ regularizacji na model.

ridge_mod2$lambda[50] # Wyświetl 50-tą wartość lambdy

## [1] 11497.57

coef(ridge_mod2)[,50] # Wyświetl współczynniki związane z 50-tą wartością lambdy

##   (Intercept)    PrivateYes          Apps        Accept        Enroll 
##  6.524147e+01  1.924300e-02  9.686567e-07  7.024066e-07 -6.112685e-07 
##     Top10perc     Top25perc   F.Undergrad   P.Undergrad      Outstate 
##  7.153334e-04  6.143038e-04 -4.144832e-07 -4.310351e-06  3.619870e-06 
##    Room.Board         Books      Personal           PhD      Terminal 
##  9.873884e-06  1.194325e-07 -1.014868e-05  4.754458e-04  5.000512e-04 
##     S.F.Ratio   perc.alumni        Expend 
## -1.970232e-03  1.009748e-03  1.899832e-06

sqrt(sum(coef(ridge_mod2)[-1,50]^2)) # Oblicz normę l2

## [1] 0.01940515

Wynik 11497.57 wskazuje na wysoki poziom regularizacji, który znacząco ogranicza wartości współczynników. Większość współczynników jest bliska zeru, co odzwierciedla silny wpływ regularizacji w redukowaniu ich wielkości. Z kolei wynik 0.01949515 dla 50-tej wartości lambda pokazuje bardzo niskie wartości współczynników, będące efektem zastosowania intensywnej regularizacji.

ridge_mod2$lambda[60] # Wyświetl 60-tą wartość lambdy

## [1] 705.4802

coef(ridge_mod2)[,60] # Wyświetl współczynniki powiązane z 60-tą wartość lambdy

##   (Intercept)    PrivateYes          Apps        Accept        Enroll 
##  6.222586e+01  2.855277e-01  1.464045e-05  1.093961e-05 -8.384396e-06 
##     Top10perc     Top25perc   F.Undergrad   P.Undergrad      Outstate 
##  1.046539e-02  9.031997e-03 -6.110823e-06 -6.531529e-05  5.334163e-05 
##    Room.Board         Books      Personal           PhD      Terminal 
##  1.453208e-04 -8.045918e-06 -1.524621e-04  6.855198e-03  7.141166e-03 
##     S.F.Ratio   perc.alumni        Expend 
## -2.809567e-02  1.502890e-02  2.710929e-05

sqrt(sum(coef(ridge_mod2)[-1,60]^2)) # Oblicz normę l2

## [1] 0.2878028

Wynik 705.482 wskazuje na niższy poziom regularizacji w porównaniu do 50-tej wartości lambda. Wartości współczynników są wyższe, co oznacza, że słabsza regularizacja pozwala na większą swobodę w ich ustalaniu. Z kolei wynik 0.2878028 jest wyższy niż dla 50-tej wartości lambda, co odzwierciedla mniejsze ograniczenie wartości współczynników przez regularizację.

#Predykcja współczynników dla zadanej wartości lambda
#Generuje współczynniki regresji `ridge` dla konkretnej wartości lambda (s = 50).
predict(ridge_mod2, s = 50, type = "coefficients")[1:18,]

##   (Intercept)    PrivateYes          Apps        Accept        Enroll 
##  4.671830e+01  1.903498e+00  1.234870e-04  1.146584e-04  3.405181e-06 
##     Top10perc     Top25perc   F.Undergrad   P.Undergrad      Outstate 
##  6.176559e-02  5.629912e-02 -3.664026e-05 -5.169350e-04  3.385363e-04 
##    Room.Board         Books      Personal           PhD      Terminal 
##  9.118980e-04 -5.822815e-04 -1.095459e-03  3.462916e-02  3.067362e-02 
##     S.F.Ratio   perc.alumni        Expend 
## -1.141984e-01  1.039575e-01  1.091441e-04

#Podział danych na zbiór treningowy i testowy 50/50
set.seed(1)

train = College %>%
  sample_frac(0.5)

test = College %>%
  setdiff(train)

#Tworzenie macierzy projektowych dla zbiorów treningowego i testowego
x_train = model.matrix(Grad.Rate~., train)[,-1]
x_test = model.matrix(Grad.Rate~., test)[,-1]

#Przygotowanie zmiennej objaśnianej
y_train = train %>%
  select(Grad.Rate) %>%
  unlist() %>%
  as.numeric()

y_test = test %>%
  select(Grad.Rate) %>%
  unlist() %>%
  as.numeric()

#Dopasowanie modelu `ridge` regression
ridge_mod2 = glmnet(x_train, y_train, alpha=0, lambda = grid, thresh = 1e-12)
#Predykcja dla zbioru testowego
ridge_pred2 = predict(ridge_mod2, s = 4, newx = x_test)
mean((ridge_pred2 - y_test)^2)

## [1] 167.8003

Utworzono model regresji grzbietowej z różnymi poziomami regularizacji. Dokonano predykcji na zbiorze testowym przy użyciu wartości lambda = 4. Obliczony błąd średniokwadratowy (MSE) dla predykcji wyniósł 167.8003. Niska wartość MSE wskazuje na dobrą jakość modelu regresji grzbietowej.

mean((mean(y_train) - y_test)^2)

## [1] 288.1504

#Porównanie wyników z różnymi wartościami lambda
ridge_pred2 = predict(ridge_mod2, s = 1e10, newx = x_test)
mean((ridge_pred2 - y_test)^2)

## [1] 288.1504

Dla wartości lambda = 10^10 (bardzo wysoka regularizacja) wynik MSE wynosi 288.1504. Tak wysoki poziom regularizacji powoduje znaczne uproszczenie modelu, co skutkuje większym błędem predykcji.

ridge_pred2 = predict(ridge_mod2, s = 0, newx = x_test)
mean((ridge_pred2 - y_test)^2)

## [1] 168.3203

#Regresja liniowa jako model bazowy
lm(Grad.Rate~., data = train)

## 
## Call:
## lm(formula = Grad.Rate ~ ., data = train)
## 
## Coefficients:
## (Intercept)   PrivateYes         Apps       Accept       Enroll    Top10perc  
##  34.3360980    2.4279951    0.0010389   -0.0013828    0.0078608    0.0371502  
##   Top25perc  F.Undergrad  P.Undergrad     Outstate   Room.Board        Books  
##   0.1233484   -0.0010036   -0.0018086    0.0009734    0.0019312   -0.0026554  
##    Personal          PhD     Terminal    S.F.Ratio  perc.alumni       Expend  
##  -0.0028311    0.1160148   -0.0463643    0.0201218    0.2333938   -0.0003622

#Predykcja współczynników dla wybranej lambda
predict(ridge_mod2, s = 0, type="coefficients")[1:18,]

##   (Intercept)    PrivateYes          Apps        Accept        Enroll 
## 34.3390554500  2.4433687719  0.0010142146 -0.0013205300  0.0076338298 
##     Top10perc     Top25perc   F.Undergrad   P.Undergrad      Outstate 
##  0.0395564696  0.1223166262 -0.0009717872 -0.0018151895  0.0009692031 
##    Room.Board         Books      Personal           PhD      Terminal 
##  0.0019333571 -0.0026600736 -0.0028335583  0.1153062666 -0.0457212412 
##     S.F.Ratio   perc.alumni        Expend 
##  0.0207052181  0.2336039536 -0.0003604071

Dla wartości lambda = 0 (brak regularizacji) wynik MSE wynosi 168.3203. Brak regularizacji umożliwia modelowi pełne dopasowanie do danych, co jednak nie zawsze prowadzi do lepszych wyników predykcji w porównaniu z modelem z optymalną wartością lambda.

#Walidacja krzyżowa dla regresji `ridge`
set.seed(1)
cv.out = cv.glmnet(x_train, y_train, alpha = 0) # Dopasuj model regresji grzbietowej na danych treningowych
bestlam_ridge = cv.out$lambda.min  # Wybierz lamdę, która minimalizuje treningowy MSE 
bestlam_ridge

## [1] 4.093486

Zatem wartość \(\lambda = 4\), która minimalizuje błąd walidacji krzyżowej, jest optymalna. Interpretacja optymalnej wartości lambda:

Optymalna wartość lambda stanowi kompromis między prostotą modelu a jego dokładnością predykcji.

Zbyt wysokie wartości lambda (nadmierna regularizacja) mogą prowadzić do zredukowania wpływu istotnych zmiennych, co skutkuje utratą ważnych informacji.
Z kolei zbyt niskie wartości lambda (lub brak regularizacji) mogą skutkować przeuczeniem modelu, gdzie dopasowuje się on zbyt dokładnie do danych uczących, tracąc zdolność do generalizacji na nowych danych.

plot(cv.out) # Narysuj wykres treningowego MSE jako funkcję lambda

Czerwona kropka wskazuje minimalny błąd MSE, odpowiadający optymalnej wartości lambda.

Szare paski reprezentują przedziały błędu dla różnych wartości lambda.

Obserwujemy, że:

Przy bardzo małych wartościach lambda, model charakteryzuje się wyższym błędem, co wynika z jego nadmiernej złożoności i ryzyka przeuczenia.
Przy bardzo dużych wartościach lambda, błąd również wzrasta, ponieważ model staje się zbyt uproszczony, tracąc zdolność do uchwycenia istotnych zależności w danych.

#Predykcja na zbiorze testowym przy użyciu optymalnej wartości lambda
ridge_pred2 = predict(ridge_mod2, s = bestlam_ridge, newx = x_test) # Użyj najlepszej lambdy do przewidywania danych testowych
mean((ridge_pred2 - y_test)^2) # Oblicz testowe MSE

## [1] 167.8427

Optymalna wartość lambda (bestlam uzyskana z walidacji krzyżowej) została zastosowana do przewidywania na danych testowych.

Wynikowy błąd średniokwadratowy (MSE) wynosi 167.8427.

Model regresji grzbietowej z optymalną regularizacją zapewnia wysoką jakość predykcji na danych testowych, skutecznie minimalizując błąd i równoważąc złożoność modelu z jego zdolnością do generalizacji.

#Dopasowanie modelu `ridge` regression do pełnego zbioru danych
out = glmnet(x1, y1, alpha = 0) # Dopasuj model regresji grzbietowej do pełnego zbioru danych
predict(out, type = "coefficients", s = bestlam_ridge)[1:18,] # Wyświetlanie współczynników przy użyciu lambda wybranego przez CV

##   (Intercept)    PrivateYes          Apps        Accept        Enroll 
## 35.7809244560  3.5091501069  0.0004522487  0.0003691766  0.0001702477 
##     Top10perc     Top25perc   F.Undergrad   P.Undergrad      Outstate 
##  0.0926780495  0.1037475794 -0.0001112942 -0.0012767754  0.0006979171 
##    Room.Board         Books      Personal           PhD      Terminal 
##  0.0017355002 -0.0021783531 -0.0018556999  0.0470636613 -0.0114667987 
##     S.F.Ratio   perc.alumni        Expend 
##  0.0319588603  0.2327598672 -0.0001406210

Współczynniki w modelu z optymalną wartością lambda są większe niż w przypadku bardzo wysokiej regularizacji, ale mniejsze niż w modelu bez regularizacji.

Taki wynik wskazuje na zrównoważony wpływ zmiennych na predykcję, co pozwala modelowi uwzględniać istotne zależności bez nadmiernego dopasowania lub uproszczenia.

Regresja Lasso

lasso_mod = glmnet(x_train, 
                   y_train, 
                   alpha = 1, 
                   lambda = grid) # Dopasuj model lasso do danych treningowych

plot(lasso_mod)    # Wykreśl współczynniki

Dopasowano model lasso regression (z alpha = 1) na zbiorze treningowym, uwzględniając różne wartości lambda.

W miarę wzrostu wartości lambda, wiele współczynników stopniowo zmniejsza się do zera, co świadczy o procesie selekcji zmiennych w modelu lasso.

Model lasso regression nie tylko regularizuje, ale także eliminuje nieistotne zmienne, co czyni go szczególnie użytecznym w analizach obejmujących dużą liczbę predyktorów, pomagając w identyfikacji najważniejszych czynników wpływających na wynik.

set.seed(1)
cv.out = cv.glmnet(x_train, y_train, alpha = 1) # Dopasuj model lasso do danych treningowych
plot(cv.out) # Narysuj wykres MSE dla próby uczącej jako funkcję lambda

bestlam_lasso = cv.out$lambda.min # Wybierz lamdę, która minimalizuje MSE w próbie uczącej
lasso_pred = predict(lasso_mod, s = bestlam_lasso, newx = x_test) # Użyj najlepszej lambdy do przewidywania danych testowych
mean((lasso_pred - y_test)^2) # Oblicz MSE w próbie testowej

## [1] 166.8336

Przeprowadzono walidację krzyżową dla modelu lasso regression (alpha = 1).

Wykres ilustruje błąd średniokwadratowy (MSE) w funkcji log(𝜆).
Czerwona kropka wskazuje wartość 𝜆, dla której błąd walidacji krzyżowej jest najmniejszy, co oznacza optymalną wartość regularizacji.
Szare paski reprezentują przedziały błędu dla różnych wartości 𝜆, a najlepsza wartość 𝜆 minimalizuje MSE i znajduje się w punkcie o najniższym błędzie.

Optymalna wartość 𝜆 dla modelu lasso została automatycznie wyznaczona na podstawie wyników walidacji krzyżowej. Model z tą wartością 𝜆 został użyty do przewidywania na zbiorze testowym.

Obliczono MSE dla predykcji na danych testowych, który wyniósł 168.8336.
Wynik pokazuje, że model lasso regression osiąga porównywalny poziom błędu z modelem `ridge` regression (MSE = 167.8427), co sugeruje podobną skuteczność obu podejść w analizowanym przypadku.

#Wyświetlenie współczynników modelu lasso dla pełnego zbioru danych

out = glmnet(x1, y1, alpha = 1, lambda = grid) # Dopasuj model lasso do pełnego zbioru danych
lasso_coef = predict(out, type = "coefficients", s = bestlam_lasso)[1:18,] # Wyświetlanie współczynników przy użyciu lambda wybranego przez CV
lasso_coef

##   (Intercept)    PrivateYes          Apps        Accept        Enroll 
## 34.9830711329  2.8429940894  0.0007166096  0.0000000000  0.0000000000 
##     Top10perc     Top25perc   F.Undergrad   P.Undergrad      Outstate 
##  0.0507860016  0.1303008938  0.0000000000 -0.0014704750  0.0009254499 
##    Room.Board         Books      Personal           PhD      Terminal 
##  0.0016861054 -0.0013031184 -0.0016745514  0.0139488663  0.0000000000 
##     S.F.Ratio   perc.alumni        Expend 
##  0.0000000000  0.2672573462 -0.0002371574

lasso_coef[lasso_coef != 0] # Wyświetlanie tylko niezerowych współczynników

##   (Intercept)    PrivateYes          Apps     Top10perc     Top25perc 
## 34.9830711329  2.8429940894  0.0007166096  0.0507860016  0.1303008938 
##   P.Undergrad      Outstate    Room.Board         Books      Personal 
## -0.0014704750  0.0009254499  0.0016861054 -0.0013031184 -0.0016745514 
##           PhD   perc.alumni        Expend 
##  0.0139488663  0.2672573462 -0.0002371574

Lasso wybrało podzbiór najważniejszych zmiennych, eliminując te, które nie miały istotnego wpływu na zmienną objaśnianą.

Wnioski

• Czy oczekiwałeś, że regresja grzbietowa będzie lepsza od lasso, czy odwrotnie? Jak wypada w stosunku do OLS? Pokaż odpowiednie raporty, miary dopasowania i krótko je omów (porównaj).

W analizie danych zawierających wiele predyktorów, często wzajemnie skorelowanych, modele regularizacyjne, takie jak regresja grzbietowa (ridge regression) i lasso, zwykle przewyższają klasyczną regresję liniową (OLS) pod względem zdolności do generalizacji na nowych danych.

Regresja grzbietowa okazuje się szczególnie skuteczna w przypadku wielu skorelowanych predyktorów, ponieważ rozkłada wagi między nimi w sposób bardziej równomierny, ograniczając problem kolinearności.

Lasso nie tylko regularizuje, ale także przeprowadza selekcję zmiennych, redukując współczynniki mniej istotnych predyktorów do zera. Dzięki temu upraszcza model i jest szczególnie użyteczne, gdy tylko kilka predyktorów ma istotny wpływ na zmienną zależną.

Na podstawie tych właściwości można oczekiwać, że:

Lasso przewyższy regresję grzbietową w sytuacjach, gdy niewielka liczba predyktorów jest istotna.
Regresja grzbietowa będzie bardziej odpowiednia, gdy wszystkie predyktory wnoszą wartość do modelu.

Porównanie modeli:

Aby zweryfikować skuteczność obu modeli, obliczamy błąd średniokwadratowy (MSE) na zbiorze testowym zarówno dla regresji grzbietowej (ridge regression), jak i regresji lasso (lasso regression). Obliczone wartości MSE pozwalają na ocenę zdolności każdego z modeli do generalizacji, czyli przewidywania wyników na nowych danych, poza zbiorem treningowym.

1. Regresja grzbietowa (ridge regression):

# Predykcja na zbiorze testowym z optymalną lambda
ridge_pred = predict(ridge_mod2, s = bestlam_ridge, newx = x_test)
# Obliczenie MSE
mse_ridge = mean((ridge_pred - y_test)^2)

MSE ridge: 167.8427

2. Regresja lasso:

# Predykcja na zbiorze testowym z optymalną lambda
lasso_pred = predict(lasso_mod, s = bestlam_lasso, newx = x_test)
# Obliczenie MSE
mse_lasso = mean((lasso_pred - y_test)^2)

MSE lasso: 168.8336

3. Regresja liniowa (OLS):

# Dopasowanie modelu OLS na zbiorze treningowym
ols_mod = lm(Grad.Rate ~ ., data = train)
# Predykcja na zbiorze testowym
ols_pred = predict(ols_mod, newdata = test)
# Obliczenie MSE
mse_ols = mean((ols_pred - y_test)^2)

MSE OLS: 168.3203

Regresja grzbietowa uzyskała najniższy MSE na zbiorze testowym, co sugeruje, że spośród trzech porównywanych modeli najlepiej przewiduje wartość Grad.Rate.

Regresja liniowa (OLS) osiągnęła MSE nieznacznie wyższe niż regresja grzbietowa, co wskazuje, że brak regularizacji nie spowodował znacznego pogorszenia wyników.
Regresja lasso miała nieco wyższe MSE niż pozostałe modele, jednak różnice były minimalne.

Wyniki wskazują, że regresja grzbietowa nieznacznie przewyższyła lasso i OLS pod względem MSE na zbiorze testowym. Oczekiwania dotyczące lepszej wydajności lasso nie potwierdziły się, co może być wynikiem charakterystyki danych, gdzie wszystkie predyktory wnoszą istotną informację do modelu.

Bliskie wyniki OLS sugerują również, że problem nadmiernego dopasowania lub wielokolinearności w tych danych nie był dominujący. Wnioski:

Wybór odpowiedniego modelu powinien zależeć nie tylko od miar dopasowania, takich jak MSE, ale również od celów analizy oraz potrzeby interpretowalności modelu. Regresja grzbietowa może być preferowana w przypadku danych o wielu skorelowanych predyktorach, podczas gdy lasso sprawdzi się lepiej w analizach wymagających selekcji zmiennych. OLS pozostaje prostym i skutecznym podejściem w sytuacjach, gdzie regularizacja nie jest konieczna.

• Które predyktory okazały się ważne w ostatecznym modelu (modelach)?

1. Regresja grzbietowa (ridge regression) W regresji grzbietowej wszystkie predyktory zachowują swoje współczynniki, ponieważ model nie eliminuje zmiennych, a jedynie regularizuje ich wartości. Ostateczne współczynniki dla optymalnej wartości𝜆 (wyznaczonej przez walidację krzyżową) są następujące:

predict(ridge_mod2, type = "coefficients", s = bestlam_ridge)[1:18,]

##   (Intercept)    PrivateYes          Apps        Accept        Enroll 
##  3.590439e+01  3.079648e+00  3.056605e-04  2.900056e-04  1.056032e-03 
##     Top10perc     Top25perc   F.Undergrad   P.Undergrad      Outstate 
##  9.179092e-02  9.254084e-02 -7.773098e-05 -1.602190e-03  6.641146e-04 
##    Room.Board         Books      Personal           PhD      Terminal 
##  1.715050e-03 -3.032236e-03 -2.722368e-03  6.887792e-02  1.412520e-02 
##     S.F.Ratio   perc.alumni        Expend 
## -7.720675e-03  2.133569e-01 -1.243617e-04

Regresja grzbietowa (ridge regression):

Wszystkie predyktory są uwzględnione w modelu, jednak ich wpływ jest regularizowany w celu zmniejszenia ryzyka przeuczenia. Wyższe wartości współczynników, takie jak PrivateYes czy Top25perc, wskazują na silniejszy wpływ tych zmiennych na wartość Grad.Rate. Regresja lasso (lasso regression):

2. Regresja lasso eliminuje mniej istotne predyktory, zmniejszając ich współczynniki do zera. Dzięki temu model automatycznie identyfikuje zmienne o największym wpływie na zmienną objaśnianą. To podejście upraszcza model i może być szczególnie przydatne w sytuacjach, gdy liczba predyktorów jest duża, a tylko niektóre z nich mają istotny wpływ.

Wybrane predyktory (współczynniki różne od zera):

predict(lasso_mod, type = "coefficients", s = bestlam_lasso)[lasso_coef != 0]

##  [1] 34.9844432060  2.0645478257  0.0004612948  0.0536851583  0.1128949262
##  [6] -0.0018452355  0.0008882663  0.0017012293 -0.0014600889 -0.0028088133
## [11]  0.0635913900  0.2347684266 -0.0001734067

Interpretacja: Lasso wyeliminowało mniej istotne predyktory, takie jak Apps, Accept, Enroll, F.Undergrad, czy P.Undergrad. Istotne zmienne: PrivateYes: Zmienna binarna, która wskazuje, czy uczelnia jest prywatna. Top10perc i Top25perc: Procent studentów z najlepszych wyników akademickich. Room.Board i Terminal: Koszty zakwaterowania oraz liczba nauczycieli z tytułem terminalnym (np. doktorem).

3. Regresja liniowa (OLS) OLS nie stosuje regularizacji ani selekcji zmiennych. Wszystkie predyktory są uwzględniane w modelu, a ich współczynniki są następujące:

Współczynniki OLS:

summary(ols_mod)$coefficients

##                  Estimate   Std. Error     t value     Pr(>|t|)
## (Intercept) 34.3360979737 7.3160760874  4.69323960 3.792381e-06
## PrivateYes   2.4279951257 2.3662117311  1.02611068 3.055098e-01
## Apps         0.0010388731 0.0006098643  1.70344958 8.932366e-02
## Accept      -0.0013827754 0.0012606191 -1.09690179 2.733977e-01
## Enroll       0.0078607547 0.0037550790  2.09336602 3.699686e-02
## Top10perc    0.0371501790 0.1060376862  0.35034883 7.262763e-01
## Top25perc    0.1233484128 0.0778745455  1.58393750 1.140623e-01
## F.Undergrad -0.0010035769 0.0006515845 -1.54021006 1.243639e-01
## P.Undergrad -0.0018085761 0.0006400996 -2.82546057 4.977280e-03
## Outstate     0.0009734411 0.0003208022  3.03439637 2.580524e-03
## Room.Board   0.0019311849 0.0008547150  2.25944891 2.443665e-02
## Books       -0.0026554051 0.0046945010 -0.56564159 5.719803e-01
## Personal    -0.0028310655 0.0010262313 -2.75870105 6.091398e-03
## PhD          0.1160148011 0.0790250027  1.46807715 1.429328e-01
## Terminal    -0.0463643398 0.0859115839 -0.53967507 5.897458e-01
## S.F.Ratio    0.0201217619 0.2415342143  0.08330812 9.336516e-01
## perc.alumni  0.2333938111 0.0706091036  3.30543512 1.041000e-03
## Expend      -0.0003622122 0.0001920999 -1.88554090 6.014007e-02

Interpretacja: Model OLS uwzględnia wszystkie zmienne, a predyktory o dużych wartościach współczynników, takie jak PrivateYes, Top10perc, i Top25perc, wydają się najbardziej istotne. Wnioski:

Lasso regression zidentyfikowało najbardziej istotne predyktory, eliminując mniej ważne zmienne, co jest szczególnie przydatne, gdy uproszczenie modelu ma kluczowe znaczenie.
`Ridge` regression uwzględniło wszystkie zmienne, regularizując ich współczynniki, co sprawdza się, gdy każda zmienna wnosi istotną informację do modelu.
Regresja liniowa (OLS) zawiera wszystkie predyktory, ale brak regularizacji może prowadzić do problemów z generalizacją w przypadku skorelowanych zmiennych.