Grafico a dispersione in R: plot(), rette di regressione e pairs()

Cosa mostra un grafico a dispersione

Un grafico a dispersione mostra la relazione tra due variabili numeriche: ogni osservazione diventa un punto, posizionato in base al primo valore lungo l'asse orizzontale e al secondo lungo l'asse verticale. Se le due variabili si muovono insieme, i punti formano un disegno; se no, ottieni una nuvola senza forma. È il primo sguardo standard prima di stimare qualsiasi modello, e R lo disegna con la stessa funzione plot() vista nella guida a plot().

Useremo mtcars, un dataset integrato con 32 automobili, e ci faremo una domanda concreta: le auto più pesanti hanno motori più potenti?

Creare il grafico

wt è il peso in migliaia di libbre, hp è la potenza in cavalli:

plot(mtcars$wt, mtcars$hp,
     main = "Horsepower vs. weight",
     xlab = "Weight (1000 lbs)",
     ylab = "Horsepower",
     pch  = 19,
     col  = "steelblue")

L'immagine: trentadue punti pieni che salgono dall'angolo basso del grafico (auto leggere, attorno a 1.5 sulla scala del peso, con una potenza tra 60 e 90 cavalli) verso le auto pesanti oltre 5, che superano i 200 cavalli. La salita è inconfondibile ma non ordinata: a ogni peso i punti si distribuiscono su una fascia ampia di potenza.

Lo stile usa i soliti strumenti dei grafici di base: pch = 19 per cerchi pieni (quello vuoto di default sparisce negli screenshot), col per il colore, cex = 1.3 se i punti devono essere più grandi. Per colorare i punti in base a una terza variabile categorica, indicizza un vettore di colori con un factor: col = c("tomato", "steelblue", "darkgreen")[factor(mtcars$cyl)] dà un colore diverso a ogni numero di cilindri.

Leggerlo: direzione, forza, forma

Tre domande, in quest'ordine, ogni volta che guardi un grafico a dispersione:

• Direzione. I punti salgono (relazione positiva) o scendono (negativa) man mano che scorri lungo l'asse orizzontale? Qui salgono: più peso significa più potenza. Se invece disegni mpg rispetto a wt, la nuvola scende: più peso significa più consumi.
• Forza. Quanto stretti stanno i punti attorno a un unico percorso? Una fascia sottile come una matita è una relazione forte; uno spruzzo sparso è una relazione debole. Questa nuvola è moderatamente stretta.
• Forma e sorprese. Il percorso è dritto o curvo? Ci sono gruppi, o punti lontani da tutto il resto? In mtcars, la Maserati Bora spicca sopra il gruppo: 335 cavalli con un peso nella media. Un punto così può spostare parecchio una retta stimata, ed è proprio per questo che guardi prima di stimare.

Aggiungere la linea di tendenza

Un grafico a dispersione mostra una relazione; una retta che lo attraversa riassume l'affermazione. Stima un modello lineare e passalo direttamente ad abline():

plot(mtcars$wt, mtcars$hp,
     pch = 19, col = "steelblue",
     xlab = "Weight (1000 lbs)", ylab = "Horsepower")

abline(lm(hp ~ wt, data = mtcars), col = "tomato", lwd = 2)

lm(hp ~ wt) stima la retta dei minimi quadrati (leggi la formula come "hp spiegato da wt", con la variabile dell'asse verticale prima di ~) e abline() la disegna sul grafico. La retta sale di circa 46 cavalli ogni mille libbre. Cosa significa quel modello, e come leggerne il riepilogo, è l'argomento della regressione lineare.

Se non vuoi presumere una linea retta, lowess() disegna una curva morbida che segue i dati ovunque vadano:

lines(lowess(mtcars$wt, mtcars$hp), col = "darkgreen", lwd = 2, lty = 2)

Quando la curva lowess e la retta più o meno coincidono, un riassunto lineare è corretto. Quando la curva si piega e se ne allontana, la relazione non è lineare e una retta la rappresenterebbe male.

Controllare i numeri con cor()

Il grafico ti dà la forma; cor() ti dà la forza in un solo numero. Questo passaggio produce solo testo, quindi eseguilo qui:

# Correlation between weight and horsepower:
print(cor(mtcars$wt, mtcars$hp))

# And the full matrix for three variables at once:
print(round(cor(mtcars[, c("mpg", "wt", "hp")]), 2))

Peso e potenza hanno una correlazione di circa 0.66: la nuvola in salita moderatamente stretta, espressa come numero. La matrice aggiunge che mpg è fortemente correlato in modo negativo con entrambi (circa −0.87 con il peso). Rispetta però l'ordine delle operazioni: prima il grafico, poi il coefficiente. Un singolo valore di r può nascondere una curva o essere gonfiato da un solo valore anomalo: vedi correlazione per i modi classici in cui inganna.

La matrice di grafici a dispersione: pairs()

Con diverse colonne numeriche, disegnare a mano ogni coppia diventa noioso. pairs() lo fa con una sola chiamata:

pairs(mtcars[, c("mpg", "wt", "hp")],
      pch = 19, col = "steelblue")

Il risultato è una griglia 3 × 3: i nomi delle variabili stanno sulla diagonale, e ogni pannello fuori dalla diagonale è il grafico a dispersione di una coppia, mpg rispetto a wt, mpg rispetto a hp, wt rispetto a hp, ognuno presente due volte con gli assi scambiati. È il modo più rapido per fare una prima selezione su un nuovo dataset: con un'occhiata vedi quali coppie sono legate, quali relazioni curvano e dove si nascondono i valori anomali. Seleziona prima le colonne, come qui: oltre sei o sette variabili i pannelli diventano illeggibili.

La versione ggplot2

In ggplot2, il grafico a dispersione con una retta stimata è fatto di due livelli:

library(ggplot2)

ggplot(mtcars, aes(x = wt, y = hp)) +
    geom_point(color = "steelblue", size = 2) +
    geom_smooth(method = "lm", color = "tomato") +
    labs(title = "Horsepower vs. weight",
         x = "Weight (1000 lbs)", y = "Horsepower")

geom_smooth(method = "lm") è abline(lm(...)) con un bonus: una fascia di confidenza ombreggiata attorno alla retta. Se non imposti method, stima invece una curva loess, l'equivalente ggplot2 di lowess(). R base vince per velocità di scrittura quando vuoi dare un'occhiata veloce; ggplot2 vince appena vuoi i punti colorati per gruppo con una legenda automatica.

Cosa ti porti a casa

• plot(x, y) con due vettori numerici è un grafico a dispersione; pch = 19 e assi con etichette lo rendono presentabile.
• Leggi direzione, forza e forma, e individua i valori anomali, prima di calcolare qualsiasi cosa.
• abline(lm(y ~ x, data = df)) aggiunge la retta di regressione; lines(lowess(x, y)) aggiunge una curva che non presume una linea retta.
• cor() quantifica quello che il grafico mostra; il grafico mantiene onesto il numero.
• pairs(df[, cols]) disegna in un colpo tutti i grafici a dispersione a coppie: la prima selezione più rapida su un nuovo dataset.

Prossimo passo: il grafico a barre, per lasciarti alle spalle le coppie numeriche e confrontare conteggi tra categorie.