Esercizio - Regressione logistica su sei punti con discesa del gradiente

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

Testo. Sei osservazioni con una feature: $x=[1,2,3,4,5,6]$ e classi $y=[0,0,1,0,1,1]$.

1. Scrivere il modello, la perdita (log-verosimiglianza negativa media) e il gradiente.
2. Con $\beta=(\beta_0,\beta_1)=(0,0)$ calcolare perdita e gradiente, ed eseguire tre passi di discesa del gradiente con $\eta=0{,}1$.
3. Dire a quale soluzione converge l'algoritmo, dove sta il bordo di decisione, e calcolare l'accuracy sul training e la probabilità prevista per $x=3$ e $x=4{,}5$.

Teoria usata: Regressione logistica e softmaxLa regressione lineare non è adatta alla classificazione (valori fuori da [0,1], retta tirata dai punti lontani). La regressione logistica passa il predittore lineare dalla sigmoide $\sigma(z)=1/(1+e^{-z})$ e interpreta $\hat y=\sigma(x^T\beta)$ come $P(y=1\mid x)$: si predice la classe 1 se $\hat y\ge0{,}5$, cioè $x^T\beta\ge0$ (bordo lineare). L'errore quadratico dà una funzione non convessa; si usa la log-verosimiglianza negativa $-\sum[y\log\hat y+(1-y)\log(1-\hat y)]$, convessa, con gradiente $X^T(\hat y-y)$ e nessuna formula chiusa (discesa del gradiente). Per più classi: one-vs-one ($C(C-1)/2$ classificatori, voto), one-vs-all ($C$ classificatori, massima probabilità), o la softmax $p_c=e^{z_c}/\sum_ke^{z_k}$ con cross-entropia. Si può regolarizzare (ridge, LASSO, Elastic Net) e la cross-validation si fa stratificata. Approfondimento: non nel programma di Telecomunicazioni.Regressione logistica e softmax → (sigmoide, log-verosimiglianza, gradiente); 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 → (discesa del gradiente); logaritmi ed esponenziali in Esponenziale e logaritmoLa funzione esponenziale a^x (base positiva diversa da 1) e la sua inversa, il logaritmo in base a, con grafici e proprietà.Esponenziale e logaritmo →; metriche in Metriche di classificazioneIn classificazione binaria ogni previsione è vero positivo (TP), vero negativo (TN), falso positivo (FP, errore di tipo I) o falso negativo (FN, errore di tipo II). Da queste quattro quantità: accuracy $=\frac{TP+TN}{TP+TN+FP+FN}$, specificità $=\frac{TN}{TN+FP}$, precision $=\frac{TP}{TP+FP}$, recall $=\frac{TP}{TP+FN}$, e la loro media armonica $F_1=\frac{2PR}{P+R}$. Con dati sbilanciati l'accuracy inganna (un modello che predice sempre la classe maggioritaria ha 99%): si usano precision, recall, F1, ROC-AUC, la cross-validation stratificata e il riequilibrio con undersampling o oversampling (non SMOTE). Cambiando la soglia sulla probabilità si ottiene la curva ROC (TPR contro FPR) e l'area AUC. Approfondimento: non nel programma di Telecomunicazioni.Metriche di classificazione →.

1. Modello, perdita, gradiente

$\hat y_i=\sigma(\beta_0+\beta_1x_i)$, con $\sigma(z)=1/(1+e^{-z})$. Con la matrice $X=[\mathbf1,x]$ ($6\times2$): $\hat y=\sigma(X\beta)$. Perdita media $$\mathcal L(\beta)=-\frac1n\sum_{i=1}^n\big[y_i\ln\hat y_i+(1-y_i)\ln(1-\hat y_i)\big],\qquad\nabla\mathcal L=\frac1nX^T(\hat y-y).$$

2. Primo passo e tre iterazioni

Con $\beta=(0,0)$. $z_i=0$, quindi $\hat y_i=\sigma(0)=0{,}5$ per tutti. Perdita: ogni termine vale $-\ln0{,}5=0{,}693$, quindi $\mathcal L=0{,}693$ (è la perdita di un modello che dice sempre «50%»).

Errore $\hat y-y=[0{,}5;\ 0{,}5;\ -0{,}5;\ 0{,}5;\ -0{,}5;\ -0{,}5]$.

• Gradiente rispetto a $\beta_0$ (prodotto con la colonna di uni): $\frac16(0{,}5+0{,}5-0{,}5+0{,}5-0{,}5-0{,}5)=\frac16\cdot0=0$.
• Rispetto a $\beta_1$ (prodotto con $x$): $\frac16(0{,}5\cdot1+0{,}5\cdot2-0{,}5\cdot3+0{,}5\cdot4-0{,}5\cdot5-0{,}5\cdot6)=\frac16(0{,}5+1-1{,}5+2-2{,}5-3)=\frac{-3{,}5}6=-0{,}5833$.

Aggiornamento: $\beta\leftarrow(0,0)-0{,}1\cdot(0;\,-0{,}5833)=(0;\ 0{,}0583)$. Nuova perdita $0{,}6655$.

La perdita cala a ogni passo (la funzione è convessa e il passo è piccolo).

3. Convergenza e soluzione

Il passo è piccolo e la discesa è lenta: dopo $5000$ iterazioni $\beta=(-4{,}247;\ 1{,}213)$, dopo $20\,000$ $\beta=(-4{,}2491;\ 1{,}2140)$, con perdita $0{,}4130$ (valore confermato dal metodo di Newton e da scikit-learn senza penalità). La perdita è scesa da $0{,}693$ (modello indeciso) a $0{,}413$: misura quanto il modello ha imparato dai dati.

Bordo di decisione. $\sigma(z)=0{,}5\iff z=0\iff\beta_0+\beta_1x=0\iff x=-\beta_0/\beta_1=4{,}249/1{,}214=3{,}5$. Quindi si predice classe $1$ per $x>3{,}5$.

Probabilità previste: $x=1,\dots,6$: $0{,}046;\ 0{,}139;\ 0{,}353;\ 0{,}647;\ 0{,}861;\ 0{,}954$. Per $x=3$: $\sigma(-4{,}249+3{,}642)=\sigma(-0{,}607)=0{,}353$; per $x=4{,}5$: $\sigma(-4{,}249+5{,}463)=\sigma(1{,}214)=0{,}771$.

Accuracy sul training. Le previsioni sono $[0,0,0,1,1,1]$ contro $y=[0,0,1,0,1,1]$: due errori (il punto $x=3$ è classe 1 ma previsto 0; il punto $x=4$ è classe 0 ma previsto 1), quindi accuracy $=4/6=0{,}667$. I dati non sono linearmente separabili in una variabile, perciò nessuna soglia può fare meglio; la regressione logistica restituisce comunque probabilità ragionevoli (vicine a $0{,}5$ nella zona di confusione). Se i dati fossero stati separabili (per esempio $y=[0,0,0,1,1,1]$) la perdita avrebbe continuato a scendere aumentando $|\beta_1|$ senza minimo finito, e servirebbe la regolarizzazione.

Verifica

import numpy as np
x = np.array([1,2,3,4,5,6.]); y = np.array([0,0,1,0,1,1.]); X = np.c_[np.ones(6), x]
sig = lambda z: 1 / (1 + np.exp(-z))
nll = lambda b: -np.mean(y * np.log(sig(X @ b)) + (1 - y) * np.log(1 - sig(X @ b)))
b = np.zeros(2)
for t in range(20000):
    b -= 0.1 * X.T @ (sig(X @ b) - y) / 6          # gradiente della log-verosimiglianza media
print(b, nll(b), -b[0] / b[1])                      # [-4.249  1.214]  0.413  3.5