Visualizzare relazioni curvilinee

Un modello statistico ha lo scopo di descrivere e prevedere uno o più fenomeni, semplificandoli, per renderli comprensibili all’uomo. Sotto questo punto di vista, i modelli statistici lineari costituiscono una soluzione ottimale per descrivere la realtà: sono semplici, facilmente interpretabili e adattabili a molte situazioni.

Nella sua formulazione più semplice, un modello lineare descrive la relazione tra due variabili attraverso una retta. Ci sono dei casi, però, in cui la relazione tra due variabili non può essere descritta attraverso una retta, ma è necessario interpolare i dati servendosi di una curva.

Descrivere una relazione di tipo curvilineo non implica necessariamente utilizzare un modello non lineare. La “linearità” di un modello, infatti, è determinata dalla forma matematica che questo assume e non dalla forma della relazione tra le variabili. Mi rendo conto che la cosa può apparire come un paradosso, ma è così: una relazione non lineare potrebbe benissimo essere descritta attraverso un modello lineare, purché questo modello sia lineare nei parametri.

In questo post vedremo un esempio di relazione non lineare e vedremo come visualizzare un modello di questo tipo utilizzando le funzioni basilari di R.

Un esempio: il compito di linea numerica


Una relazione non lineare è quella che si presenta quando un bambino piccolo, intorno ai 6/7 anni, è sottoposto al cosiddetto compito di linea numerica. Il compito consiste nel chiedere al bambino di posizionare una serie di numeri compresi entro un certo intervallo lungo un segmento. La versione più comune del compito richiede di posizionare lungo il segmento dei numeri che possono andare da 0 a 100.

Immaginiamo di presentare a un bambino trenta numeri compresi tra 0 e 100, che costituiscono gli stimoli, e di chiedergli di indicare in quale punto della linea dovrebbe essere posizionato ogni numero. Quindi, convertiamo ogni posizione indicata dal bambino nel corrispettivo numero che realmente si trova in quel punto della linea.

Registriamo nella variabile stimulus i numeri presentati e nella variable response le stime del bambino.

# Numeri presentati:
stimulus <- c(5,7,10,12,14,17,19,22,24,27,29,31,34,36,39,
            41,44,46,49,51,53,56,58,61,63,66,68,70,73,75)

# Risposte del bambino:
response <- c(38.22,42.67,40.63,45.81,44.09,44.15,45.7,52.6,
            49.51,48.13,49.68,51.83,49.4,53.48,50.06,48.76,
            53.22,51.54,53.04,56.08,53.24,50.9,52.12,53.11,
            55.95,54.81,53.96,54.62,50.92,55.61)
&#91;/code&#93;

Come ci ha insegnato Antonello Preti (<a href="">qui</a>), la primissima cosa che dovremmo fare è visualizzare la relazione tra le due variabili.

[code language="r"]
plot(stimulus, response)

Scatterplot tra stimolo e risposta nel compito di linea numerica

La relazione non sembra propriamente lineare, infatti i punti assumono un andamento leggermente curvilineo e tale curva sembra accentuata soprattutto per stimoli bassi piuttosto che elevati.

In ogni caso, se il bambino non ha risposto in maniera totalmente casuale, le sue risposte dovrebbero in qualche modo dipendere dagli stimoli presentati. Proviamo ad adattare un modello lineare:

model1 <- lm(response ~ stimulus)
summary(model1)
&#91;/code&#93;

Il modello presenta un indice R² pari a 0.6965. Non male. Vediamo allora come la retta di regressione interpola i punti:

&#91;code language="r"&#93;
plot(stimulus, response)
abline(model1, col="red")
&#91;/code&#93;

<img src="http://www.insular.it/wp-content/gallery/post-images/nlt-plot1.png" alt="Modello lineare risposta = a+b*stimolo" width="60%">

La retta non sembra interpolare benissimo i punti. Per carità, il modello non è malaccio, ma se diamo un'occhiata alla letteratura scientifica scopriremo che la relazione tra le due variabili non è lineare ma sembra essere <b>logaritmica</b>.

Quindi, a questo punto sorge un problema: come possiamo applicare un modello lineare se la relazione è non lineare (nel caso logaritmica)? Semplice: si trasforma <i>stimulus</i>, in modo da rendere lineare la sua relazione con <i>response</i>. L'assunzione è che la relazione tra il logaritmo naturale di <i>stimulus</i> e <i>response</i> sia lineare.

Adattare questo modello con R è molto semplice:

[code language="r"]
model2 <- lm(response ~ log(stimulus))
summary(model2)
&#91;/code&#93;

L'indice R² è decisamente superiore rispetto al precedente: 0.8275. La relazione potebbe davvero essere di tipo logaritmico. Come possiamo visualizzare questo modello? La funzione <b>abline</b> disegna rette, ma qui noi dobbiamo rappresentare una curva.

Per rappresentare un modello del genere si utilizza la funzione <b>lines</b>, in questo modo:

[code language="r"]
plot(stimulus, response)
lines(stimulus, predict(model2), col="blue")

Modello lineare risposta = a+b*log(stimolo)

La funzione predict estrapola i valori previsti dal modello. Alla funzione lines bisogna passare la variabile già disposta sull’asse orizzontale (stimulus) e i valori previsti dal modello, estrapolati proprio attraverso la funzione predict. Il risultato è la curva visualizzata in azzurro, che sembra interpolare i punti decisamente meglio rispetto alla retta del modello precedente.

Print Friendly, PDF & Email
Print Friendly, PDF & Email

2 Commenti per “Visualizzare relazioni curvilinee

  1. Davide MassiddaDavide Massidda Autore articolo

    Penso che quella retta derivi dal fatto che i dati non siano stati riordinati prima dell’uso del comando lines(). Provo a farti un esempio per chiarire.

    Immagina di avere queste due variabili:
    x = c(7,5,3,4,8,5,6,3,7,1,8,4,9,2,4,6)
    y = c(5.7,5.1,4.9,4.6,6.8,4.7,5.1,4.3,6.3,4.0,6.4,4.5,6.3,4.2,4.8,5.5)

    La relazione che c’è fra loro è curvilinea, come mostrato dal grafico:
    plot(x,y)

    Quindi, adatti un modello quadratico:
    model = lm(y~poly(x,2))

    Utilizzi predict() per individuare il valore di y previsto per ogni x:
    pred = predict(model)

    Quindi, provi a costruire il grafico:
    plot(x,y)
    lines(x,pred)

    Se provi a visualizzare questo grafico, noterai un gran casino, che penso sia esattamente quello che accade con i tuoi dati. Se invece utilizzi predict() non su x, ma su x riordinato in senso crescente, in questo modo:
    newdata = data.frame(x=sort(x))
    newdata$pred = predict(model, newdata)

    E provi ricostruire il grafico sfruttando i dati riordinati:
    plot(x,y)
    lines(newdata$x,newdata$pred)

    Otterrai un risultato migliore 😉

Lascia un commento