Si quieres sáltate el prologo y la introducción para ya ir al grano con el tutorial 😛

Prologo

Este tutorial lo realizo ya que hice un proyecto en el cual no encontré la información de forma fácil por internet para desarrollarlo y tuve que leer los manuales de cada cosa y eso lleva mucho tiempo, así que si en algún momento tengo que volver a hacer lo mismo ya tendré la guía en estas lineas, además, tal vez a alguno le sirva 🙂

Mis conocimientos en R aun son básicos pero suficientes para terminar el trabajo y personalizarlo a mi gusto, por lo mismo no me detendré en cada detalle a explicar, ya con tiempo podrá cada quien profundizar en los temas si así lo desea.

Bueno, la situación esta así, se me ah pedido desarrollar un “simple” mapa de color, mas propio un mapa coroplético, mi intención es hacerlo en R, tengo entendido que otros programas como QGIS, están desarrollados específicamente para estos casos, pero al no conocerlo probablemente me costaría más.

El detalle es que para ese momento de mapas en R no sabia nada, de hecho, de mapas en general. >.<

NOTA ACLARATORIA LOS DATOS SON FALSOS

Introducción

Una de las grandes ventajas de R, es la facilidad con la que se pueden crear gráficas alto nivel sin complicarse tanto, la manipulación de grandes bases de datos y su capacidad para crear simulaciones(este ultimo aun no lo manejo bien, pero lo comento).

De entra este texto esta dirigido para aquellos con conocimientos básicos de R pero intentare ser los mas explicito que pueda.

Además se utilizaran las librerías de highcharter las cuales son las necesarias para desarrollar gráficas dinámicas en R, pues es un API que generara los resultados en formato html y javascript para su visualización fácil en cualquier navegador sin la necesidad de instalar ninguna otra cosa.

La librería de tidyverse que sera la encargada de la manipulación de la información de forma sencilla, muy útil para cualquier cosa, básica para cualquier trabajo que se valla a realizar.

La librería rgdal necesaria para el manejo de los archivos shp(los mapas pué).

Y por ultimo la librería Gdal para el comando ogr2ogr con el cual, podremos realizar la transformación de nuestros archivos shp a JSON

Problema.

Mapa de colores de los casos sospechosos de Dengue en la Región VI “Selva” del estado de Chiapas, México, del año 2019.

Problemas a enfrentar:

• Búsqueda de los archivos shp para el mapa de México con su división política por municipios, ya que Highcharter solo maneja los mapas que se encuentra en su repositorio, los cuales incluye el mapa mundial, con sus divisiones por países y estados, pero hasta ahí, no viene dividido por municipios.
• Manipulación de los archivos shp para seleccionar solo el área geográfica deseada.
• Conversión de los archivos shp a JSON (que sera el formato a utilizar por la librería highcharter).
• Manipulación de la información en R
• Creación de la gráfica y mapa dinámico.

La búsqueda de los shp necesario dependerá del país que se busque, en México se cuenta con el INEGI en el cual se puede encontrar varios tipos de mapas e información abierta útil, cabe aclarar que no encontré el mapa que utilice ahí(igual no supe buscarlo bien), de hecho me costo “tanto” que perdí la fuente entre tantas pestañas que abrí en su momento 🙁 .

Comencemos con manipular el archivo shp para solo dejar en un nuevo archivo los municipios con los que trabajaremos al final, el código esta comentado paso a paso, pero al final haré una explicación para complementar.

###Libreiras a usar #####
library(rgdal) #readOGR para leer archivos shp
library(tidyverse) # para manipulacion de las bases de datos
################### MAPA GENERAL####
dirmapas <- "/home/zorro/Documentos/TEMP/mapas"    #selecionar rirectorio
setwd(dirmapas)                                    #cambiar directorio de trabajo
mexico <- readOGR("muni_2016cw.shp",layer="muni_2016cw") #abrir los archivos shp para manipularlos, creando un objeto S4, del mapa de México, con divición estatal y municipal
##### Modificar/Seleccionar partes relevantes del mapa ####
chiapas <- subset(mexico, CVE_ENT == "07" ) # Selecionar solo el estado deseado
##### Municipios a Borrar ####
mun <- c("Suchiate","Tonalá","Altamirano","Tapalapa","La Trinitaria","Villaflores",
         "Nicolás Ruíz","Chicoasén","Cintalapa","Jiquipilas","Amatenango del Valle",
         "Tenejapa","Ocotepec","Mapastepec","Chalchihuitán","Bejucal de Ocampo","Teopisca",
         "Rayón","San Andrés Duraznal","Mazapa de Madero","Berriozábal","Chiapilla","Ixtapangajoya",
         "Pijijiapan","La Independencia","Zinacantán","Totolapa","Bella Vista","Angel Albino Corzo",
         "Copainalá","Maravilla Tenejapa","El Porvenir","Mazatán","Jitotol","Chamula","Tapilula","Aldama",
         "Ostuacán","Tuxtla Chico","Ixhuatán","Montecristo de Guerrero","Huitiupán","Acala",
         "Motozintla","Santiago el Pinar","Pantelhó","Sunuapa","Chenalhó","Amatenango de la Frontera",
         "Escuintla","San Juan Cancuc","Larráinzar","La Concordia","San Cristóbal de las Casas",
         "San Fernando","San Lucas","Simojovel","Huixtán","Unión Juárez","Mitontic","Pueblo Nuevo Solistahuacá",
         "Ixtapa","Reforma","Tapachula","Juárez","Acacoyagua","Arriaga","Las Margaritas","El Bosque",
         "Ocozocoautla de Espinosa","Francisco León","Tuxtla Gutiérrez","Las Rosas","Bochil","Comitán de Domínguez",
         "Siltepec","Socoltenango","Villa Corzo","Tecpatán","Amatán","Soyaló","Metapa","Chicomuselo","Cacahoatán","Suchiapa",
         "Ixtacomitán","Coapilla","Acapetahua","Tapilula,Aldama","Tuzantán","Pantepec","Villa Comaltitlán",
         "Frontera Comalapa","Osumacinta","Chapultenango","Frontera Hidalgo","Tzimol","La Grandeza",
         "Sitalá","Chanal","Huehuetán","Huixtla","Chiapa de Corzo","Solosuchiapa","Pichucalco","Pueblo Nuevo Solistahuacán")
####formular para borar los municipios del estado seleciona con el vector previamente hecho ####
i = 0 # Declaramos un marcador
while (i != length(mun)) {
  i <- i+1 #el marcador ira aumentado con cada vuelta del ciclo
  chiapas <- chiapas[chiapas$NOM_MUN != mun[i],] #se sobre escribe la base tipo poligonos borrando cada uno de los municipios dentro del vector
}
writeOGR(chiapas, dirmapas, "Mapa.RegionVI", driver="ESRI Shapefile") #Escribir solo el estado de Chiapas en un archivo shp

0k, como podemos darnos cuenta, la primera parte es definir las librerías a usar (tidyverse y rgdal).

Hay que comentar que previamente se tiene que instalar en nuestra computadora estas librerías, si no se ah hecho antes, esto se hace con el comando install.packages(LIBRERIA) .

Segundo, se declara la variable del directorio a usar, esto se hace de esta forma pues así puede modificarse con mas facilidad en el futuro, pero también se podría a ver utilizado directamente la ubicación en el comando, es decir setwd(“/home/zorro/Documentos/TEMP/mapas”).

Se usa readOGR , para leer el archivo shp e introducirlo en una nueva base de datos llamada “mexico”.

Con subset se indica la base de datos a usar (mexico) y el filtro a usar (CVE_ENT == “07”), este filtro significa que se tomaran todas las observaciones que concuerden con el valor de “07” en la columna con nombre “CVE_ENT” (que corresponde al código del estado de Chiapas).

Además se debe aclarar que la base de datos llamada “mexico” es un objeto tipo S4(de polígonos, entre otras cosas), y para acceder de forma visual y ordenada (para manipulación de variables) a ella es de la siguiente manera: view(mexico@data)
y para observar el contenido de una columna es: view(mexico@data$CVE_ENT)

Se declara un vector llamado “mun” el cual contendrá los nombres de los municipios a borrar de nuestro mapa, para que al final solo queden los municipios que nos interesan, esto se hace de manera manual, tomando los nombres de los municipios y escribiendo cada uno de ellos(no encontré la forma de hacer automático, de seguro hay, pero no se me ocurrió 😛 )
Ahora en un bucle tipo “while” se sobrescribe la misma base de datos (“chiapas”) pero descartando de ella un municipio a la vez, por eso el tamaño del bucle es igual al numero de municipios en el vector a borrar( while (i != length(mun)) )
y para sobre escribir la misma base de datos pero descartando el municipio que corresponda al contenido del vector en la vuelta del bucle es con el comando: chiapas <- chiapas[chiapas$NOM_MUN != mun[i],] donde “NOM_MUN” corresponde al nombre del municipio.

Y al final, escribimos(writeOGR) esta nueva base de datos “chiapas” en un archivo tipo shp(driver=ESRI Shapefile) en nuestro directorio de trabajo.

Para que al final obtengamos solamente la Región VI, pero claro, en el camino podemos guardar lo obtenido de manera temporal 😉 .

0k, continuamos ahora con la conversión del archivo shp “Mapa.RegionVI.shp” que contiene la información que obtuvimos al final del código y guardamos en nuestro directorio de trabajo.

Para esto nos haremos del programa ogr2ogr, el cual lo podemos descargar directamente con el siguiente comando (probado en Debian 10): apt install gdal-bin

ogr2ogr -f "GeoJSON" JSON/Chiapas.js Mapa.RegionVI.shp

El comando utilizado en cualquier consola de nuestro linux(si tienen otro SO, con la pena, pero seguramente habrá equivalentes para esto, de hecho me parece que hay convertidores en linea) es muy simple, lo mandamos a llamar, la opción -f para especificar el formato de salida, en este caso GeoJSON, seguido de la ubicación de nuestro NUEVO archivo .js y al final el archivo de entrada (“Mapa.RegionVI.shp”), si se dan cuenta estoy ejecutando el comando en el mismo directorio donde se encuentra el archivo shp, si no fuera el caso deberíamos poner, ya sea la ruta completa o indirecta de donde se encuentre nuestro archivo shp.

Bueno ahora el código para la creación de la base de datos a usar para el mapa y el código de como crear el mapa, de igual manera esta comentado sobre el código y posterior una explicación complementaria.

library(highcharter)
library(tidyverse)
mapa <- "/home/zorro/Documentos/TEMP/mapas/Mapa.RegionVI.js" #Direcion del mapa en formato JSON
data <- jsonlite::fromJSON(mapa, simplifyVector = FALSE) #extrar la infomarción de json en una BD
casos2 <- read.csv("/home/zorro/Documentos/TEMP/Municipios.csv", head=TRUE) #se abre el archivo que contiene los casos por municipio de la region
Casos2 <- data.frame(casos %>% group_by(MUNICIPIO)%>% tally())#Base de datos con los casos por municipio
colnames(Casos2) <- c("NOM_MUN","Casos") #Renombrar las columnas de la nueva base de datos
highchart() %>% 
  hc_title(text = "Casos de Dengue en Region VI \"Selva\"",
           align = "center", style = list(color = "#000000", fontWeight = "bold")) %>% 
  hc_subtitle(text = "Fuente: HGP - Epidemiología") %>% 
  hc_add_series_map(map=data, df=Casos2, name = "Dengue",
                    value = "Casos", joinBy = "NOM_MUN",
                    dataLabels = list(enabled = TRUE,
                                      format = '{point.NOM_MUN}'),
                    borderColor = "#0EAD82", borderWidth = 0.1,
                    tooltip = list(valueSuffix = " casos")) %>% 
  hc_colorAxis(minColor = "#F57779", maxColor = "#360001")%>%
  hc_mapNavigation(enabled = TRUE)

Bueno, de entrada, de nuevo las librerías a usar.

Se prosigue con ubicar la dirección del archivo JS en una variable(“mapa”).

Se lee un archivo tipo CSV para crear una base de datos llamada “casos2”, la cual contiene lo siguiente:

Lo cual es un listado en donde cada observación es un caso sospechoso de Dengue identificado solo por el municipio de residencia del paciente, lo que hay que hacer es una nueva tabla que agrupe las repeticiones y con una nueva columna que muestre el conteo de esas repeticiones.

Por lo cual se aplica las siguientes lineas de código para la nueva base de datos llamada “Casos2” y modificando el nombre de las columnas con el comando colnames, para al final obtener:

Esta base de datos sera la que al final utilicemos para ubicar los datos en nuestro mapa.

Para la creación del mapa es a partir del highchart(), en el cual usando la tuberías(pipes) de la librería tidyverse (%>%) se irán formado las capas del mapa, dentro de los parámetros opcionales solo comentare el siguiente, que es el que realmente hace la magia

hc_add_series_map, al cual se le pasa como parámetro map la base de datos que contiene la información JSON del mapa, seguido de df para seleccionar la base de datos a usar para graficar y con el paramento value se le indica a R que columna queremos graficar de la base de datos elegida en df, en este caso la columna “Casos” de la base de datos “Casos2”.

Se le indica que busque las coincidencias entre las columnas “NOM_NUM” con el parámetro joinBy, esto quiere decir que highcart buscara entre los datos del mapa y la base de datos las observaciones que coincidan en la columna llamada NOM_MUM y al encontrar dicha observación graficará lo encontrado en la columna “Casos” de la base de datos (que previamente se le indico en el parámetro value).

De nuevo pero en español, R busca dentro de la info del mapa y la base de datos las coincidencias de la columna correspondiente al nombre del municipio (NOM_MUN), entonces, cuando en el mapa encuentre “Palenque” y observe que en la base de datos también se encuentra “Palenque” le otorgara a este el valor correspondiente que se encuentre en la columna “Casos”, en este ejemplo 392, y así sucesivamente con todos los municipios, que como en el mapa solo están los de la región VI, prácticamente compara todos, aun así hay municipios que no tienen casos.

El resultado es:(y el código del archivo también pa descargar)

El código fuente para descargar es el siguiente: (por cuestiones del servidor le agregue la extensión .txt, así qu tienes que renombrar el archivo de “Mapa.Propio.R.txt” a “Mapa.Propio.R” para poder abrirlo desde RStudio.

Bueno eso es todo de mi parte, espero te sirva de algo mi buen amigo@.

Quis custodiet ipsos custodes?
🙂