Me llamo Andy y mi dataset público preferido es el Relevamiento de Usos de Suelo.
Sin más demoras, arranquemos…
Cargo los paquetes necesarios.
library(sf)
library(dplyr)
library(here)
library(data.table)
library(ggplot2)
library(ggpointdensity)
library(viridis)
library(ggmap)
Importo el RUS 2017. Pueden descomentar el comando fread si no lo tienen descargado.
rus <- read.csv(here("data", "relevamiento-usos-del-suelo-2017.csv"),
encoding = "UTF-8")
# rus <- fread("http://cdn.buenosaires.gob.ar/datosabiertos/datasets/relevamiento-usos-del-suelo/relevamiento-usos-del-suelo-2017.csv",
# encoding = 'UTF-8',
# sep = ",",
# header = TRUE,
# showProgress = TRUE,
# data.table = FALSE)
nrow(rus)
## [1] 561722
Como vemos, tiene más de 550.000 filas.
Miro las primeras filas.
rus %>% head()
## X Y SMP CALLE NUM TIPO1_16
## 1 -58.45563 -34.59689 047-040A-001C MONTENEGRO 124 E
## 2 -58.45570 -34.59694 047-040A-002 NUEVA ZELANDIA 4820 EP
## 3 -58.45575 -34.59701 047-040A-003 NUEVA ZELANDIA 4830 E
## 4 -58.45581 -34.59707 047-040A-004 NUEVA ZELANDIA 4840 GAP
## 5 -58.45588 -34.59713 047-040A-005 NUEVA ZELANDIA 4850 GAP
## 6 -58.45588 -34.59713 047-040A-005 NUEVA ZELANDIA 4850 E
## TIPO2_16 PISOS_16 NOMBRE OBSERVACIO BARRIO COMUNA X5_DIG
## 1 VIVIENDA 1 CHACARITA 15 3
## 2 METALURGIA 1 RODOCA S.A. CHACARITA 15 28900
## 3 CASA 1 E X GAP CHACARITA 15 4
## 4 GARAGE PRIVADO 1 CHACARITA 15 63312
## 5 GARAGE PRIVADO 1 CHACARITA 15 63312
## 6 CASA 2 CHACARITA 15 4
## X4_DIG X3_DIG X2_DIG
## 1 0 0 0
## 2 2890 289 28
## 3 0 0 0
## 4 6331 633 63
## 5 6331 633 63
## 6 0 0 0
## RAMA
## 1 EDIFICIOS
## 2 FABRICACIÓN DE PRODUCTOS ELABORADOS DE METAL, EXCEPTO MAQUINARIA Y EQUIPO
## 3 EDIFICIOS
## 4 SERVICIOS COMPLEMENTARIOS AL TRANSPORTE
## 5 SERVICIOS COMPLEMENTARIOS AL TRANSPORTE
## 6 EDIFICIOS
## SUBRAMA
## 1 VIVIENDA
## 2 FABRICACIÓN DE PRODUCTOS ELABORADOS DE METAL NCP
## 3 CASA
## 4 SERVICIOS COMPLEMENTARIOS AL TRANSPORTE TERRESTRE
## 5 SERVICIOS COMPLEMENTARIOS AL TRANSPORTE TERRESTRE
## 6 CASA
## SSRAMA
## 1
## 2
## 3
## 4 ESTACIONAMIENTOS, COCHERAS Y GARAJES
## 5 ESTACIONAMIENTOS, COCHERAS Y GARAJES
## 6
Veo los principales rubros.
rus %>%
select(TIPO1_16) %>%
table() %>%
as.data.frame %>%
arrange(desc(Freq))
## . Freq
## 1 E 270677
## 2 L 116862
## 3 GAP 115083
## 4 EDU 15389
## 5 LG 13186
## 6 EP 7973
## 7 LOTE 5209
## 8 G 4788
## 9 SU 3477
## 10 U 3274
## 11 GA 2818
## 12 P 2646
## 13 ES 340
rus %>%
select(TIPO2_16) %>%
table() %>%
as.data.frame() %>%
arrange(desc(Freq)) %>%
top_n(20)
## Selecting by Freq
## . Freq
## 1 VIVIENDA 149060
## 2 GARAGE PRIVADO 115083
## 3 CASA 97541
## 4 LOCAL CERRADO 27937
## 5 VU - VIVIENDA UNICA 15738
## 6 EDIFICIO PRODUCTIVO 6622
## 7 DEPOSITO 6544
## 8 OFICINAS 6403
## 9 I FEMENINA 6098
## 10 TALLER MECANICO DE AUTOMOTORES 3523
## 11 PELUQUERIA 3187
## 12 I 2951
## 13 GARAGE 2783
## 14 SUPERMERCADO 2628
## 15 INMOBILIARIA 2387
## 16 LOTE 2311
## 17 KIOSCO 2307
## 18 EDU 2176
## 19 ABANDONADO 2146
## 20 EN OBRA 2019
table(rus$TIPO1_16, rus$TIPO2_16) %>%
as.data.frame() %>%
filter(Freq > 0) %>%
arrange(desc(Freq)) %>%
top_n(20)
## Selecting by Freq
## Var1 Var2 Freq
## 1 E VIVIENDA 149060
## 2 GAP GARAGE PRIVADO 115083
## 3 E CASA 97541
## 4 L LOCAL CERRADO 25903
## 5 E VU - VIVIENDA UNICA 15738
## 6 EP EDIFICIO PRODUCTIVO 6622
## 7 L I FEMENINA 4139
## 8 L TALLER MECANICO DE AUTOMOTORES 3050
## 9 L PELUQUERIA 3029
## 10 G DEPOSITO 2895
## 11 GA GARAGE 2783
## 12 L SUPERMERCADO 2621
## 13 L OFICINAS 2582
## 14 E OFICINAS 2466
## 15 L INMOBILIARIA 2339
## 16 LOTE LOTE 2310
## 17 L KIOSCO 2254
## 18 EDU EDU 2176
## 19 E ABANDONADO 2120
## 20 EDU DEPOSITO 2096
Acá vemos diversidad de categorías, algunas más intuitivas que otras. ¿L es local? ¿I FEMENINA es indumentaria?
Veamos las principales categorías del tipo L.
rus %>%
filter(TIPO1_16 == "L") %>%
select(TIPO2_16) %>%
table() %>%
as.data.frame() %>%
filter(Freq > 0) %>%
arrange(desc(Freq)) %>%
top_n(20)
## Selecting by Freq
## . Freq
## 1 LOCAL CERRADO 25903
## 2 I FEMENINA 4139
## 3 TALLER MECANICO DE AUTOMOTORES 3050
## 4 PELUQUERIA 3029
## 5 SUPERMERCADO 2621
## 6 OFICINAS 2582
## 7 INMOBILIARIA 2339
## 8 KIOSCO 2254
## 9 MAXIKIOSCO 1965
## 10 RESTAURANTE 1670
## 11 I 1650
## 12 VERDULERIA 1411
## 13 BAR 1364
## 14 CALZADO 1332
## 15 LAVADERO DE ROPA 1227
## 16 DEPOSITO 1224
## 17 LOTERIA 1179
## 18 FERRETERIA 1084
## 19 ALMACEN 1065
## 20 MUEBLES (VENTA) 1065
Me quedo entonces con las principales 20 categorías de locales, quitando los locales cerrados, los depósitos y las oficinas.
categorias <- rus %>%
filter(TIPO1_16 == "L" & !(TIPO2_16 %in% c("LOCAL CERRADO", "DEPOSITO", "OFICINAS"))) %>%
select(TIPO2_16) %>%
table() %>%
as.data.frame() %>%
filter(Freq > 0) %>%
arrange(desc(Freq)) %>%
top_n(20)
## Selecting by Freq
names(categorias)[1] <- "categoria"
locales <- rus %>%
filter(TIPO1_16 == "L" & TIPO2_16 %in% categorias$categoria)
locales %>% head()
## X Y SMP CALLE NUM TIPO1_16
## 1 -58.45545 -34.59703 047-040A-019 MONTENEGRO 110 L
## 2 -58.45545 -34.59703 047-040A-019 MONTENEGRO 112 L
## 3 -58.45620 -34.59678 047-040B-004 ARENAL, CONCEPCION 1826 L
## 4 -58.45645 -34.59701 047-040B-008 ARENAL, CONCEPCION 4860 L
## 5 -58.45629 -34.59713 047-040B-017 NUEVA ZELANDIA 4877 L
## 6 -58.45610 -34.59696 047-040B-020 NUEVA ZELANDIA 4851 L
## TIPO2_16 PISOS_16 NOMBRE OBSERVACIO BARRIO
## 1 TALLER MECANICO DE AUTOMOTORES 1 CHACARITA
## 2 TALLER MECANICO DE AUTOMOTORES 1 CHACARITA
## 3 TALLER MECANICO DE AUTOMOTORES 1 CHACARITA
## 4 MAXIKIOSCO 1 CHACARITA
## 5 FERRETERIA 2 CHACARITA
## 6 KIOSCO 1 CHACARITA
## COMUNA X5_DIG X4_DIG X3_DIG X2_DIG RAMA
## 1 15 50209 5020 502 50 VENTA Y REP. DE VEHICULOS AUTOMOTORES
## 2 15 50209 5020 502 50 VENTA Y REP. DE VEHICULOS AUTOMOTORES
## 3 15 50209 5020 502 50 VENTA Y REP. DE VEHICULOS AUTOMOTORES
## 4 15 52119 5211 521 52 COMERCIO
## 5 15 52363 5236 523 52 COMERCIO
## 6 15 52119 5211 521 52 COMERCIO
## SUBRAMA
## 1 MANTENIMIENTO Y REPARACION
## 2 MANTENIMIENTO Y REPARACION
## 3 MANTENIMIENTO Y REPARACION
## 4 VENTA AL POR MENOR EXEPTO LA ESPECIALIZADA
## 5 VENTA AL POR MENOR ESPECIALIZADO NCP
## 6 VENTA AL POR MENOR EXEPTO LA ESPECIALIZADA
## SSRAMA
## 1
## 2
## 3
## 4 MAXIKIOSCO
## 5 MATERIALES DE CONSTRUCCION, FERRETERIA, PINTURERIA, VIDRIOS
## 6 KIOSCO
Dato interesante: aún existen las loterías.
Veo con cuántos locales me estoy quedando.
nrow(locales)
## [1] 35436
¡Hora de mapear!
Veamos la distribución espacial de estas categorías.
Primero debo pasar el objeto locales a tipo geográfico. Y también me traigo los barrios para tener una referencia en el espacio.
(“¿te gustan los mapas?” es una gran pregunta de primera cita que podría tomar de lxs geografxs)
locales <- locales %>%
st_as_sf(coords = c("X", "Y"),
crs = 4326,
remove = FALSE
)
barrios <- st_read("http://cdn.buenosaires.gob.ar/datosabiertos/datasets/barrios/barrios.geojson")
## Reading layer `barrios_badata' from data source `http://cdn.buenosaires.gob.ar/datosabiertos/datasets/barrios/barrios.geojson' using driver `GeoJSON'
## Simple feature collection with 48 features and 4 fields
## geometry type: POLYGON
## dimension: XY
## bbox: xmin: -58.53152 ymin: -34.70529 xmax: -58.33515 ymax: -34.52649
## epsg (SRID): 4326
## proj4string: +proj=longlat +datum=WGS84 +no_defs
El parámetro remove deja las 2 columnas de latitud, longitud. Presunto que serán usadas en el futuro.
Con esto ya puedo empezar a hacer algunos mapas básicos.
ggplot() +
geom_sf(data = barrios, aes()) +
geom_sf(data = locales, aes())
Marche al Instagram de los gráficos que son horribles en el primer intento.
Va de nuevo.
ggplot() +
geom_sf(data = barrios, aes()) +
geom_sf(data = locales,
aes(color = TIPO2_16),
size = 1,
alpha = 0.5)
Mejoró rotundamente, pero existe mucha superposición, y los puntos de arriba claramente tapan a los de abajo.
Veamos creando un plano por categoría de comercio.
ggplot() +
geom_sf(data = barrios, aes()) +
geom_sf(data = locales,
aes(color = TIPO2_16),
size = 1,
alpha = 0.5) +
facet_wrap(~TIPO2_16)
Mejoró, pero la leyenda es un poco molesta y los mapas se ven muy tapados, por lo que no puedo apreciar la densidad de comercios.
ggplot() +
geom_sf(data = barrios, aes()) +
geom_sf(data = locales,
aes(color = TIPO2_16),
size = 1,
alpha = 0.2) +
facet_wrap(~TIPO2_16) +
theme(legend.position = "none")
Si se agranda la imagen ya pueden verse algunos patrones, como el predominio de ciertas categorías en el centro de la ciudad, Recoleta y Palermo. También algunas categorías poseen una distribución más uniforme a lo largo de la ciduad, como kioskos, supermercados y almacenes.
Sin embargo, aún no percibo con claridad la densidad, porque la opacidad luego de cierto nivel de superposición es engañosa.
Veamos cómo sería si usamos geom_point en lugar de geom_sf.
Y la paleta viridis, porque puedo. <3
ggplot() +
geom_sf(data = barrios, aes()) +
geom_pointdensity(data = locales,
aes(x = X, y = Y),
size = 1,
adjust = 0.1) +
scale_color_viridis() +
facet_wrap(~TIPO2_16)
# vuelvo a dejar la leyenda
# theme(legend.position = "none")
El paquete ggpointdensity tiene la virtud de que uno plotea los puntos, y se colorean por la densidad.
Eso permite ver rápidamente qué puntos de la ciudad tienen mayor concentración de cada rubro.
El parámetro adjust permite modificar qué entorno se considera para obtener la densidad. Al reducirlo pueden aparece otras modas en el espacio que habían sido suavizadas, como en los rubros de indumentaria femenina o talleres mecánicos.
Un punto interesante es que para todos los plots utiliza la misma escala de colores, por lo que cada facetado es comparable con los demás.
Otro paquete que ofrece algo interesante es ggmap. Pero primero preciso un de fondo de la bella CABA.
bbox <- barrios %>%
st_bbox() %>%
as.numeric()
# si le paso a get_stamenmap() el tipo bbox no funciona -_-
caba <- get_stamenmap(bbox = bbox,
maptype = "toner-lite",
zoom = 12)
Logré bajar las tiles, a pesar de que habían pasado cosas con una función que pedía API. Alguien tiene la gorra puesta.
Pruebo stats_density_2d, y saco viridis porque muchos valores bajos me oscurece demasiado.
Quito el geom_sf pues pruebo ggmaps, y utilizo fill en lugar de color para la densidad. Y agrego labs para darle contexto.
ggmap(caba) +
stat_density_2d(data = locales,
aes(x = X, y = Y, fill = stat(level)),
alpha = 0.6,
geom = "polygon",
size = 0.5) +
scale_fill_distiller(palette = "RdYlGn") +
facet_wrap(~TIPO2_16) +
labs(title = "Densidad de comercios por rubro",
subtitle = "Ciudad Autónoma de Buenos Aires",
x = "",
y = "",
caption = "Fuente: Elaboración propia con datos públicos y abiertos.",
fill = "Nivel")
Estaba feucho, y logré cambiar la paleta a algo más alegre.
Intento algo similar, con la paleta spectral.
ggplot() +
geom_sf(data = barrios, aes(), alpha = 0.5) +
geom_density_2d(data = locales,
aes(x = X, y = Y, color = stat(level)),
alpha = 0.8,
size = 1) +
scale_color_distiller(palette = "Spectral") +
facet_wrap(~TIPO2_16) +
labs(title = "Densidad de comercios por rubro",
subtitle = "Ciudad Autónoma de Buenos Aires",
x = "",
y = "",
caption = "Fuente: Elaboración propia con datos públicos y abiertos.",
fill = "Nivel")
También me interesa entender la interacción entre los distintos rubros de locales. Pero, suspendo por hoy, queda para el próximo episodio.
No sin antes ser agradecido.
A mi querida compa Angie, por bancarme en mis voladas y ayudarme pensando juntos. Y bajarme a tierra cuando es necesario (lo hacés menos veces de las necesarias <3). ¡Y por el feedback de mi primer Rmd! Por más almuerzos dentro y fuera del escritorio.
A mi compañero Dauny por la confianza y la oportunidad inicial.
A Mr. Tony Vazquez Brust por compartir esos análisis que se me vienen a la mente cuándo pienso ¿cómo conviene hacer esto?. Y por ser autor de mis primeras lecturas de código en datos CABA, y un faro en el ecosistema geo+R+urban. Y por ser educador de mi educadora. Angie dice que soy una groupie cholula.
A Mr. M. Alalú por fomentar este tipo de publicaciones, y que se me adelantó a publicar una punta de esta idea en este link.