METODOS MULTIVARIADOS: Análisis de correspondencia

Yeimy Maryury Montilla Montilla
Julio 15 del 2021

Introducción

En este ejercicio se estudió la existencia de relaciones entre diferentes variables del cacao que permitan identificar agrupaciones. En un primer paso se identifican los datos atípicos mediante un análisis categórico de componentes principales y posteriormente se realizó un Análisis de Correspondencia Múltiple para identificar las variables y categorías más representativas y establecer posibles agrupaciones. En particular, se trabajó con 6 variables:

  • Forma de la mazorca (fm)
  • Forma del ápice (fm)
  • Constricción basal (cb)
  • Rugosidad (rug)
  • Color de la mazorca (c_m)
  • Color de la semilla (c_s)

Datos

library(readxl)
library(Gifi)
df = read_excel('caract_cacao.xlsx', 'data')
df = data.frame(df);

dim(df)
## [1] 120   6
names(df)
## [1] "c_s" "c_m" "rug" "cb"  "fa"  "fm"

El conjunto de datos contiene la información de 6 variables para 120 registros.

Se consideran cuatro formas de mazorca (1: angoleta, 2: amelonado, 3: cundeamor, 4: calabacillo), seis formas de ápice (1: puntiagudo, 2: agudo, 3: obtuso, 4: redondeado, 5: pezón, 6: dentado), cuatro constricciones basales (1: ausente, 2: escasa, 3: intermedia, 4: clara), tres rugosidades (1: lisa, 2: intermedia, 3: áspera), dos colores de mazorca (1: verde y 2: morada) y dos colores de semilla (1: purpura obscura y 0: púrpura clara).

Análisis categórico de componentes principales (PRINCALS)

El Análisis categórico de componentes principales (PRINCALS) se uso en este caso para identificar datos atipicos dentro del conjunto de datos.

fitord <- princals(df,ordinal=TRUE,ndim = 2)  ## ordinal PCA
summary(fitord)
## 
## Loadings (cutoff = 0.1):
##     Comp1  Comp2 
## rug -0.665 -0.237
## cb  -0.942 -0.212
## fm  -0.943 -0.209
## c_s -0.259  0.860
## c_m -0.272  0.807
## fa   0.176 -0.608
## 
## Importance (Variance Accounted For):
##                  Comp1   Comp2
## Eigenvalues     2.3912  1.9057
## VAF            39.8534 31.7618
## Cumulative VAF 39.8500 71.6200
plot(fitord, "loadplot", main = "Loadings Plot ABC Data") 

plot(fitord, "biplot", main = "Biplot ABC Data")

Los componentes 1 y 2 explican el 71.6% de la varianza de los datos. En la anterior gráfica se muestra dónde se ubican los casos individuales en el espacio de los primeros dos componentes. El registro 57 se aleja totalmente de la distribución de los demás datos por lo que se puede suponer que es un dato atípico. El registro 57 es el único con forma de ápice dentado 6.

fitord <- princals(df,ordinal=TRUE,ndim = 3)  ## ordinal PCA
summary(fitord)
## 
## Loadings (cutoff = 0.1):
##     Comp1  Comp2  Comp3 
## cb  -0.904 -0.252 -0.345
## fm  -0.904 -0.253 -0.343
## rug        -0.814  0.430
## fa          0.783 -0.485
## c_s -0.506  0.437  0.632
## c_m -0.487  0.441  0.641
## 
## Importance (Variance Accounted For):
##                  Comp1   Comp2   Comp3
## Eigenvalues     2.1276  1.7880  1.4671
## VAF            35.4602 29.7996 24.4511
## Cumulative VAF 35.4600 65.2600 89.7100

Realizando el análisis con tres componentes, se explica el 89,71 % de la varianza de los datos. En la siguiente gráfica se muestra dónde se ubican los datos individuales en el espacio de los primeros tres componentes. Los registros 57 y 22 se alejan de la distribución de los demás datos por lo que se puede suponer que son datos atípicos.

library("scatterplot3d")
s3d<-scatterplot3d(fitord$objectscores[, 1], fitord$objectscores[, 2],fitord$objectscores[, 3],grid=TRUE,xlab="PC1",ylab="PC2", zlab="PC3", pch = 16,color = 'blue',angle = 45)
text(s3d$xyz.convert(fitord$objectscores[, 1:3]), labels = rownames(df),cex= 0.7, col = "black",pos=2.5)

df[57,]
##    c_s c_m rug cb fa fm
## 57   0   0   2  1  6  1
df[22,]
##    c_s c_m rug cb fa fm
## 22   1   1   1  1  2  1

Considerando lo anterior, se descartan los datos con forma de ápice 6, condición que solo cumplía un dato. por otra parte, se descartaron los datos con constitución basal ausente y rugosidad lisa, donde solo se ubicaba un dato.

df1 = df[-c(which(df$fa == 6),
            which(df$cb == 1 & df$rug == 1)),]

fitord <- princals(df1,ordinal=TRUE,ndim = 3)  ## ordinal PCA
summary(fitord)
## 
## Loadings (cutoff = 0.1):
##     Comp1  Comp2  Comp3 
## cb  -0.920  0.277  0.277
## fm  -0.920  0.277  0.277
## rug        -0.772  0.525
## fa          0.804 -0.469
## c_s  0.386  0.531  0.628
## c_m  0.407  0.510  0.637
## 
## Importance (Variance Accounted For):
##                  Comp1   Comp2   Comp3
## Eigenvalues     2.0796  1.8673  1.4476
## VAF            34.6598 31.1225 24.1270
## Cumulative VAF 34.6600 65.7800 89.9100

Se realizo el mismo procedimiento anterior, donde no se detectaron visualmente más datos atípicos.

s3d<-scatterplot3d(fitord$objectscores[, 1], fitord$objectscores[, 2],fitord$objectscores[, 3],grid=TRUE,xlab="PC1",ylab="PC2", zlab="PC3", pch = 16,color = 'blue')
text(s3d$xyz.convert(fitord$objectscores[, 1:3]), labels = rownames(df),cex= 0.7, col = "black",pos=2.5)
text(s3d$xyz.convert(fitord$objectscores[, 1:3]), labels = rownames(df1),cex= 0.7, col = "black",pos=2.5)

En la anterior imagen no se evidencian demás datos atípicos.

Análisis de Correspondencia Múltiple (ACM)

Teniendo en cuenta que los datos se espera que no sean cuantitativos continuos, deben ser llevados (convertidos) a factor.

df2 = within(df1,{
  c_s = as.factor(c_s)
  c_m = as.factor(c_m)
  rug = as.factor(rug)
  cb = as.factor(cb)
  fa = as.factor(fa)
  fm = as.factor(fm)
})

Para identificar de manera preliminar cuales de las variables pueden llegar a tener mayor representatividad, se hace un análisis de pesos en función de las respuestas por variable, en este caso sobre los datos factorizados.

library(factoextra)
## Loading required package: ggplot2
## Welcome! Want to learn more? See two factoextra-related books at https://goo.gl/ve3WBa
library(FactoMineR)
res.mca = MCA(df2)

eig.val <- get_eigenvalue(res.mca)
head(eig.val)
##       eigenvalue variance.percent cumulative.variance.percent
## Dim.1  0.4559681        19.541490                    19.54149
## Dim.2  0.4322810        18.526328                    38.06782
## Dim.3  0.3408640        14.608458                    52.67628
## Dim.4  0.2427341        10.402890                    63.07917
## Dim.5  0.2185789         9.367666                    72.44683
## Dim.6  0.1666667         7.142857                    79.58969

En el paso anterior se aprecia que las primeras tres dimensiones se explican el 52.68% de la varianza.

fviz_contrib(res.mca, choice="var", axes = c(1,2,3), fill="brown", title = "Contribución de las categorías al subespacio 1-2-3")

Conforme al gráfico de contribución se puede observar que las categorías fa_1, C_s_1, fa_5 y c_m_1, son la que menos contribuyen considerando las primeras tres dimensiones. Al revisar los datos originales, se encuentra fa_1 tiene solo un solo dato, mientras que fa_5 tiene 3 datos, por lo tanto, se consideran datos atípicos y se excluyen del análisis. Por otra parte, las categorías C_m_1 y c_s_1 tienen un número considerable de datos.

df3 = df1[-c(which(df1$fa == 1),
             which(df1$fa == 5)),]

Luego de extraer los datos considerados como atípicos, se repite el ejercicio con el propósito de identificar otros datos atípicos.

df4 = within(df3,{
  c_s = as.factor(c_s)
  c_m = as.factor(c_m)
  rug = as.factor(rug)
  cb = as.factor(cb)
  fa = as.factor(fa)
  fm = as.factor(fm)
})

res.mca = MCA(df4)

eig.val <- get_eigenvalue(res.mca)
head(eig.val)
##       eigenvalue variance.percent cumulative.variance.percent
## Dim.1  0.4564803        22.824017                    22.82402
## Dim.2  0.4302339        21.511697                    44.33571
## Dim.3  0.3410476        17.052379                    61.38809
## Dim.4  0.2417607        12.088034                    73.47613
## Dim.5  0.2133620        10.668099                    84.14423
## Dim.6  0.1233014         6.165068                    90.30929

En el paso anterior se aprecia que extrayendo los datos categorizados en fa_1 y fa_5 se logra una explicación el 61.38 % de la varianza.

fviz_contrib(res.mca, choice="var", axes = c(1,2,3), fill="brown", title = "Contribución de las categorías al subespacio 1-2-3")

El comportamiento de las categorías c_m_1 y c_s_1 es interesante. Observando las variables en los datos originales se encuentra que de 120 registros, 108 se ubican en estás dos categorías, por lo tanto, se puede considerar que estás dos variables no tienen un peso suficinente para realizar la catgorización del cacao. A explora que pasa si se omiten estas dos variables del análisis.

df5 =df3
df5$c_s <- NULL
df5$c_m <- NULL
df6 = within(df5,{
  #c_s = as.factor(c_s)
  #c_m = as.factor(c_m)
  rug = as.factor(rug)
  cb = as.factor(cb)
  fa = as.factor(fa)
  fm = as.factor(fm)
})

res.mca = MCA(df6)

eig.val <- get_eigenvalue(res.mca)
head(eig.val)
##       eigenvalue variance.percent cumulative.variance.percent
## Dim.1  0.6508525        26.034101                    26.03410
## Dim.2  0.6225073        24.900290                    50.93439
## Dim.3  0.5041268        20.165070                    71.09946
## Dim.4  0.3213547        12.854190                    83.95365
## Dim.5  0.1870061         7.480244                    91.43390
## Dim.6  0.1536657         6.146628                    97.58052

Descartando del análisis, las variables con menor contribucion en las primeras tres dimensiones se logra explicar una varianza del 71.1%

fviz_contrib(res.mca, choice="var", axes = c(1,2,3), fill="brown", title = "Contribución de las categorías al subespacio 1-2-3")

Ahora las primeras tres dimensiones explican el 71.% de todos los datos.

options(ggrepel.max.overlaps = Inf)
fviz_mca_biplot(res.mca, 
               repel = TRUE, 
               ggtheme = theme_minimal())

En la gráfica anterior se observan dos pequeños grupos, el primero compuesto por las categorías fm_1 y cb_1; el segundo compuesto por las categorias fm4 y cb4. Observando los datos originales, la primer grupo tiene 8 datos, mientras que el segundo grupo tiene 7 datos.

fviz_ellipses(res.mca,ellipse.type="t", axes=c(1,2), c("cb", "fm", "rug", "fa"), geom = "point")
## Warning in MASS::cov.trob(data[, vars]): Probable convergence failure

## Warning in MASS::cov.trob(data[, vars]): Probable convergence failure

## Warning in MASS::cov.trob(data[, vars]): Probable convergence failure

## Warning in MASS::cov.trob(data[, vars]): Probable convergence failure

## Warning in MASS::cov.trob(data[, vars]): Probable convergence failure

## Warning in MASS::cov.trob(data[, vars]): Probable convergence failure

## Warning in MASS::cov.trob(data[, vars]): Probable convergence failure

En la gráfica con elipses de concentración de las categorías sobre las nubes de datos se observa como en definitiva existe una relación estrecha entre la constricción basal (cb) y la forma de la mazorca (fm). La forma angoleta corresponde con la constitución basal ausente, la forma amelonada corresponde con la constitución basal ausente, la forma cundeamor corresponde con la constitución basal intermedia y la forma calabocillo corresponde a constitución basal clara. Solo un registro no sigue este comportamiento y corresponde un dato con forma cundeamor con constitución basal escasa.

fviz_ellipses(res.mca,ellipse.type="t", axes=c(2,3), c("cb", "fm", "rug", "fa"), geom = "point")
## Warning in MASS::cov.trob(data[, vars]): Probable convergence failure

## Warning in MASS::cov.trob(data[, vars]): Probable convergence failure
## Warning: Computation failed in `stat_ellipse()`:
## the leading minor of order 2 is not positive definite
## Warning in MASS::cov.trob(data[, vars]): Probable convergence failure

## Warning in MASS::cov.trob(data[, vars]): Probable convergence failure

## Warning in MASS::cov.trob(data[, vars]): Probable convergence failure
## Warning: Computation failed in `stat_ellipse()`:
## the leading minor of order 2 is not positive definite
## Warning in MASS::cov.trob(data[, vars]): Probable convergence failure

df5[12,]
##    rug cb fa fm
## 12   3  2  3  3

Conclusiones

A partir del anterior ejercicio se puede concluir:

La explicación del 71.1% de la varianza se logró mediante la omisión de seis datos atípicos identificados y dos variables (color de la mazorca y color del grano). De los datos atípicos omitidos dos eran únicos las categorías: 1) puntiagudo y 6) dentado y tres eran únicos en la categoría 5) pezón, las anteriores de la variable forma de ápice. La omisión de las variables color de la mazorca y color del grano se debe a que la mayoría (108 de 120) de los datos pertenecían a las categorías color de la mazorca morada y color de la semilla purpura, razón por la cual no contribuyen en la explicación de las variables. De acuerdo con el anterior hallazgo los expertos deben considerar si realmente son pertinentes estas categorías y las variables para los objetivos del estudio.

Existe una clara relación entre la constricción basal y la forma de la mazorca, lo cual permitirá generar cuatro agrupaciones Mientras que las agrupaciones considerando las variables rugosidad y forma del ápice no tienen una relación tan clara.