Se trovate utile l'impiego delle tecniche di clustering (o analisi dei gruppi) [1] come metodo per l'analisi esplorativa dei dati – ma anche se le affrontate per la prima volta – potrebbe interessarvi questo esempio.
Copiate lo script, incollatelo nella Console di R e premete ↵ Invio:
# ANALIZZARE LE DISTRIBUZIONI DEI DATI NEI CLUSTER con ggpairs
#
library(cluster) # carica il pacchetto per il clustering
library(factoextra) # carica un pacchetto addizionale per il clustering
library(GGally) # carica il pacchetto che estende le funzioni di ggplot2
#
# preparazione dei dati
#
mydata <- read.table("c:/Rdati/ema.csv", header=TRUE, sep=";", dec=".") # importa i dati dei 643 casi
newdata <- subset(mydata, select=c("gr", "hb", "mcv", "hb_a2", "ferro")) # seleziona le variabili per il clustering
#
# clustering non gerarchico con il metodo di Rousseew (k-medoids)
#
myclust <- pam(newdata, 3, metric=c("euclidean"), stand=TRUE) # genera un oggetto 'pam' che contiene i risultati del clustering
fviz_cluster(myclust, myclust$cluster, labelsize=7, main = "Grafico dei cluster - metodo di Rousseew (k-medoids)") # grafico dei cluster per le prime due componenti principali
#
Vediamo di cosa si tratta. Abbiamo nel file ema.csv una raccolta di 643 casi dei quali abbiamo raccolto i risultati di tre analisi del sangue richieste a scopo diagnostico:
→ esame emocromocitometrico – del quale ci interessano in particolare i valori di: (i) gr la concentrazione degli eritrociti nel sangue (globuli rossi) espressa in migliaia di miliardi per litro di sangue 10¹²/L o in milioni per microlitro di sangue 10⁶/µL (le due espressioni sono numericamente equivalenti); (ii) hb la concentrazione dell'emoglobina espressa in grammi per decilitro di sangue g/dL; (iii) mcv il volume medio dei globuli rossi espresso in milionesimi di miliardesimo di litro (10-¹⁵) o femtolitri (fL);
→ emoglobina A2 – hb_A2 la variante fisiologica dell'emoglobina presente nell'adulto espressa come percentuale dell'emoglobina totale;
→ ferro – ferro la concentrazione del ferro espressa in microgrammi per decilitro di siero µg/dL.
Sappiamo che i 643 casi comprendono:
→ soggetti sani, impiegati come controlli, con i risultati delle tre analisi all'interno dei valori fisiologici (valori "normali");
→ soggetti con alterazioni dell'esame emocromocitometrico e con carenza di ferro (confermata dalla diminuita concentrazione del ferro nel siero);
→ soggetti con alterazioni dell'esame emocromocitometrico e portatori del gene della beta-talassemia (confermato dall'aumento della percentuale di emoglobina A2);
→ soggetti con alterazioni dell'esame emocromocitometrico ma non classificabili ulteriormente con i dati disponibili.
Siamo interessati a una analisi esplorativa dei dati e, seguendo lo script, procediamo come segue.
1) Preparazione dei dati.
Installiamo i tre pacchetti che ci servono, cluster, factoextra e GGally, quindi ci procuriamo i dati impiegati nell'esempio (per il file di dati trovate link e istruzioni alla pagina di Download):
→ effettuando il download del file di dati ema.csv
→ salvando il file nella cartella C:\Rdati\
Con la funzione read.table() importiamo i dati nell'oggetto mydata. Quindi con la funzione subset() riportiamo nell'oggetto newdata le sole variabili che ci interessano, cioè i risultati delle analisi indicate sopra (643 record), in cinque colonne:
> newdata
gr hb mcv hb_a2 ferro
1 4.90 13.3 82.8 1.8 106
2 4.66 10.8 73.6 2.6 148
3 5.43 11.5 66.5 4.8 104
4 5.41 10.8 73.4 2.5 74
..........
642 5.28 11.0 65.4 1.9 11
643 4.03 7.2 55.0 1.9 21
2) Clustering non gerarchico (analisi dei gruppi).
A questo punto eseguiamo il clustering non gerarchico, lo facciamo impiegando il metodo di Rousseew. La prima funzione utilizzata è pam() che deriva il suo nome da "Partitioning Around Medoids", impiega l'algoritmo di clustering dei dati k-medoids (un versione robusta dell'algoritmo k-means) e prevede come argomenti:
→ i dati newdata da impiegare;
→ il numero di cluster, in questo caso 3, in cui classificare i dati;
→ la metrica da impiegare che è "euclidean" (in alternativa può essere "manhattan");
→ stand=TRUE che impone la standardizzazione dei dati, che viene eseguita automaticamente dalla funzione.
Con la successiva funzione fviz_cluster() è generato il grafico dei cluster (myclust$cluster) ottenuti con l'algoritmo k-medoids e contenuti nell'oggetto myclust.
Date le premesse riportate all'inizio abbiamo impostato 3 come numero dei gruppi (cluster) in cui classificare i dati allo scopo di separare (i) sani, (ii) carenza di ferro, (iii) beta-talassemia eterozigote. Consapevoli del fatto che data la casistica, la quantità e il tipo delle informazioni di cui disponiamo, le altre possibili e meno frequenti cause di alterazione dell'esame emocromocitometrico (e ve ne sono molte) rimarranno inestricabilmente sovrapposte tra loro e/o con qualcuno di questi tre gruppi (a causa dei vincoli imposti dalla informazione disponibile che è, come sempre, limitata).
3) Analisi delle distribuzioni dei dati nei cluster.
Eccoci arrivati al punto. Non ci basta avere, per ciascun cluster, l'indicazione dei soggetti (numerati in ordine progressivo in mydata) appartenenti a ciascun cluster. Vogliamo avere anche l'indicazione della corrispondenza tra i cluster e le distribuzioni dei risultati delle grandezze misurate. Lo facciamo in una nuova finestra grafica aperta con windows() associando a ciascun caso (record) della tabella newdata con la funzione cbind() il cluster di appartenenza myclust$clustering estratto dall'oggetto 'pam' myclust che contiene i risultati del clustering e che viene riportato aggiungendo una sesta colonna alla tabella.
#
# analisi delle distribuzioni dei dati nei cluster
#
windows() # apre e inizializza una nuova finestra grafica
newclust <- cbind(newdata, myclust$clustering) # riporta per ciascun caso il cluster di appartenenza
ggpairs(newclust, columns=1:5, aes(color=as.factor(myclust$clustering), alpha=0.5), title="Scatterplot e kernel density plot dei dati dei cluster") # scatterplot e kernel density plot dei dati nei cluster
#
Ora abbiamo la tabella che possiamo impiegare per generare il grafico della distribuzione dei dati dei cluster con la funzione ggpairs() che impiega gli argomenti:
→ newclust per indicare la tabella che contiene i dati;
→ 1:5 per indicare che i dati da rappresentare sono contenuti nelle colonne da 1 a 5 (la sesta contiene il cluster di appartenenza);
→ as.factor() in quanto le variabili numeriche 1, 2 e 3 della colonna 6 che indicano il cluster di appartenenza devono essere trasformate in variabili fattore;
→ alpha=0.5 per rappresentare i colori con una trasparenza del 50%.
Da notare che al bisogno potete anche esportare la tabella newclust, che nelle prime cinque colonne contiene i dati e nella sesta il cluster di appartenenza come qui indicato
> newclust
gr hb mcv hb_a2 ferro myclust$clustering
1 4.90 13.3 82.8 1.8 106 1
2 4.66 10.8 73.6 2.6 148 1
3 5.43 11.5 66.5 4.8 104 2
4 5.41 10.8 73.4 2.5 74 3
..........
642 5.28 11.0 65.4 1.9 11 3
643 4.03 7.2 55.0 1.9 21 3
e salvarla sotto forma di file [2].
Ed ecco il grafico risultante.
La conclusione?
1) l'analisi dei gruppi fornisce una sintesi nella quale le grandezze soggiacenti non sono più riconoscibili singolarmente in quanto le componenti principali sono una combinazione lineare di tutte le variabili originali (nel primo dei grafici la prima componente principale in ascisse è la variabile formata dalla combinazione lineare di tutte le variabili originali che spiega la maggior quantità di varianza, la seconda componente principale in ordinate è la variabile che spiega la maggior quantità di varianza in ciò che rimane una volta rimosso l'effetto della prima componente principale);
2) la distribuzione dei dati all'interno dei cluster fornisce una rappresentazione di dettaglio che consente di ricollegare questi ultimi con la distribuzione originaria delle variabili la cui combinazione ha dato loro origine, e consente di riconoscere come:
→ il cluster 1 (rosso) contiene i soggetti con emoglobina A2 e ferro normali;
→ il cluster 2 (verde) contiene i soggetti con emoglobina A2 aumentata;
→ il cluster 3 (azzurro) contiene i soggetti con ferro diminuito.
Pur con distribuzioni che in parte si sovrappongono – un fatto determinato dalla incompletezza della informazione disponibile rispetto alla complessità dei casi – l'esempio illustra come:
→ il clustering fornisce una vista sintetica dell'informazione globale contenuta nei dati;
→ possiamo aggiungere all'analisi dei cluster una vista di dettaglio delle distribuzioni dei singoli dati all'interno dei cluster;
→ le due viste degli stessi dati si integrano sinergicamente l'una con l'altra.
----------
[1] Fate click su clustering nelle Parole chiave.
