Regressione lineare semplice non parametrica

La caratteristica saliente dei metodi di regressione lineare non parametrica è costituita dalla loro robustezza. Questo significa che si tratta di metodi che risentono poco della eventuale presenza di dati apparentemente anomali, ma che non possono essere esclusi a priori dal calcolo della regressione.

Essendo n il numero delle coppie di dati/punti esistono n · (n - 1) / 2 modi di connettere due punti qualsiasi con una retta, cioè esistono n · (n - 1) / 2 possibili coefficienti angolari [1]. Nel caso più semplice il coefficiente angolare b della retta di regressione è la mediana di questi valori e l'intercetta a è la mediana degli n valori possibili calcolati mediante il coefficiente angolare trovato. Questi sono concettualmente i principi alla base della regressione lineare non parametrica, con possibili varianti da metodo a metodo.

Qui vediamo, confrontati con la regressione lineare semplice parametrica (regressione lineare ordinaria), tre modelli di regressione lineare non parametrica:

→ la regressione lineare di Thiel-Sen;

→ la regressione lineare di Siegel;

→ la regressione quantilica.

Da notare che questi tre modelli di regressione lineare non parametrica assumono la x come variabile indipendente analogamente alla regressione lineare ordinaria.

Le funzioni necessarie per il calcolo delle tre regressioni lineari non parametriche sono incluse in due pacchetti aggiuntivi, mblm e quantreg, che devono essere scaricati dal CRAN. Trovate come al solito la documentazione dei pacchetti e in particolare il loro manuale di riferimento sul repository della documentazione [2].

Copiate questo script, incollatelo nella Console di R e premete ↵ Invio:

# REGRESSIONE LINEARE ORDINARIA E NON PARAMETRICA DI THEIL-SEN, DI SIEGEL, QUANTILICA

#

library(mblm) # carica il pacchetto per le regressioni di Theil-Sen e Siegel

library(quantreg) # carica il pacchetto per la regressione quantilica

#

var_x <- c(18, 46, 104, 72, 24, 63, 23, 58, 36, 82, 94) # variabile in ascisse

var_y <- c(21, 42, 109, 64, 102, 78, 19, 66, 42, 89, 91) # variabile in ordinate

#

regpar_yx <- lm(var_y ~ var_x) # regressione lineare ordinaria x variabile indipendente

a_yx <- regpar_yx$coefficients[1] # intercetta a

b_yx <- regpar_yx$coefficients[2] # coefficiente angolare b

#

reg_TS <- mblm(var_y ~ var_x, repeated=FALSE) # regressione lineare di Theil-Sen

a_TS <- reg_TS$coefficients[1] # intercetta a

b_TS <- reg_TS$coefficients[2] # coefficiente angolare b

#

reg_S = mblm(var_y ~ var_x, repeated=TRUE) # regressione lineare di Siegel

a_S <- reg_S$coefficients[1] # intercetta a

b_S <- reg_S$coefficients[2] # coefficiente angolare b

#

reg_q <- rq(var_y ~ var_x, tau=0.5) # regressione lineare quantilica

a_q <- reg_q$coefficients[1] # intercetta a

b_q <- reg_q$coefficients[2] # coefficiente angolare b

#

# traccia il grafico dei dati con le rette di regressione e aggiunge una legenda

#

windows() # apre una nuova finestra

plot(var_x, var_y, xlim = c(0,120), ylim = c(0,120), pch=1, xlab="Variabile x", ylab="Variabile y", main="Regressioni ordinaria, Theil-Sen, Siegel, quantilica", cex.main = 0.9) # grafico dei dati

abline(a_yx, b_yx, col="blue", lty=4) # retta regressione ordinaria

abline(a_TS, b_TS, col="red", lty=4) # retta regressione di Theil-Sen

abline(a_S, b_S, col="green", lty=1) # retta regressione di Siegel

abline(a_q, b_q, col="orange", lty=1) # retta regressione quantilica

legend(60, 30, legend=c("Regressione ordinaria", "Regressione di Theil-Sen", "Regressione di Siegel", "Regressione quantilica"), col=c("blue", "red", "green", "orange"), lty=c(4,4,1,1), cex=0.8) # aggiunge al grafico la legenda

#

# crea una tabella con intercetta e coefficiente angolare delle quattro rette di regressione

#

cells <- c(a_yx,b_yx,a_TS,b_TS,a_S,b_S,a_q,b_q) # genera l'array cells con i valori numerici di a e di b

rnames <- c("Regressione ordinaria", "Regressione di Theil-Sen", "Regressione di Siegel", "Regressione quantilica") # genera l'array rnames con i nomi delle righe

cnames <- c("Intercetta a", "Coefficiente angolare b") # genera l'array cnames con i nomi delle colonne

tabris <- matrix(cells, nrow=4, ncol=2, byrow=TRUE, dimnames=list(rnames, cnames)) # costruisce la tabella con i risultati a partire dagli array cells, rnames, cnames

tabris # mostra la tabella con i risultati

#

Questo è il grafico di dispersione (xy) dei dati, con le rette calcolate:

I risultati delle quattro equazioni

y = a + b · x

sono così riportati nella Console di R:

> tabris # mostra la tabella con i risultati

                         Intercetta a Coefficiente angolare b

Regressione ordinaria       24.078677               0.7389267

Regressione di Theil-Sen     6.529412               0.9852941

Regressione di Siegel        5.873402               1.0042750

Regressione quantilica       2.581395               1.0232558

Le cose da notare sono:

→ è evidente la presenza di un dato che si discosta in modo importante da tutti gli altri;

→ i risultati dei tre metodi non parametrici sono simili tra loro;

→ i risultati delle regressione ordinaria (parametrica) si discostano da quelli delle regressioni non parametriche.

Per verificare quanto il dato, chiaramente anomalo, influenzi le regressioni, eseguiamo nuovamente lo script, sostituendo la terza e la quarta riga con queste due, nelle quali è stato eliminato il dato, il quinto, quello con le coordinate (x,y) uguali a (24, 102):

var_x <- c(18, 46, 104, 72, 63, 23, 58, 36, 82, 94) # variabile in ascisse

var_y <- c(21, 42, 109, 64, 78, 19, 66, 42, 89, 91) # variabile in ordinate

Ora i grafici delle rette sono praticamente sovrapposti (e difficilmente distinguibili tra loro)

mentre i risultati sono diventati i seguenti:

> tabris # mostra la tabella con i risultati

                         Intercetta a Coefficiente angolare b

Regressione ordinaria        1.337096                1.019512

Regressione di Theil-Sen     2.729167                1.020833

Regressione di Siegel        2.654334                1.022497

Regressione quantilica       2.581395                1.023256

Come si vede dal confronto con i risultati precedenti, con l'eliminazione del dato anomalo i risultati della regressione ordinaria sono cambiati: e questa fornisce ora risultati sovrapponibili a quelli dei tre metodi non parametrici. Poco sono cambiati, rispetti ai precedenti, i risultati della regressione di Thiel-Sen e di Siegel, mentre la regressione quantilica ha addirittura fornito identici valori nei due casi, un fatto molto interessante. I metodi non parametrici confermano la loro robustezza, mentre viene documentato quanto un singolo dato anomalo possa influire sulla regressione lineare ordinaria.

----------

[1] Per 100 coppie di dati il numero dei coefficienti angolari (sui quali calcolare la mediana) è uguale quindi a 4 950. Con 500 dati passiamo a 124 740. Con 1 000 dati passiamo a 499 500. 

[2] Available CRAN Packages By Name. URL consultato il 20/10/2018: https://goo.gl/hLC9BB