libs <- c("tidyverse", "sf", "tidygeocoder", "ggspatial", "readr")
installed_libs <- libs %in% rownames(installed.packages())
if (any(installed_libs == F)) {install.packages(libs[!installed_libs])}
invisible(lapply(libs, library, character.only = T))
rm(libs, installed_libs)
estaciones <- read_delim(file = "../data/2526020007/estaciones.csv",
delim = ";",
locale = locale(encoding = "ISO-8859-1",
decimal_mark = ","),
show_col_types = FALSE)Datos de todas las estaciones donde opera Renfe
ESTACIONES TREN
Input
Se ha obtenido el listado completo de estaciones en las que opera Renfe de la página web de dicha empresa. En esta web están disponibles los datos tanto en formato .csv como en formato .xlsx. Las observaciones disponibles son el código y nombre de la estación; la latitud y longitud de la parada; la dirección (compuesta por calle y número); el código postal, la población, provincia y país (ya que hay algunas paradas en Portugal y Francia); si son paradas de cercanías o no; si la parada pertenecía a los antiguos Ferrocarriles de Vía Estrecha, que utilizan vías más estrechas que las convencionales; y una columna con alias para la población o provincia, que está vacía para la gran mayoría de observaciones.
Si bien este dataset contiene información espacial como las direcciones y las latitudes y longitudes, sería interesante crear un dataframe de tipo sf, con geometrías que permitirían más análisis y visualizaciones. En este caso, aplicaré técnicas de limpieza de datos para facilitar la aplicación de técnicas de geocoding, como la eliminación de registros con datos faltantes que hacen imposible dicha geolocalización.
Descripción
Primero cargamos las librerías necesarias y leemos el dataset en formato .csv que hemos descargado previamente.
Ahora vamos a analizar el dataset, compuesto por 1680 registros, y los datos faltantes de este. Primero observamos que hay 280 observaciones a las que le falta dirección, pero de estas, hay 223 que sí que tienen latitud y longitud, por lo que podemos aplicarles reverse geocoding y escribir nosotros la dirección. Las 57 estaciones restantes no ofrecen dirección ni latitud y longitud, por lo que es imposible geocodificarlas con la información ofrecida por Renfe en estos momentos.
Tratamiento
Para el reverse geocoding, se ha usado la función reverse_geocode() del paquete tidygeocoder con el método “arcgis”. La dirección obtenida la almacenaremos en una columna llamada “dir”. Hay una estación a la que la función no ha podido encontrar la dirección.
Una vez hemos escrito algunas direcciones, vamos a seguir con todas aquellas observaciones que tienen dirección escrita, que en este caso son 1588, y usaremos la función geo() del mismo paquete. Antes de eso, creamos la variable “dir” para las variables que sí tenían dirección escrita, que une la dirección (calle y número) con la población, provincia y país. Además, como el formato no está unificado, vamos a “limpiar” esas direcciones para poder encontrar dobles y así hacer el proceso de geocoding más rápido al pasarle solo direcciones únicas.
estaciones1 <- estaciones %>%
filter(!is.na(DIRECION)) %>%
mutate(dir = paste0(DIRECION,", ", POBLACION, ", ", PROVINCIA, ", ", PAIS)) %>%
bind_rows(reverse)
estaciones_limpias <- estaciones1 %>%
mutate(dir_limpia = dir %>%
str_to_lower() %>% # minúsculas
str_replace_all("-", " ") %>% # quitar guiones
str_replace_all("\\s+", " ") %>% # espacios múltiples → 1
str_trim())
dirs_unicas <- estaciones_limpias %>%
distinct(dir_limpia)
arcgis <- geo(address = dirs_unicas$dir_limpia, method = "arcgis")
sum(is.na(arcgis$lat))
arcgis_sf <- st_as_sf(arcgis, coords=c("long","lat"), crs=4326, na.fail=F) Por último, volvemos a montar el dataset con la información nueva. Para ello, los pasos son los siguientes:
- Juntar las geometrías de las direcciones únicas que acabamos de encontrar con el dataset de estaciones “limpias”, es decir, aquellas estaciones de las que podemos saber sus direcciones. De esta forma, todas las paradas tienen su geometría aunque su dirección esté duplicada.
- Encontrar cuáles son las estaciones de las que no hemos podido conseguir dirección porque no tenían ni dirección escrita ni latitud y longitud y transformarlas a sf (aunque su geometría estará llena de NAs).
- Juntar en un solo sf todas las estaciones, para tener nuevamente las 1680 observaciones. 1623 de ellas tendrán su geometría y las otras 57 no.
estaciones_sf <- estaciones_limpias %>%
left_join(arcgis_sf, by = c("dir_limpia"= "address")) %>%
st_as_sf()
nrow(estaciones_sf)
estaciones_sin_geo <- estaciones %>%
anti_join(estaciones_sf %>% as.data.frame(), by = "CODIGO")
estaciones_sin_geo_sf <- estaciones_sin_geo %>%
st_as_sf(coords = c("LONGITUD", "LATITUD"), crs = 4326, na.fail = FALSE)
nrow(estaciones_sin_geo_sf)
estaciones_unificadas <- bind_rows(estaciones_sf, estaciones_sin_geo_sf)También puede ser interesante visualizar estos datos que acabamos de obtener.
ggplot() +
annotation_map_tile(type="cartolight", zoom = 6) +
geom_sf(data = estaciones_unificadas, color = rgb(0.56, 0.1, 0.46), size = 1) +
theme_void()Zoom: 6
Output
El output obtenido es el dataset original, pero con dos diferencias. La primera es que, para aquellas estaciones de las que solo se disponía de latitud y longitud se ha insertado la dirección escrita. La segunda es la columna de geometría, transformando el dataframe en un objeto espacial.
st_write(estaciones_unificadas, "estaciones_finales.gpkg", delete_dsn = TRUE)El fichero generado con este procedimiento se puede descargar de aquí.
Proyecto de Innovación Educativa Emergente (PIEE-3898312)