Componente Principales

 

EJERCICIO 4 - Práctica 4 - Parte 2

En las Tablas 2 y 3 y en el archivo P4-2-ej4-2019.txt figuran los datos de las mediciones de huesos y dientes de ratones campestres (de la especie Microtus). Las variables son:

• \(y_1\) = ancho del molar 1 superior izquierdo
• \(y_2\) = ancho del molar 2 superior izquierdo
• \(y_3\) = ancho del molar 3 superior izquierdo
• \(y_4\) = longitud de la inserción del incisivo
• \(y_5\) = longitud del hueso del paladar
• \(y_6\) = longitud del cóndilo del incisivo o longitud del cráneo
• \(y_7\) = altura del cráneo por encima de la bullae
• \(y_8\) = ancho del cráneo a través de la cara

Las variables y1 a y 5 son en mm/1000; las variables y 6 a y 8 son en mm/100. La variable grupo indica la especie de los ratones, siendo 1 la especie Microtus multiplex (Tabla2) y2 la especie Microtus subterraneus (Tabla 3). Consideremos solamente las primeras 3 variables (\(y_1\) , \(y_2\) , \(y_3\) ) que son el ancho de los molares superiores izquierdos 1, 2 y 3 respectivamente, para las 43 ratas del grupo 1.

(a) Hallar las componentes principales muestrales y los autovalores de \(S\).

(b) Hallar los porcentajes de la variabilidad total explicados por la primera y por las dos primeras componentes, e interpretarlas en función de las variables originales.

 

RESOLUCIóN

rm(list = ls())
graphics.off()

 

DATOS

## DATOS ###########################################
#ds = read.table(file.choose(), header = )
ds = read.table('/home/maurolioy/Documentos/MEL/MAESTRIA/Multivariado/Practicas/P4/II/R/datos/datos_P4_ex4_p2.txt',header = T)
# Nos quedamos con las variables Y1,Y2,Y3 del grupo 1
ds1 = ds[ds$Grupo == 1,-c(1,5:9)]

 

AUTOVALORES Y AUTOVECTORES DE S

Recordemos los estimadores de MV de \(\mathbf{\mu}\) y \(\Sigma\):

\[ \hat{\mathbf{\mu}} = \bar{ \mathbf{x} } \quad \text{y} \quad S = \frac{Q}{n} \]

n = nrow(ds1) #43
Q = (n-1) * cov(ds1)
S = Q/n

Matriz S

print(S)

          Y1        Y2       Y3
Y1 16752.435  9514.542 11839.62
Y2  9514.542  9992.946 10780.38
Y3 11839.619 10780.381 22716.02

gamma_objet = eigen(S)
autovalores = gamma_objet$values

Autovalores de S

\(\lambda_1\)	\(\lambda_2\)	\(\lambda_3\)
39040.432	7650.393	2770.575

 

COMPONENTES PRINCIPALES MUESTRALES

Sea:

\[\hat{\Sigma} = \hat{\Gamma}\hat{\Lambda}\hat{\Gamma}^T\] Entonces: \[\hat{\mathbf{v}}_i = \hat{\Gamma}^T (\mathbf{x}_i - \mathbf{\bar{x}}) \] Donde:

• \(\mathbf{x}_i\) y \(\mathbf{\bar{x}}\) \(\in R^p\) (en este caso \(p=3\)).

• \(\hat{\mathbf{v}}_i \in R^q\) es el vector de componetes principales de la medicion \(i\) (en este caso \(q=p=3\)).

• \(\hat{\Gamma}\) es la matriz de autovectores de \(S\).

• \(\hat{\Gamma} \in R^{p \times q}\) (en este caso \(p=3\) y \(q=3\)).

Gamma = gamma_objet$vectors
x_bar = matrix(colMeans(ds1), ncol = 1)

x_menos_mu = sweep(as.matrix(ds1), 2, x_bar)
v = t(Gamma)%*%t(x_menos_mu)

print(t(v))

          [,1]         [,2]         [,3]
1    59.568443   94.0583571  -28.9954830
2   297.055378   78.1333089  -39.7355281
3   205.841934  -26.2678804  -80.7416169
4    -1.842134  -26.4668684  -41.1778911
5  -252.230338   71.4937339  -62.8681576
6  -320.538704  -43.8638573  -31.2360013
7   -10.228180   95.2886095   16.3136545
8  -136.051953   60.6964004  -62.8870859
9  -113.914234  -95.2521057  -23.1380622
10  114.786196   58.6462885   48.4188899
11 -211.567605  124.0111796   47.1430237
12  135.592206  -81.3482954  -69.0316376
13   -1.337150 -145.0523394   68.4498018
14 -275.780627  -30.9217351   83.9734428
15  148.342520  143.7334199  -14.6933371
16   47.205116  -43.2032560    8.4130414
17   33.509548  126.5093019   31.0200999
18   33.987894   29.4358811   68.4035474
19   17.286013   46.9615175  -17.5900526
20 -147.343175   40.1034393  -86.4972162
21  170.804269 -175.2761259  -13.5188169
22  -58.160264  -60.0061981   74.3053613
23 -374.550795  -99.0035583  -14.1097829
24  515.937708   21.1398104   53.9622125
25   82.607084   40.4582367   41.2985288
26 -517.348002  178.1529157   12.3145738
27  143.062825  123.5922313   95.4450184
28  -81.490612   38.7692948   74.6165946
29  224.953450  -20.1423584  -17.0373947
30  -87.148744  -33.1154037   11.7709451
31  127.466137   11.8081659  -11.7093844
32  189.341256  -57.1829915   66.0575096
33  170.989584  -43.8655969  -32.2067420
34 -127.615092 -139.4962252  -24.4588147
35    2.878758 -125.6400193  -79.6798304
36 -295.836338 -146.0794312  102.7743326
37   82.539372  -26.1345961   -3.1693252
38  -90.501640 -116.9910695   17.0274438
39   55.753161   -0.2053025   -0.3108861
40  268.062205   -6.5581573  -36.4565786
41  143.047229  -33.1832280  -13.1250288
42 -220.232846   74.1256668 -115.7475813
43   53.100147  118.1388403   -1.5857862

$nbsp;

VARIABILIDAD

Las componentes principales explican una porción de la dispersión total medida a través de la traza.

• \(Tr(\Sigma) = Tr(\Gamma\Lambda\Gamma^T)=Tr(\Gamma^T\Gamma\Lambda) = Tr(\Lambda)\)

• \(\hat{\Lambda} = \text{diag}(\hat{\lambda_1},\hat{\lambda_2},\hat{\lambda_3})\)
• \(Tr(\hat{\Lambda}) = \hat{\lambda_1}+\hat{\lambda_2}+\hat{\lambda_3}\)

• \(\hat{\text{var}}(v_j)=\hat{\lambda_j}\)

Queremos calcular: \[\frac{\hat{\lambda_1}}{ \sum_{j=1}^{3} \hat{\lambda_j}}\]

y

\[\frac{\hat{\lambda_1} + \hat{\lambda_2} }{ \sum_{j=1}^{3} \hat{\lambda_j}}\]

## VARIABILIDAD ####################################
suma_p = sum(autovalores)
# variabilidad 1
pp1 = autovalores[1] / suma_p
# variabilidad 1+2
pp12 = (autovalores[1] + autovalores[2]) / suma_p

Luego:

\[\frac{\hat{\lambda_1}}{ \sum_{j=1}^{3} \hat{\lambda_j}} = 0.789\]

\[\frac{\hat{\lambda_1} + \hat{\lambda_2} }{ \sum_{j=1}^{3} \hat{\lambda_j}} = 0.943\]

 

Correlaciones

\[\text{Corr}(\mathbf{x_{j}},\mathbf{v_{l}}) = \rho_{\mathbf{x_j},\mathbf{v_l}} = \gamma_{l,j}\sqrt{\frac{\lambda_l}{\sigma_{jj}}}\]

Corr_Xjvl = matrix(0,nrow = 3,ncol=3)
## este loop me da para cada variable Xj la correlación con todas las componentes vl
for(j in 1:3){
  gammma_lj = Gamma[j,]
  correlaciones_j = gammma_lj * sqrt(autovalores / S[j,j])
  Corr_Xjvl[j,] = correlaciones_j
  
} 

data.frame(t(Corr_Xjvl),row.names = c('v1','v2','v3'))