Ciencia de datos para ciudades 2

Módulo 1: Geoprocesamiento

library(tidyverse)

## ── Attaching packages ───────────── tidyverse 1.3.0 ──

## ✓ ggplot2 3.3.0     ✓ purrr   0.3.3
## ✓ tibble  3.0.1     ✓ dplyr   0.8.5
## ✓ tidyr   1.0.2     ✓ stringr 1.4.0
## ✓ readr   1.3.1     ✓ forcats 0.5.0

## ── Conflicts ──────────────── tidyverse_conflicts() ──
## x dplyr::filter() masks stats::filter()
## x dplyr::lag()    masks stats::lag()

library(sf)

## Linking to GEOS 3.7.2, GDAL 2.4.2, PROJ 5.2.0

I. Deberán elegir una ciudad en cualquier parte del mundo que les interese y que disponga de un portal de datos abiertos que ofrece un shapefile con sus barrios. Elegimos la ciudad de Barcelona y a continuación cargaremos el dataframe que brinda la página del Ayuntamiento con los límites de sus barrios y distritos (comunas).

Barcelona_distritos <- st_read("bcn_UNITATS_ADM_POLIGONS.json")

## Reading layer `bcn_UNITATS_ADM_POLIGONS' from data source `/Users/lousil/Google Drive/Formación Académica/ECONOMIA URBANA/14. CIENCIA DE DATOS PARA CIUDADES/Ciencia de Datos para Ciudades I/Clase 1/Ciencia de Datos para ciudades 2/bcn_UNITATS_ADM_POLIGONS.json' using driver `GeoJSON'
## Simple feature collection with 1501 features and 35 fields
## geometry type:  MULTIPOLYGON
## dimension:      XY
## bbox:           xmin: 420812.5 ymin: 4574282 xmax: 435480.4 ymax: 4591066
## CRS:            25831

Vemos que el sistema de coordenadas de referencia (CRS) es 25831, por lo tanto no es Merkator. Vamos a cambiarlo:

Barcelona_distritos <- st_transform(Barcelona_distritos, crs = 4326)

names(Barcelona_distritos)

##  [1] "FID"        "ID_ANNEX"   "ANNEXDESCR" "ID_TEMA"    "TEMA_DESCR"
##  [6] "ID_CONJUNT" "CONJ_DESCR" "ID_SUBCONJ" "SCONJ_DESC" "ID_ELEMENT"
## [11] "ELEM_DESCR" "NIVELL"     "NDESCR_CA"  "NDESCR_ES"  "NDESCR_EN" 
## [16] "TERME"      "DISTRICTE"  "BARRI"      "AEB"        "SEC_CENS"  
## [21] "GRANBARRI"  "ZUA"        "AREA_I"     "LITERAL"    "PERIMETRE" 
## [26] "AREA"       "CODI_UA"    "TIPUS_UA"   "NOM"        "WEB1"      
## [31] "WEB2"       "WEB3"       "FHEX_COLOR" "Shape_Leng" "Shape_Area"
## [36] "geometry"

Como el dataset posee 36 variables, de las cual necesitamos solo unas pocas, seleccionaremos y filtraremos solo lo que necesitamos para disminuir el gran tamaño del archivo:

Barcelona_distritos <- Barcelona_distritos %>% 
  filter(!(DISTRICTE =="-")) %>% 
  group_by(DISTRICTE) %>% 
  summarise()

Del mismo portal de datos elegirán un dataset con registros geo-referenciados. Por ejemplo, las escuelas de la ciudad (o las comisarías, o las propiedades en alquiler) con sus coordenadas. Elegimos un dataset que contiene los museos y bibliotecas:

Equipamientos <- read.csv("bcn_Biblioteques_i_museus.csv")

Por alguna razón que desconocemos, aparecen varias observaciones del mismo equipamiento. Por tal motivo, agruparemos segun equipamiento para que nos quede una observación para cada uno y luego las cuente correctamente.

Equipamientos <- Equipamientos %>% 
  group_by(EQUIPAMENT, LONGITUD, LATITUD) %>% 
  summarise()

Visualicemos ambos datasets superpuestos:

ggplot()+
  geom_sf(data = Barcelona_distritos, aes(fill=DISTRICTE), color = NA) + 
  geom_point(data = Equipamientos, aes(x=LONGITUD, y=LATITUD), alpha=.8) +
         labs(title = "Equipamiento por distrito ",
         subtitle = "(museos, bibliotecas, universidades, etc)",
         fill = "Distritos",
         caption= "Fuente: Ajuntament de Barcelona | 2020",
         y="",
         x="")

Realizando un join espacial, asignar a cada registro geo-referenciado el barrio/distrito que le corresponde:

Convertimos el dataset de Bibliotecas y Museos a espacial:

Equipamientos <- Equipamientos %>%
    st_as_sf(coords = c("LONGITUD", "LATITUD"), crs = 4326)

…Y unimos al dataset que contiene los límites de barrios y distritos (comunas) de Barcelona:

Equipamientos <- st_join(Barcelona_distritos, Equipamientos)

## although coordinates are longitude/latitude, st_intersects assumes that they are planar

Utilizando ggplot():

Realizar un gráfico para mostrar los resultados de cantidad de bibliotecas y museos por distrito:

Equipamientos_distritos <- Equipamientos %>% 
      group_by(DISTRICTE) %>%
      summarise(Frecuencia = n())

ggplot(Equipamientos_distritos) +
  geom_bar(aes( x = DISTRICTE, weight = Frecuencia), fill="#69b3a2") +
  labs(title = "Cantidad de Equipamiento por distrito",
         subtitle = "(museos, bibliotecas, universidades, etc)",
         caption= "Fuente: Ajuntament de Barcelona | 2020",
         y="Cantidad",
         x="Distritos")

Ahora ponderaremos según cantidad de Museos y Bibliotecas por distrito con scale_viridis:

ggplot(Equipamientos_distritos) +
  geom_bar(aes( x = DISTRICTE, weight = Frecuencia, fill=Frecuencia)) +
  scale_fill_viridis_c() +
  labs(title = "Cantidad de Equipamiento por distrito",
       subtitle = "(museos, bibliotecas, universidades, etc)",
       fill = "Escala",
       caption= "Fuente: Ajuntament de Barcelona | 2020",
       y="Cantidad",
       x="Distritos")

Realizar un mapa con los límites de los distritos, cuyo color de relleno indique la cantidad encontrada en cada uno:

ggplot() + 
  geom_sf(data = Equipamientos_distritos, aes(fill=Frecuencia), color = NA) + 
  geom_sf_text(data = Equipamientos_distritos, aes(label = DISTRICTE), size = 2, colour = "white") +
  scale_fill_viridis_c() +
  labs(title = "Cantidad de Equipamiento por distrito",
         subtitle = "(museos, bibliotecas, universidades, etc)",
         fill = "Escala",
         caption= "Fuente: Ajuntament de Barcelona | 2020",
         y="",
         x="")

## Warning in st_point_on_surface.sfc(sf::st_zm(x)): st_point_on_surface may not
## give correct results for longitude/latitude data

Módulo 2: Acceso a información urbana georeferenciada en repositorios online

I. Descargar de OpenStreetMap la grilla de calles para la Ciudad elegida en el ejercicio 1 y mapearla por uno de sus atributos (velocidad mínima, velocidad máxima, cantidad de carriles, etc).

library(osmdata)

## Data (c) OpenStreetMap contributors, ODbL 1.0. https://www.openstreetmap.org/copyright

library(leaflet)

Barcelona <- getbb("Barcelona, España")

Barcelona

##         min       max
## x  2.052498  2.228356
## y 41.317035 41.467914

ciudad <- "Barcelona"
posicion <- 1
url <- URLencode(paste0("https://nominatim.openstreetmap.org/search.php?q=",
ciudad, "&polygon_geojson=1&format=geojson"))
destfile = paste0(tempdir(), "/limits.json")
download.file(url, destfile)
Barcelona_frontera <- st_read(destfile)[posicion,]

## Reading layer `limits' from data source `/private/var/folders/lg/vhfpqlf52vj3nv4pgd4px5pm0000gn/T/RtmpX9iYKY/limits.json' using driver `GeoJSON'
## Simple feature collection with 10 features and 9 fields
## geometry type:  GEOMETRY
## dimension:      XY
## bbox:           xmin: -80.6871 ymin: -6.069389 xmax: 124.1419 ymax: 50.3547
## CRS:            4326

leaflet(Barcelona_frontera) %>% 
  addTiles() %>% 
  addPolygons()

A continuación construiremos una consulta para que R nos de un Shapefile con las calles de la ciudad:

Barcelona_calles <- opq(Barcelona) %>% 
  add_osm_feature(key="highway")

Creamos un objeto osmdata:

Barcelona_calles <- Barcelona_calles %>% 
  osmdata_sf()

Barcelona_calles

## Object of class 'osmdata' with:
##                  $bbox : 41.3170353,2.0524977,41.4679135,2.2283555
##         $overpass_call : The call submitted to the overpass API
##                  $meta : metadata including timestamp and version numbers
##            $osm_points : 'sf' Simple Features Collection with 243987 points
##             $osm_lines : 'sf' Simple Features Collection with 45832 linestrings
##          $osm_polygons : 'sf' Simple Features Collection with 1603 polygons
##        $osm_multilines : NULL
##     $osm_multipolygons : 'sf' Simple Features Collection with 87 multipolygons

Extraemos las calles:

Barcelona_calles <- Barcelona_calles$osm_lines

Y graficamos:

ggplot() + 
  geom_sf(data = Barcelona_calles, color = "gray20")

A continuación haremos una interesección entre los límites de la ciudad y las calles obtenidas a través de OSM:

Barcelona_calles <- st_intersection(Barcelona_calles, Barcelona_frontera)

## although coordinates are longitude/latitude, st_intersection assumes that they are planar

## Warning: attribute variables are assumed to be spatially constant throughout all
## geometries

Y graficamos:

ggplot()+
  geom_sf(data=Barcelona_calles, color = "gray20")

Para mapearla por el atributo de velocidad máxima, haremos:

Barcelona_callesmax <- Barcelona_calles %>% 
  mutate(maxspeed = as.numeric(maxspeed),
         lanes = ifelse(is.na(lanes), 1, as.numeric(lanes)))

Y graficamos:

ggplot(Barcelona_callesmax) +
    geom_sf(aes(color = maxspeed), alpha = 0.5) +
    scale_color_viridis_c() +
      theme_void() +
    labs(title = "Barcelona",
         subtitle = "Vías de circulación",
         caption = "fuente: OpenStreetMap",
         color = "velocidad máxima")

Descargar de OpenStreetMap una (o más) capas de datos de tipo puntos o polígonos.

Elegimos descargar la capa de datos de las bibliotecas:

Barcelona_bibliotecas <- opq(Barcelona) %>% 
  add_osm_feature(key = "amenity", value = "library") %>%
  osmdata_sf()

Y la de museos:

Barcelona_museos <- opq(Barcelona) %>% 
  add_osm_feature(key = "tourism", value = "museum") %>%
  osmdata_sf()

Y la de universidades:

Barcelona_universidades <- opq(Barcelona) %>% 
  add_osm_feature(key = "amenity", value = "university") %>%
  osmdata_sf()

1. Proyectar los datos descargados en un mapa y comentar los resultados: ¿Cómo se distribuyen en la Ciudad?

Barcelona_bibliotecas <- Barcelona_bibliotecas$osm_points
Barcelona_museos <- Barcelona_museos$osm_points
Barcelona_universidades <- Barcelona_universidades$osm_points

Barcelona_bibliotecas <- st_intersection(Barcelona_bibliotecas, Barcelona_frontera) %>% 
  select(name, geometry) %>% 
  filter(!is.na(name))

## although coordinates are longitude/latitude, st_intersection assumes that they are planar

## Warning: attribute variables are assumed to be spatially constant throughout all
## geometries

Barcelona_museos <- st_intersection(Barcelona_museos, Barcelona_frontera) %>% 
  select(name, geometry) %>% 
  filter(!is.na(name))

## although coordinates are longitude/latitude, st_intersection assumes that they are planar

## Warning: attribute variables are assumed to be spatially constant throughout all
## geometries

Barcelona_universidades <- st_intersection(Barcelona_universidades, Barcelona_frontera) %>% 
  select(name, geometry) %>% 
  filter(!is.na(name))

## although coordinates are longitude/latitude, st_intersection assumes that they are planar

## Warning: attribute variables are assumed to be spatially constant throughout all
## geometries

ggplot()+
  geom_sf(data = Barcelona_calles, color = "gray20") +
  geom_sf(data = Barcelona_bibliotecas, color = "purple") +
  geom_sf(data = Barcelona_museos, color = "salmon") +
  geom_sf(data = Barcelona_universidades, color = "royalblue") +
   labs(title = "Barcelona",
         subtitle = "Bibliotecas, Museos y Universidades",
         caption = "fuente: OpenStreetMap")

Como podemos observar, los museos parecen estar agrupados en lugares específicos y turísticos, como ser en la antigua ciudad amurallada (Ciutat Vella) y la Rambla, así como en la montaña Montjuic. Respecto a las universidades, podemos ver un gran concentración en un extremo de la ciudad, siendo la menor distribución del lado opuesto. Sin embargo, las bibliotecas tienen una distribución más pareja en el territorio. Podemos suponer que exista alguna política pública que asegure una distribución pareja de bibliotecas por barrio.

Cargamos las siguientes librerías:

library(devtools)

## Loading required package: usethis

library(ggsflabel)

## 
## Attaching package: 'ggsflabel'

## The following objects are masked from 'package:ggplot2':
## 
##     geom_sf_label, geom_sf_text, StatSfCoordinates

Y ahora agregamos las etiquetas, sólo para el caso de las bibliotecas:

ggplot() +
  geom_sf(data = Barcelona_calles, color = "gray20") +
  geom_sf_label_repel(data = Barcelona_bibliotecas, aes(label = name), size = 2, color = "purple") +
  theme_void() +
  labs(title = "Barcelona",
       subtitle = "Bibliotecas con etiquetas",
       caption = "fuente: OpenStreetMap")

## Warning in st_point_on_surface.sfc(data$geometry): st_point_on_surface may not
## give correct results for longitude/latitude data

1. Hacer un conteo de los ítems de la capa descargada por barrio, mapearlo y compararlo con el conteo de los ítems descargados en el ejercicio anterior: ¿La distribución es similar o hay diferencias? ¿A qué se puede deber?

Primero que nada vamos a unir los shapes de bibliotecas, museos y universidades en un solo:

Barcelona_bibliotecas_museos_universidades <- rbind(Barcelona_bibliotecas, Barcelona_museos, Barcelona_universidades)

Haremos un join espacial entre el shape de bibliotecas, universidades y museos unido anteriormente y el shape que contiene los distritos:

osm_bmu_distritos <- st_join(Barcelona_bibliotecas_museos_universidades, Barcelona_distritos)

## although coordinates are longitude/latitude, st_intersects assumes that they are planar
## although coordinates are longitude/latitude, st_intersects assumes that they are planar

Hacemos un conteo agrupando por distrito y sacándole la geometría:

osm_bmu_distritos <- osm_bmu_distritos %>%
group_by(DISTRICTE) %>%
summarise(Frecuencia = n()) %>%
st_set_geometry(NULL)

Con left_join() unimos el conteo anterior con el shape que tiene los límites de los distritos:

osm_bmu_distritos <- left_join(Barcelona_distritos, osm_bmu_distritos, by="DISTRICTE")

Y graficamos:

ggplot() +
  geom_sf(data = osm_bmu_distritos, aes(fill = Frecuencia), color=NA) +
  geom_sf_text(data = osm_bmu_distritos, aes(label = Frecuencia), size = 2, colour = "white") +
  scale_fill_viridis_c() +
  labs(title = "Barcelona",
       subtitle = "Bibliotecas, Museos y Universidades por Distritos",
       caption = "fuente: OpenStreetMap",
        y="",
        x="")

## Warning in st_point_on_surface.sfc(data$geometry): st_point_on_surface may not
## give correct results for longitude/latitude data

Vemos que la frecuencia máxima es de 31 para el Distrito Nº1, coloreado en amarillo.

Comparemos ahora con el gráfico hecho anteriormente:

ggplot() + 
  geom_sf(data = Equipamientos_distritos, aes(fill=Frecuencia), color = NA) + 
  geom_sf_text(data = Equipamientos_distritos, aes(label = Frecuencia), size = 2, colour = "white") +
  scale_fill_viridis_c() +
  labs(title = "Barcelona",
         subtitle = "Bibliotecas, Museos y Universidades por Distritos",
         fill = "Escala",
         caption= "Fuente: Ajuntament de Barcelona | 2020",
         y="",
         x="")

## Warning in st_point_on_surface.sfc(data$geometry): st_point_on_surface may not
## give correct results for longitude/latitude data

En este mapa, hecho en el módulo anterior con data abierta que brinda el Ayuntamiento de Barcelona, la Frecuencia máxima llega a los 70. Además, la cantidad de observaciones en cada Distrito es mayor a la observada con los datos de OSM. Por tal motivo entendemos que los datos brindados por el ayuntamiento incluyen otros equipamientos, además de bibliotecas museos y universidades, y esto hace que difiera la frecuencia en cada distrito.

Comprobémemoslo ahora numéricamente:

summary(Equipamientos_distritos)

##    DISTRICTE   Frecuencia             geometry
##  01     :1   Min.   : 7.00   MULTIPOLYGON :2  
##  02     :1   1st Qu.:10.75   POLYGON      :8  
##  03     :1   Median :20.50   epsg:4326    :0  
##  04     :1   Mean   :26.20   +proj=long...:0  
##  05     :1   3rd Qu.:26.00                    
##  06     :1   Max.   :70.00                    
##  (Other):4

En este dataset hecho en el módulo Nº1, la frecuencia máxima alcanza las 70 observaciones en el Distrito 1 y además el promedio es 26.20.

summary(osm_bmu_distritos)

##    DISTRICTE          geometry   Frecuencia   
##  01     :1   MULTIPOLYGON :2   Min.   : 2.00  
##  02     :1   POLYGON      :8   1st Qu.: 3.50  
##  03     :1   epsg:4326    :0   Median : 7.00  
##  04     :1   +proj=long...:0   Mean   :10.60  
##  05     :1                     3rd Qu.:17.25  
##  06     :1                     Max.   :31.00  
##  (Other):4

En cambio, en este dataset hecho con información de OSM vemos que la frecuencia máxima es 31 y el promedio es 10,60.

Ciencia de datos para ciudades 2

Cuviello_Silano

03/05/2020

Módulo 1: Geoprocesamiento

Módulo 2: Acceso a información urbana georeferenciada en repositorios online