Spesso vogliamo spiegare la relazione tra una variabile di interesse \(Y\) ed un vettore di variabili esplicative \(X\), oppure vogliamo fare una previsione di \(Y\) sulla base dei valori noti di \(X\).

Shmueli spiega meglio di me la distinzione tra i due approcci: to explain or to predict? (Reference: https://www.stat.berkeley.edu/~aldous/157/Papers/shmueli.pdf).

Io provo a sintetizzarla in breve.

Se esiste una vera relazione (teoria sottostante) del tipo: \[Y = f(X) + \epsilon\]

dove \(\epsilon\) è il termine di errore imprevedibile e \(f\) è il segnale, a seconda dell’obiettivo (predictive vs explanatory) cambia la stima della funzione \(\widehat{f}\).

Modello esplicativo: si cerca il modello migliore minimizzando la distorsione \(B(\widehat{f})=E(f - \widehat{f})\), per avere una rappresentazione più accurata della teoria sottostante.

Modello previsionale: si cerca il modello migliore minimizzando congiuntamente la distorsione e la varianza \(V(\widehat{f})=E\left[\left(\widehat{f}-E(\widehat{f})\right)^2\right]\) del modello.

La conseguenza è che il miglior modello previsionale, se è quello che produce le migliori previsioni per nuovi valori di \(Y\) sulla base dei valori noti \(X\), potrebbe anche essere teoreticamente sbagliato.

Quando analizziamo dati reali, dobbiamo fronteggiare tutti i possibili problemi: il dataset iniziale non è quasi mai pulito, ci sono valori mancanti, valori anomali (e bisogna vedere di che tipo… rimuoverli alla cieca è molto pericoloso), errori di misurazione, variabili multicollineari, variabile target molto asimmetrica, relazioni tra variabili non lineari e difficili da modellare a priori, matrice dei dati sparsa, \(p \ge n\), …

..la lista è lunga.

I modelli non parametrici basati sugli alberi decisionali (Random Forest, Bagging, Gradient Boosting, …) hanno il vantaggio di essere più robusti e - nel caso di dati tabulari - più performanti in termini di capacità previsionali.

In particolare XGBoost (Extreme Gradient Boosting) è un algoritmo che ha vinto diverse competizioni e suscitano interesse per i risultati robusti anche in presenza di molti dei problemi citati sopra (Reference: https://www.sciencedirect.com/science/article/pii/S0169207021001874)

Il prezzo da pagare per non imporre a priori restrizioni alla possibile funzione \(f\) è quello di un costo computazionale più alto.

Un singolo albero decisionale cerca di prevedere il valore di \(Y\) (sia essa numerica o no) partizionando lo spazio degli input in modo da ottenere sotto-gruppi di unità omogenee nella variabile target.

XGBoost parte da un singolo albero e - ad ogni iterazione - lo migliora applicando alberi successivi (https://xgboost.readthedocs.io/en/stable/tutorials/model.html)

Esempio 1: Multicollinearità e outliers

\[Y = X_1 + 2X_2 + \pi X_3 + \epsilon, \quad \epsilon \sim N(0,1), \quad Cor(X_1,X_2)\ge0.9 \]

## 
## Call:
## lm(formula = y ~ ., data = data)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -19.131  -6.812  -1.937   2.144 141.560 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept)  48.8826    21.8655   2.236   0.0277 *  
## x1          -14.2533     5.5860  -2.552   0.0123 *  
## x2            0.7731     5.7826   0.134   0.8939    
## x3            2.5929     0.2073  12.507   <2e-16 ***
## x4            2.1198     1.7963   1.180   0.2409    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 18.77 on 95 degrees of freedom
## Multiple R-squared:  0.6467, Adjusted R-squared:  0.6318 
## F-statistic: 43.47 on 4 and 95 DF,  p-value: < 2.2e-16

Dalla tabella di summary del modello, vediamo che la regressione non riesce ad attribuire correttamente gli effetti delle variabili. Il vero vettore dei parametri è \(\beta = (0,1,2,\pi,0)\).

Anche guardando alle previsioni puntuali e ai grafici dei residui, si vede chiaramente quale modello riesce meglio a prevedere \(Y\) in presenza di multicollinearità e solo \(3\) osservazioni anomale. Le metriche di mean squared error sono:

## Esempio 1: Collinearità e outliers
## OLS MSE: 334.5319
## XGBoost MSE: 1.661577e-06

Esempio 2: Relazioni non lineari

\[ Y = 25\sin^2(X_1) + 10X_2 + 3\sqrt{X_2}+ \epsilon, \epsilon \sim N(0,1)\]

## 
## Call:
## lm(formula = y ~ ., data = data)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -14.0123  -7.3771  -0.3097   7.6633  15.0453 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept)  27.3857    15.8264   1.730   0.0868 .  
## x1           -0.1846     0.1539  -1.199   0.2335    
## x2           10.0036     0.5001  20.001   <2e-16 ***
## x3           -0.3012     0.2074  -1.452   0.1498    
## x4            0.4661     0.2665   1.749   0.0835 .  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 8.553 on 95 degrees of freedom
## Multiple R-squared:  0.8249, Adjusted R-squared:  0.8175 
## F-statistic: 111.9 on 4 and 95 DF,  p-value: < 2.2e-16

Dalla tabella di summary del modello, vediamo che la regressione non riesce ad attribuire correttamente gli effetti delle variabili. L’unico parametro stimato correttamente è il coefficiente di \(X_2\).

Anche guardando alle previsioni puntuali e ai grafici dei residui, si vede chiaramente quale modello riesce meglio a prevedere \(Y\) in presenza di relazioni non lineari. Le metriche di mean squared error sono:

## Esempio 2: Relazioni non lineari
## OLS MSE: 69.50064
## XGBoost MSE: 1.187448e-06