Al momento de escribir esto, en México estamos muy próximos a las elecciones para elegir Presidente del país. Entre los muchos temas de interés para la ciudadanía que forman parte de la agenda de los candidatos a la presidencia, uno muy importante es la seguridad.
Datos recientes revelan que hasta 76% de los mexicanos se sienten inseguros en sus ciudades de residencia,
lo cual es alarmante. Seguramente una razón por la que los ciudadanos
se sienten más inseguros es por el incremento en el número de crímenes
cometidos en años recientes, en particular los más violentos.
Lo anterior me llevó a preguntarme ¿dónde se encuentran situados los crímenes violentos en México? ¿Se encuentran concentrados en alguna región o se distribuyen de manera similar en todo el territorio nacional?
Aprovecharemos esta curiosidad personal para revisar cómo crear mapas temáticos coropléticos, es decir, coloreados de acuerdo a un criterio cuantitativo, con R.
Exploraremos dónde han ocurrido los homicidios en México durante el año 2017.
Comencemos preparando nuestro entorno de trabajo.
Preparación del entorno de trabajo
Estos son los paquetes que necesitamos:
- tidyverse: Un conjunto de paquetes para importar, manipular, exportar y visualizar datos.
- data.table: Un paquete extremadamente eficiente para leer y manipular tablas de datos.
- tmap: Paquete para crear mapas temáticos.
- tmaptools: Funciones auxiliares para tmap.
- Hmisc: De este paquete usaremos sólo la función
cut2(), pero lo cargamos igual.
libary(tidyverse)
library(data.table)
library(tmap)
library(tmaptools)
library(Hmisc)Si no tienes algún paquete, puedes instalarlo como es usual con install.packages().
En caso de tener problemas para instalar tmap, generalmente debidos a la configuración de tus sitema operativo, puedes instalar primero rgdal y rgeos con el parámetro type = "binary" de install.packages(), en una sesión de R fuera de RStudio.
install.packages(rgdal, type = "binary")
install.packages(rgeos, type = "binary")Descarga de datos
Lo primero que necesitamos es un mapa de México con división política por municipio. En https://datos.gob.mx encontramos un paquete de mapas proporcionado por la Comisión Nacional para el Conocimiento y Uso de la Biodiversidad (CONABIO) que cumple con este requisito.
Los datos de incidencia delictiva los encontramos en la página del Secretariado Ejecutivo del Sistema Nacional de Seguridad Pública (SESNSP).
Descargamos específicamente los datos que corresponden a Cifras de Incidencia Delictiva Municipal, 2015 — abril 2018.
Además, descargamos datos de población por municipio, que nos serán de utilidad más adelante, desde el sitio del Consejo Nacional de Población (CONAPO).
Específicamente, descargamos este archivo, que nos hará la vida más fácil.
Podemos hacer lo anterior directamente desde R.
download.file("http://www.conapo.gob.mx/work/models/CONAPO/Proyecciones/Datos/Bases_de_Datos/Proyecciones_Municipios/R/baseprymunMX_r.zip", mode = "wb")Mapa de México dividido por municipio
El conjunto de mapas que hemos obtenido de CONABIO se encuentra contenido en un archivo .zip, así que usaremos unzip() para extraer su contenido.
unzip("mapa_division_politica.zip")Usaremos el archivo tipo shapefile con extensión .shp para importar el mapa que deseamos a nuestro entorno de trabajo.
Importamos el mapa con la función read_shape() de tmaptools, con el argumento as.sf = TRUE.
mex <-
read_shape("mapa_division_politica/division_politica_municipal/Municipios_2010_5A.shp", as.sf = TRUE)Vale la pena mencionar que es necesario un archivo con extensión .shx para poder leer un archivo .shp, por lo que debes asegurarte de contar con ambos. En nuestro caso este archivo es Municipios_2010_5A.shx.
Nuestro resultado es un objeto que contiene información geográfica y además una tabla de datos que podemos manipular directamente, como si se tratara de un data frame.
mex## Simple feature collection with 2456 features and 4 fields
## geometry type: MULTIPOLYGON
## dimension: XY
## bbox: xmin: 907821.8 ymin: 319149.1 xmax: 4082997 ymax: 2349615
## epsg (SRID): NA
## proj4string: +proj=lcc +lat_1=17.5 +lat_2=29.5 +lat_0=12 +lon_0=-102 +x_0=2500000 +y_0=0 +ellps=GRS80 +units=m +no_defs
## First 10 features:
## CVE_ENT CVE_MUN NOM_MUN OID
## 0 09 012 Tlalpan 1
## 1 09 013 Xochimilco 2
## 2 09 008 La Magdalena Contreras 3
## 3 09 002 Azcapotzalco 4
## 4 09 014 Benito Juárez 5
## 5 09 015 Cuauhtémoc 6
## 6 09 010 Álvaro Obregón 7
## 7 09 005 Gustavo A. Madero 8
## 8 09 004 Cuajimalpa de Morelos 9
## 9 09 016 Miguel Hidalgo 10
## geometry
## 0 MULTIPOLYGON (((2793795 814...
## 1 MULTIPOLYGON (((2800674 804...
## 2 MULTIPOLYGON (((2788802 811...
## 3 MULTIPOLYGON (((2797499 836...
## 4 MULTIPOLYGON (((2799251 820...
## 5 MULTIPOLYGON (((2795912 825...
## 6 MULTIPOLYGON (((2794396 824...
## 7 MULTIPOLYGON (((2801416 846...
## 8 MULTIPOLYGON (((2787230 825...
## 9 MULTIPOLYGON (((2796917 831...Cambiamos los nombres a minúscula para simplificar su manejo (y la verdad, para que sea más cómodo escribir código).
names(mex) <-
tolower(names(mex))Unimos la clave de entidad y la de municipio en una sola, de cinco caracteres de largo. Esto nos será conveniente en el futuro.
mex <-
mutate(mex, cve = paste0(cve_ent, cve_mun))Pasemos a los datos sobre homicidios.
Incidencia delictiva
Dado que están guardados en un archivo xlsx, el procedimiento más intuitivo sería usar el paquete readxl para leerlos. Sin embargo, todos mis intentos con este procedimiento terminaron con un crash de R y ocasionalmente de mi sistema operativo.
Por fortuna, es posible abrir el archivo xlsx en cualquier programa para manejar Hojas de Cálculo y exportarlo como un documento .csv sin mayor complicación. El resultado de esta operación será un archivo Municipal_Delitos_2015_2018_abr18.csv.
Ahora sí, leemos estos datos con la función fread() de data.table. Aprovechamos para de una vez omitir algunas columnas que no son de nuestro interés con el argumento drop y los nombres de las columnas que no conservaremos.
crimen <-
fread("Municipal_Delitos_2015_2018_abr18.csv",
drop = c("Entidad","Municipio", "Bien jurídico afectado", "Modalidad"))Los nombres de las columnas de esta tabla de datos nos pueden dar problemas más adelante por contener espacios, así que les haremos algunos cambios con tolower() y gsub().
# Los nombres originales
names(crimen)## [1] "Año" "Clave_Ent" "Cve. Municipio"
## [4] "Tipo de delito" "Subtipo de delito" "Enero"
## [7] "Febrero" "Marzo" "Abril"
## [10] "Mayo" "Junio" "Julio"
## [13] "Agosto" "Septiembre" "Octubre"
## [16] "Noviembre" "Diciembre"# Cambios
names(crimen) <-
names(crimen) %>%
tolower() %>%
gsub("\\W+", "_", .) %>%
gsub("año", "periodo", .)# Los nombres modificados
names(crimen)## [1] "periodo" "clave_ent" "cve_municipio"
## [4] "tipo_de_delito" "subtipo_de_delito" "enero"
## [7] "febrero" "marzo" "abril"
## [10] "mayo" "junio" "julio"
## [13] "agosto" "septiembre" "octubre"
## [16] "noviembre" "diciembre"Homicidios dolosos
Los datos de incidencia delictiva que hemos obtenido incluyen todo tipo de crímenes. Para este análisis, nos concentraremos en los homicidios dolosos, es decir, aquellos en los que una persona tiene la intención de dañar a otra, a diferencias los homicidios culposos, en los que la muerte ocurre por negligencia o descuido.
Haremos el filtrado usando la sintaxis de data.table. Elegiremos los homicidios dolosos del 2017.
homicidio <-
crimen[periodo == 2017 & subtipo_de_delito == "Homicidio doloso"]Se nos proporcionan los datos por mes del año. Para obtener los homicidios de todo el año, agruparemos por entidad y municipio, para obtener el total anual. Usaremos las funciones tbl_df(), group_by(), select(), summarise() y ungroup() de dplyr, y la función gather() de tidyr.
homicidio <-
homicidio %>%
tbl_df() %>%
group_by(clave_ent, cve_municipio) %>%
select(enero:diciembre) %>%
gather(key = "mes", value = "numero", enero:diciembre) %>%
summarise(homicidios = sum(numero)) %>%
ungroup()## Adding missing grouping variables: `clave_ent`, `cve_municipio`Combinamos las claves de entidad y de municipio en una sola, de cinco caracteres de largo. Este es el mismo formato que tenemos en nuestro mapa y nos servirá para unir las diferentes fuentes de datos que tenemos. También cambiamos los nombres de las columnas, para facilitar su manejo.
Para esta tarea, usamos las funciones rename(), mutate() y mutate_at() de dplyr.
homicidio <-
homicidio %>%
rename(cve_ent = clave_ent, cve_mun = cve_municipio) %>%
mutate_at(c("cve_ent", "cve_mun"), as.character) %>%
mutate(cve_ent = ifelse(nchar(cve_ent) == 1, paste0("0", cve_ent), cve_ent),
cve_mun = ifelse(nchar(cve_mun) == 4, substr(cve_mun, 2, 4), substr(cve_mun, 3, 5)),
cve = paste0(cve_ent, cve_mun))Crearemos un objeto en el que sólo tengamos el número de homicidios agrupados por la clave que acabamos de general. Aunque podríamos sobrescribir el objeto homicidio, prefiero conservarlo, en caso de que quiera volver a él.
homicidio_cve <-
homicidio %>%
group_by(cve) %>%
summarise(homicidios = sum(homicidios)) %>%
ungroup()Veamos como se distribuyen los homicidios con un histograma.
hist(homicidio_cve$homicidios)
Usemos la función cut2() de Hmisc para crear una variable discreta, h_cat, a partir de lo anterior. De esta manera tendremos dos opciones para visualizar nuestra información, el conteo y las categorías.
homicidio_cve <-
homicidio_cve %>%
mutate(h_cat = cut2(homicidios, cuts = c(0, 1, 5, 10, 50, 100)))Tenemos lo homicidios dolosos, obtengamos ahora la población.
Población por municipio
Como obtuvimos un conjunto de datos con los habitantes por municipio compatible con R, podemos usar la función load() para importar fácilmente su contenido, después de extraer el contenido del archivo que descargamos con unzip()
unzip("baseprymunMX_r.zip")
load("baseprymunMX_r.rdata")Lo anterior agregará a nuestro espacio de trabajo el objeto baseprymunMX, que contiene la población por municipio por grupo de edad. Usaremos la variable cvegeo, que es una clave de cinco caracteres para identificar por entidad y municipio cada dato, para obtener la población total estimado del 2017.
pob_cve <-
baseprymunMX %>%
filter(año == 2017) %>%
group_by(cvegeo) %>%
summarise(poblacion = sum(pob))¡Gracias a CONAPO por facilitar todo!
Generación de mapas temáticos
Hecho todo lo anterior, crear un mapa temático es relativamente sencillo. Basta añadir los datos que nos interesan al objeto que contiene el mapa.
Añadir datos a un mapa
Usamos la función append_data() de tmaptools para agregar datos a un mapa.
Lo único que necesitamos es indicar con los argumentos key.shp y key.data cuales son las columnas en común entre los datos de ambos objetos para poder hacer la unión. En nuestro caso, esta columna será la clave de cinco caracteres con el identificador de entidad y municipio.
Empezamos por anexar los datos de homicidio_cve al mapa mex.
mex <-
append_data(shp = mex, data = homicidio_cve,
key.shp = "cve",
key.data = "cve")## Under coverage: 535 out of 2456 shape features did not get appended data. Run under_coverage() to get the corresponding feature id numbers and key values.## Over coverage: 10 out of 1931 data records were not appended. Run over_coverage() to get the corresponding data row numbers and key values.Podemos ver que se ha agregado una columna a los datos del objeto mex.
mex## Simple feature collection with 2456 features and 7 fields
## geometry type: MULTIPOLYGON
## dimension: XY
## bbox: xmin: 907821.8 ymin: 319149.1 xmax: 4082997 ymax: 2349615
## epsg (SRID): NA
## proj4string: +proj=lcc +lat_1=17.5 +lat_2=29.5 +lat_0=12 +lon_0=-102 +x_0=2500000 +y_0=0 +ellps=GRS80 +units=m +no_defs
## First 10 features:
## cve_ent cve_mun nom_mun oid cve homicidios h_cat
## 1 09 012 Tlalpan 1 09012 68 [ 50, 100)
## 2 09 013 Xochimilco 2 09013 34 [ 10, 50)
## 3 09 008 La Magdalena Contreras 3 09008 13 [ 10, 50)
## 4 09 002 Azcapotzalco 4 09002 34 [ 10, 50)
## 5 09 014 Benito Juárez 5 09014 21 [ 10, 50)
## 6 09 015 Cuauhtémoc 6 09015 110 [ 100,1615]
## 7 09 010 Álvaro Obregón 7 09010 96 [ 50, 100)
## 8 09 005 Gustavo A. Madero 8 09005 196 [ 100,1615]
## 9 09 004 Cuajimalpa de Morelos 9 09004 7 [ 5, 10)
## 10 09 016 Miguel Hidalgo 10 09016 42 [ 10, 50)
## geometry
## 1 MULTIPOLYGON (((2793795 814...
## 2 MULTIPOLYGON (((2800674 804...
## 3 MULTIPOLYGON (((2788802 811...
## 4 MULTIPOLYGON (((2797499 836...
## 5 MULTIPOLYGON (((2799251 820...
## 6 MULTIPOLYGON (((2795912 825...
## 7 MULTIPOLYGON (((2794396 824...
## 8 MULTIPOLYGON (((2801416 846...
## 9 MULTIPOLYGON (((2787230 825...
## 10 MULTIPOLYGON (((2796917 831...Haremos el mismo procedimiento para los datos de habitantes por municipio.
mex <-
append_data(shp = mex, data = pob_cve,
key.shp = "cve",
key.data = "cvegeo")## Over coverage: 1 out of 2457 data records were not appended. Run over_coverage() to get the corresponding data row numbers and key values.Estamos listos para generar nuestro primer mapa temático.
Creación de un mapa temático coroplético con tmap
Generamos un mapa coroplético usando las funciones tm_shape() y tm_fill() de tmap. La primera de estas funciones indica qué mapa será usado, la segunda, cuales serán los datos que se usarán para rellenar sus formas.
Para nuestro mapa, usaremos los datos en la columna “homicidios”.
tm_shape(mex) + tm_fill(col = "homicidios")
Pues… es un mapa, sin duda, pero puede mejorar. No podemos hacer mucho para recuperar los datos perdidos en Oaxaca, marcados con gris, pero si podemos hacer más atractivo nuestro resultado.
Mejorando la presentación de nuestros mapas
Crearemos una función para producir paletas de colores graduales. En este caso, una que genere una paleta que vaya del blanco al rojo, la cual usaremos para colorear nuestros mapas.
paleta_rojo <- colorRampPalette(c("white", "red"))Crearemos el mismo mapa en la ocasión anterior, pero ahora agregaremos los argumentos breaks y palette a la función tm_fill() para controlar mejor el resultado que obtendremos. Con breaks definimos los puntos de corte en los datos y con palette los colores de relleno.
tm_shape(mex) +
tm_fill(col = "homicidios",
breaks = c(0, 1, 10, 50, 100, 500, 1000),
palette = paleta_rojo(6))
Mucho mejor! Empezamos a darnos una idea de dónde se concentran los homicidios dolosos en México. Los municipios con 0 homicidios ahora aparecen en blanco, lo cual facilita la lectura y permite distinguir mejor tendencias.
Desde luego, podemos usar la variable homicidio_cat que hemos creado para obtener un resultado un poco más sencillo para su interpretación.
tm_shape(mex) +
tm_fill(col = "h_cat",
palette = paleta_rojo(6))
Sin embargo, obtener el número de homicidios por municipio nos da un resultado sesgado. No tiene las implicaciones para la seguridad pública tener 100 homicidios en una municipio de 3000 habitantes que en uno de 100000. Necesitamos una medida proporcional.
Tasa de homicidios por 1000 habitantes
Mencionamos antes que podemos manipular los datos de un objeto de mapa creado con tmtools directamente, así que con ayuda de mutate() obtendremos la tasa de homicidios por cada 1000 habitantes de cada municipio.
mex <-
mex %>%
mutate(tasa_hom = (homicidios / (poblacion / 1000)) )Veamos cómo se distribuye el resultado de la operación anterior.
hist(mex$tasa_hom)
Nos conviene crear una variable discreta con la tasa de homicidios para una mejor visualización de la información, así que la anexamos a nuestro mapa, apoyándonos de cut2() y mutate().
mex <-
mex %>%
mutate(tasa_cat = cut2(tasa_hom, cuts = c(0, 0.001, .5, 1, 1.5, 2),
digits = 2))Veamos el mapa temático con estos datos.
tm_shape(mex) +
tm_fill(col = "tasa_cat",
palette = paleta_rojo(5))
Aun podemos mejorar la presentación de nuestro mapa usando diferentes argumentos dentro de la función tm_fill(), pero por el momento, este es un resultado satisfactorio.
Conclusiones
Creo que hemos obtenido una primera respuesta a mis preguntas iniciales sobre la ubicación de los crímenes violentos en México. Hay regiones de México con una mayor proporción de homicidios dolosos que el resto del país. Esto, desde luego, no es ninguna revelación para quienes siguen de cerca el tema del crimen y la violencia en México, pero ha resultado interesante explorarlo nosotros mismos.
Es posible analizar otros crímenes de manera similar a la que hemos trabajado con los homicidios dolosos, incluso podríamos hacer un análisis de tendencias a través del tiempo con los conjuntos de datos con los que contamos.
El procedimiento que hemos seguido es relativamente sencillo y tiene la flexibilidad suficiente para explorar este u otros temas. Pero eso, será tema para otras ocasiones.
***
Consultas, dudas, comentarios y correcciones son bienvenidas:
El código y los datos usados en este documento se encuentran en Github:
Esta entrada fue publicada originalmente en junio de 2018 en:
Deja una respuesta