Mapeando Ingresos en las Ciudades de Uruguay

· 2019/04/13 · 5 minute read

En un post anterior hice un mapa de Montevideo combinando capas de Open Street View y Shapefiles publicados por el INE. Hace un tiempo que estoy trabajando sobre un curso de Charlotte Wickham sobre mapas y pensé en replicar el proyecto final del curso con datos de Uruguay.

La principal pregunta que me interesa contestar es ¿cómo es la distribución del ingreso en las ciudades de Uruguay? En Montevideo sabemos que la zona costera es donde se concentran los hogares de mayores ingresos (Carrasco, Punta Gorda, Pocitos, etc.).

Para visualizar esta estructura, linkeo los datos de la Encuesta Continua de Hogares con un mapa de los segmentos censales del país.

## OGR data source with driver: ESRI Shapefile 
## Source: "/Users/rafa/blog_data/mapas_vectoriales_ine", layer: "ine_seg_11"
## with 4313 features
## It has 13 fields

El mapa

readOGR lee los archivos publicados por el INE como un objecto SpatialPolygonsDataFrame.

# Cargar mapa
ine_seg_11 <- readOGR(dsn = '~/blog_data/mapas_vectoriales_ine/', 
                      layer = 'ine_seg_11')
## OGR data source with driver: ESRI Shapefile 
## Source: "/Users/rafa/blog_data/mapas_vectoriales_ine", layer: "ine_seg_11"
## with 4313 features
## It has 13 fields
class(ine_seg_11)
## [1] "SpatialPolygonsDataFrame"
## attr(,"package")
## [1] "sp"

Este objeto es un data.frame que además de los tipo de datos usuales (AREA, PERIMETER, departamento, localidad) tiene una columna con la forma y la ubicación de cada segmento censal.

head(ine_seg_11@data)
##        AREA PERIMETER DEPTO SECCION SEGMENTO LOCALIDAD CODSEC CODSEG
## 0  14476.50 1785.1361     1       0        0        20    100 100000
## 1 122200.01 1563.4714     1       1        1        20    101 101001
## 2 151608.67 1676.3479     1       1        2        20    101 101002
## 3  99377.28 1456.9389     1       1        3        20    101 101003
## 4  54395.51 1040.6595     1       1      104        20    101 101104
## 5  41055.69  838.3599     1       1      105        20    101 101105
##   CODLOC  NOMBDEPTO    NOMBLOC CDEPTO_ISO CLOC_ISO
## 0   1020 MONTEVIDEO MONTEVIDEO       UYMO  UYMOMON
## 1   1020 MONTEVIDEO MONTEVIDEO       UYMO  UYMOMON
## 2   1020 MONTEVIDEO MONTEVIDEO       UYMO  UYMOMON
## 3   1020 MONTEVIDEO MONTEVIDEO       UYMO  UYMOMON
## 4   1020 MONTEVIDEO MONTEVIDEO       UYMO  UYMOMON
## 5   1020 MONTEVIDEO MONTEVIDEO       UYMO  UYMOMON

Para saber más sobre estos objetos, recomiendo el capítulo 2 de Geocomputation with R.

La ECH tiene los ingresos de los hogares identificados por segmento censal, por lo que podemos calcular la mediana de ingresos por segmento con group_by y summarize de dplyr:

Sys.setlocale("LC_ALL", "UTF-8")
## [1] "en_US.UTF-8/UTF-8/en_US.UTF-8/C/en_US.UTF-8/en_US.UTF-8"
ingresos_ech <- read_tsv("~/blog_data/H_2018_TERCEROS.dat") %>% 
  transmute(id       =as.character(numero),
            depto_id = dpto, 
            depto = nomdpto,
            localidad_id = as.numeric(locagr),
            localidad = nom_locagr,
            seccion_id = as.numeric(secc),
            segmento_id = as.numeric(segm),
            ingreso = YSVL)

medianas <- ingresos_ech %>% 
    group_by(depto_id, seccion_id, segmento_id) %>% 
    summarize(mediana = median(ingreso))

Elegir una localidad

Una vez calculada esta estimación, elijo una localidad para subsetear el mapa.

# Salto es la localidad 1120
cod_loc <- "15120"
mapa_loc <- ine_seg_11[ine_seg_11$CODLOC == cod_loc,]

Pegar las estimaciones

Para visualizar las medianas por segmento, le agrego la columna con la mediana de cada segmento usando depto_id, seccion_id y segmento_id como identificadores:

mapa_loc@data <- mapa_loc@data %>% 
    left_join(medianas, by = c("DEPTO" = "depto_id",
                               "SECCION" = "seccion_id",
                               "SEGMENTO" = "segmento_id"))

En el post anterior usé leaflet para hacer los mapas. Como Charlotte usa tmap voy a usar esa librería para hacer el choropleth:

tm_shape(mapa_loc) +
    tm_fill(col = "mediana", 
            style="quantile", 
            textNA = "S/D",
            colorNA = "gray80",
            labels = c("0-20%", "20%-40%", "40%-60%", "60%-80%", "80%-100%"),
            title = "Quintiles de Ingreso") +
    tm_layout("Salto", legend.position = c("right", "top"),
              legend.text.size = 0.6,
              legend.title.size = 0.8, 
              title.size = 1)

Una función

Para poder hacer el mapa para varias localidades, lo puse todo en una función, a la que hay que pasarle el nombre y el código de la localidad para que haga el mapa:

mapa_ingresos <- function(cod_loc, nom_loc, legend_pos = c("right","top")) {
  
  # segmento censal
  ine_seg_11 <- readOGR(dsn = '~/blog_data/mapas_vectoriales_ine/', layer = 'ine_seg_11')
  crswgs84 <- CRS("+proj=longlat +ellps=WGS84 +datum=WGS84 +no_defs")
  # proyeccion ine
  prj <- CRS("+proj=utm +zone=21 +south +ellps=WGS84 +datum=WGS84 +units=m +no_defs")
  proj4string(ine_seg_11) <- prj
  #ine_seg_11 <- spTransform(states, crswgs84)
  
  Sys.setlocale("LC_ALL", "UTF-8")
  ingresos_ech <- read_tsv("~/blog_data/H_2018_TERCEROS.dat") %>% 
    transmute(id       =as.character(numero),
              depto_id = dpto, 
              depto = nomdpto,
              localidad_id = as.numeric(locagr),
              localidad = nom_locagr,
              seccion_id = as.numeric(secc),
              segmento_id = as.numeric(segm),
              ingreso = YSVL)
  
  medianas <- ingresos_ech %>% 
    group_by(depto_id, seccion_id, segmento_id) %>% 
    summarize(mediana = median(ingreso))
  
  mapa_loc <- ine_seg_11[ine_seg_11$CODLOC == cod_loc,]
  
  mapa_loc@data <- mapa_loc@data %>% 
    left_join(medianas, by = c("DEPTO" = "depto_id",
                               "SECCION" = "seccion_id",
                               "SEGMENTO" = "segmento_id"))
  
  tm_shape(mapa_loc) +
    tm_fill(col = "mediana", 
            style="quantile", 
            textNA = "S/D",
            colorNA = "gray80",
            labels = c("0-20%", "20%-40%", "40%-60%", "60%-80%", "80%-100%"),
            title = "Quintiles de Ingreso") +
    tm_layout(nom_loc, legend.position = legend_pos,
              legend.text.size = 0.6,
              legend.title.size = 0.8, 
              title.size = 1)
}

Para sacar varios mapas a la vez podemos usar tmap_arrange:

tmap_arrange(
            mapa_ingresos("14320", "Rocha", c("left", "bottom")),
            mapa_ingresos("6220", "Durazno", c("left", "bottom")),
            mapa_ingresos("18220", "Tacuarembó", c("left", "bottom")),
            mapa_ingresos("9220", "Minas", c("left", "bottom")))
## OGR data source with driver: ESRI Shapefile 
## Source: "/Users/rafa/blog_data/mapas_vectoriales_ine", layer: "ine_seg_11"
## with 4313 features
## It has 13 fields
## OGR data source with driver: ESRI Shapefile 
## Source: "/Users/rafa/blog_data/mapas_vectoriales_ine", layer: "ine_seg_11"
## with 4313 features
## It has 13 fields
## OGR data source with driver: ESRI Shapefile 
## Source: "/Users/rafa/blog_data/mapas_vectoriales_ine", layer: "ine_seg_11"
## with 4313 features
## It has 13 fields
## OGR data source with driver: ESRI Shapefile 
## Source: "/Users/rafa/blog_data/mapas_vectoriales_ine", layer: "ine_seg_11"
## with 4313 features
## It has 13 fields

Lo que más me gustaría mejorar de estos mapas es poder agregar referencias como calles importantes o cuerpos de agua que permitan ubicarse mejor en la geografía de la ciudad. Dejo eso para otro post.

La estimación de las medianas probablemente deja mucho que desear desde el punto de vista estadístico también.