Práctica 3

Parte 1

Aplique lo aprendido sobre Clúster jerárquico a la base de datos denominada Poblacion2010_2022. Esta base de datos contiene la población por municipios en Puerto Rico para los años 2010 y 2022.

Base de datos

Para esto importaremos la base base de datos denominada Poblacion2010_2022. Luego procederemos a escalar los datos.

library(readxl)
datos<-read_excel("Poblacion2010_2022.xlsx")
head(datos)

## # A tibble: 6 × 3
##   Municipios   Pob_2010 Pob_2022
##   <chr>           <dbl>    <dbl>
## 1 Adjuntas        19483    17905
## 2 Aguada          41959    37666
## 3 Aguadilla       60949    53931
## 4 Aguas Buenas    28659    23538
## 5 Aibonito        25900    24555
## 6 Añasco          29261    25026

data.scaledpob <- cbind(datos[,1], scale(x = datos[,-1],
                     center = TRUE,
                     scale = TRUE))

head(data.scaledpob)

##     Municipios   Pob_2010    Pob_2022
## 1     Adjuntas -0.5227021 -0.50868616
## 2       Aguada -0.1073876 -0.07910698
## 3    Aguadilla  0.2435122  0.27447357
## 4 Aguas Buenas -0.3531467 -0.38623186
## 5     Aibonito -0.4041279 -0.36412356
## 6       Añasco -0.3420229 -0.35388462

Ya tenemos la base de datos importada y escalada, lo que significa que cada variable ahora tendrá una media cero y una desviación estándar uno.

Matriz de distancia

Calcularemos la matriz de distancia y la visualizaremos.

dist      <- dist(data.scaledpob,method = "euclidean")
# matriz de distancias
dist_mat  <- as.matrix(round(dist,3))
as.dist(dist_mat[1:6,1:6])

##       1     2     3     4     5
## 2 0.732                        
## 3 1.342 0.610                  
## 4 0.256 0.482 1.090            
## 5 0.229 0.504 1.114 0.068      
## 6 0.291 0.443 1.052 0.042 0.077

Visualizamos los primeros 6 datos y no se observan números claramente lejanos.

Número óptimo de clusters

Calcularemos el número óptimo de clúster. Para esto emplearemos distintos métodos, ya programados de la librería factoextra.

library(factoextra)
fviz_nbclust(data.scaledpob, FUN = hcut, method = "silhouette")

Según el método silhouette, el número óptimo de clusters es \(k=2\).

fviz_nbclust(data.scaledpob, FUN = hcut, method = "wss")

Según el método wss, el número óptimo de clusters podriamos interpretar que es \(k=3\).

fviz_nbclust(data.scaledpob[,-1], FUN = hcut, method = "gap_stat")

Según el método gap statistics, el número óptimo de clusters es \(k=1\).

Tomaremos 2 clusters de acuerdo al resultado del método de silhouette, y porque las otras 2 metodologías llegan a resultados cercanos.

Visualización de dendogramas

Ahora calcularemos el modelo

modelo2 <- hclust(dist, method = "ward.D")

Procedemos a visualizar los dendogramas

fviz_dend(modelo2, cex =.8, k = 2,
                 rect = TRUE,
                 k_colors = "jco",
                 rect_border = "jco",
                 rect_fill = TRUE,
                 type = "circular")

fviz_dend(modelo2, cex = 0.8, k = 2,
                  rect = TRUE,
                  k_colors = "jco",
                  rect_border = "jco",
                  rect_fill = TRUE,
                  type = "phylogenic")

Podemos observar que al realizar 2 clusters de la población con el metodo ward.D, se estratifica o dividen los datos de una forma más clara. Aunque no llega a ser un resultado satisfactoio.

Nueva columna

Agregamos la nueva columna de clasificación a la base de datos.

grupos<-cutree(modelo2, k=2) 
table(grupos)

## grupos
##  1  2 
## 67 11

library(DT)
grupos_pueblos <- cbind(datos,grupos)
datatable(grupos_pueblos, options = list(order = list(3, "desc")), class = "display")

Podemos notar que el grupo 2 se compone de los 11 pueblos con mayor cantidad de habitantes en Puerto Rico.

Parte 1: Utilizando res.nbclust

Presionar

Aunque utilizando la función dada en los “Métodos de particionamiento” se vería así:

Número óptimo de clusters

library(NbClust)
set.seed(1234)
res.nbclust <- NbClust(data.scaledpob[-1], distance = "euclidean",
                  min.nc = 2, max.nc = 10, 
                  method = "complete", index ="all")

## *** : The Hubert index is a graphical method of determining the number of clusters.
##                 In the plot of Hubert index, we seek a significant knee that corresponds to a 
##                 significant increase of the value of the measure i.e the significant peak in Hubert
##                 index second differences plot. 
##

## *** : The D index is a graphical method of determining the number of clusters. 
##                 In the plot of D index, we seek a significant knee (the significant peak in Dindex
##                 second differences plot) that corresponds to a significant increase of the value of
##                 the measure. 
##  
## ******************************************************************* 
## * Among all indices:                                                
## * 5 proposed 2 as the best number of clusters 
## * 10 proposed 3 as the best number of clusters 
## * 4 proposed 6 as the best number of clusters 
## * 1 proposed 8 as the best number of clusters 
## * 4 proposed 10 as the best number of clusters 
## 
##                    ***** Conclusion *****                            
##  
## * According to the majority rule, the best number of clusters is  3 
##  
##  
## *******************************************************************

Tomando la recomendación, estaríamos utilizando 3 clusters.

Visualización de dendogramas

Ahora calcularemos el modelo

modelo2 <- hclust(dist, method = "ward.D")

Procedemos a visualizar los dendogramas

fviz_dend(modelo2, cex = 0.3, k=3, 
          rect = TRUE,  
          k_colors = "jco",
          rect_border = "jco", 
          rect_fill = TRUE, 
          horiz = TRUE)

fviz_dend(modelo2, cex =.8, k = 3,
                 rect = TRUE,
                 k_colors = "jco",
                 rect_border = "jco",
                 rect_fill = TRUE,
                 type = "circular")

fviz_dend(modelo2, cex = 0.8, k = 3,
                  rect = TRUE,
                  k_colors = "jco",
                  rect_border = "jco",
                  rect_fill = TRUE,
                  type = "phylogenic")

Podemos observar que al realizar 3 clusters de la población con el metodo ward.D, se estratifica o dividen los datos de una forma más clara. Resulta en un dendograma filogénico con mayor sentido.

Nueva columna

Agregamos la nueva columna de clasificación a la base de datos.

grupos<-cutree(modelo2, k=3) 
table(grupos)

## grupos
##  1  2  3 
## 33 34 11

library(DT)
grupos_pueblos <- cbind(datos,grupos)
datatable(grupos_pueblos, options = list(order = list(3, "desc")), class = "display")

Podemos notar que como en el ejercicio realizado con solo 2 clusters, el grupo 3 se compone de los 11 pueblos con mayor cantidad de habitantes en Puerto Rico. Esta partición divide el grupo 1 del ejercicio anterior en 2 de acuerdo a la población. Lo que podemos traducir a

Grupo 3: Pueblos Altamente Poblados
Grupo 2: Pueblos Moderadamente Poblados
Grupo 1: Pueblos Escasamente Poblados

Parte 2

Aplique lo aprendido sobre métodos de particionamiento a la base de datos Kanga, de la librería Faraway. Esta base de datos contiene diferentes medidas de los cráneos de 148 ejemplares de canguros de 3 especies. Considera solo las medidas asociadas a su mandíbula (columnas 17 a 19), y como hay valores NA en estas variables, usa solo las filas que no contengan ningún NA.

Base de datos

Para esto importaremos la base base de datos denominada Kanga, de la librería Faraway. Luego procederemos a eliminar las filas que contienen NA en las columnas 17:19 de Kanga, por último vamos a escalar los datos.

library(faraway)
library(dplyr)
nrow(kanga)

## [1] 148

data.scaledkan <- na.omit(cbind(kanga[1:2],
                             scale(kanga[17:19],
                                   center = TRUE,
                                   scale = TRUE)))
head(data.scaledkan)

##     species  sex mandible.length mandible.width mandible.depth
## 1 giganteus Male      -1.1098270     -0.5768333    -0.74474030
## 2 giganteus Male      -0.6138240      0.6530926    -0.65684432
## 3 giganteus Male      -0.7998251     -1.6620620    -1.18422023
## 4 giganteus Male      -1.0822713     -0.5044847    -0.30526038
## 5 giganteus Male      -0.4967122     -0.5768333     0.04632356
## 6 giganteus Male      -2.3842792     -2.7472907    -2.54660799

nrow(data.scaledkan)

## [1] 136

Verificamos cuántas filas contiene la base de datos kanga 148 y tras omitir los NA 136.

Número óptimo de particiones

NbClust

Trabajaremos con la función NbClust() que prueba 30 índices diferentes para encontrar el número óptimo de particiones

library(NbClust)
set.seed(1234)
res.nbclustkan <- NbClust(data.scaledkan[3:5], distance = "euclidean",
                  min.nc = 2, max.nc = 10, 
                  method = "complete", index ="all")

## *** : The Hubert index is a graphical method of determining the number of clusters.
##                 In the plot of Hubert index, we seek a significant knee that corresponds to a 
##                 significant increase of the value of the measure i.e the significant peak in Hubert
##                 index second differences plot. 
##

## *** : The D index is a graphical method of determining the number of clusters. 
##                 In the plot of D index, we seek a significant knee (the significant peak in Dindex
##                 second differences plot) that corresponds to a significant increase of the value of
##                 the measure. 
##  
## ******************************************************************* 
## * Among all indices:                                                
## * 3 proposed 2 as the best number of clusters 
## * 7 proposed 3 as the best number of clusters 
## * 5 proposed 4 as the best number of clusters 
## * 1 proposed 5 as the best number of clusters 
## * 5 proposed 8 as the best number of clusters 
## * 2 proposed 9 as the best number of clusters 
## 
##                    ***** Conclusion *****                            
##  
## * According to the majority rule, the best number of clusters is  3 
##  
##  
## *******************************************************************

Tomando la recomendación, estaremos utilizando 3 clusters.

dist      <- dist(data.scaledkan[3:5],method = "euclidean")
modelo    <- hclust(dist, method = "complete")

fviz_dend(modelo, cex = 0.5, k=3, 
          rect = TRUE,  
          k_colors = "jco",
          rect_border = "jco", 
          rect_fill = TRUE,
          horiz = TRUE,
          ggtheme = theme_bw())

clValid

library(clValid)
validcluskan  <- clValid(data.scaledkan[3:5], nClust = 2:5, 
                clMethods = c("hierarchical","kmeans","diana","fanny","pam","clara","agnes"),
                validation = "internal")
summary(validcluskan)

## 
## Clustering Methods:
##  hierarchical kmeans diana fanny pam clara agnes 
## 
## Cluster sizes:
##  2 3 4 5 
## 
## Validation Measures:
##                                  2       3       4       5
##                                                           
## hierarchical Connectivity   4.2782 11.2738 14.2028 23.5849
##              Dunn           0.1635  0.1363  0.1363  0.1082
##              Silhouette     0.5384  0.4287  0.3925  0.3654
## kmeans       Connectivity  21.8536 23.5095 35.5675 51.6306
##              Dunn           0.0433  0.0512  0.0726  0.0860
##              Silhouette     0.4208  0.3917  0.3823  0.3342
## diana        Connectivity  20.5925 29.4718 42.2274 48.0603
##              Dunn           0.0623  0.0685  0.0984  0.1105
##              Silhouette     0.4215  0.3736  0.3440  0.3319
## fanny        Connectivity  22.0718 40.5500 47.7984 59.0040
##              Dunn           0.0438  0.0360  0.0674  0.0516
##              Silhouette     0.4207  0.3335  0.2611  0.2297
## pam          Connectivity  26.4821 23.8437 41.3623 43.0528
##              Dunn           0.0390  0.0518  0.0531  0.0678
##              Silhouette     0.4153  0.3932  0.2980  0.3341
## clara        Connectivity  22.8349 23.2115 48.3885 50.9563
##              Dunn           0.0270  0.0568  0.0445  0.0632
##              Silhouette     0.4191  0.3751  0.2896  0.3119
## agnes        Connectivity   4.2782 11.2738 14.2028 23.5849
##              Dunn           0.1635  0.1363  0.1363  0.1082
##              Silhouette     0.5384  0.4287  0.3925  0.3654
## 
## Optimal Scores:
## 
##              Score  Method       Clusters
## Connectivity 4.2782 hierarchical 2       
## Dunn         0.1635 hierarchical 2       
## Silhouette   0.5384 hierarchical 2

Tomando la recomendación, estaremos utilizando 2 clusters.

dist      <- dist(data.scaledkan[3:5],method = "euclidean")
modelo    <- hclust(dist, method = "complete")

fviz_dend(modelo, 
          cex = 0.5, 
          k=2,
          k_colors = "RdBu",
          rect_border = "RdBu",
          horiz = TRUE,
          ggtheme = theme_minimal())

Conclusión

Visualmente tiene más sentido el agrupar la información en 3 clusters sin embargo, según la comparación realizada con clValid podemos determinar que lo correcto es utilizar solo 2 clusters.