#paquetes <- c("cluster", "dbscan", "e1071", "factoextra", "fpc", "GGally", "kernlab", "mclust", "mlbench", "scatterpie", "seriation", "tidyverse")
#install.packages(paquetes, dependencies = TRUE)
¿Para qué sirven estos paquetes?
Abrimos Tidyverse
library(tidyverse)
## ── Attaching core tidyverse packages ──────────────────────── tidyverse 2.0.0 ──
## ✔ dplyr 1.1.4 ✔ readr 2.1.4
## ✔ forcats 1.0.0 ✔ stringr 1.5.0
## ✔ ggplot2 3.4.2 ✔ tibble 3.2.1
## ✔ lubridate 1.9.2 ✔ tidyr 1.3.0
## ✔ purrr 1.0.1
## ── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
## ✖ dplyr::filter() masks stats::filter()
## ✖ dplyr::lag() masks stats::lag()
## ℹ Use the ]8;;http://conflicted.r-lib.org/conflicted package]8;; to force all conflicts to become errors
Abrimos la base de datos de la población penitenciaria
library(readxl)
# Lee la hoja por defecto del archivo Excel
datos <- read_excel("Denuncias.xlsx")
summary (datos)
## Ubigeo Hecho Departamento Provincia Distrito
## Min. : 10101 Length:1435 Length:1435 Length:1435
## 1st Qu.: 50621 Class :character Class :character Class :character
## Median :110210 Mode :character Mode :character Mode :character
## Mean :112629
## 3rd Qu.:151030
## Max. :250401
##
## ESTAFA EXTORSION HURTO ROBO
## Min. : 1.00 Min. : 1.0 Min. : 1.0 Min. : 1.00
## 1st Qu.: 2.00 1st Qu.: 1.0 1st Qu.: 2.0 1st Qu.: 2.00
## Median : 5.00 Median : 4.5 Median : 8.0 Median : 4.00
## Mean : 40.31 Mean : 31.2 Mean : 143.4 Mean : 95.77
## 3rd Qu.: 20.00 3rd Qu.: 19.0 3rd Qu.: 42.0 3rd Qu.: 24.00
## Max. :1278.00 Max. :1092.0 Max. :9260.0 Max. :4428.00
## NA's :646 NA's :869 NA's :123 NA's :389
## HOMICIDIO
## Min. : 1.000
## 1st Qu.: 1.000
## Median : 2.000
## Mean : 5.265
## 3rd Qu.: 4.000
## Max. :131.000
## NA's :1050
head (datos)
## # A tibble: 6 × 9
## `Ubigeo Hecho` Departamento Provincia Distrito ESTAFA EXTORSION HURTO ROBO
## <dbl> <chr> <chr> <chr> <dbl> <dbl> <dbl> <dbl>
## 1 10101 AMAZONAS CHACHAPOYAS CHACHAPO… 127 34 172 31
## 2 10103 AMAZONAS CHACHAPOYAS BALSAS NA 1 3 1
## 3 10105 AMAZONAS CHACHAPOYAS CHILIQUIN NA NA NA 1
## 4 10106 AMAZONAS CHACHAPOYAS CHUQUIBA… NA NA 2 NA
## 5 10109 AMAZONAS CHACHAPOYAS LA JALCA 3 NA 1 1
## 6 10110 AMAZONAS CHACHAPOYAS LEIMEBAM… 1 NA 6 NA
## # ℹ 1 more variable: HOMICIDIO <dbl>
tail (datos)
## # A tibble: 6 × 9
## `Ubigeo Hecho` Departamento Provincia Distrito ESTAFA EXTORSION HURTO ROBO
## <dbl> <chr> <chr> <chr> <dbl> <dbl> <dbl> <dbl>
## 1 250301 UCAYALI PADRE ABAD PADRE ABAD 16 7 85 14
## 2 250302 UCAYALI PADRE ABAD IRAZOLA 6 2 27 11
## 3 250303 UCAYALI PADRE ABAD CURIMANA 1 NA 27 3
## 4 250304 UCAYALI PADRE ABAD NESHUYA 3 4 36 7
## 5 250305 UCAYALI PADRE ABAD ALEXANDER… 2 NA 12 4
## 6 250401 UCAYALI PURÚS PURÚS 2 NA 9 NA
## # ℹ 1 more variable: HOMICIDIO <dbl>
library(dplyr)
# Eliminar filas con NA en cualquier columna
data_limpia <- datos %>% drop_na()
# Seleccionamos las variables relevantes (EDAD y Nº DE INGRESOS)
variables_numericas <- data_limpia[, c("ESTAFA", "EXTORSION", "HURTO", "ROBO", "HOMICIDIO")]
# Escalamos las variables
variables_escaladas <- scale(variables_numericas)
# Aplicamos k-means con 4 conglomerados y 10 inicios aleatorios
set.seed(123)
km <- kmeans(variables_numericas, centers = 3, nstart = 10)
# Mostramos los resultados
km
## K-means clustering with 3 clusters of sizes 13, 215, 43
##
## Cluster means:
## ESTAFA EXTORSION HURTO ROBO HOMICIDIO
## 1 724.76923 405.00000 4830.3077 2800.23077 47.61538
## 2 38.48837 29.20465 187.8791 90.91163 3.60000
## 3 241.46512 110.81395 1513.4884 868.93023 11.32558
##
## Clustering vector:
## [1] 2 2 2 2 2 2 2 2 3 2 2 2 2 2 2 2 3 2 2 3 2 2 2 2 3 2 2 2 2 3 2 2 2 2 2 3 2
## [38] 2 2 2 2 3 2 2 2 2 2 3 2 2 1 2 2 2 2 3 2 2 2 3 2 2 2 2 2 2 2 2 2 2 2 2 2 2
## [75] 2 2 3 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 3 3 2 2 2 2 2 2 2 2 2 2 2 2 2 1 2 2
## [112] 2 3 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 1 3 3 2 2 2 2 2 2 2 2 2 2 2
## [149] 2 2 2 1 2 1 3 3 3 2 3 2 1 3 3 3 3 1 3 1 3 2 2 1 2 2 3 3 3 3 1 3 1 3 3 2 1
## [186] 3 3 3 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2
## [223] 1 3 2 2 2 2 3 2 2 2 2 2 2 2 3 2 2 2 3 2 3 2 2 2 2 2 2 2 2 2 2 2 2 3 2 2 2
## [260] 2 2 2 2 2 2 2 3 2 2 2 2
##
## Within cluster sum of squares by cluster:
## [1] 58241224 11614616 23279782
## (between_SS / total_SS = 81.8 %)
##
## Available components:
##
## [1] "cluster" "centers" "totss" "withinss" "tot.withinss"
## [6] "betweenss" "size" "iter" "ifault"
# Creamos un data frame escalado con la asignación de clusters
clustered_data <- as.data.frame(variables_numericas) %>%
add_column(cluster = factor(km$cluster))
# Visualizamos las primeras filas
head(clustered_data)
## ESTAFA EXTORSION HURTO ROBO HOMICIDIO cluster
## 1 127 34 172 31 2 2
## 2 2 1 8 5 4 2
## 3 86 16 552 203 9 2
## 4 17 1 69 14 2 2
## 5 1 1 8 1 1 2
## 6 43 18 189 56 1 2
library(ggplot2)
# Visualizamos los clusters
ggplot(clustered_data, aes(x = HURTO, y = ROBO, color = cluster)) +
geom_point(size = 2) +
labs(title = "Distribución de Conglomerados k-means", x = "Hurto", y = "Robo") +
theme_minimal()
# Visualizamos con 8 clusters
library(factoextra)
## Welcome! Want to learn more? See two factoextra-related books at https://goo.gl/ve3WBa
fviz_cluster(kmeans(variables_numericas, centers = 3, nstart = 10),
data = variables_numericas, centroids = TRUE, geom = "point", ellipse.type = "norm") +
labs(title = "Distribución con 4 Clusters")
#3.2. Análisis Jerárquico
# Seleccionamos una muestra de 10,000 registros (ajustar según la memoria disponible)
set.seed(123) # Semilla para reproducibilidad
# Creamos la matriz de distancias y continuar con el análisis
d <- dist(variables_escaladas)
# Agrupación jerárquica aglomerativa con enlace completo
hc <- hclust(d, method = "complete")
# Visualizamos el dendrograma
plot(hc, main = "Dendrograma - Método Complete", xlab = "", sub = "", cex = 0.8)
library(factoextra)
# Visualizamos el dendrograma cortado en 4 conglomerados
fviz_dend(hc, k = 3, rect = TRUE, show_labels = FALSE, cex = 0.6)
## Warning: The `<scale>` argument of `guides()` cannot be `FALSE`. Use "none" instead as
## of ggplot2 3.3.4.
## ℹ The deprecated feature was likely used in the factoextra package.
## Please report the issue at <]8;;https://github.com/kassambara/factoextra/issueshttps://github.com/kassambara/factoextra/issues]8;;>.
## This warning is displayed once every 8 hours.
## Call `lifecycle::last_lifecycle_warnings()` to see where this warning was
## generated.
# Extraemos los clusters
clusters <- cutree(hc, k = 3)
# Creamos un data frame con las asignaciones
cluster_complete <- as.data.frame(variables_escaladas) %>%
add_column(cluster = factor(clusters))
# Vemos los primeros registros con clusters asignados
head(cluster_complete)
## ESTAFA EXTORSION HURTO ROBO HOMICIDIO cluster
## 1 0.13158691 -0.2052215 -0.38528446 -0.4494547 -0.3489628 1
## 2 -0.57182267 -0.4638960 -0.52603857 -0.4867504 -0.2076042 1
## 3 -0.09913143 -0.3463167 -0.05914688 -0.2027289 0.1457924 1
## 4 -0.48741352 -0.4638960 -0.47368490 -0.4738404 -0.3489628 1
## 5 -0.57744995 -0.4638960 -0.52603857 -0.4924883 -0.4196421 1
## 6 -0.34110433 -0.3306394 -0.37069409 -0.4135934 -0.4196421 1
# Visualización de los conglomerados
ggplot(cluster_complete, aes(x = HURTO, y = ROBO, color = cluster)) +
geom_point(size = 2) +
labs(title = "Conglomerados jerárquicos (3 clusters)", x = "Hurto (escalada)", y = "Robo (escalado)") +
theme_minimal()
# Visualización con 3 clusters
fviz_cluster(list(data = variables_escaladas, cluster = cutree(hc, k = 3)), geom = "point", ellipse.type = "convex") +
labs(title = "Distribución con 3 Clusters")
# Agrupación jerárquica con enlace simple
hc_single <- hclust(d, method = "single")
# Visualizamos el dendrograma
fviz_dend(hc_single, k = 3, rect = TRUE, show_labels = FALSE, cex = 0.6)
# Visualización de clusters con enlace simple
fviz_cluster(list(data = variables_escaladas, cluster = cutree(hc_single, k = 3)), geom = "point") +
labs(title = "Clusters - Método Single")
#3.3. DBSCAN
DBSCAN son las siglas de Density-Based Spatial Clustering of Applications with Noise (agrupación espacial de aplicaciones con ruido basada en la densidad). Agrupa los puntos que están muy juntos y trata los puntos en regiones de baja densidad como valores atípicos.
library(dbscan)
##
## Attaching package: 'dbscan'
## The following object is masked from 'package:stats':
##
## as.dendrogram
Realizamos el gráfico llamado “distancias kNN”. El eje x del gráfico representa los puntos de datos ordenados por su distancia al k-ésimo vecino más cercano (kNN), y el eje y muestra la distancia al k-ésimo vecino más cercano. En el contexto de DBSCAN, este tipo de gráfico puede ser útil para identificar un valor de epsilon (ε), que es un parámetro clave en DBSCAN que define la distancia máxima entre dos puntos para que uno sea considerado vecino del otro. Con el “abline” agregamos una línea horizontal al gráfico en la posición y = 0.32. Este valor está relacionado con la distancia epsilon (ε) que estamos considerando como umbral para definir vecinos. La línea roja podría ayudarte a visualizar y seleccionar un valor de epsilon que sea adecuado para tu análisis.
# Crear el gráfico de distancias kNN
kNNdistplot(variables_escaladas, k = 3) # k = minPts - 1
abline(h = 0.32, col = "red") # Valor sugerido para epsilon
Usamos los siguientes parámetros:
Valor de ε (epsilon): 0.32 Número mínimo de puntos (minPts): 4
minPts define cuántos puntos de la vecindad épsilon son necesarios para que un punto sea un punto central. A menudo se elige como parámetro el suavizado. Aquí utilizamos minPts = 4.
Para decidir el épsilon, se utiliza el gráfico de distancia kNN. Si la distancia entre dos puntos es menor o igual a 0.32, entonces estos puntos se consideran vecinos cercanos y se agruparán juntos en el mismo clúster si ambos tienen suficientes vecinos (determinado por el parámetro minPts). Si la distancia entre dos puntos es mayor que 0.32, entonces estos puntos no se consideran vecinos cercanos y pueden estar en diferentes clústeres, o si no tienen suficientes vecinos dentro de este radio, se considerarán puntos de ruido (outliers)
# Ejecutar DBSCAN
db <- dbscan(variables_escaladas, eps = 0.32, minPts = 4)
# Ver los resultados
db
## DBSCAN clustering for 271 objects.
## Parameters: eps = 0.32, minPts = 4
## Using euclidean distances and borderpoints = TRUE
## The clustering contains 1 cluster(s) and 77 noise points.
##
## 0 1
## 77 194
##
## Available fields: cluster, eps, minPts, metric, borderPoints
Vemos la estructura de lo obtenido
# Agregamos los clusters a las variables escaladas
clustered_data <- as.data.frame(variables_escaladas) %>%
add_column(cluster = factor(db$cluster))
# Visualización básica de los clusters
ggplot(clustered_data, aes(x = HURTO, y = ROBO, color = cluster)) +
geom_point() +
labs(title = "Resultados DBSCAN", x = "Hurto (escalada)", y = "Robo (escalado)") +
theme_minimal()
library(factoextra)
# Visualizar los clusters con factoextra
fviz_cluster(db, variables_escaladas, geom = "point") +
labs(title = "Visualización de Clusters DBSCAN")
# Seleccionamos una muestra aleatoria de 10,000 registros (ajustar según la memoria disponible)
set.seed(123)
#datos_muestra <- datos[sample(1:nrow(datos), 10000), ]
# Preparamos las variables mixtas (numéricas y categóricas)
data_limpia$Departamento <- as.factor(data_limpia$Departamento)
# Calculamos la matriz de distancias Gower
library(cluster)
d <- daisy(data_limpia[, c("ESTAFA", "EXTORSION", "HURTO", "ROBO", "HOMICIDIO", "Departamento")], metric = "gower")
# Aplicamos PAM
p <- pam(d, k = 3)
# Revisamos los resultados
p
## Medoids:
## ID
## [1,] 197 197
## [2,] 116 116
## [3,] 176 176
## Clustering vector:
## [1] 1 1 2 1 1 1 1 1 3 1 1 1 2 1 1 1 3 2 2 3 1 1 2 2 2 2 1 2 1 2 1 1 1 1 1 2 1
## [38] 1 1 1 1 2 1 1 1 1 1 2 1 1 3 2 2 3 1 2 1 2 2 2 1 1 1 1 1 1 1 1 2 1 1 2 1 1
## [75] 1 1 3 1 1 1 1 1 2 1 1 1 1 1 2 1 1 1 1 3 3 1 1 1 1 1 1 1 1 1 1 1 1 1 2 2 2
## [112] 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 3 3 2 1 2 1 1 1 1 1 1 2 1 1
## [149] 1 1 1 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3
## [186] 3 3 3 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 2 1 1 1 1 2 2 1 2 1 1 1 1
## [223] 3 3 1 1 1 1 2 1 2 1 1 1 1 2 3 1 2 1 3 1 3 1 2 1 1 1 1 1 2 1 1 1 1 2 1 1 2
## [260] 2 1 1 1 1 1 1 3 2 2 1 1
## Objective function:
## build swap
## 0.1605865 0.1520545
##
## Available components:
## [1] "medoids" "id.med" "clustering" "objective" "isolation"
## [6] "clusinfo" "silinfo" "diss" "call"
# Añadimos los clusters al subconjunto
data_limpia$cluster <- factor(p$clustering)
# Vemos las primeras filas del subconjunto con los clusters
head(data_limpia)
## # A tibble: 6 × 10
## `Ubigeo Hecho` Departamento Provincia Distrito ESTAFA EXTORSION HURTO ROBO
## <dbl> <fct> <chr> <chr> <dbl> <dbl> <dbl> <dbl>
## 1 10101 AMAZONAS CHACHAPOYAS CHACHAPO… 127 34 172 31
## 2 10203 AMAZONAS BAGUA COPALLIN 2 1 8 5
## 3 10701 AMAZONAS UTCUBAMBA BAGUA GR… 86 16 552 203
## 4 10702 AMAZONAS UTCUBAMBA CAJARURO 17 1 69 14
## 5 20507 ÁNCASH BOLOGNESI COLQUIOC 1 1 8 1
## 6 20801 ÁNCASH CASMA CASMA 43 18 189 56
## # ℹ 2 more variables: HOMICIDIO <dbl>, cluster <fct>
library(ggplot2)
# Visualizamos en un gráfico de dispersión
ggplot(data_limpia, aes(x = HURTO, y = ROBO, color = cluster)) +
geom_point(size = 2) +
labs(title = "Clusters PAM con Variables Mixtas", x = "Hurto", y = "Robo") +
theme_minimal()