El análisis de componentes principales parte de la matriz de correlaciones de las variables. El problema radica cuando, en el conjunto de variables, se contienen variables ordinales pues es necesario calcular las correlaciones entre dichas variables y correlaciones de las variables continuas con ellas.
A continuación se presentan alternativas para calcular dichas correlaciones. Para el caso de dos variables ordinales usaremos el Polychoric Correlation Coefficient y para calcular la correlación entre una variable ordinal y una continua se utilizará el Polyserial Correlation Coefficient.
Estos métodos están implementados en el paquete en el paquete psych de R.
La correlación polychoric es una técnica que permite estimar la correlación entre pares variables latentes (normalmente distribuida) que correspondan a dos variables ordinales. En el caso particular de que las variables seas dicotómicas entonces se le conoce con el nombre de correlación tetrachoric.
A continuación se muestra el cálculo para el caso particular de una variable dicotómica. En la siguiente elipse se muestra la distribución conjunta de un evento desde la perspectiva de dos jueces ($Y_{1}$ y $Y_{2}$) y las tolerancias $t_{1}, t_{2}$.
Los valores de $a, b, c$ y $d$ corresponden a la proporción de casos que caen en cada una de los casos mostrados en la elipse.
Por lo que sólo es necesario estimar el modelo representado por la elipse para encontrar los puntos de $t_{1}$ y $t_{2}$ así como el parámetro $\rho$ (rho) que determina la forma de la elipse. La forma en que se calcula dicha correlación es encontrando la cambinación de $t_{1}, t_{2}$ y $\rho$ que hacen que las proporciones de $a, b, c$ y $d$ observadas sean similares a las de la elipse.
El cálculo para el caso general no es tan sencillo, inlcuso en la historia de la estadística este tipo de correlación era conocidad por la dificulatad que implicaba cálculo. Para simplificar el cálculo Pearson creó tablas que ayudaban a calcular las series, sin embargo hoy en día se obtiene utilizando softwares especializados.
Suponiendo que la tabla de contigencia tiene la estructura anterior, entonces el coeficiente de correlación se puede calcular utilizando la regla de Cole (1949):
$$ \rho = \cos \left(\frac{180°}{1 + \sqrt{\frac{bc}{ad}}} \right) $$
El coeficiente de correlación polyserial es la correlación latente entre una variable continua X y una variable ordinal Y pero asumida para representar una variable continua con distribución normal.
Supuestos:
Principalmente se utilizan dos métodos para estimar la correlación:
knitr::opts_chunk$set(fig.dpi = 300, fig.height = 7, fig.width = 7) library(psych)
Este paquete cuenta un diversas funciones, entre ellas tenemos las siguientes:
Asi se veria una correlación tetrachorica fuerte
draw.tetra(.9,2,0.5)
Asi se veria una correlación tetrachorica debil
draw.tetra(.1,0,0)
library(ggplot2)
Se utilizará la base de datos diamonds del paquete ggplot2. Esta base cuenta con el precio y otros 9 atributos de aproximadamente 53,940 diamantes. Las variables son:
Primero cargarmos la base a un data frame.
Las variables cut, color y clarity son ordinales. Por lo tanto las convertimos en un ordered factor
dfOrd <- data.frame(diamonds) head(dfOrd)
Este conjunto de datos es de particular interés en R ya que permite realizar predicciones respecto al precio de un diamante en función de variables continuas y ordinales con más de 5 categorias. Para entender un poco el conjunto de datos incluimos la variación del precio en función de la raiz cubica de los kilates de la piedra dependiendo de su claridad, corte y su color.
library(scales) cuberoot_trans = function() trans_new('cuberoot', transform = function(x) x^(1/3), inverse = function(x) x^3)
ggplot(aes(x = carat, y = price), data = diamonds) + geom_point(alpha = 0.5, size = 1, position = 'jitter',aes(color=clarity)) + scale_color_brewer(type = 'div', guide = guide_legend(title = 'Clarity', reverse = T, override.aes = list(alpha = 1, size = 2))) + scale_x_continuous(trans = cuberoot_trans(), limits = c(0.2, 3), breaks = c(0.2, 0.5, 1, 2, 3)) + scale_y_continuous(trans = log10_trans(), limits = c(350, 15000), breaks = c(350, 1000, 5000, 10000, 15000)) + ggtitle('Price (log10) by Cube-Root of Carat and Clarity') ggplot(aes(x = carat, y = price), data = diamonds) + geom_point(alpha = 0.5, size = 1, position = 'jitter',aes(color=cut)) + scale_color_brewer(type = 'div', guide = guide_legend(title = 'Cut', reverse = T, override.aes = list(alpha = 1, size = 2))) + scale_x_continuous(trans = cuberoot_trans(), limits = c(0.2, 3), breaks = c(0.2, 0.5, 1, 2, 3)) + scale_y_continuous(trans = log10_trans(), limits = c(350, 15000), breaks = c(350, 1000, 5000, 10000, 15000)) + ggtitle('Price (log10) by Cube-Root of Carat and Cut') ggplot(aes(x = carat, y = price), data = diamonds) + geom_point(alpha = 0.5, size = 1, position = 'jitter',aes(color=color)) + scale_color_brewer(type = 'div', guide = guide_legend(title = 'Color', reverse = F, override.aes = list(alpha = 1, size = 2))) + scale_x_continuous(trans = cuberoot_trans(), limits = c(0.2, 3), breaks = c(0.2, 0.5, 1, 2, 3)) + scale_y_continuous(trans = log10_trans(), limits = c(350, 15000), breaks = c(350, 1000, 5000, 10000, 15000)) + ggtitle('Price (log10) by Cube-Root of Carat and Color')
Posteriormente cambiamos las categorías de las variables ordinales a números consecutivos comenzando desde 1.
ordNum <- data.matrix(dfOrd) head(ordNum)
Podemos comparar los resultados de aplicar la correlación de Pearson como si los datos ordinales fueran en realidad de intervalo con los resultados de buscar la $\rho$ que mejor describe la elipse de los datos.Los resultados son los siguientes:
Por lo general la correlación de Pearson tiende a subestimar la correlación entre variables ordinales.
cor(ordNum[,2],ordNum[,3]) prop.table(table(ordNum[,2],ordNum[,3])) polychoric(ordNum[,c(2,3)])
Calculamos la matriz de correlaciones de las variables utilizando la función mixed.cor.
pc <- mixed.cor(x=ordNum,smooth=TRUE, correct=TRUE)
pc
Calculamos los primeros 5 componentes principales para la matriz anterior.
faPC <- principal(r=pc$rho, nfactors=5, rotate="none",scores=TRUE)
faPC
Para encontrar la representación de alguna observación en las nuevas coordenadas podemos utilizar la función predict.
predict(faPC,ordNum[1,],ordNum) fa.plot(faPC) fa.diagram(faPC)
Para determinar que variables contribuyen a cada una de las componentes principales podemos utilizar un Biplot. Les proponemos nombrar cada una de las componentes en base a estas variables.
faPC$scores <- psych::factor.scores(ordNum,faPC) biplot.psych(faPC,choose=c(1,2),main = "Componente 1 vs Componente 2") biplot.psych(faPC,choose=c(1,3),main = "Componente 1 vs Componente 3") biplot.psych(faPC,choose=c(2,3),main = "Componente 2 vs Componente 3") biplot.psych(faPC,choose=c(4,5),main = "Componente 4 vs Componente 5")
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.