Cap # 10, FIABILIDAD - Oscar Padilla
Las filas del conjunto de datos contienen los 50 estados, en orden alfabetico.
states =row.names(USArrests )
states [1] "Alabama" "Alaska" "Arizona" "Arkansas" "California"
[6] "Colorado" "Connecticut" "Delaware" "Florida" "Georgia"
[11] "Hawaii" "Idaho" "Illinois" "Indiana" "Iowa"
[16] "Kansas" "Kentucky" "Louisiana" "Maine" "Maryland"
[21] "Massachusetts" "Michigan" "Minnesota" "Mississippi" "Missouri"
[26] "Montana" "Nebraska" "Nevada" "New Hampshire" "New Jersey"
[31] "New Mexico" "New York" "North Carolina" "North Dakota" "Ohio"
[36] "Oklahoma" "Oregon" "Pennsylvania" "Rhode Island" "South Carolina"
[41] "South Dakota" "Tennessee" "Texas" "Utah" "Vermont"
[46] "Virginia" "Washington" "West Virginia" "Wisconsin" "Wyoming" El conjunto de datos contienen las cuatro variables.
names(USArrests )[1] "Murder" "Assault" "UrbanPop" "Rape" Se usará la función apply() que permite aplicar una función a cada columna o fila del dataset, en este caso usaremos la función mean().
apply(USArrests , 2, mean) Murder Assault UrbanPop Rape
7.788 170.760 65.540 21.232 Tenga en cuenta que la funcion apply() permite aplicar la funcion mean(), a cada fila o columna del conjunto de datos. La segunda entrada indica si se desea calcular la media de las filas o las columnas. -Vemos que hay en promedio tres veces mas violaciones que asesinatos, y mas de ocho veces mas asaltos que violaciones.
apply(USArrests , 2, var) Murder Assault UrbanPop Rape
18.97047 6945.16571 209.51878 87.72916 La variable UrbanPop mide el porcentaje de la poblacion en cada estado que vive en un area urbana, que no es un numero comparable al numero de violaciones en cada estado por cada 100,000 personas. Si no logramos escalar las variables antes de realizar PCA entonces la mayoria de los componentes principales que observamos estaria impulsada por la variable Assault, ya que tiene, por mucho, la media y la varianza mas grande. Por lo tanto, es importante estandarizar las variables para que tengan una media de cero y una desviacion estandar antes de realizar la PCA.
pr.out =prcomp (USArrests , scale =TRUE)La funcion prcomp() centra las variables para que tengan media cero. Al usar la opcion scale = TRUE, escalamos las variables para tener una desviacion estandar de uno. La salida de prcomp() contiene cantidades utiles.
names(pr.out )[1] "sdev" "rotation" "center" "scale" "x" El centro y los componentes de escala corresponden a las medias y desviaciones estandar de las variables quese utilizaron para escalar antes de implementar PCA.
pr.out$center Murder Assault UrbanPop Rape
7.788 170.760 65.540 21.232 pr.out$scale Murder Assault UrbanPop Rape
4.355510 83.337661 14.474763 9.366385 La matriz de rotacion proporciona las cargas principales de los componentes; Cada columna de rotacion pr.out$ contiene el vector de carga del componente principal correspondiente.
pr.out$rotation PC1 PC2 PC3 PC4
Murder -0.5358995 0.4181809 -0.3412327 0.64922780
Assault -0.5831836 0.1879856 -0.2681484 -0.74340748
UrbanPop -0.2781909 -0.8728062 -0.3780158 0.13387773
Rape -0.5434321 -0.1673186 0.8177779 0.08902432En general hay min(n - 1, p) componentes principales informativos en un conjunto de datos con n observaciones y p variables. Usando la funcion prcomp(), no necesitamos multiplicar explicitamente los datos por el componente principal de carga de vectores con el fin de obtener el Componente principal en los vectores de puntuacion. Mas bien, la matriz x de 50×4 tiene como columnas los vectores de puntuacion del componente principal.
dim(pr.out$x )[1] 50 4Podemos trazar los dos primeros componentes principales de la siguiente manera:
biplot (pr.out , scale =0)
Observe que esta figura es una imagen reflejada de la Figura 10.1. Recordar que los componentes principales solo son unicos hasta un cambio de signo, por lo que podemos reproducir la Figura 10.1 realizando algunos pequeños cambios:
pr.out$rotation=-pr.out$rotation
pr.out$x=-pr.out$x
biplot (pr.out , scale =0)
La funcion prcomp() tambien genera la desviacion estandar de cada componente principal. Por ejemplo, en el conjunto de datos de USArrests, podemos acceder a estas desviaciones estandar de la siguiente manera:
pr.out$sdev[1] 1.5748783 0.9948694 0.5971291 0.4164494La varianza explicada por cada componente principal se obtiene de esta forma:
pr.var =pr.out$sdev ^2
pr.var[1] 2.4802416 0.9897652 0.3565632 0.1734301Para calcular la proporcion de varianza explicada por cada componente principal, simplemente dividimos la varianza explicada por cada componente principal por la varianza total explicada por los cuatro componentes principales
pve=pr.var/sum(pr.var )
pve[1] 0.62006039 0.24744129 0.08914080 0.04335752Vemos que el primer componente principal explica el 62.0% de la varianza en los datos, el siguiente componente principal explica el 24.7% de la varianza, y asi sucesivamente. Podemos trazar el PVE explicado por cada componente, asi como el PVE acumulado, de la siguiente manera:
plot(pve, xlab = "Principal Component", ylab = "Proportion of Variance Explained", ylim = c(0,1), type = 'b')
plot(cumsum(pve), xlab = "Principal component", ylab = "Cumulative Proportion of Variance Explained", ylim = c(0,1), type = 'b')
La función cumsum() calcula la suma acumulativa de los elementos de un vector numerico.
Por ejemplo:
a=c(1,2,8,-3)
cumsum (a)[1] 1 3 11 8