Usare i gRafici in R

Tavolozza

Predicare ai convertiti. Così dicono gli Americani.

Noi diremmo “sfondare una porta aperta”.

L’argomento è stato già trattato più volte. Comunque, repetita iuvant.

Uno dei vantaggi di R, e non dei minori, è la sua estrema ricchezza e flessibilità sul piano grafico. Nel suo formato base R ha un sofisticato sistema di rappresentazione grafica dei dati, che consente di avere una rappresentazione dei dati utile alla comprensione delle relazioni sottese.

Questa disponibilità grafica è ulteriormente arricchita dagli sviluppatori, che non mancano di dedicare degli script alla rappresentazione grafica dei dati secondo l’analisi proposta nel pacchetto cui contribuiscono.

Alcuni pacchetti sono stati sviluppati in modo specifico per rendere ancora più versatili le capacità grafiche di R. Alcuni di questi pacchetti gestiscono le forme geometriche, altri sono dedicati alle palettes, le disposizioni dei colori in un grafico, una utilità che chi deve creare grafici per lavoro è in grado di apprezzare con gratitudine.

L’importanza della rappresentazione grafica dei dati è testimoniata dalla frequenza con la quale sono scaricati i pacchetti specificamente dedicati.

Tra i pacchetti di R maggiormente scaricati figurano: ggplot2, colorspace, RColorBrewer, Scales, lattice.

Un passo indietro.

La rappresentazione grafica dei dati è importante, e permette di raggiungere una intuitiva e più immediata consapevolezza delle relazioni statistiche tra le variabili.

Usare R facilmente conduce alla scoperta di questa ovvietà.

Appena raggiunto il satori la prima reazione è quella di incazzarsi come una iena del Serengeti (copyleft, Paolo Villaggio, per i più distratti).

Se il programma statistico su cui avete studiato avesse avuto una migliore interfaccia grafica – vi verrà spontaneo pensare – avreste compreso meglio e più in fretta le basi della statistica. Ma la maggior parte dei programmi statistici commerciali concentra il focus sull’ottenimento della preziosa p < .05, e molte cose sono ad essa sacrificate.

Di seguito, facendo riferimento a dataset disponibili in R, riassumo qualche giochino utile; potrete copiare e incollare i comandi direttamente in R.

Caricate il pacchetto (già disponibile in R) e usate attach per rendere il dataset immediatamente accessibile. Dopo aver usato attach, potrete chiamare le variabili direttamente con il loro nome, risparmiandovi la seccatura di dovere specificare di continuo iris$variabile.

data(iris)
attach(iris)

Ora ottenete una descrizione del dataset:

dim(iris)
head(iris)

La funzione pairs vi consente di ottenere una matrice di correlazione delle variabili del dataset. Ogni variabile è posta in correlazione con le altre: intuitiva la comprensione di come la matrice deve essere letta. Sconsigliabile usarla con dataset con più di 10 variabili; nel caso, scomponete il dataset in vari insiemi. Come? Ripassate la lezione introduttiva su R… è possibile “forzare” pairs a fare cose parecchio interessanti date un’occhiata al manuale online.

pairs(iris)

pairs

Volete sapere se esiste una relazione tra lunghezza dei petali e lunghezza dei sepali in un fiore? Usate plot, che utilizza il formato plot(x,y): prima la variabile che va sull’asse delle x, poi quella che va sull’asse delle y.

plot(Petal.Length, Sepal.Length)

Un metodo alternativo per specificare l’ordine delle variabili in un grafico è usare il formato plot(y~x): la variabile che precede la tilde (il simbolo ~) è la variabile dipendente, che verrà posta sull’asse y, mentre la variabile inserita dopo la tilde è la variabile indipendente, che verrà posta sull’asse delle x.

plot(Sepal.Length ~ Petal.Length, data = iris)

Sembra proprio che esista una relazione lineare, con varianza spiegata pari al 75.8%. Volete esserne sicuri? Provate una regressione lineare:

fit  <- lm(Sepal.Length ~ Petal.Length, data = iris)
summary(fit)

Aggiungete la retta di regressione al grafico:

abline(fit)

Aggiungete la curva della local regression, una regressione non-parametrica. Cosa è lowess? Studiatevelo, vi tornerà utile!

fit2 <- lowess(Sepal.Length ~ Petal.Length)
lines(fit2, col = "red")

plot

Non male, la curva della local regression non si discosta troppo dalla retta della regressione lineare.

Volete l’intervallo di confidenza della retta della regressione lineare?

IntConf <- predict(fit, interval="confidence")

Nella matrice risultante, nelle colonne 2 e 3 si trovano i valori rispettivamente lower e upper dell’intervallo di confidenza della predizione. Le due linee dell’intervallo di confidenza vanno centrate sull’asse delle x (dove è il Petal.Length).

plot(Sepal.Length ~ Petal.Length)
abline(fit)
lines(Petal.Length, IntConf[,2], lty=3, col="blue")
lines(Petal.Length, IntConf[,3], lty=3, col="blue")

lines

Volete essere sicuri che la vostra regressione sia valida? Controllate le diagnostiche del modello. Sono quattro: nelle prime due (in alto) i residui dovrebbero avere una distribuzione casuale e la curva che li attraversa dovrebbe essere più o meno orizzontale; nel QQplot i vostri residui dovrebbero cadere sulla retta o non discostarsene troppo; infine, il quarto plot mostra i casi che violano la distanza di Cook, che dovrebbero essere pochi (sono numerati).

par(mfcol = c (2,2))
plot(fit)

La funzione par(mfcol = c(2,2)) serve a ripartire la lavagna grafica in quattro riquadri; se vi servono più riquadri appaiati, semplicemente variate gli input in c(num righe, num colonne).

plot(fit)

Nel complesso, non male.

Una migliore descrizione della violazione della distanza di Cook si ottiene con un plot del pacchetto car (da scaricare, se non lo avete), dove car sta per Companion to Applied Regression; creiamo un influence plot con punti proporzionali alla distanza di Cook (casi rilevanti identificati in rosso)

library(car)
layout(1)
influencePlot(fit, id.n=3,id.col="red")

Influence plot

Dato che prima, con il comando par(), avevamo ripartito la lavagna grafica in quattro riquadri, con layout(1) la riportiamo a un solo riquadro.

Sull’argomento delle regression diagnostics, utile riassunto qui, qui e qui. John Fox è il masterchef dell’argomento.

Volete avere una distribuzione della lunghezza dei petali nelle tre specie di fiori del dataframe? Perché non usare un bel density plot, sempre dal pacchetto car? La variabile di raggruppamento deve essere trattata come fattore e Species è già di classe “factor”. Non ci credete? Chiedetelo:

is.factor(Species)

La risposta è TRUE, mentre per is.numeric(Species) la risposta è FALSE. Se la vostra variabile non è un fattore, rendetela tale così (meglio assegnare un nuovo nome alla variabile, non si sa mai):

nuovo_nome_vostra_variabile <- as.factor(nome_vostra_variabile)

Ora possiamo andare sicuri a costruire i density plot:

par(mfcol = c(1,2))
densityPlot(Petal.Length~Species,data=iris)
densityPlot(Sepal.Length~Species,data=iris)

density plot

Bella differenza, eh?

Altra rappresentazione, ottenibile con un altro pacchetto, graphics (già installato e attivo nella distribuzione base di R) è il conditioning plot:

coplot(Sepal.Length ~ Petal.Length | Species, data = iris, panel = panel.smooth, rows = 1)

coplot

Il risultato suggerisce che per le tre specie esiste una relazione tra lunghezza dei petali e dei sepali, meno marcata per la specie setosa.

E adesso provate voi: caricare il dataframe mtcars:

data(mtcars)

Rapido check alle variabili:

head(mtcars)

E ora create un plot della relazione tra la variabile mpg (Miles/(US) gallon, miglia con un litro) e quella wt (weight, peso della macchina):

plot(mpg,wt)

Aha! messaggio di errore:

Errore in plot(disp, wt) : oggetto “mpg” non trovato

Capito perché è utile usare attach? Con attach potete nominare direttamente le variabili:

attach(mtcars)
plot(mpg,wt)

Ora funziona. Bene, proseguite con l’esercizio oppure variate a piacimento. Sbizzarritevi chiamando la lista dei dataset disponibili in R, digitando:

library(help = "datasets")

Per caricare il dataset, basta chiamarlo con: data(nome_del_dataset).

That’s all, folks!

Antonello Preti

Stay Tuned for our next episode!

Print Friendly

Lascia un Commento