Ricerca Sociale con R

Una wiki per l'analisi dei dati con R

Strumenti Utente

Strumenti Sito


Il Laboratorio di Analisi dei Dati con R, dell'Università di Teramo su piattaforma Meet, inizia il 9 aprile 2021 - Iscrizione - email
r:analisi_bivariata:regressione_lineare_bivariata

Regressione lineare bivariata

Per studiare l’esistenza, la forza e la direzione della relazione fra due variabili cardinali, possono essere utilizzati diversi tests statistici, in riferimento tanto all’andamento delle due variabili, quanto alle ipotesi della ricerca.

Rappresentazione grafica

Per controllare che la relazione fra le due variabili abbia un andamento almeno grosso modo lineare, è necessario partire dalla loro rappresentazione grafica: in questo caso, il grafico di dispersione o Scatterplot (grafico a dispersione).

Esempio

Utilizziamo il dataset di esempio cars incluso nel pacchetto datasets, e che raccoglie i dati relativi a velocità (speed) e spazio di frenata (dist) di 50 automobili. Ipotizziamo che lo spazio di frenata aumenti con la velocità: $Y=f(X)$, quindi X (indipendente) = speed e Y (dipendente) = dist

Eseguiamo le seguenti righe di comando:

# costruisco il grafico delle due variabili
# prima la variabile X, poi la variabile Y
plot(cars$speed, cars$dist, ylab = "dist", xlab="speed")
 
# aggiungo la retta dei minimi quadrati
abline(lm(cars$dist ~ cars$speed), col = "red")
&nolink |Scatterplot

con ggplot (vedi I grafici con ggplot2)

ggplot(cars, aes(speed, dist)) +    # sistema di riferimento
  geom_point() +                    # punti
  geom_smooth(method = "lm")        # geom_smooth: interpolazione
&nolink |Scatterplot

Il grafico prodotto conferma la linearità della relazione: lo spazio di frenata aumenta in maniera uniforme all'aumentare della velocità.

Osserviamo però anche che i punti non sono perfettamente allineati e si discostano dalla retta. Nel secondo grafico, l'area in grigio intorno alla retta rappresenta l'area dell'errore standard (*se*; vedi oltre).

Il modello

Il modello di regressione utilizza la funzione della retta per rappresentare la forma della relazione fra due variabili. Più la relazione è lineare (direttamente o indirettamente proporzionale), più i valori di Y potranno essere calcolati (predetti) in base a quelli di X, grazie alla funzione della retta $$Y = a+bX+e$$, dove:

  • $b$ = coefficiente di regressione, e indica l’inclinazione della retta (coefficiente angolare);
  • $a$ = intercetta della retta sulle ascisse;
  • $e$ = errore, ovvero la distanza fra valori attesi e valori reali.

La funzione utilizzata per lo studio dei modelli lineari è lm(), il cui argomento principale è la formula che rappresenta il modello da utilizzare, nella forma lm(dist ~ speed).

Eseguiamo dunque le seguenti righe di comando:

res <- lm(dist ~ speed, data=cars)
summary(res)

Analizziamo i risultati. In primo luogo, troveremo il comando che ha definito il modello stesso:

##  Call:
##  lm(formula = dist ~ speed, data = cars)

I residui

Seguono le statistiche descrittive della distribuzione dei residui, ovvero delle differenze fra valori osservati e valori attesi sulla base del modello:

##  Residuals:
##      Min      1Q  Median      3Q     Max 
##  -29.069  -9.525  -2.272   9.215  43.201 

Ci aspettiamo che i residui siano casuali e abbiano un andamento normale. In questo caso, i risultati indicano che i residui sono solo leggermente asimmetrici, come si nota dalla differenza interquartile e dal valore della mediana (la media dei residui è pari a 0).

L'errore standard dei residui è una misura dello scostamento dei valori attesi rispetto a quelli osservati, ed è pari a $\sqrt{\text{residui}^2/df}$.

sqrt(sum(res$residuals^2)/48)
##  [1] 15.37959

Possiamo ottenere questo valore con la funzione

sigma(res)
##  [1] 15.37959
##  Residual standard error: 15.38 on 48 degrees of freedom

I parametri (coefficienti)

Passiamo a questo punto a considerare i parametri del modello, ovvero i coefficienti a (intercetta, intercept) e b (indicato con il nome della variabile esplicativa, speed) della retta di regressione:

##  Coefficients:
##              Estimate Std. Error t value Pr(>|t|)    
##  (Intercept) -17.5791     6.7584  -2.601   0.0123 *  
##  speed         3.9324     0.4155   9.464 1.49e-12 ***
##  ---
##  Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Per ciascun parametro della retta, vengono forniti: il valore stimato (Estimate), l’errore standard (Std. Error), il valore t (t value), e il valore Pr(>|t|) che è il valore di probabilità del test t a due code.

L'intercetta

L’intercetta (-17,58) indica lo spazio di frenata quando la velocità è pari a 0, ovvero quando l’auto è ferma: valore che naturalmente non ha alcun senso in questi termini, ma che dipende essenzialmente dalle relative scale e non ha un significato sostantivo.

Il coefficiente di regressione

Il parametro che interessa è però il coefficiente di regressione è pari a 3,9324, e ci informa sulla forma della retta che interpola i punti, ma nulla sulla forza della relazione in quanto tale.

Se le due variabili covariassero in maniera lineare, la retta avrebbe questa forma.

I valori attesi

Le funzioni fitted() e predict() restituiscono dunque i valori attesi in base al modello:

$$\hat {\text{speed}} = -17,58+3,9 \text{dist}$$

-17.5791 + (3.9324 * cars$speed)

restituisce — a meno di approssimazioni — gli stessi valori di

fitted(lm.res)

La funzione predict() può però essere applicata ad un diverso dataset (appunto per predire i valori attesi in base al modello).

Intervallo di confidenza dei coefficienti

Per calcolare gli intervalli di confidenza dei parametri del modello (di questo e di altri), possiamo utilizzare la funzione confint():

confint(res)
 
##                  2.5 %    97.5 %
## (Intercept) -31.167850 -3.990340
## speed         3.096964  4.767853

Il coefficiente di determinazione R^2

Ma è il coefficiente di determinazione R2 la misura della qualità della stima prodotta dal modello, in termini di varianza spiegata: Nel nostro caso, quindi, R2 è pari a 0,65, a dire che la varianza di Y (dist) spiegata dal modello lineare è pari al 65% della sua varianza complessiva; con i tre outliers in meno, il valore di R2 sale a 0,72.

##  Multiple R-squared:  0.6511,	Adjusted R-squared:  0.6438 
##  F-statistic: 89.57 on 1 and 48 DF,  p-value: 1.49e-12 

L'output include anche il valore dell’R2 corretto, utile in caso di regressione multivariati.

Per valutare la significatività statistica di R2, si utilizza invece il test F: il numero di gradi di libertà del modello è pari al numero dei coefficienti (in caso di regressione bivariata, 1); il numero dei gradi di libertà della distribuzione congiunta delle due variabili è pari al numero dei casi, meno il numero delle medie (48).

Scomposizione dei quadrati (analisi della devianza) e test F

Coefficienti di regressione, determinazione e correlazione

In caso di regressione bivariata, il coefficiente di regressione può essere calcolato come rapporto fra covarianza e varianza della variabile assunta come indipendente:

$$b_{yx} = \frac{covXY}{varX} = \frac{\sigma_{xy}}{\sigma^2_x}$$

$$b_{xy} = \frac{covXY}{varY} = \frac{\sigma_{xy}}{\sigma^2_y}$$

cov(cars$speed, cars$dist) / var(cars$speed)
##  [1] 3.932409
 
cov(cars$speed, cars$dist) / var(cars$dist)
##  [1] 0.1655676

Il prodotto dei due coefficienti di regressione è uguale al coefficiente di determinazione $R^2$:

$$b_{yx} b_{xy} = R^2 = \frac{\text {devianza spiegata}}{\text {devianza totale}}$$

3.932409 * 0.1655676
##  [1] 0.6510795

La radice quadrata di $R^2$, e dunque del prodotto dei due coefficienti di regressione è il coefficiente di correlazione $r$.

$$ r = \sqrt{b_{yx} b_{xy}} $$

sqrt(3.932409 * 0.1655676)
##  [1] 0.806895

Considerando le relazioni fra queste grandezze, in caso di relazione bivariata, e una volta controllati gli assunti della correlazione, è sufficiente in fase esplorativa svolgere il test di correlazione (per avere una misura di significatività), e valutare l’adeguatezza del modello.

Script di esempio

E' possibile scaricare ed eseguire lo script dell'esempio:

Es-lm-bivariata.R
# grafico della retta dei minimi quadrati
plot(cars$speed, cars$dist, ylab = "dist", xlab="speed")
abline(lm(dist ~ speed), col = "red")
 
# con ggplot2
library(ggplot2)
 
ggplot(cars, aes(speed, dist)) +
  geom_point() +
  geom_smooth(method = "lm")
 
# lm  
res <- lm(dist ~ speed, data = cars)
summary(res)
 
# intervalli di confidenza
confint(res)
 
# errore standard dei residui
sigma(res)
 
# regressione e correlazione
cov(cars$speed, cars$dist) / var(cars$speed)
cov(cars$speed, cars$dist) / var(cars$dist)
3.932409 * 0.1655676
sqrt(3.932409 * 0.1655676)
r/analisi_bivariata/regressione_lineare_bivariata.txt · Ultima modifica: 25/09/2021 12:49 da admin