Regressione lineare multipla

La regressione lineare semplice (regressione bivariata) costruita sul piano cartesiano

y = a + b·x

può essere estesa a uno spazio n-dimensionale: abbiamo così la regressione lineare multipla nella quale la variabile dipendente y viene fatta dipendere non più da una, bensì da k variabili indipendenti, indicate come x₁, x₂, x₃ … xₖ secondo l'equazione

y = a + b₁·x₁ + b₂·x₂ + b₃·x₃ + ... + bₖ·xₖ

Il risultato è rappresentato quindi da una intercetta a e da tanti coefficienti angolari b (indicati con b₁, b₂, b₃ ... bₖ) quante sono le variabili indipendenti (x₁, x₂, x₃ ... xₖ) che contribuiscono a determinare il valore della variabile dipendente y.

Come esempio impieghiamo il set di dati ais nel quale, tra i dati ematologici rilevati in un gruppo di atleti australiani, abbiamo tre grandezze:

→ la concentrazione degli eritrociti (globuli rossi) espressa in milioni per microlitro di sangue (10⁶/µL) e contenuta nella variabile rcc;

→ il valore ematòcrito cioè la frazione del volume del sangue occupata dagli eritrociti espressa in percentuale (%) e contenuta nella variabile hc;

→ la concentrazione dell'emoglobina espressa in grammi per decilitro di sangue (g/dL) e contenuta nella variabile hg.

Le tre grandezze sono strettamente legate tra loro dal punto di vista biologico, in quanto gli eritrociti hanno ciascuno un certo volume e contengono ciascuno una data quantità di emoglobina: pertanto all'aumentare del numero di eritrociti aumenta la frazione percentuale del volume del sangue occupata dagli eritrociti (il valore ematòcrito) e aumenta la quantità di emoglobina. Mentre la diminuzione dell'ematòcrito e la diminuzione dell'emoglobina sono indice della diminuzione del numero degli eritrociti.

L'esempio è stato scomposto in una serie di script indipendenti - onde facilitare il riutilizzo del codice per analizzare i propri dati - che impiegano i pacchetti aggiuntivi gvlma, car e rgl, che devono essere scaricati dal CRAN. Accertatevi di avere installato inoltre il pacchetto DAAG che contiene i dati impiegati nello script, o in alternativa  procedete come indicato nel post Il set di dati ais. 

Iniziamo con una analisi esplorativa preliminare dei dati, copiate questo primo script, incollatelo nella Console di R e premete ↵ Invio:

# REGRESSIONE LINEARE MULTIPLA - analisi esplorativa preliminare dei dati

#

library(DAAG) # carica il pacchetto che include il set di dati ais

library(car) # carica il pacchetto per la grafica

mydata <- ais[,c(1,3,4)] # le colonne 1, 3, 4 contengono le variabili rcc, hc, hg

windows() # apre e inizializza una nuova finestra grafica

#

scatterplotMatrix(~rcc+hc+hg, col="black", pch=20, regLine=list(method=lm, lty=1, lwd=2, col="chartreuse3"), smooth=FALSE, diagonal=list(method="histogram", breaks="FD"), main="Istogrammi e grafici xy con rette di regressione", data=mydata) # istogrammi e grafici xy con rette di regressione

#

Dopo avere caricato i pacchetti necessari, copiato nell'oggetto mydata le tre variabili che ci interessano, aperto e inizializzato una nuova finestra grafica, l'analisi preliminare dei dati viene effettuata con la funzione scatterplotMatrix() nella quale:

→ ~rcc+hc+hg specifica le variabili da rappresentare, che sono la concentrazione degli eritrociti rcc, il valore ematòcrito hc e la concentrazone di emoglobia hg;

→ col="..." specifica il colore dei grafici; 

→ pch=... specifica il tipo di simbolo da impiegare (un piccolo cerchio pieno);

→ per la retta di regressione regLine sono impiegati il modello lineare (lm), una linea continua (lty=1) con uno spessore doppio (lwd=2) e un colore (col="chartreuse3"), mentre l'argomento smooth=FALSE fa si che non venga aggiunta ai grafici la rappresentazione della funzione non lineare prevista di default nel pacchetto car;

→ l'argomento diagonal=list() indica la rappresentazione della variabile riportata nella diagonale, che qui è l'istogramma (method="histogram") con il numero di classi (breaks="FD") calcolato mediante la regola di Freedman-Diaconis. 

Da notare che le espressioni che è possibile impiegare per l'argomento  diagonal=... sono:

→ list(method="boxplot");

→ list(method="density", bw="nrd0", adjust=1, kernel="gaussian", na.rm=TRUE);

→ list(method="histogram");

→ list(method="oned");

→ list(method="qqplot").

Vista la relazione lineare molto buona che intercorre tra le variabili sembra ragionevole continuare ed effettuare il calcolo della regressione lineare multipla, copiate questo secondo script, incollatelo nella Console di R e premete ↵ Invio:

# REGRESSIONE LINEARE MULTIPLA - calcolo della regressione lineare multipla

#  

library(DAAG) # carica il pacchetto che include il set di dati ais 

mydata <- ais[,c(1,3,4)] # le colonne 1, 3, 4 contengono le variabili rcc, hc, hg

reglin <- lm(rcc~hg+hc, data=mydata) # calcolo della regressione 

#

coefficients(reglin) # intercetta e coefficienti angolari 

confint(reglin, level=0.95) # intervalli di confidenza al 95% 

#  

La funzione lm(...) viene impiegata per calcolare la regressione lineare multipla specificando con l'espressione rcc~hg+hc la variabile dipendente (rcc), le due variabili indipendenti (hg e hc) e l'equazione 

rcc = a + b₁·hg + b₂·hc

quindi possiamo ricavare con la funzione coefficients() i valori dell'intercetta (a), del coefficiente angolare di hg (b₁) e del coefficiente angolare di hc (b₂):

(Intercept)          hg          hc 

-0.24269563  0.03283747  0.10403395 

e riportare l'equazione della regressione lineare multipla come

rcc = -0.24269563 + 0.03283747·hg + 0.10403395·hc

Successivamente con la funzione confint() sono riportati gli intervalli di confidenza al 95% dell'intercetta e dei due coefficienti angolari:

                  2.5 %     97.5 %

(Intercept) -0.53163317 0.04624192

hg          -0.02459817 0.09027311

hc           0.08267073 0.12539718

Da notare che gli intervalli di confidenza forniscono un test di significatività per intercetta e coefficienti angolari:

→ dato un intervallo di confidenza al 95% che va da -0.53163317 a 0.04624192 (che quindi include il valore 0) l'intercetta non è significativamente diversa da 0, cosa molto logica in quanto al diminuire fino ad azzerarsi dei valori di emoglobina ed ematòcrito deve necessariamente diminuire e azzerarsi il numero degli eritrociti;

→ con un intervallo di confidenza al 95% che va da -0.02459817 a 0.09027311 (che quindi include il valore 0) il coefficiente angolare della concentrazione di emoglobina hg non è significativamente diversa da 0, cosa meno logica in quanto non solo dal punto di vista biologico deve esistere una proporzionalità tra numero degli eritrociti e concentrazione dell'emoglobina, ma anche nell'analisi esplorativa preliminare dei dati riportata con lo script precedente risulta evidente la proporzionalità tra numero degli eritrociti e concentrazione dell'emoglobina che è una dato di fatto dal punto di vista biologico;

→ con un intervallo di confidenza al 95% che va da 0.08267073 a 0.12539718 (che quindi non include il valore 0) il coefficiente angolare del valore ematòcrito hc risulta significativamente diverso da 0.

Come detto all'inizio, oltre ad effettuare una analisi esplorativa preliminare dei dati è opportuno sottoporre i risultati ad una valutazione critica a posteriori, che effettuiamo ora sotto forma di una diagnostica numerica della regressione lineare multipla, copiate il terzo script, incollatelo nella Console di R e premete ↵ Invio:

# REGRESSIONE LINEARE MULTIPLA - diagnostica numerica della regressione

#

library(gvlma) # carica il pacchetto per la diagnostica numerica 

library(car) # carica il pacchetto per la grafica 

library(DAAG) # carica il pacchetto che include il set di dati ais

mydata <- ais[,c(1,3,4)] # le colonne 1, 3, 4 contengono le variabili rcc, hc, hg

reglin <- lm(rcc~hg+hc, data=mydata) # calcolo della retta di regressione lineare multipla

#

summary(gvlma(reglin)) #  test globale per l'assunto di linearità 

outlierTest(reglin) # valore p di Bonferroni per la presenza di dati aberranti (outliers) 

#

Ecco la prima parte del riepilogo dei risultati fornito dalla funzione summary():

Call:

lm(formula = rcc ~ hg + hc, data = mydata)

Residuals:

     Min       1Q   Median       3Q      Max 

-0.55791 -0.11284 -0.00699  0.09950  0.53595 

Coefficients:

            Estimate Std. Error t value Pr(>|t|)    

(Intercept) -0.24270    0.14652  -1.656   0.0992 .  

hg           0.03284    0.02913   1.127   0.2609    

hc           0.10403    0.01083   9.603   <2e-16 ***

---

Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 0.1744 on 199 degrees of freedom

Multiple R-squared:  0.8565,    Adjusted R-squared:  0.855 

F-statistic: 593.8 on 2 and 199 DF,  p-value: < 2.2e-16

Il valore di probabilità Pr(>|t|) ottenuto con il test t di Student conferma (ovviamente) i risultati della significatività ricavati dagli intervalli di confidenza riportati sopra, per i quali valgono le considerazioni già fatte. E se la significatività statistica non è coerente con una realtà biologica assodata, dovrà ovviamente essere quest'ultima a prevalere: per cui assumiamo che l'equazione calcolata sia, fino a prova contraria, una descrizione sufficientemente adeguata della relazione di funzione che intercorre tra le tre variabili esaminate.

Che questo sia vero ce lo conferma la funzione gvlma(), acronimo di global validation of linear models assumptions, che fornisce questi risultati:

→ i dati sono distribuiti in modo lineare (Global Stat), assunto accettabile;

→ i dati non presentano asimmetria (Skewness), assunto accettabile;

→ i dati non presentano curtosi (Skewness), assunto accettabile;

→ nell'ambito dei valori osservati la varianza è costante (Heteroscedasticity), assunto accettabile.

ASSESSMENT OF THE LINEAR MODEL ASSUMPTIONS

USING THE GLOBAL TEST ON 4 DEGREES-OF-FREEDOM:

Level of Significance =  0.05 

Call:

 gvlma(x = reglin) 

                    Value p-value                   Decision

Global Stat        6.5801 0.15981    Assumptions acceptable.

Skewness           0.5149 0.47304    Assumptions acceptable.

Kurtosis           0.7559 0.38461    Assumptions acceptable.

Link Function      3.9844 0.04592 Assumptions NOT satisfied!

Heteroscedasticity 1.3249 0.24972    Assumptions acceptable.

L'assenza nei dati di asimmetria e di curtosi, e la loro omoschedasticità, sono condizioni sufficienti per garantire la loro adeguata descrizione mediante una regressione lineare. Il primo test (Global Stat) di valutazione globale degli assunti di linearità, che sostanzialmente comprende tutti e tre i test precedenti in forma compatta, li conferma. In questo contesto il test rimanente e l'unico discordante (Link Function) può essere ignorato.

Con la funzione outlierTest() viene  effettuato il test di Bonferroni, che non evidenzia dati anomali (outliers) ma segnala il dato 78 come quello che si discosta maggiormente dagli altri:

No Studentized residuals with Bonferroni p < 0.05

Largest |rstudent|:

   rstudent unadjusted p-value Bonferroni p

78 -3.28855          0.0011917      0.24073

In tema di valutazione critica a posteriori, aggiungiamo ora alla diagnostica numerica anche una diagnostica grafica della regressione lineare multipla, per effettuarla copiate questo quarto script, incollatelo nella Console di R e premete ↵ Invio:

# REGRESSIONE LINEARE MULTIPLA - diagnostica grafica della regressione

#

library(car) # carica il pacchetto per la grafica

library(DAAG) # carica il pacchetto che include il set di dati ais

mydata <- ais[,c(1,3,4)] # le colonne 1, 3, 4 contengono le variabili rcc, hc, hg

reglin <- lm(rcc~hg+hc, data=mydata) # calcolo della retta di regressione lineare multipla

#

windows() # apre e inizializza una nuova finestra grafica

qqPlot(reglin, main="Grafico dei quantili per i residui", xlab="t-quantili", ylab="Residui studentizzati") # grafico globale per linearità e dati aberranti

#

windows() # apre e inizializza una nuova finestra grafica

ceresPlots(reglin, ask=FALSE, main="Grafico di Ceres", ylab="Residui di Ceres (rcc)") # grafico di Ceres per la valutazione separata della linearità di hg e hc vs. rcc

#

windows() # apre e inizializza una nuova finestra grafica

plot(mydata$rcc, reglin$fitted.values, xlim=c(4,7), ylim=c(4,7)) # grafico valori originari vs valori stimati

abline(0, 1, lty=2) # linea di indentità y = x

#

La funzione qqPlot() riporta una distribuzione dei residui lineare, evidenziando nel contempo i dati 78 e 161 come quelli che dei quali sarebbe utile controllare la validità, ovvero che sono situati in intervalli di valori per i quali potrebbe essere opportuno acquisire più dati. 

La funzione ceresPlot() genera due grafici di Ceres separati per emoglobina ed ematòcrito, e consente di  visualizzare una relazione lineare tra ematòcrito ed eritrociti migliore di quella che intercorre tra emoglobina ed eritrociti, confermando in un certo senso i valori di significatività riportati sopra.

La diagnostica grafica della regressione lineare multipla si conclude con questo grafico minimalista, ma che entra nell'essenza dei risultati.

Il grafico, generato con la funzione plot(), riporta in ascisse i valori originari della concentrazione degli eritrociti (mydata$rcc) e in ordinate i corrispondenti valori reglin$fitted.values ricalcolati da hg e hc impiegando l'equazione della regressione lineare multipla ove

reglin$fitted.values = -0.24269563 + 0.03283747·hg + 0.10403395·hc

La linea di identità y = x che viene riportata evidenzia, senza la necessità di ulteriori statistiche, che l'informazione acquisita misurando la concentrazione dell'emoglobina e l'ematòcrito può essere compressa in un unico valore, che fornisce una stima mediamente adeguata della concentrazione degli eritrociti.

Nota: con "mediamente adeguata" si intende che, data la misura diretta della concentrazione degli eritrociti riportata in ascisse, la sua misura indiretta - ottenuta misurando emoglobina ed ematòcrito ed applicando la regressione lineare multipla - comporta una perdita di informazione [aumento dell'incertezza della misura] che si riflette nella "dispersione media" dei valori sull'asse delle ordinate. Se non ci fosse perdita di informazione i punti sarebbero perfettamente allineati sulla retta y = x.

Ora, visto che con tre variabili la cosa è ancora fattibile, non resta che aggiungere un grafico tridimensionale (3D).

Copiate quest'ultimo script, incollatelo nella Console di R e premete ↵ Invio:

# GRAFICO 3D ANIMATO con il pacchetto rgl

#  

library(rgl) # carica il pacchetto per la grafica 3D

library(DAAG) # carica il pacchetto che include il set di dati ais 

mydata <- ais[,c(1,3,4)] # le colonne 1, 3, 4 contengono le variabili rcc, hc, hg

windows() # apre e inizializza una nuova finestra grafica

#  

# grafico 3D animato

#  

plot3d(mydata$rcc, y=mydata$hg, z=mydata$hc, type="p", col="red", size=3) # realizza il grafico

play3d(spin3d(axis=c(0,0,1), rpm=10), duration=10) # anima il grafico

#  

Dopo avere caricato il pacchetto necessario con il comando library(rgl) il grafico viene realizzato con la funzione plot3d().

Quindi con la funzione play3d() viene realizzata una animazione della durata di 10 secondi (duration=10) con i parametri specificati dalla funzione spin3d() e cioè l'asse sul quale fare ruotare l'animazione, nel nostro caso  l'asse delle z (0,0,1), e il numero di rotazioni al minuto (rpm=10).

Quando l'animazione che compare nella finestra "RGL device" aperta nella Console di R dallo script si ferma, se “afferrate” il grafico 3D facendo click con il tasto sinistro del mouse e tenendolo premuto senza rilasciarlo potete ruotarlo a vostro piacimento muovendo il mouse. 

Anche se non aggiunge nulla di nuovo, il grafico [1] ci conforta documentando la linearità della relazione tra le tre variabili.

Conclusione: mettere in relazione tra di loro troppe variabili contemporaneamente può portare a risultati di difficile comprensione. Tuttavia una analisi esplorativa preliminare dei dati e una attenta valutazione critica a posteriori dei risultati possono aiutare a identificare i casi nei quali l'impiego della regressione lineare multipla può essere appropriato.

----------

[1] Trovate altre modalità di rappresentazione sotto forma di grafici tridimensionali nel post Grafici 3D.

Grafici 3D

(l'animazione del grafico 3D è stata realizzata con l'ultimo script riportato in questo post,

per riavviarla riaprire la pagina web dal browser)

Con R è possibile realizzare grafici tridimensionali (grafici 3D) impiegando una sola riga di codice.

Per vedere come sia possibile farlo è necessario:

→ scaricare i quattro pacchetti aggiuntivi scatterplot3d, plot3D, rgl e Rcmdr [1];

→ installare il pacchetto DAAG che contiene nella tabella ais i dati impiegati come esempio [2].

Le tre variabili x, y, z utilizzate nelle rappresentazioni grafiche che seguono sono rispettivamente la concentrazione nel sangue degli eritrociti (rcc), espressa in milioni per microlitro di sangue (10⁶/µL), la concentrazione nel sangue dell'emoglobina (hg), espressa in grammi per decilitro di sangue (g/dL) e il valore ematòcrito (hc), cioè la quota di sangue occupata dalla parte corpuscolata, espresso in percentuale (%) [3].

Ho riportato per i quattro pacchetti altrettanti script, limitati all'essenziale, che hanno il solo scopo di interessare a R per questo tipo di rappresentazione. Per i dettagli si rimanda ai manuali di riferimento dei pacchetti, che contengono la documentazione completa delle funzioni e dei rispettivi argomenti [4]. Ho poi aggiunto un quinto script che realizza un grafico 3D animato (quello riportato in apertura del post) per consentire di valutare anche questo aspetto di R.

Negli script le prime due righe di codice servono per caricare con library(DAAG) il pacchetto che contiene i dati (prima riga) e, di volta in volta, caricare il pacchetto di grafica (seconda riga) necessario per lo specifico grafico 3D. Questo viene infine realizzato con la funzione impiegata nell'ultima riga di codice (nel caso del grafico 3D animato nelle ultime due).

Poiché le variabili da rappresentare (rcc, hg, hc) sono contenute nella tabella ais (del pacchetto DAAG) il nome completo e univoco delle variabili è rappresentato dal prefisso ais più ($) il nome della variabile, pertanto in tre script trovate le variabili indicate come ais$rcc, ais$hg, ais$c.

In altri due script è stata invece impiegata la funzione attach() che, una volta specificato come argomento il nome della tabella che contiene le variabili con attach(ais), consente di impiegare i nomi originari delle variabili (quindi i nomi senza il prefisso ais$), che pertanto trovate riportati come rcc, hg, hc. 

Copiate e incollate nella Console di R il primo script e premete ↵ Invio:

# GRAFICO 3D SEMPLICE impiega il pacchetto scatterplot3d

#  

library(DAAG) # carica il pacchetto che include il set di dati/tabella ais 

library(scatterplot3d) # carica il pacchetto per la rappresentazione grafica 3D

#  

# grafico 3D semplice 

#  

scatterplot3d(ais$rcc, y=ais$hg, z=ais$hc, color="red", type="p", main="Grafico 3d semplice", xlab="Eritrociti", ylab="Emoglobina", zlab="Ematòcrito")  

#  

Dopo avere caricato i due pacchetti necessari, il primo contenente i dati e il secondo contenente le funzioni grafiche, il grafico viene realizzato con la funzione scatterplot3d() specificando come argoment:

→ le variabili della tabella ais da rappresentare (ais$rcc, ais$hg, ais$hc);

→ il tipo (type="p") di rappresentazione dei dati ("p" = punti, ma sono possibili le alternative "h" = punti con linee verticali che congiungono ciascun punto al piano x-y, e "l" = linee che congiungono i punti tra loro);

→ il colore dei punti (color="red");

→ il titolo del grafico (main="Grafico 3d semplice") e le etichette degli assi. 

Questo è il grafico risultante:

Prima di eseguire ognuno dei successivi script uscite da R e poi rientrate (procedura raccomandata), oppure in alternativa chiudete la finestra contenente il grafico, quindi digitate nella Console di R

rm(list=ls(all=TRUE))

per effettuare la pulizia dell'area di lavoro [5].

Ora copiate e incollate nella Console di R questo secondo script e premete ↵ Invio:

# GRAFICO 3D CON VISTA DA ANGOLI VARIABILI impiega il pacchetto plot3D

#  

library(DAAG) # carica il pacchetto che include il set di dati/tabella ais 

library(plot3D) # carica il pacchetto per la rappresentazione grafica 3D

#  

# grafico 3D con vista da angoli variabili

#  

scatter3D(ais$rcc, y=ais$hg, z=ais$hc, theta=20, phi=0, ticktype="detailed", xlab="Eritrociti", ylab="Emoglobina", zlab="Ematòcrito")  

# 

Dopo avere caricato i due pacchetti necessari il grafico viene realizzo con la funzione scatter3D() specificando come argoment:

→ le variabili della tabella ais da rappresentare (ais$rcc, ais$hg, ais$hc);

→ l'angolo di rotazione del sistema di coordinate sull'asse verticale (theta=20);

→ l'angolo di rotazione del sistema di coordinate sull'asse orizzontale (phi=0);

→ la presenza sugli assi delle scale numeriche (ticktype="detailed");

→ le etichette degli assi.

I punti sono riportati con una scala di colori che indica la loro distanza dal piano x-y, misurata lungo l'asse z. La scala colore/distanza è riportata nella legenda che compare sulla destra del grafico.  Questo è il grafico:

Copiate e incollate nella Console di R il terzo script e premete ↵ Invio:

# GRAFICO 3D INTERATTIVO CHE PUO' ESSERE RUOTATO CON IL MOUSE impiega il pacchetto rgl

#  

library(DAAG) # carica il pacchetto che include il set di dati/tabella ais 

library(rgl) # carica il pacchetto per la rappresentazione grafica 3D

#

attach(ais) # consente di impiegare direttamente i nomi delle variabili

#  

# grafico 3D che può essere ruotato con il mouse

#  

plot3d(rcc, y=hg, z=hc, type="p", col="red", size=3) 

#  

detach(ais) # termina l'impiego diretto dei nomi delle variabili 

#

Dopo avere caricato i due pacchetti necessari, la funzione attach(ais) consente questa volta di realizzare il grafico nella riga successiva impiegando direttamente i nomi delle variabili contenute nella tabella ais senza specificare per ciascuna di esse la tabella nella quale sono incluse (notare quindi la differenza con i due script precedenti).

Il grafico è realizzato con la funzione plot3d() specificando come argomenti:

→ le variabili della tabella ais da rappresentare (rcc, hg, hc);

→ type="p" che consente di rappresentare i dati come punti;

→ col="red" che specifica il colore dei punti;

→ size=3 che specifica la dimensione dei punti [6].

La funzione detach() tecnicamente, come riporta la documentazione "Detach a database, i.e., remove it from the search() path of available R objects. Usually this is either a data.frame which has been attached or a package which was attached by library", praticamente ripristina le condizioni precedenti al comando attach().

Notare anche che, in questo e nei due script successivi, non essendo stati specificati gli argomenti xlab, ylab e zlab, R riporta automaticamente come etichette degli assi i nomi delle variabili.

 

Questa volta il grafico 3D è interattivo. Infatti se “afferrate” il grafico risultante facendo click su di esso e tenendo premuto il tasto sinistro del mouse senza rilasciarlo, potete ruotare il grafico a vostro piacimento muovendo il mouse. In questo modo sono state realizzate le due immagini che, qui mostrate affiancate, forniscono due delle possibili viste [della distribuzione] dei dati. 

Copiate e incollate nella Console di R questo quarto script e premete ↵ Invio:

# GRAFICO 3D INTERATTIVO CHE PUO' ESSERE RUOTATO CON IL MOUSE impiega il pacchetto Rcmdr

#  

library(DAAG) # carica il pacchetto che include il set di dati/tabella ais 

library(Rcmdr) # carica il pacchetto per la rappresentazione grafica 3D

#

attach(ais) # consente di impiegare direttamente i nomi delle variabili

#  

# grafico 3D che può essere ruotato con il mouse

#  

scatter3d(rcc, y=hg, z=hc) 

#  

detach(ais) # termina l'impiego diretto dei nomi delle variabili

#

Dopo avere caricato i due pacchetti necessari, la funzione attach(ais) consente anche questa volta di realizzare il grafico impiegando direttamente i nomi delle variabili senza specificare per ciascuna di esse la tabella (ais) nella quale sono incluse (notare nuovamente la differenza con i primi due script).

Il grafico viene realizzato con la funzione scatter3d() specificando come argomento solamente le variabili della tabella ais da rappresentare (rcc, hg, hc), mentre per tutti gli altri argomenti, colori inclusi, sono stati lasciati i valori di default.

Come nel caso del grafico precedente si tratta di un grafico interattivo, e se lo “afferrate” [7] facendo click su di esso con il tasto sinistro del mouse e tenendolo premuto senza rilasciarlo, potete ruotarlo a vostro piacimento spostando il mouse. In questo modo sono state realizzate le immagini che, qui mostrate affiancate, forniscono due delle viste [della distribuzione] dei dati che è possibile ottenere. 

Concludiamo con un grafico 3D animato, molto semplice, ma che apre in R un mondo, al quale si rimandano gli interessati [8].

 

Copiate e incollate nella Console di R il quinto e ultimo script e premete ↵ Invio:

# GRAFICO 3D ANIMATO impiega il pacchetto rgl

#  

library(DAAG) # carica il pacchetto che include il set di dati/tabella ais 

library(rgl) # carica il pacchetto per la rappresentazione grafica 3D

#  

# grafico 3D animato

#  

plot3d(ais$rcc, y=ais$hg, z=ais$hc, type="p", col="red", size=3) 

play3d(spin3d(axis=c(0,0,1), rpm=10), duration=30) 

#  

Dopo avere caricato i due pacchetti necessari, il primo contenente i dati e il secondo contenente le funzioni grafiche, il grafico 3D viene realizzato con la funzione plot3d() del pacchetto rgl già impiegata nel terzo degli script riportati sopra (vedi).

Nella riga successiva con la funzione play3d() viene realizzata una animazione del grafico:

→ con i parametri specificati dalla funzione spin3d() che sono l'asse sul quale fare ruotare l'animazione, nel nostro caso  l'asse delle z (0,0,1), e il numero di rotazioni al minuto (rpm=10). 

→ delladurata di 30 secondi (duration=30)

Il grafico 3D animato risultante, che è quello riportato anche in apertura del post, è interessante non tanto per l'aspetto coreografico, ma piuttosto in quanto offre una serie di visuali che consentono di apprezzare appieno la relazione di linearità che intercorre tra le tre variabili rappresentate.

Nonostante gli script riportati sopra siano stati limitati all'essenziale [3], potete immediatamente riutilizzarli per una rappresentazione non banale dei vostri dati semplicemente sostituendo ai nomi delle variabili qui impiegati i nomi delle vostre variabili x, y e z.

Per una guida rapida all'importazione dei dati in R potete consultare i link:

→ importazione dei dati da un file .csv

→ importazione dei dati da un file .xls o .xlsx

→ gestione dei dati mancanti

Infine ricordo che nel post Salvare i grafici di R in un file è riportato uno script che consente di trasformare i grafici in immagini e salvarli sotto forma di file .bmp, .jpeg, .png, .pdf, .ps per poterli stampare, archiviare, inserire in una pubblicazione, in un post o in un sito web. 

-----------

[1] Pacchetti aggiuntivi in quanto vanno ad aggiungersi al pacchetto base di R. Dopo avere scaricato il pacchetto Rcmdr chiudere e riavviare R, quindi digitare nella Console di R library(Rcmdr) e premere ↵ Invio per scaricare dal CRAN gli ulteriori pacchetti aggiuntivi necessari al suo funzionamento, infine nuovamente chiudere e riavviare R.

[2] Se non volete installare il pacchetto potete in alternativa fare il download di un file contenente gli stessi dati, come riportato nel post Il set di dati ais e alla pagina Dati.

[3] La parte corpuscolata del sangue include eritrociti, leucociti e piastrine, ma dato che gli ultimi due in condizioni normali occupano un volume irrilevante, la quota di sangue occupata dalla parte corpuscolata (il valore ematòcrito) è di fatto quella occupata dagli eritrociti (globuli rossi).

[4] Per la documentazione completa delle funzioni e dei rispettivi argomenti vedere i manuali di riferimento dei pacchetti scatterplot3d, plot3D, rgl e Rcmdr sul repository della documentazione: Available CRAN Packages By Name. URL consultato il 10/11/2021: https://bit.ly/3zLwHWH

[5] Vedere il post Pulizia dell'area di lavoro e uscita dal programma

[6] Altre cose interessanti le trovate su Impressive package for 3D and 4D graph - R software and data visualization. URL consultato il 10/11/2021: https://bit.ly/3mKEeC1

[7] Attenzione perché il grafico potrebbe essere nascosto in parte o completamente dalla finestra di R Commander che si apre, al bisogno minimizzatela.

[8] Vedere ad esempio Animated 3d chart with R in: The R Graph Gallery. URL consultato il 10/11/2021: https://bit.ly/3BXySI4