Quante componenti estrarre? La parallel analysis

L’analisi in componenti principali (ACP) è una celebre tecnica di analisi multivariata che, dato un gruppo di variabili quantitative misurate sulle medesime unità statistiche, consente di identificare un insieme ridotto di nuove variabili in grado di rappresentare i dati in maniera più sintetica. Date n variabili osservate, fra loro più o meno correlate, la ACP le trasforma in altrettante variabili fra loro linearmente indipendenti. Queste nuove variabili sono dette “componenti” e vengono estratte in ordine d’importanza: da quella che spiega la quota maggiore di varianza dei dati osservati a quella che ne spiega la quota minore.

Lo scopo della ACP è ridurre i dati trasformandoli, cercando di veicolare la medesima quantità d’informazione utilizzando un numero inferiore di variabili. Ne abbiamo già parlato in un precedente post, al quale si rimanda per ulteriori approfondimenti.

Quello che vorrei affrontare in questo articolo è un vecchissimo problema che ancora sembra affliggere molte ricerche, ovvero la determinazione del numero di componenti da estrarre. Una volta individuate le componenti, bisogna scegliere quante tenerne: se le tenessimo tutte non avremmo certo ridotto i dati, ma se ne tenessimo troppo poche perderemmo preziose informazioni. Si tratta di un problema non da poco, anche perché spesso è attraverso la ACP che viene definita la dimensionalità di strumenti di misura come i test psicologici utilizzati nella prassi psicodiagnostica.

Prima di addentrarci nel problema, vediamo di procurarci dei dati da utilizzare per svolgere qualche analisi. Il pacchetto Flury contiene tutti i dataset utilizzati in un volume del 1997 di B. Flury dedicato alla statistica multivariata; la libreria va installata perché non è distribuita con la versione standard di R. Noi utilizzeremo il dataset microtus, che contiene la registrazione di otto variabili biometriche su 288 roditori appartenenti appunto al genere “Microtus”.

install.pakcages("Flury") # Installa la libreria Flury
library(Flury)            # Carica la libreria Flury
data(microtus)            # Attiva il dataset

head(microtus)            # Visualizza l'anteprima del dataset
      Group M1Left M2Left M3Left Foramen Pbone Length Height Rostrum
1 multiplex   2078   1649   1708    3868  5463   2355    805     475
2 multiplex   1929   1551   1550    3825  4741   2305    760     450
3 multiplex   1888   1613   1674    4440  4807   2388    775     460
4 multiplex   2020   1670   1829    3800  4974   2370    766     460
5 multiplex   2223   1814   1933    4222  5460   2470    815     475
6 multiplex   2190   1800   2066    4662  4860   2535    838     521

La variabile Group (prima colonna) classifica i roditori in tre gruppi: “multiplex”, “subterraneus” e “unknown”. Noi ci interesseremo del gruppo “subterraneus”, per cui isoliamo i dati di questi esemplari in un nuovo dataset:


subterraneus <- microtus[microtus$Group=="subterraneus", -1] [/code]

R offre diversi comandi per realizzare una ACP. Qui utilizzeremo la funzione principal del pacchetto psych, anche questo da installare manualmente. Utilizzare la funzione principal è molto semplice: oltre alla matrice di dati osservati, dobbiamo specificare anche il numero di componenti da estrarre. Dato che ancora non abbiamo deciso quante componenti utilizzare (non abbiamo nessun aggancio teorico), ne estraiamo tante quante sono le variabili osservate.


install.packages(“psych”) # Installa la libreria psych
library(psych) # Carica la libreria psych
k <- ncol(subterraneus) # Numero di componenti (= num. variabili) PCA <- principal(subterraneus, nfactors=k, rotate="none") [/code]

Il problema che sorge a questo punto è capire quante sono effettivamente le componenti che andrebbero estratte. Come già visto, una cosa che possiamo fare è visionare gli autovalori:

> PCA$values
[1] 4.02 1.21 0.85 0.68 0.46 0.38 0.27 0.12

Gli autovalori rappresentano la varianza delle singole componenti, e il rapporto tra ogni autovalore e la somma degli autovalori di tutte le componenti fornisce la varianza spiegata da ogni componente. Del 100% della variabilità originaria, quanto è spiegato da ogni componente? Lo si può capire così:

> PCA$values / sum(PCA$values)
[1] 0.50 0.15 0.11 0.08 0.06 0.05 0.03 0.02

Una regola molto utilizzata per stabilire il numero di componenti è quella di estrarne tante quanti sono gli autovalori maggiori di 1. Nel nostro caso dovremmo tenere solo le prime due componenti, che complessivamente spiegano il 65% della varianza (0.50 + 0.15). Siamo però davvero sicuri che tenere le prime due componenti sia la soluzione migliore? L’autovalore della seconda componente vale 1.21, un valore molto vicino alla fatidica soglia 1: non è che magari vi è così vicino da dover essere escluso?

Ragionando in ottica di analisi in componenti principali potremmo anche tenere la seconda componente: avrà pure un autovalore troppo vicino a 1, ma spiega pur sempre un 15% di varianza. Molto spesso, però, l’analisi degli autovalori calcolati dalla ACP è utilizzata per stabilire il numero di componenti in un’analisi fattoriale, la quale oltre che alla percentuale di varianza spiegata deve guardare soprattutto all’economia della soluzione trovata, alla generalizzabilità del risultato e alla congruenza con assunzioni teoriche più o meno fondate.

C’è un tipo di analisi molto raffinata che aiuta nella scelta del numero di componenti da estrarre: la parallel analysis. Non si tratta di una tecnica innovativa, anzi; essa risale agli anni ’60 ed è basata su un’osservazione: talvolta può capitare che alcuni campioni di dati producano degli autovalori maggiori di 1 per puro effetto del caso. La differenza tra un autovalore di poco superiore a 1 e un altro di poco inferiore può essere molto sottile, così sottile che può essere stato il caso a determinarla.

La parallel analysis consiste nel confrontare gli autovalori calcolati sui dati osservati non solo con il valore soglia 1, ma anche con degli autovalori calcolati su matrici di dati casuali. L’idea è che una componente, perché sia considerata significativa, non solo debba presentare un autovalore maggiore di 1, ma il suo autovalore deve anche essere superiore al corrispettivo autovalore prodotto da dei dati casuali. In sostanza, la componente deve dimostrare di essere più forte del caso.

La tecnica prevede di realizzare, parallelamente all’analisi in componenti principali, delle altre analisi su matrici di dati generati casualmente. Definito un numero n di iterazioni, si producono n matrici di dati casuali dalle caratteristiche simili a quelle dei dati osservati: stesso numero di unità campionarie e stesso numero di variabili. Su ogni matrice si esegue un’analisi in componenti principali e si calcolano gli autovalori; come risultato si avrà, per ogni componente, una serie di autovalori di numerosità n.

Per esempio, nel dataset subterraneus abbiamo 8 variabili misurate su 46 unità. Possiamo definire un n = 1000 e quindi costruire 1000 matrici di dati casuali con 46 righe e 8 colonne. Su ogni matrice si esegue un’analisi in componenti principali e si calcolano gli autovalori, ottenendo così 1000 autovalori per ognuna delle 8 componenti possibili.

Al termine dell’operazione, per ognuna delle 8 componenti, sui 1000 autovalori si calcola un indice riassuntivo che può essere la media o la mediana. Avremo come risultato 8 medie (o mediane), che rappresentano gli autovalori generati dai dati casuali; questi indici rispondono alla domanda: «Quali autovalori avremmo ottenuto se i dati osservati fossero stati determinati completamente dal caso?».

A questo punto si confrontano autovalori “reali” e autovalori “casuali”: il numero di componenti da estrarre corrisponde al numero di autovalori reali che risultano più elevati degli autovalori casuali.

Realizzare una parallel analysis in R è abbastanza semplice. La libreria psych contiene l’ottima funzione fa.parallel, ma per diversi motivi (che magari approfondiremo in un altro post) io preferisco utilizzare uno script realizzato da me artigianalmente. Nel mio account GitHub ho caricato uno script R che consente di realizzare in maniera completamente automatica la procedura descritta sopra. Per scaricare lo script è sufficiente andare in basso a destra nella pagina e cliccare su “Donwload ZIP” oppure cliccare direttamente su questo link. Lo script su chiama parallelPCA.R e va estratto dall’archivio zip.

Vediamo ora di applicare la parallel analysis sul dataset microtus. Lo script deve preventivamente caricato in R; il caricamento può essere eseguito direttamente da menù: File > Sorgente codice R… per renderlo disponibile nell’area di lavoro. Quindi, utilizzando il comando parallelPCA si esegue la parallel analysis.


parpca <- parallelPCA(subterraneus) [/code]

L’output è prevalentemente di tipo grafico, e mostra uno scree-graph al quale sono sovrapposti, in grigio, i quantili delle distribuzioni degli autovalori dei dati casuali: i puntini rappresentano le mediane mentre le barre rappresentano il 2.5° e 97.5° percentile.

Gli stessi dati sono riportati nell’output testuale: autovalori reali, medie, intervalli di confidenza e quantili degli autovalori calcolati sui dati casuali.

> parpca
Parallel Analysis by Random Data 
Correlation matrix: pearson 

Observed values:
        c1   c2   c3   c4   c5   c6   c7   c8
eigen 4.02 1.21 0.85 0.68 0.46 0.38 0.27 0.12

Parallel Confidence Intervals:
         c1   c2   c3   c4   c5   c6   c7   c8
CI Sup 1.93 1.61 1.38 1.19 1.02 0.85 0.69 0.52
Mean   1.68 1.40 1.20 1.03 0.89 0.74 0.60 0.45
CI Inf 1.43 1.19 1.02 0.88 0.76 0.63 0.51 0.39

Parallel Quantiles:
        c1   c2   c3   c4   c5   c6   c7   c8
97.5% 2.00 1.59 1.34 1.16 1.01 0.87 0.73 0.59
50%   1.67 1.40 1.20 1.03 0.89 0.74 0.60 0.45
2.5%  1.45 1.23 1.06 0.91 0.76 0.61 0.47 0.32

A questo punto abbiamo un’arma in più per decidere il numero di componenti da estrarre. Dei due autovalori maggiori di 1, solo uno si dimostra effettivamente “più forte” del caso: già il secondo autovalore prodotto dai dati reali (1.21) è inferiore alla mediana degli autovalori prodotti dai dati casuali (1.40), e gli altri di conseguenza a seguire.

Quindi, quante componenti estrarre? Premesso che certe decisioni dovrebbe essere la teoria a guidarle, la parallel analysis ci indica come potenzialmente rischiosa l’assunzione di due componenti. Se sotto analisi ci fosse un questionario psicometrico e questa analisi fosse il preludio di un’analisi fattoriale, con le componenti da estrarre che rappresentano dei costrutti latenti, io sceglierei un unico fattore.

Bibliografia essenziale

Horn, J. L. (1965). A rationale and test for the number of factors in factor analysis. Psychometrika, 30(2), 179-185.

Zwick, W. R., & Velicer, W. F. (1986). Comparison of five rules for determining the number of components to retain. Psychological Bulletin, 99, 432-442.

Franklin, S. B., Gibson, D. J., Robertson, P. A., Pohlmann, J. T., & Fralish, J. S. (1995). Parallel analysis: a method for determining significant principal components. Journal of Vegetation Science, 6(1), 99-106.

Crawford, A. V., Green, S. B., Levy, R., Lo, W. J., Scott, L., Svetina, D., & Thompson, M. S. (2010). Evaluation of parallel analysis methods for determining the number of factors. Educational and Psychological Measurement, 70(6), 885–901.

Print Friendly

Lascia un Commento