1 Analisi delle componenti principali

Vengono caricati i packages necessari per realizzare questo documento

1.1 Primo esempio con tre sole variabili

source

view(dfSummary(dati),method = "render")

Data Frame Summary

dati

Dimensions: 1427 x 3
Duplicates: 103
No Variable Stats / Values Freqs (% of Valid) Graph Valid Missing
1 Altezza [numeric] Mean (sd) : 151.9 (10.1) min < med < max: 127 < 151 < 183 IQR (CV) : 15 (0.1) 54 distinct values 1427 (100%) 0 (0%)
2 Span [numeric] Mean (sd) : 153.6 (11.2) min < med < max: 123 < 153 < 184 IQR (CV) : 16 (0.1) 60 distinct values 1427 (100%) 0 (0%)
3 Torace [numeric] Mean (sd) : 75.6 (7.8) min < med < max: 57 < 74 < 104 IQR (CV) : 10 (0.1) 44 distinct values 1427 (100%) 0 (0%)

Generated by summarytools 0.9.6 (R version 4.0.2)
2020-09-12

1.2 Variabili osservate.

matrice di grafici per le sole variabili Altezza, Span, Torace rilevate su 1427 righe

source

round(      cor(dati),3   )

output

        Altezza  Span Torace
Altezza   1.000 0.892  0.587
Span      0.892 1.000  0.578
Torace    0.587 0.578  1.000

source

MLA.explor.pairs(dati)

plot

1.3 Principal Component Analysis

Primo esempio con tre sole variabili:

  • Come misurare la correlazione globalmente presente fra le 3 variabili?

  • Possiamo trovare una nuova variabile combinazione lineare delle variabili originali che in qualche modo le riassuma meglio

  • cerchiamo la retta che minimizza la somma dei quadrati degli scarti dei punti dalla retta stessa (misurata ortogonalmente) (retta rossa nella figura)

  • Vedremo dopo che questa nuova variabile corrispondente a questa retta (ossia le nuove coordinate ottenute dalle proiezioni dei punti) è quella che ha la massima varianza

source

rid =scale(dati[,c(1:3)])
n   =nrow(rid)
round(cor(rid),3)

output

        Altezza  Span Torace
Altezza   1.000 0.892  0.587
Span      0.892 1.000  0.578
Torace    0.587 0.578  1.000

source

     autov=eigen(cor(rid))
     round(autov$values,3)

output

[1] 2.383 0.509 0.108

source

     sum(autov$values)

output

[1] 3

source

     v=autov$vectors
     round(v,3)

output

       [,1]   [,2]   [,3]
[1,] -0.609 -0.352  0.711
[2,] -0.606 -0.372 -0.703
[3,] -0.512  0.859 -0.013

source

open3d()

output

glX 
  1

source

MLA.pca3d(rid,segment=FALSE, princ.axes=TRUE,ellipse = FALSE)

output

Importance of components:
                          Comp.1    Comp.2     Comp.3
Standard deviation     1.5430776 0.7134693 0.32828463
Proportion of Variance 0.7942527 0.1697985 0.03594879
Cumulative Proportion  0.7942527 0.9640512 1.00000000

You must enable Javascript to view this page properly.

1.4 Retta di regressione principale

Voglio trovare la retta che meglio passa attraverso i dati.

Sottolineo la differenza con la regressione semplice e multipla: adesso voglio minimizzare le distanze dai punti ortogonalmente rispetto a questa retta, che chiameremo Retta di regressione principale

1.5 Rotazione assi

Dopo aver fatto questa operazione una prima volta, possiamo ripeterla cercando un’altra retta che minimizzi la sommadei quadrati delle distanze dei punti da questa retta, ma limitando la ricerca alle retteortogonali rispetto alla prima retta di regressione principale

L’operazione si può ripetere fino alla \(k-\)esima retta…

Analisi delle componenti principali: con tre variabili non vi è alcuna perdita di informazione ovviamente; la trasformazione corrisponde ad una rotazione degli assi

2 Esempio numerico:

source

     pca1=princomp(rid)
     print(str(pca1))

output

List of 7
 $ sdev    : Named num [1:3] 1.543 0.713 0.328
  ..- attr(*, "names")= chr [1:3] "Comp.1" "Comp.2" "Comp.3"
 $ loadings: 'loadings' num [1:3, 1:3] 0.609 0.606 0.512 0.352 0.372 ...
  ..- attr(*, "dimnames")=List of 2
  .. ..$ : chr [1:3] "Altezza" "Span" "Torace"
  .. ..$ : chr [1:3] "Comp.1" "Comp.2" "Comp.3"
 $ center  : Named num [1:3] 8.90e-16 -8.87e-16 7.67e-16
  ..- attr(*, "names")= chr [1:3] "Altezza" "Span" "Torace"
 $ scale   : Named num [1:3] 1 1 1
  ..- attr(*, "names")= chr [1:3] "Altezza" "Span" "Torace"
 $ n.obs   : int 1427
 $ scores  : num [1:1427, 1:3] -1.95 -1.36 -1.83 -2.02 -1.26 ...
  ..- attr(*, "dimnames")=List of 2
  .. ..$ : NULL
  .. ..$ : chr [1:3] "Comp.1" "Comp.2" "Comp.3"
 $ call    : language princomp(x = rid)
 - attr(*, "class")= chr "princomp"
NULL

source

     MLA.explor.pairs(pca1$scores)

plot

source

     kable(cov(pca1$scores))
Comp.1 Comp.2 Comp.3
Comp.1 2.382758 0.0000000 0.0000000
Comp.2 0.000000 0.5093955 0.0000000
Comp.3 0.000000 0.0000000 0.1078464

source

     tot  =sum(pca1$sdev^2) #pca1$sdev sono le dev. standard delle nuove variabili
     perc =pca1$sdev^2/tot
     print(pca1$loadings[,1],cutoff=0)

output

  Altezza      Span    Torace 
0.6085715 0.6063977 0.5117838

source

     print(sum(pca1$loadings[,1]^2))

output

[1] 1

source

     print(pca1$loadings,cutoff=0)

output


Loadings:
        Comp.1 Comp.2 Comp.3
Altezza  0.609  0.352  0.711
Span     0.606  0.372 -0.703
Torace   0.512 -0.859 -0.013

               Comp.1 Comp.2 Comp.3
SS loadings     1.000  1.000  1.000
Proportion Var  0.333  0.333  0.333
Cumulative Var  0.333  0.667  1.000

source

     print(round(t(pca1$loadings)%*%pca1$loadings,4))

output

       Comp.1 Comp.2 Comp.3
Comp.1      1      0      0
Comp.2      0      1      0
Comp.3      0      0      1

2.1 Confronto fra spazio originale e trasformato:

Landing spot

Mettiamo a fianco le due figure precedenti (in 3d):

La prima variabile principale, lungo il primo asse principale, in rosso, spiega il 79.43% della variabilità complessiva.

Di contro l’ultima variabile principale, la terza, asse blu, spiega soltanto il 3.59% della variabilità.

output

glX 
  5

output

Importance of components:
                          Comp.1    Comp.2     Comp.3
Standard deviation     1.5430776 0.7134693 0.32828463
Proportion of Variance 0.7942527 0.1697985 0.03594879
Cumulative Proportion  0.7942527 0.9640512 1.00000000

output

Importance of components:
                          Comp.1    Comp.2     Comp.3
Standard deviation     1.5430776 0.7134693 0.32828463
Proportion of Variance 0.7942527 0.1697985 0.03594879
Cumulative Proportion  0.7942527 0.9640512 1.00000000

You must enable Javascript to view this page properly.

2.2 Analisi in Componenti principali per 7 variabili

Utilizziamo ancora lo stesso data frame antropometric ma prendiamo stavolta 7 variabili, di cui riportiamo ancora delle statistiche e grafici descrittivi, insieme con la matrice dei grafici a coppia:

Data Frame Summary

dati

Dimensions: 1427 x 7
Duplicates: 0
No Variable Stats / Values Freqs (% of Valid) Graph Valid Missing
1 Altezza [numeric] Mean (sd) : 151.9 (10.1) min < med < max: 127 < 151 < 183 IQR (CV) : 15 (0.1) 54 distinct values 1427 (100%) 0 (0%)
2 Peso [numeric] Mean (sd) : 45 (10.7) min < med < max: 21 < 43 < 100 IQR (CV) : 14 (0.2) 65 distinct values 1427 (100%) 0 (0%)
3 Torace [numeric] Mean (sd) : 75.6 (7.8) min < med < max: 57 < 74 < 104 IQR (CV) : 10 (0.1) 44 distinct values 1427 (100%) 0 (0%)
4 Cranio [numeric] Mean (sd) : 54.8 (1.6) min < med < max: 50 < 55 < 60 IQR (CV) : 2 (0) 11 distinct values 1427 (100%) 0 (0%)
5 Bisacrom [numeric] Mean (sd) : 34.5 (3) min < med < max: 23 < 34 < 46 IQR (CV) : 4 (0.1) 21 distinct values 1427 (100%) 0 (0%)
6 Bitrocan [numeric] Mean (sd) : 26.3 (2.8) min < med < max: 20 < 26 < 38 IQR (CV) : 4 (0.1) 18 distinct values 1427 (100%) 0 (0%)
7 Span [numeric] Mean (sd) : 153.6 (11.2) min < med < max: 123 < 153 < 184 IQR (CV) : 16 (0.1) 60 distinct values 1427 (100%) 0 (0%)

Generated by summarytools 0.9.6 (R version 4.0.2)
2020-09-12

plot

Altezza Peso Torace Cranio Bisacrom Bitrocan Span
Altezza 1.00 0.72 0.59 0.48 0.75 0.74 0.89
Peso 0.72 1.00 0.91 0.54 0.75 0.84 0.69
Torace 0.59 0.91 1.00 0.48 0.69 0.77 0.58
Cranio 0.48 0.54 0.48 1.00 0.50 0.49 0.47
Bisacrom 0.75 0.75 0.69 0.50 1.00 0.76 0.78
Bitrocan 0.74 0.84 0.77 0.49 0.76 1.00 0.71
Span 0.89 0.69 0.58 0.47 0.78 0.71 1.00

Le variabili sono a due a due correlate, ma ovviamente stavolta non possiamo rappresentarle graficamente tutte.

L’analisi delle componenti principali ci permette però di ricavare le combinazioni di maggior varianza delle variabili originarie. Questo ci consente di rappresentare i dati in uno spazio di due o tre dimensioni ma anche di misurare in qualche modo una sorta di correlazione multipla fra le variabili.

Potremo anche vedere il grado di collinearità fra le variabili.

Calcoliamo i 7 autovalori, e riportiamo anche le deviazioni standard (ricavate dalla teoria asintotica, per la quale \(Var(\hat{\lambda_j})=\frac{2\lambda_j^2}{n}\), nel caso di campioni provenienti da distribuzioni normali multivariate: ci danno comunque un’ordine di grandezza dell’eventuale errore campionario)

output

[1] "Autovalori"

output

Comp.1 Comp.2 Comp.3 Comp.4 Comp.5 Comp.6 Comp.7 
 5.082  0.666  0.624  0.246  0.207  0.105  0.065

output

[1] "Errori standard approssimati"

output

Comp.1 Comp.2 Comp.3 Comp.4 Comp.5 Comp.6 Comp.7 
 0.190  0.025  0.023  0.009  0.008  0.004  0.002

In sintesi possiamo vedere che la variabilità globale delle sette variabili del dataset dati è spiegata per una percentuale di 72.65% dalla prima componente, mentre le prime tre componenti, rappresentate nel grafico 3d successivo, spiegano una percentuale di 91.09%

plot

output

Importance of components:
                         Comp.1    Comp.2     Comp.3
Standard deviation     2.254309 0.8160813 0.78985197
Proportion of Variance 0.797567 0.1045219 0.09791109
Cumulative Proportion  0.797567 0.9020889 1.00000000

You must enable Javascript to view this page properly.

2.3 I loadings (autovettori) e i primi due assi principali

output

       Comp.1 Comp.2 Comp.3 Comp.4 Comp.5 Comp.6 Comp.7
Comp.1      1      0      0      0      0      0      0
Comp.2      0      1      0      0      0      0      0
Comp.3      0      0      1      0      0      0      0
Comp.4      0      0      0      1      0      0      0
Comp.5      0      0      0      0      1      0      0
Comp.6      0      0      0      0      0      1      0
Comp.7      0      0      0      0      0      0      1

output

         Comp.1 Comp.2 Comp.3
Altezza   0.388  0.416  0.245
Peso      0.409 -0.194 -0.325
Torace    0.376 -0.333 -0.488
Cranio    0.283 -0.674  0.681
Bisacrom  0.393  0.161  0.023
Bitrocan  0.399 -0.001 -0.226
Span      0.384  0.446  0.285

plot

2.4 Complementi. Il biplot: relazioni fra osservazioni e variabili

plot

2.5 Parte teorica

Vedremo che se \(\Sigma\) è la matrice di varianza e covarianza dei dati, le varianze delle componenti principali sono gli autovalori di \(\Sigma\), mentre i coefficienti delle nuove variabili (loadings), sono i corrispondenti autovettori normalizzati.

Indichiamo con \(\boldsymbol{ \gamma} _{j}\) un autovettore di \(\mathrm{V}\left[\mathbf{X}\right]\), (normalizzato, ossia con \(\boldsymbol{ \gamma} _{j}^{\mathsf{T}}\boldsymbol{ \gamma} _{j}=1\)) e con \(\lambda_j\) il corrispondente autovalore; allora si ha:

\[ \boldsymbol{ \gamma} _{j}^{\mathsf{T}} \mathrm{V}\left[\mathbf{X}\right] \boldsymbol{ \gamma} _{j}=\boldsymbol{ \gamma} _{j}^{\mathsf{T}} \lambda_j \boldsymbol{ \gamma} _{j}= \lambda_j \] La prima delle precedenti eguaglianze deriva dalla definizione degli autovettori e degli autovalori; la seconda eguaglianza deriva dalla condizione di normalizzazione degli autovettori \(\boldsymbol{ \gamma} _{j}^{\mathsf{T}}\boldsymbol{ \gamma} _{j}=1.\)

Deriva dall’equazione fondamentale: \[ \mathrm{V}\left[\mathbf{X}\right]\boldsymbol{ \gamma} _{j}=\lambda_j \boldsymbol{ \gamma} _{j} \] e con la convenzione che gli autovalori siano ordinati in senso decrescente: \(\lambda_1 \ge \lambda_2 \ge \dots \ge \lambda_p\) (in effetti per alcune delle scomposizioni fatte occorrerebbe anche ipotizzare che siano tutti distinti, ma per ora non è necessario precisare altro)

dal momento che una matrice di varianze e covarianze è sempre semidefinita positiva

2.6 Estensioni e generalizzazioni.

L’analisi delle componenti principali è la più semplice delle cosiddette tecniche fattoriali. Appartiene alla famiglia dell’analisi delle correlazioni canoniche (massimizzare la correlazione fra due gruppi di variabili), come l’analisi delle corrispondenze (per variabili categoriali); sono possibili estensioni anche per dati suddivisi in gruppi (analisi delle funzioni discriminanti, recentemente ribattezzate come tecniche di SVM, Support Vector Machine).