L’analisi di regressione multipla, spiegata semplice.

I fenomeni cui assistiamo, e che vogliamo studiare per approfondirne la comprensione, raramente si presentano in maniera così semplice da potersi definire attraverso due sole variabili, di cui una predittiva (indipendente) e una responso (dipendente).

Per questo, se pure l’analisi di regressione lineare presenta una fondamentale importanza teorica, nella pratica fornisce poca informazione in più rispetto allo studio attraverso il semplice coefficiente di correlazione.

E’ per questo che in letteratura scientifica è la regressione multipla a svolgere un ruolo preponderante, potendo fornire un set completo di strumenti in grado di spiegare la variazione della nostra variabile dipendente per ogni diversa variabile predittiva presente nel modello, e per il complesso delle interazioni delle variabili indipendenti.

L’equazione della regressione multipla

Partiamo dall’equazione di regressione multipla, che ovviamente si presenta come una “espansione” di quella della regressione lineare semplice, ed ha questa forma generale:

\( \\ y = a_1 x_1 + a_2 x_2 + … + a_i x_i + b \\ \)

dove
y è la variabile responso (n.b. : è unica)
a1,a2…ai sono i coefficienti di regressione per le variabili predittive
x1,x2…xi sono le variabili predittive
b è l’intercetta (ed è anch’essa unica)

Quali informazioni posso ricavare?

Posso e devo, in prima istanza, capire quanto le mie variabili predittive, combinate, siano relate in maniera significativa alla variabile dipendente. E quanta parte del risultato è spiegata dalla combinazione delle variabili predittive usate nel modello.

Poi devo capire quanto ciascuna variabile predittiva sia legata alla variabile dipendente, mentre controllo l’altra variabile indipendente (suppongo per semplicità che le variabili predittive siano solo due).

Devo poi capire quale dei miei due predittori sia quello più “forte” nello stimare la variabile dipendente.

Queste e molte altre cose possono essere esaminate con gli strumenti messi a disposizione dall’analisi di regressione multipla e questi primi accenni mostrano già la potenza e l’utilità pratica di questa tecnica.

In pratica, come procedere?

Il mio suggerimento per impratichirsi è quello di procedere per passaggi. Io suggerisco questi:

  1. Tracciamo i grafici di dispersione (scatterplot) per ogni variabile predittiva rispetto alla variabile responso, per valutare l’esistenza o meno di correlazione.
  2. Calcoliamo l’equazione di regressione multipla.
  3. Verifichiamo che siamo rispettate le assunzioni di normalità e omoschedasticità.
  4. Analizziamo il coefficiente di determinazione.

Per ognuno di questi passaggi, R si dimostra come di consueto un compagno fidatissimo, capace di fornirci le elaborazioni richieste con pochi comandi.

Mettiamoci all’opera!

Scaldiamo i motori di R e vediamo un esempio pratico, sia pur semplificando al massimo. Ovviamente, sfioriamo solo la superficie dell’argomento. Lo scopo infatti è quello di gettare le prime basi per iniziare la navigazione in un argomento sterminato. Starà alla passione di chi legge continuare ad approfondire…

Basta con i preamboli: scegliamo un dataset di esempio.
R ci fornisce tra gli altri il Longley’s Economic Regression Data, che userò in questo articolo.
Di cosa si tratta? Di un data frame che contiene 7 variabili di tipo economico, dall’anno 1947 al 1962.

Apriamo R Studio, carichiamo il dataset e diamogli una prima occhiata:

data(longley)
dim(longley)
head(longley)

Vediamo le varie metriche. Per il nostro esempio, andremo a testare come variabile responso il numero degli occupati, e come variabili predittive indipendenti il GNP (Gross National Product, il prodotto interno lordo) e la popolazione.

Ok, vediamo come si presentano gli scatterplot:

plot(longley$Employed, longley$GNP)
impiegati gnp
impiegati /prodotto interno lordo
plot(longley$Employed, longley$Population)
Residui regressione multipla
impiegati / popolazione

La correlazione crescente di tipo lineare appare evidente. Procediamo con la regressione multipla. Il comando lo conosciamo già: è lm(). La sintassi per più variabili predittive è la seguente:

regressione = lm(Employed ~ GNP + Population, data=longley)
summary(regressione)

Il risultato è una miniera di informazioni!

Call:
lm(formula = Employed ~ GNP + Population, data = longley)

Residuals:
     Min       1Q   Median       3Q      Max 
-0.80899 -0.33282 -0.02329  0.25895  1.08800 

Coefficients:
            Estimate Std. Error t value Pr(>|t|)    
(Intercept) 88.93880   13.78503   6.452 2.16e-05 ***
GNP          0.06317    0.01065   5.933 4.96e-05 ***
Population  -0.40974    0.15214  -2.693   0.0184 *  
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 0.5459 on 13 degrees of freedom
Multiple R-squared:  0.9791,	Adjusted R-squared:  0.9758 
F-statistic: 303.9 on 2 and 13 DF,  p-value: 1.221e-11

Notiamo per prima cosa il valore p della Statistica F. E’ un valore molto piccolo (1.211e-11), altamente significante. Quindi, almeno una delle nostre variabili predittive è relata in maniera statisticamente significativa alla variabile risultato. Dunque, procediamo.

Tracciamo anche il grafico dei residui:

plot(longley$Year, summary(regressione)$res, type='o')
abline(h=0, lty=2, col="red", lwd=2)
res

Idealmente, la somma dei residui, cioè le differenze tra i valori reali e quelli predetti, dovrebbe tendere a zero, o quantomeno essere la più bassa possibile. Uno sguardo a questa sezione dell’output del comando summary ci dice che in questo caso la condizione è rispettata:

residuals
la mediana dei residui è vicina allo 0, vero?

Controlliamo graficamente l’assunzione di normalità dei residui:

hist(resid(regressione),main='Istogramma dei residui',xlab='Residui Standardizzati',ylab='Frequenza')

Il grafico che otteniamo ci mostra che l’assunzione di normalità è rispettata:

istogramma dei residui
i residui rispettano l’assunzione di normalità

E’ giunto il momento di andare a vedere in dettaglio cosa ci dice l’output del comando summary relativo alla nostra regressione.

La mediana è vicina allo 0, i residui sono distribuiti in maniera approssimativamente normale. Quindi possiamo procedere.

I valori dei coefficienti angolari per le nostre due variabili predittive e quello dell’intercetta li troviamo nel nostro output:

coefficienti
R ci fornisce…la pappa pronta!

La retta di regressione è dunque:

\( \\ y = 0.06317x_1 – 0.40974 x_2 + 88.9388 \\ \)

Prestiamo attenzione ai valori p per quanto riguarda i coefficienti angolari. L’output ci dice che le nostre variabili indipendenti sono significative per predire il valore della variabile dipendente. I p-value infatti sono inferiori al valore standard 0.05. Il prodotto nazionale lordo ha addirittura un valore inferiore al livello dello 0.001. Possiamo dunque rigettare l’ipotesi nulla (le variabili predittive non sono significative) e ritenere valido il nostro modello, notando come il GNP sia l’elemento più “fedele” per la stima rispetto alla popolazione. Molto utile la rappresentazione con gli asterischi restituita dal nostro output, capace di fornirci subito il colpo d’occhio del risultato.

pvalue 1
Occhio agli asterischi!

Possiamo anche chiedere a R il calcolo dell’intervallo di confidenza per il nostro modello, con il comando:

confint(regressione)

Nel nostro esempio otterremo questo output:

                  2.5 %       97.5 %
(Intercept) 59.15805807 118.71953854
GNP          0.04017053   0.08617434
Population  -0.73841447  -0.08107138

Quanto è valido il mio modello?

A questo punto entra in scena un protagonista che abbiamo già incontrato nei precedenti articoli, il coefficiente di determinazione, o r2

E’ proprio lui a dirci una cosa fondamentale, cioè quanto vicini sono i dati alla nostra retta di regressione. In altri termini più pratici, quale percentuale dei “movimenti” della variabile dipendente sono spiegati dalle nostre variabili predittive. Il valore è standardizzato tra 0 e 1, ed è superfluo notare come il nostro modello sia tanto più utile quanto più questo valore si avvicina a 1.

R, come al solito formidabile, ci ha già dato in output questo utilissimo valore. Eccolo:

coefficiente di determinazione
quando due variabili spiegano tutto…o quasi!

Ok, ma cos’è il valore Adjusted R-squared? E’ il valore da considerare, per superare un paradosso del valore r2, che ha un valore sempre crescente al crescere del numero delle variabili (anche se quelle variabili non sono affatto significative). Il valore r quadro adjusted corregge questa anomalia e restituisce un valore (sempre inferiore a r quadro) perfettamente utilizzabile.

Sintesi finale

Sarò brevissimo: un valore elevato di r quadrato e un valore molto basso, tendente a zero, dei residui indicano che il modello è buono.

Non è tutto, anzi è solo l’inizio. Ma è il primo passo per iniziare a padroneggiare uno strumento, come l’analisi di regressione lineare multipla, che ha grandissima utilità pratica. Vi lascio, spero, la curiosità di approfondire e andare oltre. Buon lavoro!




Lascia un commento

Il tuo indirizzo email non sarà pubblicato. I campi obbligatori sono contrassegnati *