Objetivo

Se presentan algunos ejercicios básicos de la cartografía digital con R y haciendo uso de algunas paqueterías para integrar herramientas geográficas y estadísticas. Así bien, para visualizar y modelar datos geográficos con el software de código abierto y tratar de entender principalmente el manejo de Sistemas de Información Geográfica (SIG) y visualizar datos e información que requiera de una localización específica.

Sistema de Información Geográfica

El Instituto Nacional de Estadística, Geografía e Informática (INEGI) tiene puestos a disposición del público los shapefiles con los distintos niveles de información geográfica. Estos archivos digitales se encuentran disponibles en el Marco Geoestadístico Nacional (MGN), el cual es un sistema único y se presenta la división del territorio nacional en diferentes niveles de desagregación para referir geográficamente la información estadística de los censos y encuestas. Se integra al Sistema Nacional de Información Estadística y Geográfica (SNIEG).

Este producto integra información vectorial, tablas de atributos y catálogos. Muestra la división geoestadística del territorio nacional en sucesivos niveles del territorio. Esta división está dada por los llamados límites estadísticos, que pueden coincidir con los límites político-administrativos oficiales, los cuales tienen sustento legal.

Tabla: Subdivisiones territoriales básicas del INEGI.
División	Unidades territoriales	Datos vectoriales
AGEE
Estatales	32	Polígonos
AGEM
Municipales	2 469	Polígonos
Localidades
Localidades rurales amanzanadas	45 397	Polígonos
Localidades urbanas	4 911	Polígonos
Localidades rurales amanzanadas y no amanzanadas	295 779	Puntos
AGEB
Áreas geoestadísticas básicas rurales	17 469	Polígonos
Áreas geoestadísticas básicas urbanas	63 982	Polígonos
MZA
Manzanas urbanas y rurales (incluyendo caserío disperso)	2 513 853	Polígonos
TI
Territorio insular	350	Polígonos

Especificación del sistema de coordenadas

Ocurre a menudo que las capas de información vectorial se han obtenido de fuentes con sistemas de coordenadas distintas; por lo que se debe transformar la información para homogenizarla y representarla en el proyecto.

La Norma Técnica para el Sistema Geodésico Nacional define las disposiciones con el fin de que el marco sea homogéneo, compatible y comparable; y establece que el Marco de Referencia oficial para los Estados Unidos Mexicanos es el International Terrestrial Reference Frame 2008, con datos de época 2010.0 (ITRF08 ) asociado al elipsoide de referencia definido en el GRS80.

A continuación se enlista los datos del sistema de coordenadas para la estandarización de la información:

Tabla: Características técnicas de la información
Parámetros
Cónica Conforme de Lambert (`CCL`)
Dátum:	ITRF2008
Elipsoide	GRS80
Meridiano origen	102° 00’ 00’ W.
Latitud origen	12° 00’ 00’ N.
Primer paralelo estándar	17° 30’ N
Segundo paralelo estándar	29° 30’N
Falso Este	2 500 000
Falso Norte	0

Geospatial Data Abstraction Library (`GDAL`)

En este documento se va a hacer uso del paquete rgdal el cual permite importar, exportar o bien manipular los datos especiales:

readOGR() y writeOGR(): para datos vectoriales
readGDAL() y writeGDAL(): para rasters

También el paquete gdalUtils permite reproyectar, transformar, reclasificar, etc.

La función ogrDrivers permite visualizar los tipos de formatos que el paquete rgdal permite trabajar:

name	long_name	write	copy	isVector
AeronavFAA	Aeronav FAA	FALSE	FALSE	TRUE
AmigoCloud	AmigoCloud	TRUE	FALSE	TRUE
ARCGEN	Arc/Info Generate	FALSE	FALSE	TRUE
AVCBin	Arc/Info Binary Coverage	FALSE	FALSE	TRUE
AVCE00	Arc/Info E00 (ASCII) Coverage	FALSE	FALSE	TRUE
BAG	Bathymetry Attributed Grid	TRUE	TRUE	TRUE
BNA	Atlas BNA	TRUE	FALSE	TRUE
CAD	AutoCAD Driver	FALSE	FALSE	TRUE
Carto	Carto	TRUE	FALSE	TRUE
Cloudant	Cloudant / CouchDB	TRUE	FALSE	TRUE
CouchDB	CouchDB / GeoCouch	TRUE	FALSE	TRUE
CSV	Comma Separated Value (.csv)	TRUE	FALSE	TRUE
CSW	OGC CSW (Catalog Service for the Web)	FALSE	FALSE	TRUE
DGN	Microstation DGN	TRUE	FALSE	TRUE
DXF	AutoCAD DXF	TRUE	FALSE	TRUE
EDIGEO	French EDIGEO exchange format	FALSE	FALSE	TRUE
EEDA	Earth Engine Data API	FALSE	FALSE	TRUE
Elasticsearch	Elastic Search	TRUE	FALSE	TRUE
ESRI Shapefile	ESRI Shapefile	TRUE	FALSE	TRUE
ESRIC	Esri Compact Cache	FALSE	FALSE	TRUE
ESRIJSON	ESRIJSON	FALSE	FALSE	TRUE
FITS	Flexible Image Transport System	TRUE	FALSE	TRUE
FlatGeobuf	FlatGeobuf	TRUE	FALSE	TRUE
Geoconcept	Geoconcept	TRUE	FALSE	TRUE
GeoJSON	GeoJSON	TRUE	FALSE	TRUE
GeoJSONSeq	GeoJSON Sequence	TRUE	FALSE	TRUE
Geomedia	Geomedia .mdb	FALSE	FALSE	TRUE
GeoRSS	GeoRSS	TRUE	FALSE	TRUE
GML	Geography Markup Language (GML)	TRUE	FALSE	TRUE
GPKG	GeoPackage	TRUE	TRUE	TRUE
GPSBabel	GPSBabel	TRUE	FALSE	TRUE
GPSTrackMaker	GPSTrackMaker	TRUE	FALSE	TRUE
GPX	GPX	TRUE	FALSE	TRUE
HTF	Hydrographic Transfer Vector	FALSE	FALSE	TRUE
HTTP	HTTP Fetching Wrapper	FALSE	FALSE	TRUE
Idrisi	Idrisi Vector (.vct)	FALSE	FALSE	TRUE
JML	OpenJUMP JML	TRUE	FALSE	TRUE
JP2OpenJPEG	JPEG-2000 driver based on OpenJPEG library	FALSE	TRUE	TRUE
JPEG2000	JPEG-2000 part 1 (ISO/IEC 15444-1), based on Jasper library	FALSE	TRUE	TRUE
KML	Keyhole Markup Language (KML)	TRUE	FALSE	TRUE
LVBAG	Kadaster LV BAG Extract 2.0	FALSE	FALSE	TRUE
MapInfo File	MapInfo File	TRUE	FALSE	TRUE
MapML	MapML	TRUE	FALSE	TRUE
MBTiles	MBTiles	TRUE	TRUE	TRUE
Memory	Memory	TRUE	FALSE	TRUE
MSSQLSpatial	Microsoft SQL Server Spatial Database	TRUE	FALSE	TRUE
MVT	Mapbox Vector Tiles	TRUE	FALSE	TRUE
MySQL	MySQL	TRUE	FALSE	TRUE
netCDF	Network Common Data Format	TRUE	TRUE	TRUE
NGW	NextGIS Web	TRUE	TRUE	TRUE
OAPIF	OGC API - Features	FALSE	FALSE	TRUE
ODBC	ODBC	TRUE	FALSE	TRUE
ODS	Open Document/ LibreOffice / OpenOffice Spreadsheet	TRUE	FALSE	TRUE
OGCAPI	OGCAPI	FALSE	FALSE	TRUE
OGR_GMT	GMT ASCII Vectors (.gmt)	TRUE	FALSE	TRUE
OGR_PDS	Planetary Data Systems TABLE	FALSE	FALSE	TRUE
OGR_SDTS	SDTS	FALSE	FALSE	TRUE
OGR_VRT	VRT - Virtual Datasource	FALSE	FALSE	TRUE
OpenAir	OpenAir	FALSE	FALSE	TRUE
OpenFileGDB	ESRI FileGDB	FALSE	FALSE	TRUE
OSM	OpenStreetMap XML and PBF	FALSE	FALSE	TRUE
PCIDSK	PCIDSK Database File	TRUE	FALSE	TRUE
PDF	Geospatial PDF	TRUE	TRUE	TRUE
PDS4	NASA Planetary Data System 4	TRUE	TRUE	TRUE
PGDUMP	PostgreSQL SQL dump	TRUE	FALSE	TRUE
PGeo	ESRI Personal GeoDatabase	FALSE	FALSE	TRUE
PLSCENES	Planet Labs Scenes API	FALSE	FALSE	TRUE
PostgreSQL	PostgreSQL/PostGIS	TRUE	FALSE	TRUE
REC	EPIInfo .REC	FALSE	FALSE	TRUE
S57	IHO S-57 (ENC)	TRUE	FALSE	TRUE
SEGUKOOA	SEG-P1 / UKOOA P1/90	FALSE	FALSE	TRUE
SEGY	SEG-Y	FALSE	FALSE	TRUE
Selafin	Selafin	TRUE	FALSE	TRUE
SQLite	SQLite / Spatialite	TRUE	FALSE	TRUE
SUA	Tim Newport-Peace’s Special Use Airspace Format	FALSE	FALSE	TRUE
SVG	Scalable Vector Graphics	FALSE	FALSE	TRUE
SXF	Storage and eXchange Format	FALSE	FALSE	TRUE
TIGER	U.S. Census TIGER/Line	TRUE	FALSE	TRUE
TopoJSON	TopoJSON	FALSE	FALSE	TRUE
UK .NTF	UK .NTF	FALSE	FALSE	TRUE
VDV	VDV-451/VDV-452/INTREST Data Format	TRUE	FALSE	TRUE
VFK	Czech Cadastral Exchange Data Format	FALSE	FALSE	TRUE
VICAR	MIPL VICAR file	TRUE	TRUE	TRUE
Walk	Walk	FALSE	FALSE	TRUE
WAsP	WAsP .map format	TRUE	FALSE	TRUE
WFS	OGC WFS (Web Feature Service)	FALSE	FALSE	TRUE
XLS	MS Excel format	FALSE	FALSE	TRUE
XLSX	MS Office Open XML spreadsheet	TRUE	FALSE	TRUE
XPlane	X-Plane/Flightgear aeronautical data	FALSE	FALSE	TRUE

Temas cubiertos en este documento.
Tema	Métodos	Librería y función
Importación de datos	Archivo shapefiles .shp	`rgdal::readOGR()`
Mapas coropléticos	Función para crear polígonos	`ggplot2:: geom_polygon()`
	`fortify()` + datos	`ggplot() + geom_map()`
	.shp + datos	`ggspatial::layer_spatial()`
	.shp + datos	`ggplot2::geom_sf()`
	.shp + datos	`tmap::shape()`
Etiquetado de mapas	Como texto	`ggplot2::geom_label`
	Ubicación en centroides	`sp::coordinates()`
	Etiquetado parcial	`ifelse() + ggreppel::geom_label_repel()`

Base de datos

Índice de marginación

Para conocer las desigualdades territoriales, desde 1993 la Secretaría General del Consejo Nacional de Población (CONAPO) realiza un ejercicio que permite construir el Índice de marginación con base en los resultados censales, una medida que permite identificar las zonas y regiones con más carencias dentro del país. Este índice resume y permite diferenciar los estados, municipios, localidades y áreas geoestadísticas básicas del país según el impacto global de las carencias que padece la población como resultado de la falta de acceso a la educación, la residencia en viviendas inadecuadas, la percepción de ingresos monetarios insuficientes y las relacionadas con la residencia en localidades pequeñas.

Hoy este índice se ha convertido en una de las principales herramientas analíticas y operativas para la definición y focalización de políticas públicas enfocadas al abatimiento de las carencias socioeconómicas de la población mexicana.

Para este ejercicio se hace uso de las bases de datos que contienen las estimaciones del índice de marginación en los diferentes niveles de desagregación geográfica censal. Permitiendo de esta manera hacer de manera visual los territorios de desigualdad que enfrenta el país.

La base de datos del índice de marginación en los diferentes niveles de desagregación se encuentra disponible en la página de datos abiertos.

Estructura de un mapa

Entidad Federativa

.shp

La función readOGR del paquete rgdal, extrae automáticamente la información utilizada por otros paquetes SIG de código abierto como QGIS y permite a R manejar una gama más amplia de formatos de datos espaciales. Esta función lee datos OGR y datos vectoriales, pero solamente permite manejar capas con características geométricas (no mezcla puntos, líneas o polígonos en una sola capa) y a su vez establecerá un sistema de referencia espacial si la capa tiene dichos metadatos.
Para leer un archivo shapefile, se establecen los siguientes argumentos, como dsn, en donde se indica el directorio que contiene los shapes y layer que es el nombre explícito de la capa a trabajar y dichas capas deben de ir sin la extensión .shp.

A continuación se lee el archivo .shp que contiene de manera integrada la división de el área geoestadística estatal agee.

require(rgdal)
shape_estados <- readOGR(dsn ="MGN 2020/conjunto_de_datos", 
                          layer = "00ent",
                           encoding = "UTF-8",
                            use_iconv = TRUE)

La función rename() del paquete dplyr permite cambiar el nombre de la columna de la clave geoestadística a nivel estatal dentro de la base de datos del shape.

shape_estados@data <- shape_estados@data %>%
                       rename("CVE_GEO" = "CVEGEO")

Un objeto SpatialPolygons se puede combinar con una data.frame para crear lo que se llama SpatialPolygonsDataFrame. La diferencia está en los atributos asociados con los polígonos. Un objeto SpatialPolygonsDataFrame tiene información adicional asociada a los polígonos (ejemplo, nombre de la entidad, población, etc.) mientras que SpatialPolygons contiene solo la información espacial (vértices) sobre el polígono.

class(shape_estados)
## [1] "SpatialPolygonsDataFrame"
## attr(,"package")
## [1] "sp"

Si el archivo no tuviera Sitema de Referencia de Coordenadas CRS asignado se debería utilizal la función proj4string() para asignar el CRS adecuado. En cambio, si se quisiera cambiar el sistema de referencia (transormar) debe utilizarse la función spTransform() del paquete sp.

La función proj4string()representa el sistema de referencia de coordenadas utilizado en los datos.

proj4string(shape_estados)
## [1] "+proj=lcc +lat_0=12 +lon_0=-102 +lat_1=17.5 +lat_2=29.5 +x_0=2500000 +y_0=0 +ellps=GRS80 +units=m +no_defs"

IME_2020

La base de datos del índice de marginación por estados se encuentra disponible en la página oficial de CONAPO o bien se puede consultar en la página de Datos Abiertos y se presenta en formato .xlsx Consulta.

Para poder trabajar con la base de datos del índice de marginación a nivel estatal, se elimina la fila que contiene los datos a nivel nacional con la función ::filter() del paquete dplyr y por otro lado se cambia el nombre de la columna CVE_ENT que hace referencia a la clave geográfica del estado por CVE_GEO para fines prácticos.

IME_2020 <- read_xlsx("Data/IME_2020.xlsx", sheet = "IME_2020") %>% 
             filter(NOM_ENT != "Nacional") %>% 
              mutate(CVE_GEO = CVE_ENT) %>% 
               as.data.frame()

Ahora bien, para la realización de los mapas cloropléticos que muestran áreas geográficas divididas o regiones que están coloreadas en relación con una variable de interés. Para este documento se decide hacer uso del grado de marginación que permite englobar a las unidades geográficas que presentan las mismas condiciones de exclusión social en cinco categorías.

Tabla: Entidades según el grado de marginación.
Grado de marginación	Estados
Muy alto	3
Alto	9
Medio	8
Bajo	8
Muy bajo	4

Layers

Existen dos maneras de trabajar los datos de los mapas:

Es trabajar con los objetos SpatialPolygonsDataFrame o bien spdf y tratar de unir los datos que se desean analizar, llamando algunas de las funciones de join_ o merge del paquete dplyr.
O bien utilizar la función fortify() del paquete ggplot2, que permite convertir un objeto (ej.spdf) a un data.frame. Al hacer este tipo de conversión permite aplicar un join para unir ambas estructuras de datos al nuevo formato del data.frame.

\[SpatialPolygons \Rightarrow data.frame\]

Cuando es necesario cambiar los ID’s de las características en los objetos SpatialLines o SpatialPolygons, se puede usar la función spChFIDs() de la paquetería sp. Donde los nuevos ID’s deben ser un vector de caracteres únicos de la longitud correcta.

La función fortify() convierte un objeto S3 para convertir objetos diversos a un data.frame para ggplot2. El data.frame resultante de la conversión de un objeto spdf es sumamente largo, tiene una fila por cada cambio en un polígono y un columna de grupo para separar cada polígono.

capas_estados <- shape_estados %>%
                  sp::spChFIDs(., str_pad(shape_estados@data$CVE_GEO, 2, "left", pad = "0")) %>%
                   fortify(., id = "CVE_GEO") %>%
                    left_join(., IME_2020 %>% dplyr::select(CVE_GEO, NOM_ENT, GM_2020),
                               by = c("id" = "CVE_GEO"))

\[SpatialPolygons \Rightarrow SpatialPolygons + Datos\]

layer_estados <- merge(shape_estados,
                        IME_2020 %>%
                         dplyr::select(c(-CVE_ENT, -NOM_ENT)) %>%
                          mutate(GM_2020 = fct_relevel(.$GM_2020,"Muy alto", "Alto", "Medio", "Bajo", "Muy bajo")),
                           by = "CVE_GEO")

\[Data.frame \Rightarrow SpatialPolygonsDataFrame\]

Para poder emparejar polígonos y puntos, los objetos espaciales deben tener el mismo CRS. Para después poder calcular el área de cada polígono o bien se puede reproyectar los datos espaciales usando la función spTransform:

tabla <- SpatialPointsDataFrame(data = capas_estados, 
                                 coords = capas_estados %>% 
                                             dplyr::select(long, lat), 
                                  match.ID = FALSE,
                                   proj4string = CRS("+proj=lcc +lat_0=12 +lon_0=-102 +lat_1=17.5 +lat_2=29.5 +x_0=2500000 +y_0=0 +ellps=GRS80 +units=m +no_defs"))

Municipio

.shp

A continuación se lee el archivo .shp que contiene de manera integrada la división de el área geoestadística municipal agem.

shape_municipios <- readOGR(dsn ="MGN 2020/conjunto_de_datos", 
                             layer = "00mun",
                              encoding = "UTF-8",
                               use_iconv = TRUE)

La función rename() del paquete dplyr permite cambiar el nombre de la columna de la clave geoestadística a nivel estatal dentro de la base de datos del shape.

shape_municipios@data <- shape_municipios@data %>%
                          rename("CVE_GEO" = "CVEGEO")

IMM_2020

La base de datos del índice de marginación por municipios se encuentra disponible en la página oficial de CONAPO o bien se puede consultar en la página de Datos Abiertos y se presenta en formato .xlsx Consulta.

Para poder trabajar con la base de datos del índice de marginación a nivel municipal, se elimina la fila que contiene los datos a nivel nacional con la función ::filter() del paquete dplyr y por otro lado se cambia el nombre de la columna CVE_MUN que hace referencia a la clave geográfica del municipio por CVE_GEO para fines prácticos.

IMM_2020 <- read_xlsx("Data/IMM_2020.xlsx", sheet = "IMM_2020") %>%
             filter(NOM_ENT != "Nacional") %>%
               mutate(CVE_GEO = CVE_MUN) %>%
                as.data.frame()

Tabla: Municipios según el grado de marginación.
Grado de marginación	Municipios
Muy alto	204
Alto	586
Medio	494
Bajo	530
Muy bajo	655

Layers

Existen dos maneras de trabajar los datos de los mapas:

Es trabajar con los objetos SpatialPolygonsDataFrame o bien spdf y tratar de unir los datos que se desean analizar, llamando algunas de las funciones de join_ o merge del paquete dplyr.
O bien utilizar la función fortify() del paquete ggplot2, que permite convertir un objeto (ej.spdf) a un data.frame. Al hacer este tipo de conversión permite aplicar un join para unir ambas estructuras de datos al nuevo formato del data.frame.

\[SpatialPolygons \Rightarrow data.frame\]

capas_municipio <-  shape_municipios %>%
                     sp::spChFIDs(., str_pad(shape_municipios@data$CVE_GEO, 5, "left", pad = "0")) %>%
                      fortify(., id = "CVE_GEO") %>%
                       left_join(., IMM_2020 %>% dplyr::select(CVE_GEO, CVE_ENT, NOM_ENT, NOM_MUN, GM_2020), by = c("id" = "CVE_GEO"))

\[SpatialPolygons \Rightarrow SpatialPolygons + Datos\]

layer_municipios <- merge(shape_municipios,
                           IMM_2020 %>% dplyr::select(c(-CVE_ENT, -CVE_MUN)) %>%
                                         mutate(GM_2020 = fct_relevel(.$GM_2020,"Muy alto", "Alto", "Medio", "Bajo", "Muy bajo")),
                            by = "CVE_GEO")

Zonas Metropolitanas

ZM_2015 <- read_xlsx("Data/ZM_2015.xlsx", sheet = "ZM_2015")

Base: Zonas Metropolitanas 2015.
CVE_MUN	CVE_ENT	NOM_ENT	MUN	NOM_MUN	ORDEN_ZM	NOM_ZM	CVE_ZM	SUN	TIPO_CD
01001	01	Aguascalientes	001	Aguascalientes	1	Aguascalientes	01.01	1	1
01005	01	Aguascalientes	005	Jesús María	1	Aguascalientes	01.01	1	1
01011	01	Aguascalientes	011	San Francisco de los Romo	1	Aguascalientes	01.01	1	1
02001	02	Baja California	001	Ensenada	2	Ensenada	02.01	68	1
02002	02	Baja California	002	Mexicali	3	Mexicali	02.02	3	1
02003	02	Baja California	003	Tecate	4	Tijuana	02.03	2	1
02004	02	Baja California	004	Tijuana	4	Tijuana	02.03	2	1
02005	02	Baja California	005	Playas de Rosarito	4	Tijuana	02.03	2	1
03003	03	Baja California Sur	003	La Paz	5	La Paz	03.01	69	1
04002	04	Campeche	002	Campeche	6	Campeche	04.01	70	1
Fuentes: Estimaciones del CONAPO con base en el INEGI, Censo de Población y Vivienda 2020

Localidad

.shp

A continuación se lee el archivo .shp que contiene de manera integrada la división de el área geoestadística por localidad.

shape_localidad <- readOGR(dsn ="MGN 2020/conjunto_de_datos", 
                            layer = "00l",
                             encoding = "UTF-8",
                              use_iconv = TRUE)

La función rename() del paquete dplyr permite cambiar el nombre de la columna de la clave geoestadística a nivel estatal dentro de la base de datos del shape.

shape_localidad@data <- shape_localidad@data %>%
                         rename("CVE_GEO" = "CVEGEO")

IML_2020

La base de datos del índice de marginación por localidades se encuentra disponible en la página oficial de CONAPO o bien se puede consultar en la página de Datos Abiertos y se presenta en formato .xlsx Consulta.

Para poder trabajar con la base de datos del índice de marginación a nivel localidad, se cambia el nombre de la columna CVE_LOC que hace referencia a la clave geográfica de la localidad por CVE_GEO para fines prácticos.

IML_2020 <- read_xlsx("Data/IML_2020.xlsx", sheet = "IML_2020") %>%
              mutate(CVE_GEO = CVE_LOC) %>% 
               as.data.frame()

Tabla: Localidades según el grado de marginación.
Grado de marginación	Localidades
Muy alto	13256
Alto	15985
Medio	25294
Bajo	31615
Muy bajo	24089

Layers

Existen dos maneras de trabajar los datos de los mapas:

Es trabajar con los objetos SpatialPolygonsDataFrame o bien spdf y tratar de unir los datos que se desean analizar, llamando algunas de las funciones de join_ o merge del paquete dplyr.
O bien utilizar la función fortify() del paquete ggplot2, que permite convertir un objeto (ej.spdf) a un data.frame. Al hacer este tipo de conversión permite aplicar un join para unir ambas estructuras de datos al nuevo formato del data.frame.

\[SpatialPolygons \Rightarrow data.frame\]

capas_localidad <- shape_localidad %>%
                    sp::spChFIDs(., str_pad(shape_localidad@data$CVE_GEO, 9, "left", pad = "0")) %>%
                     fortify(., id = c("CVE_GEO")) %>%
                      right_join(., IML_2020 %>% dplyr::select(CVE_GEO, ENT, NOM_ENT, MUN, NOM_MUN, NOM_LOC, GM_2020), 
                                  by = c("id" = "CVE_GEO")) 
length(unique(capas_localidad$id))
## [1] 110239

\[SpatialPolygons \Rightarrow SpatialPolygons + Datos\]

layer_localidad <- merge(shape_localidad,
                          IML_2020 %>%
                              dplyr::select(c(-CVE_LOC)) %>%
                                mutate(GM_2020 = fct_relevel(.$GM_2020,"Muy alto", "Alto", "Medio", "Bajo", "Muy bajo")),
                           by = "CVE_GEO")

AGEB

.shp

A continuación se lee el archivo .shp que contiene de manera integrada la división de el área geoestadística municipal ageb.

shape_ageb <- readOGR(dsn ="MGN 2020/conjunto_de_datos", 
                       layer = "00a",
                        encoding = "UTF-8",
                         use_iconv = TRUE)

La función rename() del paquete dplyr permite cambiar el nombre de la columna de la clave geoestadística a nivel estatal dentro de la base de datos del shape.

shape_ageb@data <- shape_ageb@data %>%
                    rename("CVE_GEO" = "CVEGEO")

IMU_2020

La base de datos del índice de marginación por ageb se encuentra disponible en la página oficial de CONAPO o bien se puede consultar en la página de Datos Abiertos y se presenta en formato .xlsx Consulta.

Para poder trabajar con la base de datos del índice de marginación a nivel ageb, se cambia el nombre de la columna `CVE_AGEB que hace referencia a la clave geográfica de la ageb por CVE_GEO para fines prácticos.

IMU_2020 <- read_xlsx("Data/IMU_2020.xlsx", sheet = "IMU_2020") %>%
             mutate(CVE_GEO = CVE_AGEB) %>% 
              as.data.frame()

Tabla: Agebs según el grado de marginación.
Grado de marginación	Ageb
Muy alto	5910
Alto	8274
Medio	13209
Bajo	14105
Muy bajo	9292

Layers

Existen dos maneras de trabajar los datos de los mapas:

Es trabajar con los objetos SpatialPolygonsDataFrame o bien spdf y tratar de unir los datos que se desean analizar, llamando algunas de las funciones de join_ o merge del paquete dplyr.
O bien utilizar la función fortify() del paquete ggplot2, que permite convertir un objeto (ej.spdf) a un data.frame. Al hacer este tipo de conversión permite aplicar un join para unir ambas estructuras de datos al nuevo formato del data.frame.

\[SpatialPolygons \Rightarrow data.frame\]

capas_ageb <- shape_ageb %>%
               sp::spChFIDs(., str_pad(shape_ageb@data$CVE_GEO, 13, "left", pad = "0")) %>%
                fortify(., id = "CVE_GEO") %>%
                 right_join(., IMU_2020 %>% dplyr::select(CVE_GEO, ENT, NOM_ENT, MUN, NOM_MUN, LOC, NOM_LOC, GM_2020), 
                             by = c("id" = "CVE_GEO"))

\[SpatialPolygons \Rightarrow SpatialPolygons + Datos\]

layer_ageb <- merge(shape_ageb,
                     IMU_2020 %>%
                        dplyr::select(c(-CVE_AGEB)) %>%
                         mutate(GM_2020 = fct_relevel(.$GM_2020,"Muy alto", "Alto", "Medio", "Bajo", "Muy bajo")),
                    by = "CVE_GEO")

Zona Metropolitana del Valle de México

capas_zm <- shape_ageb %>%
             sp::spChFIDs(., str_pad(shape_ageb@data$CVE_GEO, 13, "left", pad = "0")) %>%
              fortify(., id = "CVE_GEO") %>%
               right_join(., IMU_2020 %>% dplyr::select(CVE_GEO, ENT, NOM_ENT, MUN, NOM_MUN, LOC, NOM_LOC, GM_2020), 
                           by = c("id" = "CVE_GEO")) %>% 
                 mutate(CVE_MUN = paste0(ENT, MUN)) %>%
                  right_join(., ZM_2015 %>% dplyr::select(CVE_MUN, CVE_ZM, NOM_ZM), 
                              by = "CVE_MUN")

Mapas por desagregación geográfica

Estructura básica de un mapa

La forma más básica de mapear datos en R es crear un objeto ggplot() regular y mapear los argumentos x = longuitude y la y = latitude. Se puede usar esta técnica para crear mapas de áreas geográficas, como estados o países, y para gráficar ubicaciones como puntos, líneas y otras formas. El paquete ggplot2 incluye algunos conjuntos de funciones geom que permiten graficar con información geográfica. Teniendo en cuenta que el data.frame incluye las columnas con ubicación (long y lat). También incluye una columna que describe el orden en que se deben conectar estos puntos para formar un polígono (order), el nombre del estado (NOM_ENT) y una columna de group que separa los puntos en polígonos únicos que se deben trazar.

p <- ggplot() + 
      geom_polygon(data = capas_estados,
                    mapping = aes(x = long, 
                                   y = lat, 
                                    group = group),
                     fill = "#ECECEC", 
                     color = "dark grey",
                      size = 0.5) +
        expand_limits(x = capas_estados$long,  
                      y = capas_estados$lat) + 
          theme(plot.title = element_text(size = 20, hjust = 0.15, family = "montserrat", face = "bold"),
                 plot.caption = element_text(size = 9, hjust = 0.2, vjust = 1, family = "montserrat")) +
       labs(title = "Estructura básica de un mapa",
              caption = expression(paste("Fuentes: Estimaciones del CONAPO con base en el INEGI, Censo de Población y Vivienda 2020")))
p

Entidad Federativa

`geom_map`

La función ggmap() devuelve un objeto ggplot. Puede usar este objeto ggplot resultante como lo haría con cualquier otro objeto ggplot (ejemplo, agregar geoms, cambiar tema). GitHub: fresques/ggmap

Una de las bondades del paquete ggmap, permite hacer mediante el uso de herramientas de Google Maps directamente desde R. Este paquete permite extraer un mapa base de Google Maps y algunos otros servidores de mapas, para que luego se puedan trazar y agregar puntos, polígonos, líneas y otros elementos usando funciones ggplot2. Este paquete usa la API de Google Maps. Otras funciones en el paquete ggmap, necesita una clave API de Google Maps para que funcione.

# Paquete
require(ggmap) 

p <- IME_2020 %>% 
      mutate(GM_2020 = fct_relevel(.$GM_2020,"Muy alto","Alto","Medio","Bajo","Muy bajo")) %>%
       ggplot() + 
        geom_map(map = capas_estados, 
                  aes(map_id = CVE_GEO, fill = GM_2020,  color = "GM_2020")) +
          expand_limits(x = capas_estados$long, 
                        y = capas_estados$lat) + 
           theme_transparent() + 
            theme(plot.title = element_text(size = 20, hjust = 0.15, family = "montserrat", face = "bold"),
                   plot.caption = element_text(size = 9, hjust = 0.2, vjust = 1, family = "montserrat"),
                    legend.key.size = unit(0.5, "cm"),
                     legend.text = element_text(size = 8, family = "montserrat"), 
                      legend.title = element_text(size = 10, hjust = 0.5, family = "montserrat", face = "bold"),
                       legend.position = c(0.8, 0.7)) + 
              scale_fill_viridis_d(option = "A", begin = 0.3, end = 0.8) + 
               scale_color_manual(values = c("#BDBDBD")) + 
                guides(color = 'none') +
       labs(title = "Índice de marginación a nivel estatal, 2020",
             fill = stringr::str_wrap("Grado de marginación", 10),
              caption = expression(paste("Fuentes: Estimaciones del CONAPO con base en el INEGI, Censo de Población y Vivienda 2020")))

#path = "Images/IME_2020_gmap.pdf"
#ggexport (p, width = 8, height = 5, res = 400, filename = path)

p

`geom_polygon`

La función geom_polygon() de la paquetería ggplot2, devuelve un objeto ggplot en el cual dibuja líneas entre puntos y los cierra ( es decir, dibuja una línea desde el último punto hasta el primer punto). En el cual es importante tener en cuenta que el data.frame incluye columnas con la ubicación de las áreas geográficas y se agregan como argumentos en la función aes(x = long, y = lat). También se incluye una columna que describe el orden en que se deben conectar estos puntos para formar un polígono (order), el nombre del estado (id) y una columna de grupo que separa los puntos en polígonos únicos que se deben trazar aes (group = group).

# Paleta de colores
require(unikn) 

p <- ggplot() + 
      geom_polygon(data = capas_estados,
                    mapping = aes(x = long, 
                                   y = lat, 
                                    group = group,
                                     fill = fct_relevel(GM_2020,"Muy alto","Alto","Medio","Bajo","Muy bajo")),
                     color = "dark grey",
                      size = 0.5) +
        expand_limits(x = capas_estados$long,  
                      y = capas_estados$lat) + 
         theme_transparent() + 
          theme(plot.title = element_text(size = 20, hjust = 0.15, family = "montserrat", face = "bold"),
                 plot.caption = element_text(size = 9, hjust = 0.2, vjust = 1, family = "montserrat"), 
                  legend.key.size = unit(0.5, "cm"),
                   legend.text = element_text(size = 8, family = "montserrat"), 
                    legend.title = element_text(size = 10, hjust = 0.5, family = "montserrat", face = "bold"),
                     legend.position = c(0.8, 0.7)) + 
           scale_fill_manual(values = rev(usecol(pal = pal_petrol, n = 5))) + 
            scale_color_manual(values = c("#BDBDBD")) +  
             guides(color = guide_legend(override.aes = list(fill = usecol(pal = pal_petrol, n = 5)))) +
       labs(title = "Índice de marginación a nivel estatal, 2020",
             fill = stringr::str_wrap("Grado de marginación", 10),
              caption = expression(paste("Fuentes: Estimaciones del CONAPO con base en el INEGI, Censo de Población y Vivienda 2020")))

#path = "Images/IME_2020_geom_polygon.pdf"
#ggexport (p, width = 8, height = 5, res = 400, filename = path)

p

`tmap`

Con el paquete tmap se pueden generar mapas temáticos con gran flexibilidad. La sintaxis para crear gráficos es similar a la de ggplot2, pero adaptada a los mapas. Además, ofrece un sistema de trazado coherente para mapas temáticos que se basan en la estructura de las capas en los gráficos. Los mapas temáticos se crean apilando capas, donde cada capa y cada dato se pueden asignar a una o más aes (). También es posible escalarlos y agregando posibles atributos del mapa, como una barra de escala o una brújula.

Una de las ventajas a diferencia de ggplot(), es más rápido de generar los gráficos y permite trabajar con múltiples capas a la vez. Además, permite hacerlo de manera interactiva como un widget HTML.

La función tm_shape, crea un elemento que especifica un objeto de SpatialPolygonsDataFrame, sin la necesidad de convertir los datos a un data.frame. GitHub: r-tmap/tmap

# Paquetes 
require(tmap)
require(tmaptools)

p <- tm_shape(layer_estados) + 
      tm_polygons("GM_2020", 
                   title = stringr::str_wrap("Grado de marginación", 10),
                    #palette = hcl.colors(18, palette = "Inferno")[3:9], 
                     palette = "BuPu") + 
       tm_borders(col = "grey20", lwd = 0.5) + 
        tm_layout(title.snap.to.legend = TRUE,
                   title = "Índice de marginación a nivel estatal, 2020",
                    title.fontfamily = "montserrat",
                     title.fontface = "bold",
                      title.size = 1.9,
                       main.title.fontface = "bold",
                        main.title.position = c("center", "top"),
                         main.title.size = 7, 
                          main.title.fontfamily = "montserrat",
                           legend.outside.position = "right",
                            legend.title.fontfamily = "montserrat",
                             legend.title.fontface = "bold", 
                              legend.text.fontfamily = "montserrat",
                               legend.title.size = 1.2, 
                                legend.text.size = 1, 
                                 inner.margins = c(t = 0.05, r = 0.02, b = 0.11, l = 0.15),
                                  fontfamily = "montserrat",
                                   frame = FALSE) + 
         tm_legend(legend.position = c("right", "center")) +
          tm_credits(expression(paste("Fuentes: Estimaciones del CONAPO con base en el INEGI, Censo de Población y Vivienda 2020")),
                      position=c("left", "bottom"), 
                       size = 1.1)
     
#path = "Images/IME_2020_tmap.pdf"
#ggexport (p, width = 8, height = 5, res = 400, filename = path)

p

`ggspatial`

Cualquier capa espacial se puede agregar a un ggplot() usando la función layer_spatial() del paquete ggspatial (bueno, cualquier objeto de los paquetes sf, sp o raster…). Estas capas entrenarán las escalas y a diferencia de las funciones geom_ o stat_, layer_spatial() siempre toma sus datos primero sin necesidad de especificar los límites de las escalas. Estructura ggspatial.

# Paquete 
require(ggspatial)

p <- ggplot() + 
      layer_spatial(layer_estados, aes(fill = GM_2020)) + 
       theme_transparent() + 
        theme(plot.title = element_text(size = 20, hjust = 0.15, family = "montserrat", face = "bold"),
               plot.caption = element_text(size = 9, hjust = 0.2, vjust = 1, family = "montserrat"), 
                legend.key.size = unit(0.5, "cm"),
                 legend.text = element_text(size = 10, family = "montserrat"), 
                  legend.title = element_text(size = 12, hjust = 0.5, family = "montserrat", face = "bold"),
                   legend.position = c(0.8, 0.7)) + 
          scale_fill_brewer(palette = rev("BuPu")) +
           scale_color_manual(values = c("#BDBDBD")) + 
            guides(color = guide_legend(override.aes = list(fill = usecol(pal = pal_petrol, n = 5)))) +
     labs(title = "Índice de marginación a nivel estatal, 2020",
           fill = stringr::str_wrap("Grado de marginación", 10), 
            caption = expression(paste("Fuentes: Estimaciones del CONAPO con base en el INEGI, Censo de Población y Vivienda 2020")))

#path = "Images/IME_2020_ggspatial.pdf"
#ggexport (p, width = 8, height = 5, res = 400, filename = path)

p

`geom_sf`

Conversión de datos: st_as_sf, convierte un objeto externo en un spatial object.

Los objetos espaciales creados por el paquete sp (que hasta hace poco era la forma estándar de manejar datos espaciales en R) se pueden formatear en un formato de características simples usando la función st_as_sf():

#Paquetes
require(sf)
estados_sf <- st_as_sf(capas_estados, 
                        coords = c("long","lat"),
                         crs = "+proj=lcc +lat_0=12 +lon_0=-102 +lat_1=17.5 +lat_2=29.5 +x_0=2500000 +y_0=0 +ellps=GRS80 +units=m +no_defs")

A continuación, se muestra un ejemplo de cómo visualizar objetos sf. Este tipo de función es un poco más amigable de usar.

La función as_tibble de la paquetería tibble, genera una base de datos genérica de S3 con métodos más eficientes para matrices y data.frame. Su envoltorio delgado permite trabajar de manera más rápida que un data.frame.

# Paquetes 
require(ggplot2)
require(tibble)
require(sf)
require(unikn)
ENT_tbl <- as_tibble(st_as_sf(shape_estados)) %>%
             merge(., IME_2020 %>% dplyr::select(CVE_GEO, GM_2020))

p <- ggplot(data = ENT_tbl) +  
      geom_sf(mapping = aes(geometry = geometry, 
                             fill = fct_relevel(GM_2020,"Muy alto","Alto","Medio","Bajo","Muy bajo")),
               color = "#BDBDBD", 
                size = 0.1)  + 
       theme_transparent() + 
        theme(plot.title = element_text(size = 16, hjust = 0.15, family = "montserrat", face = "bold"),
               plot.caption = element_text(size = 9, hjust = 0.2, vjust = 1, family = "montserrat"),  
                legend.key.size = unit(0.5, "cm"),
                 legend.text = element_text(size = 8, family = "montserrat"),   
                  legend.title = element_text(size = 10, hjust = 0.5, family = "montserrat", face = "bold"),
                   legend.position = c(0.8, 0.7)) + 
          scale_fill_manual(values = rev(usecol(pal = pal_petrol, n = 5))) + 
           scale_color_manual(values = rev(usecol(pal = pal_petrol, n = 5))) + 
            guides(color = guide_legend(override.aes = list(fill = rev(usecol(pal = pal_petrol, n = 5))))) +
      labs(title = "Índice de marginación a nivel estatal, 2020",
            fill = stringr::str_wrap("Grado de marginación", 10),
             caption = expression(paste("Fuentes: Estimaciones del CONAPO con base en el INEGI, Censo de Población y Vivienda 2020")))

#path = "Images/IME_2020_geom_sf.pdf"
#ggexport (p, width = 8, height = 5, res = 400, filename = path)

p

Municipio

`geom_map`

p <- IMM_2020 %>% 
      mutate(GM_2020 = fct_relevel(.$GM_2020,"Muy alto", "Alto", "Medio", "Bajo", "Muy bajo")) %>%
       ggplot() + 
        geom_map(map = capas_municipio, 
                  aes(map_id = CVE_GEO, fill = GM_2020), 
                   color = 'transparent')  +
         expand_limits(x = capas_municipio$long, y = capas_municipio$lat) + 
          theme_transparent() + 
           theme(plot.title = element_text(size = 20, hjust = 0.15, family = "montserrat", face = "bold"), 
                  plot.caption = element_text(size = 9, hjust = 0.2, vjust = 1, family = "montserrat"), 
                   legend.key.size = unit(0.5, "cm"),
                    legend.text = element_text(size = 8, family = "montserrat"), 
                     legend.title = element_text(size = 10, hjust = 0.5, family = "montserrat", face = "bold"),
                      legend.position = c(0.8, 0.7)) + 
            scale_fill_viridis_d(option = "A", begin = 0.3, end = 0.8) + 
             scale_color_manual(values = c("#BDBDBD")) + 
              guides(color = 'none') +
      labs(title = "Índice de marginación a nivel municipal, 2020",
            fill = stringr::str_wrap("Grado de marginación", 10),
             caption = expression(paste("Fuentes: Estimaciones del CONAPO con base en el INEGI, Censo de Población y Vivienda 2020")))

#path = "Images/IMM_2020_gmap.pdf"
#ggexport (p, width = 8, height = 5, res = 400, filename = path)

p

`geom_polygon`

La función geom_polygon() de la paquetería ggplot2, devuelve un objeto ggplot en el cual dibuja líneas entre puntos y los cierra ( es decir, dibuja una línea desde el último punto hasta el primer punto). En el cual es importante tener en cuenta que el data.frame incluye columnas con la ubicación de las áreas geográficas y se agregan como argumentos en la función aes(x = long, y = lat). También se incluye una columna que describe el orden en que se deben conectar estos puntos para formar un polígono (order), el nombre del municipio (id) y una columna de grupo que separa los puntos en polígonos únicos que se deben trazar aes (group = group).

# Paleta de colores 
require(RColorBrewer)

p <- ggplot() + 
      geom_polygon(capas_municipio,
                    mapping = aes(x = long, 
                                  y = lat, 
                                  group = group, 
                                  fill = fct_relevel(GM_2020, "Muy alto", "Alto", "Medio", "Bajo", "Muy bajo")),
                     color = "transparent",  
                      size = .15,
                       alpha = 0.9) +
      geom_polygon(capas_estados,
                    mapping = aes(x = long, 
                                  y = lat, 
                                  group = group),
                     color = "dark grey",
                      fill = "transparent",
                       size = 0.6) +
       coord_equal() +
        theme_transparent() + 
         theme(plot.title = element_text(size = 22, hjust = 0.15, family = "montserrat", face = "bold"),     
                plot.caption = element_text(size = 9, hjust = 0.2, vjust = 1, family = "montserrat"),
                 legend.key.size = unit(0.5, "cm"), 
                  legend.text = element_text(size = 10, family = "montserrat"), 
                   legend.title = element_text(size = 12, hjust = 0.5, family = "montserrat", face = "bold"),
                    legend.position = c(0.8, 0.7)) + 
          scale_fill_brewer(palette = "YlGnBu") +
           guides(color = 'none') +
     labs(title = expression(paste("Índice de marginación a nivel municipal, 2020")),
           fill = stringr::str_wrap("Grado de marginación", 10),
            caption = expression(paste("Fuentes: Estimaciones del CONAPO con base en el INEGI, Censo de Población y Vivienda 2020")))

#path = "Images/IMM_2020_geom_polygon.pdf"
#ggexport (p, width = 8, height = 5, res = 400, filename = path)

p

`tmap`

La función tm_shape, crea un elemento que especifica un objeto de SpatialPolygonsDataFrame, sin la necesidad de convertir los datos a un data.frame. GitHub: r-tmap/tmap

# Paquetes
require(tmap)
require(tmaptools)

p <- tm_shape(layer_municipios, id = "CVE_GEO") +
      tm_polygons("GM_2020", 
                   title = stringr::str_wrap("Grado de marginación", 10),
                    palette = "BuPu", 
                     border.col = 'transparent') + 
     tm_shape(layer_estados) + 
       tm_borders("grey20", lwd = 0.5) + 
        tm_layout(title.snap.to.legend = TRUE,
                   title = "Índice de marginación a nivel municippal, 2020",
                    title.fontfamily = "montserrat",
                     title.fontface = "bold",
                      title.size = 1.8,
                       main.title.fontface = "bold",
                        main.title.position = c("center", "top"),
                         main.title.size = 7, 
                          main.title.fontfamily = "montserrat",
                           legend.outside.position = "right",
                            legend.title.fontfamily = "montserrat",
                             legend.title.fontface = "bold", 
                              legend.text.fontfamily = "montserrat",
                               legend.title.size = 1.2, 
                                legend.text.size = 1, 
                                 inner.margins = c(t = 0.06, r = 0.04, b = 0.12, l = 0.15),
                                  fontfamily = "montserrat",
                                   frame = FALSE) + 
         tm_legend(legend.position = c("right", "center")) +
          tm_credits(expression(paste("Fuentes: Estimaciones del CONAPO con base en el INEGI, Censo de Población y Vivienda 2020")),
                      position=c("left", "bottom"), 
                       size = 1.4)
     
#path = "Images/IMM_2020_tmap.pdf"
#ggexport (p, width = 8, height = 5, res = 400, filename = path)

p

`geom_sf`

A continuación, se muestra un ejemplo de cómo visualizar objetos sf. Este tipo de función es un poco más amigable de usar.

require(ggplot2)
require(tibble)
require(sf)
MUN_tbl <- as.tibble(st_as_sf(shape_municipios)) %>%
            merge(., IMM_2020 %>% dplyr::select(CVE_GEO, GM_2020))
ENT_tbl <- as_tibble(st_as_sf(shape_estados))


p <- ggplot() +
      geom_sf(data = MUN_tbl,
               mapping = aes(geometry = geometry, 
                             fill = fct_relevel(GM_2020, "Muy alto","Alto","Medio","Bajo","Muy bajo"),
                             color = fct_relevel(GM_2020, "Muy alto","Alto","Medio","Bajo","Muy bajo")),
                  size = 0.1) + 
      geom_sf(data = ENT_tbl,
               mapping = aes(geometry = geometry), 
                fill = 'transparent',
                 color = "black") +
        coord_sf() +
        theme_transparent() + 
         theme(plot.title = element_text(size = 22, hjust = 0.15, family = "montserrat", face = "bold"),
                plot.caption = element_text(size = 10, hjust = 0.2, vjust = 1, family = "montserrat"),  
                 legend.key.size = unit(0.5, "cm"),
                  legend.text = element_text(size = 10, family = "montserrat"),   
                   legend.title = element_text(size = 12, hjust = 0.5, family = "montserrat", face = "bold"),
                    legend.position = c(0.8, 0.7)) + 
          scale_fill_manual(values = usecol(pal_bordeaux, n = 5)) + 
           scale_color_manual(values = usecol(pal_bordeaux, n = 5)) + 
      labs(title = "Índice de marginación a nivel municipal, 2020",
            fill = stringr::str_wrap("Grado de marginación", 10),
             color = stringr::str_wrap("Grado de marginación", 10),
             caption = expression(paste("Fuentes: Estimaciones del CONAPO con base en el INEGI, Censo de Población y Vivienda 2020")))

#path = "Images/IMM_2020_geom_sf.pdf"
#ggexport (p, width = 8, height = 5, res = 400, filename = path)

p

Localidad

`geom_polygon`

La función geom_polygon() de la paquetería ggplot2, devuelve un objeto ggplot en el cual dibuja líneas entre puntos y los cierra ( es decir, dibuja una línea desde el último punto hasta el primer punto). En el cual es importante tener en cuenta que el data.frame incluye columnas con la ubicación de las áreas geográficas y se agregan como argumentos en la función aes(x = long, y = lat). También se incluye una columna que describe el orden en que se deben conectar estos puntos para formar un polígono (order), el nombre de la localidad (id) y una columna de grupo que separa los puntos en polígonos únicos que se deben trazar aes (group = group).

p <- ggplot() + 
      geom_polygon(capas_estados,
                    mapping = aes(x = long, 
                                  y = lat,  
                                  group = group), 
                     fill = 'transparent',
                      color = "black", 
                       size = .15,
                        alpha = 0.9) +
       geom_polygon(capas_localidad,
                     mapping = aes(x = long, 
                                   y = lat, 
                                   group = group, 
                                   fill = fct_relevel(GM_2020, "Muy alto","Alto","Medio","Bajo","Muy bajo"), 
                                   color = fct_relevel(GM_2020, "Muy alto", "Alto", "Medio", "Bajo", "Muy bajo")),
                      size = .15,
                       alpha = 0.9) +
        coord_equal() +
         theme_transparent() + 
          theme(plot.title = element_text(size = 21, hjust = 0.15, family = "montserrat", face = "bold"),
                 plot.caption = element_text(size = 9, hjust = 0.2, vjust = 1, family = "montserrat"),  
                  legend.key.size = unit(0.5, "cm"),
                   legend.text = element_text(size = 10, family = "montserrat"), 
                    legend.title = element_text(size = 12, hjust = 0.5, family = "montserrat", face = "bold"),
                     legend.position = c(0.8, 0.7)) + 
           scale_fill_viridis_d(option = "A", begin = 0.3, end = 0.8) + 
            scale_color_viridis_d(option = "A", begin = 0.3, end = 0.8) + 
             guides(color = 'none') +
      labs(title = "Índice de marginación a nivel localidad, 2020",
            fill = stringr::str_wrap("Grado de marginación", 10),
             caption = expression(paste("Fuentes: Estimaciones del CONAPO con base en el INEGI, Censo de Población y Vivienda 2020")))

#path = "Images/IML_2020_geom_polygon.pdf"
#ggexport (p, width = 8, height = 5, res = 400, filename = path)

p

Por estados

# Clave única de estados 
estados <- unique(IML_2020$ENT)
nom_ent <- unique(IML_2020$NOM_ENT)

p <- NULL 
q <- NULL
#for(i in 1:32){
for(i in 22){
  
municipio <- IML_2020 %>% 
              mutate(CVE_MUN = paste0(ENT, MUN)) %>%
               filter(ENT == estados[i]) %>% 
                dplyr::select(CVE_MUN) %>%
                 unique()

nom_mun <- IML_2020 %>% 
            filter(ENT == estados[i]) %>% 
             dplyr::select(NOM_MUN) %>%
               unique()

p[[i]] <- ggplot() + 
           geom_polygon(capas_estados %>% filter(id == estados[i]),
                         mapping = aes(x = long, 
                                       y = lat, 
                                       group = group), 
                          fill = 'transparent',
                           color = "black", 
                            size = .5,
                             alpha = 1) +
           geom_polygon(capas_municipio %>% filter(CVE_ENT == estados[i]),
                         mapping = aes(x = long, 
                                       y = lat, 
                                       group = group), 
                          fill = 'transparent',
                           color = "#7A7A7A", 
                            size = .2,
                             alpha = 0.9) +
           geom_polygon(capas_localidad %>% filter(ENT == estados[[i]]),
                         mapping = aes(x = long, 
                                        y = lat, 
                                         group = group, 
                                          fill = fct_relevel(GM_2020, "Muy alto", "Alto", "Medio", "Bajo", "Muy bajo"),
                                           color = fct_relevel(GM_2020, "Muy alto", "Alto", "Medio", "Bajo", "Muy bajo")),
                          size = 0.5,
                           alpha = 0.6) +
             coord_equal() +
              theme_transparent() + 
               theme(plot.title = element_text(size = 15, hjust = 0.5, family = "montserrat", face = "bold"),
                      plot.caption = element_text(size = 9, hjust = 0, family = "montserrat"), 
                       legend.key.size = unit(0.3, "cm"),
                        legend.text = element_text(size = 10, family = "montserrat"), 
                         legend.title = element_text(size = 12, hjust = 0.5, family = "montserrat", face = "bold"),
                          legend.position = "bottom") + 
                scale_fill_viridis_d(option = "A", begin = 0.3, end = 0.8) + 
                 scale_color_viridis_d(option = "A", begin = 0.3, end = 0.8) + 
                  guides(color = 'none', 
                         fill = guide_legend(override.aes = list(alpha = 1))) +
          labs(title = stringr::str_wrap(paste(nom_ent[i], ": índice de marginación a nivel localidad, 2020"), 80),
                fill = stringr::str_wrap("Grado de marginación", 10),
                 #caption = expression(paste("Fuentes: Estimaciones del CONAPO con base en el INEGI, Censo de Población y Vivienda 2020"))
               )

for(j in 1:nrow(municipio)){
  
q[[j]] <- ggplot() + 
           geom_polygon(data = capas_municipio %>% filter(id == municipio$CVE_MUN[j]), 
                         mapping = aes(x = long, 
                                       y = lat, 
                                       group = group), 
                          fill = 'transparent',
                           color="black", 
                            size = .15,
                             alpha = 0.9) +
           geom_polygon(data = capas_localidad %>% 
                                 mutate(CVE_MUN = paste0(ENT, MUN)) %>% 
                                  filter(CVE_MUN == municipio$CVE_MUN[j]),
                         mapping = aes(x = long, 
                                       y = lat, 
                                       group = group, 
                                       fill = fct_relevel(GM_2020, "Muy alto", "Alto", "Medio", "Bajo", "Muy bajo"),
                                       color = fct_relevel(GM_2020, "Muy alto", "Alto", "Medio", "Bajo", "Muy bajo")),
                          size = .15,
                           alpha = 0.9) +
             coord_equal() +
              theme_transparent() + 
               theme(plot.title = element_text(size = 18, hjust = 0.5, family = "montserrat", face = "bold"),
                      plot.caption = element_text(size = 9, hjust = 0, family = "montserrat"), 
                       legend.key.size = unit(0.3, "cm"),
                        legend.text = element_text(size = 10, family = "montserrat"), 
                         legend.title = element_text(size = 12, hjust = 0.5, family = "montserrat", face = "bold"),
                          legend.position = "right") + 
                scale_fill_viridis_d(option = "A", begin = 0.3, end = 0.8) + 
                 scale_color_viridis_d(option = "A", begin = 0.3, end = 0.8) + 
                  guides(color = 'none', fill = "none") + 
           labs(title = str_wrap(paste(nom_mun$NOM_MUN[j]),50),
                 fill = stringr::str_wrap("Grado de marginación", 10))
  }
}

# Se genera un gráfico con su estado y municipio 

### Se extrae la estructura de  la legenda del gráfico de entidad 

legend <- get_legend(p[[22]] + theme(legend.key.size = unit(0.5, "cm"),
                                      legend.text = element_text(size = 12, family = "montserrat"), 
                                       legend.title = element_text(size = 10, hjust = 0.5, family = "montserrat", face = "bold"),                                legend.position = "right"))

tabla <- lapply(1:nrow(municipio), function(x){
                                      ggarrange(p[[22]] + guides(fill = 'none') + labs(title = ""),
                                                 q[[x]], 
                                                  legend , 
                                                    ncol = 3, nrow = 1,
                                                     widths = c(1, 1, 0.5),
                                                      heights = c(1, 0.3, 1))
})

### Se generan los títulos y los pie de página del gráfico 

tabla <- lapply(1:nrow(municipio), function(x){
                                     annotate_figure(tabla[[x]], 
                                                     top = text_grob(stringr::str_wrap(paste(nom_ent[22], ": índice de marginación a nivel localidad, 2020"), 50),
                                                                     family = "montserrat", 
                                                                     face = "bold", 
                                                                     size = 16), 
                                                     bottom = text_grob(expression(paste("Fuentes: Estimaciones del CONAPO con base en el INEGI, Censo de Población y Vivienda 2020")),
                                                                     family = "montserrat", 
                                                                     face = "plain", 
                                                                     size = 9,
                                                                     hjust = 0.5))
})

#path = "Images/IML_2020_por estado_geom_polygon.pdf"
#ggexport (list = tabla, width = 8, height = 5, res = 400, filename = path)

tabla[[1]]

`geom_sf`

A continuación, se muestra un ejemplo de cómo visualizar objetos sf. Este tipo de función es un poco más amigable de usar.

require(ggplot2)
require(tibble)
require(sf)
LOC_tbl <- as.tibble(st_as_sf(shape_localidad)) %>%
             right_join(., IML_2020 %>% dplyr::select(CVE_GEO, GM_2020), by = "CVE_GEO")
MUN_tbl <- as.tibble(st_as_sf(shape_municipios))
ENT_tbl <- as_tibble(st_as_sf(shape_estados))

p <- list()
#for(i in 1:32){
for(i in 7){
  limites <- capas_estados %>%
               filter(id == estados[i]) %>%
                 summarise(min.x = min(.$long),
                           max.x = max(.$long),
                           min.y = min(.$lat),
                           max.y = max(.$lat))

p[[i]] <-  ggplot() + 
            geom_sf(data = LOC_tbl,
                     mapping = aes(geometry = geometry, 
                                   fill = fct_relevel(GM_2020, "Muy alto","Alto","Medio","Bajo","Muy bajo"),
                                   color = fct_relevel(GM_2020, "Muy alto","Alto","Medio","Bajo","Muy bajo"))) +
            geom_sf(data = MUN_tbl,
                     mapping = aes(geometry = geometry), 
                     fill = 'transparent',
                     color = "#D4D4D4",
                     size = 0.3) + 
            geom_sf(data = ENT_tbl,
                     mapping = aes(geometry = geometry), 
                      fill = 'transparent',
                       color = "black") +
             coord_sf(xlim = c(limites$min.x,limites$max.x),
                      ylim = c(limites$min.y,limites$max.y)) +
              theme_transparent() + 
               theme(plot.title = element_text(size = 22, hjust = 0, family = "montserrat", face = "bold"),
                      plot.caption = element_text(size = 9, hjust = 0.2, vjust = 1, family = "montserrat"),  
                       plot.margin = unit(c(t = 0.5, r = 0, b = 0, l = 0), "cm"),
                        legend.key.size = unit(0.5, "cm"),
                         legend.text = element_text(size = 10, family = "montserrat"),   
                          legend.title = element_text(size = 12, hjust = 0.5, family = "montserrat", face = "bold"),
                           legend.position = "right") + 
                scale_fill_manual(values = usecol(pal_bordeaux, n = 5)) + 
                 scale_color_manual(values = usecol(pal_bordeaux, n = 5)) + 
            labs(title = paste(nom_ent[i]),
                  fill = stringr::str_wrap("Grado de marginación", 10),
                   color = stringr::str_wrap("Grado de marginación", 10),
                   caption = expression(paste("Fuentes: Estimaciones del CONAPO con base en el INEGI, Censo de Población y Vivienda 2020")))
}

#path = "Images/IML_2020_geom_sf.pdf"
#ggexport (list = p, width = 8, height = 5, res = 400, filename = path)

p[[7]]

AGEB

`geom_polygon`

La función geom_polygon() de la paquetería ggplot2, devuelve un objeto ggplot en el cual dibuja líneas entre puntos y los cierra ( es decir, dibuja una línea desde el último punto hasta el primer punto). En el cual es importante tener en cuenta que el data.frame incluye columnas con la ubicación de las áreas geográficas y se agregan como argumentos en la función aes(x = long, y = lat). También se incluye una columna que describe el orden en que se deben conectar estos puntos para formar un polígono (order), el nombre de la ageb (id) y una columna de grupo que separa los puntos en polígonos únicos que se deben trazar aes (group = group).

p <- ggplot() + 
      geom_polygon(capas_estados,
                    mapping = aes(x = long, 
                                  y = lat, 
                                  group = group), 
                     fill = 'transparent',
                      color="black", 
                       size = .15,
                        alpha = 0.9) +
      geom_polygon(capas_ageb,
                    mapping = aes(x = long, 
                                  y = lat, 
                                  group = group, 
                                  fill = fct_relevel(GM_2020, "Muy alto", "Alto", "Medio", "Bajo", "Muy bajo"), 
                                  color = fct_relevel(GM_2020, "Muy alto", "Alto", "Medio", "Bajo", "Muy bajo")),
                      size = .15,
                       alpha = 0.9) +
        coord_equal() +
         theme_transparent() + 
          theme(plot.title = element_text(size = 22, hjust = 0.15, family = "montserrat", face = "bold"), 
                 plot.caption = element_text(size = 9, hjust = 0.2, vjust = 1, family = "montserrat"), 
                  legend.key.size = unit(0.5, "cm"),
                   legend.text = element_text(size = 10, family = "montserrat"), 
                    legend.title = element_text(size = 12, hjust = 0.5, family = "montserrat", face = "bold"),
                     legend.position = c(0.8, 0.7)) + 
            scale_fill_viridis_d(option = "A", begin = 0.3, end = 0.8) + 
             scale_color_viridis_d(option = "A", begin = 0.3, end = 0.8) +  
              guides(color = 'none', 
                      fill = guide_legend(override.aes = list(alpha = 1))) +
       labs(title = "Índice de marginación a nivel AGEB, 2020",
             fill = stringr::str_wrap("Grado de marginación", 10),
              caption = expression(paste("Fuentes: Estimaciones del CONAPO con base en el INEGI, Censo de Población y Vivienda 2020")))


#path = "Images/IMU_2020_geom_polygon.pdf"
#ggexport (p, width = 8, height = 5, res = 400, filename = path)

p

Zona Metropolitana del Valle de México (ZMVM)

tabla <- ZM_2015 %>%
          filter(CVE_ZM == "09.01")

p <-  ggplot() + 
       geom_polygon(capas_estados %>% filter(id == c("09", "13", "15")),
                     mapping = aes(x = long, 
                                   y = lat, 
                                   group = group), 
                      fill = 'transparent',
                       color = "black", 
                        size = .5,
                         alpha = 0.9) +
       geom_polygon(capas_municipio %>% filter(id == tabla$CVE_MUN),
                     mapping = aes(x = long, 
                                   y = lat, 
                                   group = group), 
                      fill = 'transparent',
                       color = "#C7C7C7", 
                        size = .3,
                         alpha = 0.9) +
       geom_polygon(capas_zm %>% filter(CVE_ZM == "09.01"),
                     mapping = aes(x = long, 
                                   y = lat, 
                                   group = group, 
                                   fill = fct_relevel(GM_2020, "Muy alto", "Alto", "Medio", "Bajo", "Muy bajo"),
                                   color = fct_relevel(GM_2020, "Muy alto", "Alto", "Medio", "Bajo", "Muy bajo")),
                      size = .15,
                       alpha = 0.7) +
        coord_equal() +
         theme_transparent() + 
          theme(plot.title = element_text(size = 20, hjust = 0.15, family = "montserrat", face = "bold"), 
                 plot.caption = element_text(size = 9, hjust = 0.2, vjust = 1, family = "montserrat"), 
                  legend.key.size = unit(0.5, "cm"),
                   legend.text = element_text(size = 10, family = "montserrat"), 
                    legend.title = element_text(size = 12, hjust = 0.5, family = "montserrat", face = "bold"),
                     legend.position = "right") + 
            scale_fill_viridis_d(option = "A", begin = 0.3, end = 0.8) + 
             scale_color_viridis_d(option = "A", begin = 0.3, end = 0.8) +  
              guides(color = 'none', 
                      fill = guide_legend(override.aes = list(alpha = 1))) +
       labs(title = "ZMVM: Índice de marginación a nivel AGEB, 2020",
             fill = stringr::str_wrap("Grado de marginación", 10),
              caption = expression(paste("Fuentes: Estimaciones del CONAPO con base en el INEGI, Censo de Población y Vivienda 2020")))

#path = "Images/IMU_2020_ZMVM_geom_polygon.pdf"
#ggexport (p, width = 8, height = 5, res = 400, filename = path)

p

`geom_sf`

A continuación, se muestra un ejemplo de cómo visualizar objetos sf. Este tipo de función es un poco más amigable de usar.

# Paquetes
require(ggplot2)
require(tibble)
require(sf)
AGEB_tbl <- as.tibble(st_as_sf(shape_ageb)) %>%
              right_join(., IMU_2020 %>% dplyr::select(CVE_GEO, GM_2020), by = "CVE_GEO")
LOC_tbl <- as.tibble(st_as_sf(shape_localidad)) 
MUN_tbl <- as.tibble(st_as_sf(shape_municipios))
ENT_tbl <- as_tibble(st_as_sf(shape_estados))

q <- NULL
p <- list()
#for(i in 1:32){
for(i in 6){
  limites <- capas_estados %>%
               filter(id == estados[i]) %>%
                 summarise(min.x = min(.$long),
                           max.x = max(.$long),
                           min.y = min(.$lat),
                           max.y = max(.$lat))
q[[i]] <-  ggplot() + 
            geom_sf(data = AGEB_tbl %>% filter(CVE_ENT == estados[i]),
                     mapping = aes(geometry = geometry, 
                                   fill = fct_relevel(GM_2020, "Muy alto","Alto","Medio","Bajo","Muy bajo"),
                                   color = fct_relevel(GM_2020, "Muy alto","Alto","Medio","Bajo","Muy bajo"),
                                   linetype = "Line1"),
                      size = 0.7) +
            geom_sf(data = LOC_tbl %>% filter(CVE_ENT == estados[i]),
                     mapping = aes(geometry = geometry, linetype = "Line2"), 
                     fill = 'transparent',
                     color = "#C70E8E", 
                      size = 0.03) + 
            geom_sf(data = MUN_tbl %>% filter(CVE_ENT == estados[i]),
                     mapping = aes(geometry = geometry), 
                     fill = 'transparent',
                     color = "#D4D4D4",
                     size = 0.3) + 
            geom_sf(data = ENT_tbl %>% filter(CVE_ENT == estados[i]),
                     mapping = aes(geometry = geometry), 
                      fill = 'transparent',
                       color = "black") +
             coord_sf(xlim = c(limites$min.x,limites$max.x),
                      ylim = c(limites$min.y,limites$max.y)) +
              theme_transparent() + 
               theme(plot.title = element_text(size = 22, hjust = 0.15, family = "montserrat", face = "bold"),
                      plot.caption = element_text(size = 9, hjust = 0.2, vjust = 1, family = "montserrat"),  
                       legend.key.size = unit(0.5, "cm"),
                        legend.text = element_text(size = 10, family = "montserrat"),   
                         legend.title = element_text(size = 12, hjust = 0.5, family = "montserrat", face = "bold"),
                          legend.position = "right") + 
                scale_fill_manual(values = usecol(pal_bordeaux, n = 5)) + 
                 scale_color_manual(values = usecol(pal_bordeaux, n = 5)) +  
                  scale_linetype_manual("Nivel", labels = c("AGEB", "Localidad"), values = c("Line1" = "solid", "Line2" = "solid")) +
                   guides(linetype = guide_legend(override.aes = list(linetype = "solid", fill = 'transparent', color = c("#771434", "#C70E8E"), size = c(0.7, 0.03)))) +
            labs(title = paste(nom_ent[i]),
                  fill = stringr::str_wrap("Grado de marginación", 10),
                   color = stringr::str_wrap("Grado de marginación", 10),
                   caption = expression(paste("Fuentes: Estimaciones del CONAPO con base en el INEGI, Censo de Población y Vivienda 2020")))
if(i == 6){
  p[[i]] <- q[[i]] + 
             coord_sf(xlim = c(limites$min.x + 1060000,limites$max.x),
                      ylim = c(limites$min.y,limites$max.y)) 
} else {
  p[[i]] <- q[[i]]
}
}

#path = "Images/IMU_2020_geom_sf.pdf"
#ggexport (list = p, width = 8, height = 5, res = 400, filename = path)

p[[6]]

Complementos

Entidad

Etiquetas

Se utiliza la función coordinates() del paquete sp que permite localizar los centroides de cada estado. Permitiendo de esta manera encontrar la localización para las etiquetas de cada área geográfica y utilizando la función geom_label_repel() mediante los atributos mapping = aes(x = Longitude, y = Latitude) y de esta manera no permitir que las etiquetas del mapa se sobrelapen entre sí y se alejen de los puntos de datos generados por la función coordinates().

#Paquete
require(ggrepel)
require(viridisLite)

# "coordinates()" extrae los centroides de los polígonos, en el orden indicado en shape@data
centroids.df <- as.data.frame(coordinates(shape_estados)) %>%
                  rename("Longitude" = "V1","Latitude" = "V2")

p <- IME_2020 %>% 
      cbind(., centroids.df) %>%
       ggplot() + 
        geom_polygon(data = capas_estados,
                      mapping = aes(x = long, 
                                    y = lat, 
                                    group = group,
                                    fill = fct_relevel(GM_2020, "Muy alto", "Alto", "Medio", "Bajo", "Muy bajo")),
                       color = "dark grey",
                        size = 0.5) +
         geom_label_repel(aes(label = stringr::str_wrap(NOM_ENT, 15),
                               x = Longitude,
                                y = Latitude,
                                 color = fct_relevel(GM_2020,"Muy alto", "Alto", "Medio", "Bajo", "Muy bajo")), 
                           force = 0.6,
                            box.padding = 0.25,
                             min.segment.length = 0.5,
                              size = 3.5, 
                               fill = "transparent",
                                point.padding = 0.5,
                                 segment.size = 0.7,
                                  segment.colour= "dark grey",
                                   label.size = 0,
                                    alpha = 0.8) +
           expand_limits(x = capas_estados$long,  
                         y = capas_estados$lat) + 
            theme_transparent() + 
             theme(plot.title = element_text(size = 20, hjust = 0.15, family = "montserrat", face = "bold"),
                    plot.caption = element_text(size = 9, hjust = 0.2, vjust = 1, family = "montserrat"), 
                     legend.key.size = unit(0.5, "cm"),
                      legend.text = element_text(size = 8, family = "montserrat"), 
                       legend.title = element_text(size = 10, hjust = 0.5, family = "montserrat", face = "bold"),
                        legend.position = c(0.9, 0.7)) + 
               scale_fill_viridis_d(option = "A", begin = 0.3, end = 0.8) + 
                scale_color_viridis_d(option = "A", begin = 0.3, end = 0.8) + 
                 guides(fill = guide_legend(override.aes = list(color = viridis(n = 5, option = "A", begin = 0.3, end = 0.8))),
                         color = 'none') +
         labs(title = "Índice de marginación a nivel estatal, 2020",
               fill = stringr::str_wrap("Grado de marginación", 10),
                caption = expression(paste("Fuentes: Estimaciones del CONAPO con base en el INEGI, Censo de Población y Vivienda 2020")))

#path = "Images/IME_2020_geom_polygon_etiquetas.pdf"
#ggexport (p, width = 8, height = 5, res = 400, filename = path)

p

`Annotation_scale`

Usando las funciones annotation_scale() y annotation_north_arrow() del paquete ggspatial se puede agregar una barra de escala o bien una flecha de norte. El sistema de posicionamiento de estas características está controlado por el argumento de ubicación (ejemplo, location = “tl” ~ topleft). Usando los argumentos pad_x y pad_y se pueden hacer ajustes menores al posicionamiento de los iconos.

#Paquetes
require(ggspatial)

#Ejemplos
#style = north_arrow_orienteering
#style = north_arrow_fancy_orienteering
#style = north_arrow_minimal
#style = north_arrow_nautical

p <- IME_2020 %>% 
       ggplot() + 
        geom_polygon(data = capas_estados,
                      mapping = aes(x = long, 
                                    y = lat, 
                                    group = group,
                                    fill = fct_relevel(GM_2020, "Muy alto", "Alto", "Medio", "Bajo", "Muy bajo")),
                       color = "dark grey",
                        size = 0.5) +
         annotation_north_arrow(location = "tl", 
                                 style = north_arrow_fancy_orienteering, 
                                  pad_x = unit(0, "cm"), 
                                   height = unit(1, "cm"),
                                    width = unit(1, "cm")) +
          annotation_scale(location = "bl", 
                            plot_unit = "m" , 
                             pad_y = unit(0, "cm"), 
                              text_family = "montserrat") +
           expand_limits(x = capas_estados$long,  
                         y = capas_estados$lat) + 
            theme_transparent() + 
             theme(plot.title = element_text(size = 20, hjust = 0.15, family = "montserrat", face = "bold"),
                    plot.caption = element_text(size = 9, hjust = 0.2, vjust = 1, family = "montserrat"), 
                     legend.key.size = unit(0.5, "cm"),
                      legend.text = element_text(size = 10, family = "montserrat"), 
                       legend.title = element_text(size = 12, hjust = 0.5, family = "montserrat", face = "bold"),
                        legend.position = c(0.8, 0.7)) + 
              scale_fill_viridis_d(option = "A", begin = 0.3, end = 0.8) + 
               scale_color_viridis_d(option = "A", begin = 0.3, end = 0.8) + 
                guides(fill = guide_legend(override.aes = list(color = viridis(n = 5, option = "A", begin = 0.3, end = 0.8))),
                        color = 'none') +
       labs(title = "Índice de marginación a nivel estatal, 2020",
             fill = stringr::str_wrap("Grado de marginación", 10),
               caption = expression(paste("Fuentes: Estimaciones del CONAPO con base en el INEGI, Censo de Población y Vivienda 2020")))

#path = "Images/IME_2020_annotation.pdf"
#ggexport (p, width = 8, height = 5, res = 400, filename = path)

p

Municipios

Etiquetas

Se utiliza la función coordinates() del paquete sp que permite localizar los centroides de cada muncipio. Permitiendo de esta manera encontrar la localización para las etiquetas de cada área geográfica y utilizando la función geom_label_repel() mediante los atributos mapping = aes(x = Longitude, y = Latitude) y de esta manera no permitir que las etiquetas del mapa se sobrelapen entre sí y se alejen de los puntos de datos generados por la función coordinates().

capas_zm <- left_join(capas_municipio, 
                       ZM_2015 %>% dplyr::select(-NOM_ENT, -NOM_MUN, -CVE_ENT), 
                        by = c("id" = "CVE_MUN")) %>%
             filter(CVE_ENT %in% c("09","13","15")) %>%
              mutate(condicional = ifelse(CVE_ZM == "09.01", .$GM_2020, NA))

# "coordinates()" extrae los centroides de los polígonos, en el orden indicado en shape@data
centroids.df <- as.data.frame(coordinates(shape_estados)) %>%
                  rename("Longitude" = "V1","Latitude" = "V2")

tabla <- ZM_2015 %>%
          filter(CVE_ZM == "09.01") %>%
           dplyr::select(CVE_MUN) %>%
            unique() %>% 
             as.vector()

p <- ggplot() + 
      geom_polygon(data = capas_zm,
                    mapping = aes(x = long, 
                                  y = lat, 
                                  group = group,
                                  fill = fct_relevel(condicional, c("Muy alto","Alto","Medio","Bajo","Muy bajo"))),
                     color = "dark grey",
                      size = 0.3) +
      geom_polygon(data = capas_estados %>% filter(id %in% c("09","13","15")),
                    mapping = aes(x = long, 
                                  y = lat, 
                                  group = group),
                     fill = 'transparent',
                      color = "black",
                       size = 0.8) +
      geom_label_repel(data = data.frame(IME_2020, centroids.df) %>%
                                filter(CVE_ENT %in% c("09","13","15")),
                        mapping = aes(label = stringr::str_wrap(NOM_ENT, 15),
                                      x = Longitude,
                                      y = Latitude), 
                         box.padding = 0.25,
                          point.padding = 3,
                           size = 5, 
                            fill = "#EDEDED",
                             color = "#000036",
                              segment.colour= "dark grey",
                               alpha = 0.8) +
         coord_equal() + 
          theme_transparent() + 
           theme(plot.title = element_text(size = 14, hjust = 0.15, family = "montserrat", face = "bold"),
                  plot.caption = element_text(size = 10, hjust = 0.2, vjust = 1, family = "montserrat"), 
                   legend.key.size = unit(0.5, "cm"),
                    legend.text = element_text(size = 8, family = "montserrat"), 
                     legend.title = element_text(size = 10, hjust = 0.2, family = "montserrat", face = "bold")) +
             scale_fill_manual(labels = c("Medio","Bajo","Muy bajo"),
                               values = rev(usecol(pal = pal_petrol, n = 3)),
                               na.translate = F) + 
               guides(color = 'none') +
       labs(title = "ZMVM: Zona Metropolitana del Valle de México, 2020",
             fill = stringr::str_wrap("Grado de marginación", 10),
              caption = expression(paste("Fuentes: Estimaciones del CONAPO con base en el INEGI, Censo de Población y Vivienda 2020")))
 
#path = "Images/IMM_2020_geom_polygon_etiquetas.pdf"
#ggexport (p, width = 6, height = 5, res = 400, filename = path)

p

# "coordinates()" extrae los centroides de los polígonos, en el orden indicado en shape@data
centroids.df <- as.data.frame(coordinates(shape_municipios)) %>%
                 rename("Longitude" = "V1","Latitude" = "V2") %>%
                  mutate(CVE_ENT = shape_municipios@data$CVE_ENT,
                         CVE_MUN = shape_municipios@data$CVE_GEO,
                         NOM_MUN = iconv(shape_municipios@data$NOMGEO, from = "utf8", to = "latin9"))
p <- NULL
#for(i in 1:32){
for(i in 7){
# Se generan los límites de cada estados, en caso de que se requiera modificar alguno
limites <- capas_municipio %>%
            subset(CVE_ENT == estados[i]) %>%
             summarise(min.x = min(.$long),
                       max.x = max(.$long),
                       min.y = min(.$lat),
                       max.y = max(.$lat))

# Se genera una muestra aleatoria   
sample <- centroids.df %>%
            filter(CVE_ENT == estados[i]) %>%
             sample_n(., min(30, nrow(.)))

q <- ggplot() + 
      geom_polygon(data = capas_municipio %>% filter(CVE_ENT == estados[i]),
                    mapping = aes(x = long, 
                                  y = lat, 
                                  group = group,
                                  fill = fct_relevel(GM_2020, c("Muy alto", "Alto", "Medio", "Bajo", "Muy bajo"))),
                     color = "dark grey",
                      size = 0.3) +
       geom_polygon(data = capas_estados %>% filter(id == estados[i]),
                     mapping = aes(x = long, 
                                   y = lat, 
                                   group = group),
                      fill = 'transparent',
                       color = "black",
                        size = 0.8) +
        geom_label_repel(data = sample,
                          mapping = aes(label = stringr::str_wrap(NOM_MUN, 15),
                                        x = Longitude,
                                        y = Latitude), 
                           size = 3.5, 
                            fill = "#EDEDED",
                             color = "#000036",
                              alpha = 0.8,
                               force = 0.2,
                                force_pull = 0.5,
                                 direction = c("both"),
                                  box.padding = 0.25,
                                   point.padding = 0,
                                    segment.colour= "dark grey",
                                     segment.size = 1,
                                      min.segment.length = 0.5) +
           coord_equal() +
           theme_transparent() + 
            theme(plot.title = element_text(size = 14, hjust = 0.15, family = "montserrat", face = "bold"),
                   plot.caption = element_text(size = 10, hjust = 0.2, vjust = 1, family = "montserrat"), 
                    legend.key.size = unit(0.5, "cm"),
                     legend.text = element_text(size = 8, family = "montserrat"), 
                      legend.title = element_text(size = 10, hjust = 0.2, family = "montserrat", face = "bold")) +
             scale_fill_manual(values = brewer.pal(n = 9, name = "Spectral"),
                                na.translate = T) + 
              guides(color = 'none') +
        labs(title = paste(nom_ent [i]),
              fill = stringr::str_wrap("Grado de marginación", 10),
               caption = expression(paste("Fuentes: Estimaciones del CONAPO con base en el INEGI, Censo de Población y Vivienda 2020")))

if(i %in% c(1:5,7:32)){
  p[[i]] <- q  
} else {
  p[[i]] <- q + coord_cartesian(xlim = c(limites$min.x + 1060000, limites$max.x),
                                ylim = c(limites$min.y, limites$max.y)) 
  }
}

#path = "Images/IMM_2020_etiquetas por estado.pdf"
#ggexport (list = p, width = 8, height = 5, res = 400, filename = path)

p[[7]]

2010 - 2020

Son nueve indicadores simples, los que conforman el índice de marginación a nivel municipal. De los cuales se enfocan en cuatro dimensiones de exclusión social.

Indicadores <- c('Porcentaje de población de 15 años o más analfabeta', 
                 'Porcentaje de población de 15 años o más sin educación básica',
                 'Porcentaje de ocupantes en viviendas sin drenaje ni excusado',
                 'Porcentaje de ocupantes en viviendas sin energía eléctrica', 
                 'Porcentaje de ocupantes en viviendas sin agua entubada', 
                 'Porcentaje de ocupantes en viviendas con piso de tierra',
                 'Porcentaje de población en localidades con menos de 5 000 habitantes', 
                 'Porcentaje de viviendas con algún nivel de hacinamiento',
                 'Porcentaje de población ocupada con ingresos de hasta 2 salarios mínimos')

Se cargan las bases de datos de los indicadores simples del índice de marginación 2010-2020

tablas <- c("2010","2015","2020")  

for(i in tablas){
  assign(paste0("Ind_", i), read_xlsx(paste0("Data/Indicadores_", i, ".xlsx"), sheet = "Indicadores") %>% 
                              rename("ANIO" = "AÑO"))
}

### Se crea un data.frame (se incluyen todo los años)
Ind <- do.call(dplyr::bind_rows,list(Ind_2010, Ind_2015, Ind_2020))

for(i in tablas){
  rm(list = paste0("Ind_",i))
}

Utilizando la función cut() y quantile se dividen en cuantiles cada una de las variables simples de todos los años y de esta manera poder visualizar la evolución de los indicadores de rezago.

Ind <- Ind %>% 
        mutate(Q_ANALF = cut(ANALF, breaks = quantile(unique(ANALF), probs = seq(0, 1, 0.2)), include.lowest = TRUE, dig.lab = 2),
               Q_SBASC = cut(SBASC, breaks = quantile(unique(SBASC), probs = seq(0, 1, 0.2)), include.lowest = TRUE, dig.lab = 2),
               Q_OVSDE = cut(OVSDE, breaks = quantile(unique(OVSDE), probs = seq(0, 1, 0.2)), include.lowest = TRUE, dig.lab = 2),
               Q_OVSEE = cut(OVSEE, breaks = quantile(unique(OVSEE), probs = seq(0, 1, 0.2)), include.lowest = TRUE, dig.lab = 2),
               Q_OVSAE = cut(OVSAE, breaks = quantile(unique(OVSAE), probs = seq(0, 1, 0.2)), include.lowest = TRUE, dig.lab = 3),
               Q_OVPT = cut(OVPT, breaks = quantile(unique(OVPT), probs = seq(0, 1, 0.2)), include.lowest = TRUE, dig.lab = 2),
               Q_PL500 = cut(PL.5000, breaks = quantile(unique(PL.5000), probs = seq(0, 1, 0.2)), include.lowest = TRUE, dig.lab = 3),
               Q_VHAC = cut(VHAC, breaks = quantile(unique(VHAC), probs = seq(0, 1, 0.2)), include.lowest = TRUE, dig.lab = 2),
               Q_PO2SM = cut(PO2SM, breaks = quantile(unique(PO2SM), probs = seq(0, 1, 0.2)), include.lowest = TRUE, dig.lab = 3)
               )

Para poder visualizar el mapa a nivel municipal dentro de conjunto de datos del Marco Geoestadístico Nacional [MGN], donde este conjunto abarca varias capas de datos vectoriales.
Se toma como referencia el Marco Geoestadístico, 2010 Versión 4.3, Marco Geoestadístico, junio 2016 y Marco Geoestadístico, 2020, que es producto de la actualización de censos, conteos, encuestas y registros administrativos.

Shape de municipios 2010

shape_estados_2010 <- readOGR(dsn = "MGN 2010/32_Entidades_Federativas",
                               layer = "ESTADOS",
                              #  encoding = "UTF-8",
                               #   use_iconv = TRUE
                              )

shape_municipios_2010 <- readOGR(dsn = "MGN 2010/32_Entidades_Federativas",
                                 layer = "MUNICIPIOS",
                               #  encoding = "UTF-8",
                                #   use_iconv = TRUE
                                )

# Se crea un ID único o una Clave Geoestadística 
shape_municipios_2010@data <- shape_municipios_2010@data %>%
                                mutate(CVE_GEO = paste0(CVE_ENT, CVE_MUN))

Shape municipio 2015

shape_estados_2015 <- readOGR(dsn = "MGN Junio 2016/conjunto_de_datos",
                               layer = "areas_geoestadisticas_estatales"
                              #  encoding = "UTF-8",
                                #  use_iconv = TRUE
                              )

shape_municipios_2015 <- readOGR(dsn = "MGN Junio 2016/conjunto_de_datos",
                                  layer = "areas_geoestadisticas_municipales"
                                #  encoding = "UTF-8",
                                 #   use_iconv = TRUE
                                )

shape_municipios_2015@data <- shape_municipios_2015@data %>%
                                mutate(CVE_GEO = paste0(CVE_ENT, CVE_MUN))

Se elimina al municipio 23011 Puerto Morelos del estado de Quintana Roo, debido a que no tiene calculado el IM 2015.

require(spdplyr)
shape_municipios_2015 <- shape_municipios_2015 %>%
                           spdplyr:::filter.Spatial(CVE_GEO != 23011)

Se unen las bases de los indicadores simples a los shapes

capas_estados_2010 <- fortify(shape_estados_2010, id = "CVE_ENT") 
capas_municipios_2010 <- dplyr::right_join(fortify(shape_municipios_2010, id = "CVE_GEO"), 
                                            Ind %>% dplyr::select(-CVE_ENT, -NOM_ENT, -NOM_MUN) %>%
                                             filter(ANIO == "2010") %>%
                                              arrange(.$CVE_MUN) %>%
                                               mutate(id = as.character(seq(0, to = (nrow(.) - 1)))),
                                                by =  c("id"))

capas_estados_2015 <- fortify(shape_estados_2015, id = "CVE_ENT") 
capas_municipios_2015 <- dplyr::right_join(fortify(shape_municipios_2015, id = "CVE_GEO"), 
                                            Ind %>% dplyr::select(-CVE_ENT, -NOM_ENT, -NOM_MUN) %>%
                                             filter(ANIO == "2015") %>%
                                              arrange(.$CVE_MUN) %>%
                                               mutate(id = as.character(seq(0, to = (nrow(.) - 1)))),
                                                by =  c("id"))

capas_estados_2020 <- fortify(shape_estados, id = "CVE_ENT") 
capas_municipios_2020 <- dplyr::right_join(fortify(shape_municipios, id = "CVE_GEO"), 
                                            Ind %>% dplyr::select(-CVE_ENT, -NOM_ENT, -NOM_MUN) %>%
                                             filter(ANIO == "2020") %>%
                                              arrange(.$CVE_MUN) %>%
                                               mutate(id = as.character(seq(0, to = (nrow(.) - 1)))),
                                                by =  c("id"))

Utilizando dos capas

En este tipo de gráficos se utiliza la capa de municipios y la capa de estados. Para cada caso se usa la función geom_polygon() del paquete ggplot2:: y así poder encontrar las delimitaciones de los municipios en su respectiva entidad.

Variables <- colnames(capas_municipios_2010)

p <- NULL
for(i in tablas){
#p[[i]] <- lapply(c(20:28), function(x){ 
p[[i]] <- lapply(c(20), function(x){ 
            ggplot() + 
              geom_polygon(get(paste0("capas_municipios_", i)) ,
                            mapping=aes(x = long, y = lat, group = group, fill = get(Variables[x])),
                             color = "transparent", 
                              size = 0.5,
                               alpha = 0.9) +
               geom_polygon(get(paste0("capas_estados_", i)),
                             mapping = aes(x = long, y = lat, group = group),
                              color = "dark grey",
                               fill = "transparent",
                                size = 0.5) +
                  expand_limits(x = get(paste0("capas_estados_", i))$long, 
                                  y = get(paste0("capas_estados_", i))$lat) + 
                    coord_equal() +
                     theme_void() + 
                      theme(plot.title = element_text(size = 15, family = "montserrat"),
                             plot.subtitle = element_text(size = 15, hjust = 0.08 , family = "montserrat"),
                              plot.margin = margin(1, 1, 1, 1, "cm"),
                               plot.caption = element_text(size = 10, hjust = 0.1, vjust = -1.5, family = "montserrat"), 
                                axis.title.x = element_text(size = 8, family = "montserrat"),
                                 axis.text.x = element_blank(),
                                  axis.text.y = element_blank(), 
                                   legend.text = element_text(size = 10, family = "montserrat"),
                                    legend.key.size = unit(0.6, "cm"),
                                     legend.position = c(0.8, 0.7)) + 
                           scale_fill_manual(values = usecol(pal_petrol, n = 5)) +
                             guides(color = FALSE) +
            labs(subtitle = paste(i),
                  x = "",
                   y = "",
                    fill = paste(Variables[x-9]))
    }
  )   
}

Mapas de los indicadores simples 2010-2020

# Titulo de los gráficos
#for(i in 1:9){
for(i in 1){
legend <- get_legend(p[["2010"]][[i]] + theme(legend.position = c(0.5, 0.5),
                                               legend.key.size = unit(1, "cm"),
                                                legend.title = element_blank()))

q <- ggarrange(p[["2010"]][[i]] + theme(plot.margin = margin(0.5, -1, -1, 0, "cm"),
                                         legend.position = 'none'), 
               p[["2015"]][[i]] + theme(plot.margin = margin(0.5, 0, -1, -1, "cm"),
                                        legend.position = 'none'), 
               p[["2020"]][[i]] + theme(plot.margin = margin(0.5, -1, 0, 0, "cm"),
                                         legend.position = 'none'), 
               legend, 
               ncol = 2, nrow = 2, widths = c(0.5, 0.5, 3, 0.6), heights = c(0.5, 0.5, 3, 0.6))
        
#path = paste0("Images/Mapa_", colnames(Ind)[6 + i], ".pdf")
#ggexport (annotate_figure(q, bottom =  text_grob(paste(stringr::str_wrap(Indicadores[i], 25, indent = 1, exdent = 1),"\n"),
 #                                                family = "montserrat", 
  #                                               face = "plain", 
   #                                              size = 18
    #                                             #hjust = 0,
     #                                            #vjust = 0
      #                                           )),
       #   width = 10, height = 8, res = 400, filename = path)
}

annotate_figure(q, bottom =  text_grob(paste(stringr::str_wrap(Indicadores[i], 25, indent = 1, exdent = 1),"\n"),
                                                 family = "montserrat", 
                                                 face = "plain", 
                                                 size = 20))

Librerías

Librerías que se usaron en el trabajo

##  [1] "sf"           "ggspatial"    "doMC"         "iterators"    "foreach"     
##  [6] "viridisLite"  "unikn"        "tmaptools"    "tmap"         "stringr"     
## [11] "purrr"        "readr"        "tidyr"        "tibble"       "tidyverse"   
## [16] "spdplyr"      "rgdal"        "sp"           "RColorBrewer" "readxl"      
## [21] "openxlsx"     "lubridate"    "kableExtra"   "knitr"        "gridExtra"   
## [26] "ggpubr"       "ggmap"        "ggrepel"      "ggthemes"     "ggplot2"     
## [31] "forcats"      "dplyr"        "showtext"     "showtextdb"   "sysfonts"

Referencias

Datos Abiertos de México - Índice de marginación (carencias poblacionales) por localidad, municipio y entidad. (2021). Retrieved February 13, 2022, from https://datos.gob.mx/busca/dataset/indice-de-marginacion-carencias-poblacionales-por-localidad-municipio-y-entidad

Wickham, H. (2022). A Grammar of Data Manipulation • dplyr. Retrieved February 12, 2022, from https://dplyr.tidyverse.org/

Garnett, R. (2018). cheatsheets/sf.pdf at main · rstudio/cheatsheets · GitHub. Retrieved February 12, 2022, from https://github.com/rstudio/cheatsheets/blob/main/sf.pdf

Lovelace, R., Nowosad, J., & Muenchow, J. (2019). Geocomputation with R. CRC Press. https://geocompr.robinlovelace.net/

This work by Diana Villasana Ocampo is licensed under a Creative Commons Attribution 4.0 International License.

Mapas por desagregaciones geográficas de México

Diana Villasana Ocampo

17/02/2022

Objetivo

Sistema de Información Geográfica

Especificación del sistema de coordenadas

Geospatial Data Abstraction Library (GDAL)

Base de datos

Índice de marginación

Estructura de un mapa

Entidad Federativa

.shp

IME_2020

Layers

\[SpatialPolygons \Rightarrow data.frame\]

\[SpatialPolygons \Rightarrow SpatialPolygons + Datos\]

\[Data.frame \Rightarrow SpatialPolygonsDataFrame\]

Municipio

.shp

IMM_2020

Layers

\[SpatialPolygons \Rightarrow data.frame\]

\[SpatialPolygons \Rightarrow SpatialPolygons + Datos\]

Zonas Metropolitanas

Localidad

.shp

IML_2020

Layers

\[SpatialPolygons \Rightarrow data.frame\]

\[SpatialPolygons \Rightarrow SpatialPolygons + Datos\]

AGEB

.shp

IMU_2020

Layers

\[SpatialPolygons \Rightarrow data.frame\]

\[SpatialPolygons \Rightarrow SpatialPolygons + Datos\]

Zona Metropolitana del Valle de México

Mapas por desagregación geográfica

Estructura básica de un mapa

Entidad Federativa

geom_map

geom_polygon

tmap

ggspatial

geom_sf

Municipio

geom_map

geom_polygon

tmap

geom_sf

Localidad

geom_polygon

Por estados

geom_sf

AGEB

geom_polygon

Zona Metropolitana del Valle de México (ZMVM)

geom_sf

Complementos

Entidad

Etiquetas

Annotation_scale

Municipios

Etiquetas

2010 - 2020

Librerías

Referencias

Geospatial Data Abstraction Library (`GDAL`)

`geom_map`

`geom_polygon`

`tmap`

`ggspatial`

`geom_sf`

`geom_map`

`geom_polygon`

`tmap`

`geom_sf`

`geom_polygon`

`geom_sf`

`geom_polygon`

`geom_sf`

`Annotation_scale`