Esercizio - Cross-validation annidata per il LASSO sul dataset prostate

Questa pagina non ha ancora la versione ripasso: qui sotto c'è il testo completo.

Testo (laboratorio LAB3, parte 3). Dataset prostate: 96 osservazioni, 8 feature cliniche (lcavol, lweight, age, lbph, svi, lcp, gleason, pgg45) e target lpsa (logaritmo del PSA). Si vuole addestrare un LASSO con discesa del gradiente e scegliere gli iperparametri $\eta$ (passo) e $\lambda$ con la cross-validation annidata:

• passi $\eta\in\{10^{-5},10^{-4},10^{-3},10^{-2}\}$ e $\lambda\in\{10^{-4},10^{-3},10^{-2},0{,}1,1,10,100\}$ (28 combinazioni);
• ciclo esterno $K=8$ fold, ciclo interno $D=4$ fold, massimo 1000 iterazioni di discesa del gradiente.

Stimare MSE e $R^2$ del procedimento e commentare la scelta degli iperparametri.

Teoria usata: LASSO e discesa del gradienteIl LASSO è la regressione regolarizzata con penalità $L_1$: minimizza $\sum_i(y_i-x_i^T\beta)^2+\lambda\sum_{j\ge1}|\beta_j|$. A differenza della ridge, porta alcuni coefficienti esattamente a zero (soluzione sparsa, selezione delle feature): geometricamente le curve di livello dell'errore toccano il vincolo $\sum|\beta_j|\le s$ (un rombo) in uno spigolo. Non ha formula chiusa, quindi si minimizza con la discesa del gradiente $W\leftarrow W-\eta,\nabla J(W)$, usando il subgradiente $\operatorname{sign}(\beta_j)$ per il valore assoluto (nel punto 0 qualunque valore in $[-1,1]$). Il passo $\eta$ è critico: per l'errore quadratico converge se $\eta<1/\mu_{\max}(X^TX)$; si può usare $\eta_t=\eta_0/(1+\gamma t)$. L'Elastic Net combina le penalità $L_1$ e $L_2$ con $\lambda_1=\alpha\lambda$, $\lambda_2=(1-\alpha)\lambda$.LASSO e discesa del gradiente → (subgradiente, passo), Overfitting, ridge regression e cross-validationUna buona prestazione sul training non basta: serve stimare quella su dati nuovi. La cross-validation (K-fold: $k$ parti, ciascuna a turno come test, errore medio; Monte Carlo: $k$ divisioni casuali con quota di test $q$; leave-one-out se $k=n$) evita di dipendere da una sola divisione casuale. L'errore atteso si scompone in $\text{bias}^2+\text{varianza}+\sigma^2$: i modelli semplici fanno underfitting (bias alto), quelli complessi overfitting (varianza alta). La regolarizzazione aggiunge alla perdita una penalità: la ridge regression minimizza $|y-X\beta|^2+\lambda\sum_{j\ge1}\beta_j^2$ e ha soluzione $\beta=(X^TX+\lambda\tilde I)^{-1}X^Ty$ (l'intercetta non si penalizza, le feature si standardizzano): riduce i coefficienti, rende l'inversa stabile con feature collineari, e $\lambda$ è un iperparametro scelto con la validazione (cross-validation annidata per non contaminare il test).Overfitting, ridge regression e cross-validation → (cross-validation annidata, standardizzazione nei fold).

Procedura (per ogni fold esterno $i=1,\dots,8$)

1. Fold esterno: 12 campioni di test, 84 di training (e train+validation).
2. Ciclo interno: per ciascuna delle 28 combinazioni $(\eta,\lambda)$, si dividono gli 84 campioni in 4 fold da 21; per ogni fold interno: i 63 restanti sono il training interno, si standardizzano le feature con media e deviazione standard di quel training interno (e si applicano al fold di validazione), si addestra il LASSO con la discesa del gradiente ($\beta$ iniziale casuale $\sim N(0,1)$, aggiornamento $\beta\leftarrow\beta-\eta\,[2X^T(X\beta-y)+\lambda\operatorname{sign}(\beta)]$, intercetta non penalizzata) e si calcola l'MSE sul fold di validazione.
3. Per ogni combinazione si fa la media dei 4 MSE di validazione; si sceglie la combinazione con la media minore (le combinazioni che divergono danno inf o nan e non possono essere scelte).
4. Con la combinazione scelta si riaddestra il modello su tutti gli 84 campioni (standardizzati con le loro statistiche) e lo si valuta sui 12 campioni del fold esterno: MSE e $R^2=1-\mathrm{MSE}/\operatorname{Var}(y_{\text{test}})$.
5. Finiti gli 8 fold si riassumono MSE e $R^2$ (media e deviazione standard) e si riaddestra un modello finale su tutto il dataset con la combinazione più frequente.

Risultati (seme 0)

Media sui fold esterni: $\mathrm{MSE}=0{,}512$ (deviazione standard $0{,}202$) e $R^2=0{,}509$. Il modello spiega circa la metà della varianza del target, e l'errore dipende molto dal fold (da $0{,}23$ a $0{,}88$): con 12 campioni di test la stima di un singolo fold è molto rumorosa, ed è per questo che si fa la media.

Commento

• Passo. La scelta cade quasi sempre su $\eta=10^{-4}$ (7 fold su 8), una volta su $10^{-3}$. Il passo $10^{-2}$ diverge per ogni $\lambda$: per il training interno la matrice $X^TX$ (con la colonna di uni) ha autovalore massimo $\approx242$, quindi il limite di convergenza è $\eta<1/242\approx0{,}004$ (LASSO e discesa del gradienteIl LASSO è la regressione regolarizzata con penalità $L_1$: minimizza $\sum_i(y_i-x_i^T\beta)^2+\lambda\sum_{j\ge1}|\beta_j|$. A differenza della ridge, porta alcuni coefficienti esattamente a zero (soluzione sparsa, selezione delle feature): geometricamente le curve di livello dell'errore toccano il vincolo $\sum|\beta_j|\le s$ (un rombo) in uno spigolo. Non ha formula chiusa, quindi si minimizza con la discesa del gradiente $W\leftarrow W-\eta,\nabla J(W)$, usando il subgradiente $\operatorname{sign}(\beta_j)$ per il valore assoluto (nel punto 0 qualunque valore in $[-1,1]$). Il passo $\eta$ è critico: per l'errore quadratico converge se $\eta<1/\mu_{\max}(X^TX)$; si può usare $\eta_t=\eta_0/(1+\gamma t)$. L'Elastic Net combina le penalità $L_1$ e $L_2$ con $\lambda_1=\alpha\lambda$, $\lambda_2=(1-\alpha)\lambda$.LASSO e discesa del gradiente →); $10^{-2}$ lo supera e i valori esplodono (nan), $10^{-5}$ è così piccolo che in 1000 iterazioni non si arriva alla soluzione. Questo è il compromesso «passo troppo grande: diverge; troppo piccolo: lento».
• $\lambda$. Le scelte sono piccole ($10^{-4}$ ... $1$): con 8 feature standardizzate e 63-84 campioni l'OLS ha già poca varianza, e la penalità conta poco per i valori scelti. Le combinazioni con $\lambda$ grande non vengono mai scelte (l'MSE di validazione è peggiore). Nei modelli finali i coefficienti standardizzati più grandi sono quelli di lcavol ($\approx0{,}6$-$0{,}7$), poi svi ($0{,}24$-$0{,}34$) e lweight ($0{,}2$-$0{,}3$), con intercetta $\approx2{,}4$-$2{,}5$; nessun coefficiente è esattamente zero (il subgradiente non dà zeri esatti, vedi LASSO e discesa del gradienteIl LASSO è la regressione regolarizzata con penalità $L_1$: minimizza $\sum_i(y_i-x_i^T\beta)^2+\lambda\sum_{j\ge1}|\beta_j|$. A differenza della ridge, porta alcuni coefficienti esattamente a zero (soluzione sparsa, selezione delle feature): geometricamente le curve di livello dell'errore toccano il vincolo $\sum|\beta_j|\le s$ (un rombo) in uno spigolo. Non ha formula chiusa, quindi si minimizza con la discesa del gradiente $W\leftarrow W-\eta,\nabla J(W)$, usando il subgradiente $\operatorname{sign}(\beta_j)$ per il valore assoluto (nel punto 0 qualunque valore in $[-1,1]$). Il passo $\eta$ è critico: per l'errore quadratico converge se $\eta<1/\mu_{\max}(X^TX)$; si può usare $\eta_t=\eta_0/(1+\gamma t)$. L'Elastic Net combina le penalità $L_1$ e $L_2$ con $\lambda_1=\alpha\lambda$, $\lambda_2=(1-\alpha)\lambda$.LASSO e discesa del gradiente →).
• Perché annidata. Gli MSE di validazione (circa $0{,}46$-$0{,}56$) sono ottimistici rispetto a quelli del test (la media è $0{,}51$ ma con grande varianza): la combinazione scelta dipende dai fold interni, e il ciclo esterno dà una stima non influenzata dalla scelta.
• Modello finale. Si sceglie la combinazione più frequente (qui $(10^{-4},\ 0{,}1)$ e $(10^{-4},\ 1)$, due volte ciascuna) e si riaddestra su tutti i 96 campioni.

Codice (struttura)

for i in range(K):                                         # ciclo esterno
    X_test, y_test = fold i;  X_tr, y_tr = gli altri
    best, best_mse = None, np.inf
    for lr, lam in combos:                                 # 28 combinazioni
        mses = []
        for j in range(D):                                 # ciclo interno
            X_in, y_in, X_val, y_val = split(X_tr, y_tr, j)
            mu, sd = X_in.mean(0), X_in.std(0)             # statistiche del solo training interno
            b = train_lasso_gd((X_in - mu) / sd, y_in, lam, lr, max_iter=1000)
            mses.append(np.mean((y_val - (b[0] + ((X_val - mu) / sd) @ b[1:])) ** 2))
        if np.isfinite(np.mean(mses)) and np.mean(mses) < best_mse:
            best, best_mse = (lr, lam), np.mean(mses)
    # riaddestra su X_tr con `best` e valuta su X_test -> MSE, R2

Il numero totale di addestramenti è $8\times28\times4=896$ per la ricerca più $8$ per la valutazione. Il codice completo (con train_lasso_gd e lasso_gradient) è in LASSO e discesa del gradienteIl LASSO è la regressione regolarizzata con penalità $L_1$: minimizza $\sum_i(y_i-x_i^T\beta)^2+\lambda\sum_{j\ge1}|\beta_j|$. A differenza della ridge, porta alcuni coefficienti esattamente a zero (soluzione sparsa, selezione delle feature): geometricamente le curve di livello dell'errore toccano il vincolo $\sum|\beta_j|\le s$ (un rombo) in uno spigolo. Non ha formula chiusa, quindi si minimizza con la discesa del gradiente $W\leftarrow W-\eta,\nabla J(W)$, usando il subgradiente $\operatorname{sign}(\beta_j)$ per il valore assoluto (nel punto 0 qualunque valore in $[-1,1]$). Il passo $\eta$ è critico: per l'errore quadratico converge se $\eta<1/\mu_{\max}(X^TX)$; si può usare $\eta_t=\eta_0/(1+\gamma t)$. L'Elastic Net combina le penalità $L_1$ e $L_2$ con $\lambda_1=\alpha\lambda$, $\lambda_2=(1-\alpha)\lambda$.LASSO e discesa del gradiente →.