An interactive map is a great tool to support spatial analysis. This
project shows how to create an interactive map using the
leaflet package (library(leaflet)) in R. The
leaflet package works with spatial data, such as
coordinates to add markers or polygons to delimit areas. Here, I present
maps that shows the distribution of veterinary services in Mexico
City.
I use data from the National Statistical Directory of Economic Units (DENUE in spanish) of the National Institute of Statistics and Geography (INEGI in spanish). This directory contains information of all the establishments in Mexico.
To download the data, you can visit the DENUE website. The data can be filter by Economic Activity, Size of Establishments and/or Geographic Area.
For this project, I focus on veterinary services in the Mexico City. I use the exact locations provided by latitude and longitude, as well as the name of establishments.
Read the data from DENUE. This file contains all establishments in Mexico City.
denue_inegi <- read_csv("denue_inegi_09_.csv", show_col_types = F, locale = locale(encoding = "LATIN1"))
head(denue_inegi[c("nom_estab", "nombre_act","entidad", "municipio", "latitud", "longitud")])## # A tibble: 6 × 6
## nom_estab nombre_act entidad municipio latitud longitud
## <chr> <chr> <chr> <chr> <dbl> <dbl>
## 1 AJOLOTARIO APATLACO Piscicultura e… Ciudad… Xochimil… 19.3 -99.1
## 2 AJOLOTARIO CARRIZAL Otra acuicultu… Ciudad… Xochimil… 19.3 -99.1
## 3 AJOLOTARIO CHINAMPA ONKALI Otra acuicultu… Ciudad… Xochimil… 19.3 -99.1
## 4 AJOLOTARIO TLAZOCAMA TONANA Otra acuicultu… Ciudad… Xochimil… 19.3 -99.1
## 5 AJOLOTERIA APANTLI Otra acuicultu… Ciudad… Xochimil… 19.3 -99.1
## 6 AMILLI TLAPALLI Piscicultura e… Ciudad… Xochimil… 19.3 -99.1Filter by veterinary services. Based on the North American Industry Classification System, Mexico (SCIAN 2023) used by INEGI, the codes for these activities are 541941, 541942, 541943, and 541944.
data_vet <- denue_inegi %>%
filter(codigo_act >= 541941 & codigo_act <= 541944)
head(data_vet[c("nom_estab", "nombre_act","entidad", "municipio", "latitud", "longitud")])## # A tibble: 6 × 6
## nom_estab nombre_act entidad municipio latitud longitud
## <chr> <chr> <chr> <chr> <dbl> <dbl>
## 1 A PATA DE PERRO Servicios vete… Ciudad… Venustia… 19.4 -99.1
## 2 ACUA MASCOTA Servicios vete… Ciudad… Iztapala… 19.3 -99.0
## 3 ACUAMASCOTAS VETERINARIA Servicios vete… Ciudad… Azcapotz… 19.5 -99.2
## 4 ADOPTANDO UN CORAZON CANINO Servicios vete… Ciudad… Miguel H… 19.4 -99.2
## 5 ALFA PETS Servicios vete… Ciudad… Gustavo … 19.5 -99.1
## 6 AMALTEA Servicios vete… Ciudad… Coyoacán 19.3 -99.1For a map by geographic area, I aggregate the data by municipality.
data_vet$vet <- 1
data_vet_mun <- data_vet %>%
group_by(municipio) %>%
summarise(total_vet = sum(vet))
head(data_vet_mun)## # A tibble: 6 × 2
## municipio total_vet
## <chr> <dbl>
## 1 Azcapotzalco 107
## 2 Benito Juárez 199
## 3 Coyoacán 170
## 4 Cuajimalpa de Morelos 40
## 5 Cuauhtémoc 130
## 6 Gustavo A. Madero 286Read the shapefile with the municipality polygons. You can download the polygons from Data Portal of CDMX.
Now, we merge the data with the polygons by municipality.
Get the centroids; this step helps with the elements on the map.
centroides <- st_centroid(mun_sf)
coords <- st_coordinates(centroides)
mun_sf$lon <- coords[,1]
mun_sf$lat <- coords[,2]Finally, the data by municipality is as follows.
## Simple feature collection with 6 features and 7 fields
## Geometry type: POLYGON
## Dimension: XY
## Bounding box: xmin: -99.36492 ymin: 19.23257 xmax: -99.09908 ymax: 19.51514
## Geodetic CRS: WGS 84
## NOMGEO CVEGEO CVE_ENT CVE_MUN total_vet
## 1 Álvaro Obregón 09010 09 010 155
## 2 Azcapotzalco 09002 09 002 107
## 3 Benito Juárez 09014 09 014 199
## 4 Coyoacán 09003 09 003 170
## 5 Cuajimalpa de Morelos 09004 09 004 40
## 6 Cuauhtémoc 09015 09 015 130
## geometry lon lat
## 1 POLYGON ((-99.18906 19.3955... -99.24683 19.33617
## 2 POLYGON ((-99.18231 19.5074... -99.18211 19.48533
## 3 POLYGON ((-99.14762 19.4040... -99.16113 19.38064
## 4 POLYGON ((-99.13427 19.3565... -99.15038 19.32667
## 5 POLYGON ((-99.25738 19.4011... -99.31094 19.32431
## 6 POLYGON ((-99.12951 19.4626... -99.14906 19.43137The first map shows the total number of veterinary services by municipality in Mexico City. Below, I explain how the code works:
leaflet package uses %>% as pipe to build the map
step by step.
addProviderTiles(providers$CartoDB.Positron adds the
basemap.
addPolygons() draws the municipality areas and sets
their appearance. popup = ~ NOMGEO creates a pop-up label
that shows the municipality name.
addCircles() adds the circles to the map. Each
circle is centered on a municipality using the previously calculated
centroid, and its radius represents the total number of veterinary
services in the municipality.
label = ~paste0(NOMGEO, ": ", total_vet) adds a label that
shows the municipality name and the total number of veterinary
services.
addControl() places the map title in the top-right
corner.
map1 <- leaflet(mun_sf) %>%
addProviderTiles(providers$CartoDB.Positron) %>%
addPolygons( color = "darkblue", weight = 2, opacity = 0.5,
fillColor = "white", fillOpacity = 0.5,
popup = ~NOMGEO) %>%
addCircles(lng = ~lon, lat = ~lat,radius = ~total_vet*10,
color = "white",fillColor = "blue",fillOpacity = 0.8,
label = ~paste0(NOMGEO, ": ", total_vet) ) %>%
addControl(html = "<div style='text-align: center; font-size: 16px;'> <b> Total of veterinary Services in Mexico City </b> </div>", position = "topright" )
map1 htmlwidgets::saveWidget(map1, "html/map_01.html",
title = "G1.Total of veterinary Services in Mexico City", selfcontained = TRUE)This map shows the distribution of all veterinary establishments in Mexico city. Below, I explain how the code works. Some functions were already used in previous maps.
makeIcon() upload a custom icon to use as marker on
the map.
addMarkers places markers on the map at the
specified coordinates (latitude and longitude). Each marker represents a
veterinary establishment.
setView() sets the initial view of the map, defining
the center coordinates and zoom level.
icon_gato <- makeIcon(iconUrl = "pets.png", iconWidth = 15, iconHeight = 15)
map2 <- leaflet(data = data_vet) %>%
addProviderTiles(providers$CartoDB.Positron) %>%
addMarkers(lng = ~longitud, lat = ~latitud, icon = icon_gato, label = ~nom_estab) %>%
addControl(html = "<div style='text-align: center; font-size: 16px;'> <b> veterinary Services in Mexico City </b> </div>", position = "topright" ) %>%
setView(lng = -99.140487, lat = 19.430072, zoom = 13)
map2 htmlwidgets::saveWidget(map2, "html/map_02.html",
title = "G2.veterinary Services in Mexico City", selfcontained = TRUE)Finally, I present a cluster map. The clusters are created based on the proximity of establishments. Note that the clusters are adjusted as you zoom in the map or click on a cluster.
In this map, I used the same functions as before, but added
clusterOptions = markerClusterOptions(), which creates the
clusters on the map.
map3 <- leaflet(data = data_vet) %>%
addProviderTiles(providers$CartoDB.Positron) %>%
addPolygons( data = mun_sf, color = "green", weight = 2, opacity = 0.5,
fillColor = "white", fillOpacity = 0.5,
popup = ~ NOMGEO) %>%
addMarkers(lng = ~longitud, lat = ~latitud, icon = icon_gato, label = ~nom_estab,
clusterOptions = markerClusterOptions(showCoverageOnHover = F)) %>%
addControl(html = "<div style='text-align: center; font-size: 16px;'> <b> veterinary Services in Mexico City </b> </div>", position = "topright" )
map3