In this gist we show with a reproducible example how to create an animation of public transport networks using GTFS data in R. We use a few packages to do this. One of the core packages here is the new gtfs2gps. The gtfs2gps package converts public transport data in GTFS format to GPS-like records in a data.frame/data.table, which we will be using to create a .gif with the gganimate package.
The first step is to process a GTFS.zip file. The function gtfs2gps{gtfs2gps} interpolates the space-time position of each vehicle in each trip considering the network distance and average speed between stops. The output is a data.table where each row represents the timestamp of each vehicle at a given spatial resolution. In this example, we use a sample of the public transport network of Sao Paulo (Brazil) mapped every 30 meters.
# load libraries
library(gtfs2gps)
library(geobr)
library(ggplot2)
library(gganimate)
library(data.table)
library(ggthemes)
library(sf)
library(viridis)
library(sfheaders)
###### 1. process public transport data ------------------
# get local GTFS.zip file
gtfs_zip <- system.file("extdata/saopaulo.zip", package="gtfs2gps" )
# read gtfs
gtfs_dt <- gtfs2gps::read_gtfs(gtfs_zip)
# filter transport services on weekdays
gtfs_dt <- gtfs2gps::filter_week_days(gtfs_dt)
# get transport network as sf object
shapes_sf <- gtfs2gps::gtfs_shapes_as_sf(gtfs_dt)
# Convert GTFS data into a data.table with GPS-like records
gps_dt <- gtfs2gps::gtfs2gps(gtfs_dt, spatial_resolution = 30, progress = T, parallel = T)
head(gps_dt)
# subset time interval. Get only vehicle positions between 7am and 7:30am
gps_dt2 <- gps_dt[ between(departure_time, as.ITime("07:00:"), as.ITime("07:31"))]
# remove observations with 'too high speed' (in km/h)
gps_dt2 <- gps_dt2[speed < 20 , ]
# Convert "GPS" points into sf
gps_sf <- sfheaders::sf_point(gps_dt2, x = "shape_pt_lon" , y = "shape_pt_lat", keep = T)
sf::st_crs(gps_sf) <- 4326
head(gps_sf)
Once you have created the GPS-like dataset of the public transport network, creating an animation gif is simple as this:
###### 2. gif with blank background [fully reproducible] ------------------
# download Sao Paulo border for spatial context
sao_paulo <- geobr::read_municipality(code_muni = 3550308)
anim <- ggplot() +
geom_sf(data = sao_paulo, color='gray50', fill=NA) +
geom_sf(data = shapes_sf, color='gray90', size=0.01) +
geom_point(data = gps_dt2, aes(x = shape_pt_lon, y=shape_pt_lat, colour = trip_id), size=1.5, alpha = 0.5, show.legend = FALSE) +
scale_colour_viridis( discrete = T ) +
coord_sf(ylim = c(min(gps_dt2$shape_pt_lat), max(gps_dt2$shape_pt_lat)) ) +
theme_map() +
# gganimate specificatons
labs(title = 'Time: {frame_time}') +
transition_time(as.POSIXct(departure_time) + 10800) + # need to fix issue with time zone
shadow_wake(wake_length = 0.015, alpha = FALSE) +
ease_aes('linear')
# save gif
anim_save(animation = anim, "./tests_rafa/spo_gtfs_sample_20fps__GIST.gif", fps = 20)
In case you'd like to add more context info to the gif, you can use rasters as spatial tiles. I've added here some background info using a mapbox tile created by Kaue Braga for the Access to Opportunities Project. Thanks Kaue.
###### 3. gif tile background ------------------
# read map tile for spatial context
map_tiles <- readr::read_rds("L:/Proj_acess_oport/data/maptiles_crop/2019/mapbox/maptile_crop_mapbox_spo_2019.rds")
# convert gps data to UTM
gps_sf_utm <- sf::st_transform(gps_sf, 3857)
gps_dt_utm <- sfheaders::sf_to_df( gps_sf_utm, fill = TRUE )
anim2 <- ggplot() +
geom_raster(data = map_tiles, aes(x, y, fill = hex), alpha = 1) +
scale_fill_identity() +
geom_sf(data= st_transform(sao_paulo, 3857), color='gray50', fill=NA) +
geom_point(data = gps_dt_utm, aes(x = x, y=y, colour = trip_id), size=1.5, alpha = 0.6, show.legend = FALSE) +
scale_colour_viridis(discrete = T) +
coord_sf(ylim = c(min(gps_dt_utm$y), max(gps_dt_utm$y)) ) +
theme_map() +
# gganimate specificatons
labs(title = 'Time: {frame_time}') +
transition_time(as.POSIXct(departure_time) + 10800) + # need to fix issue with time zone
shadow_wake(wake_length = 0.015, alpha = FALSE) +
ease_aes('linear')
# save gif
anim_save(animation = anim2, "./tests_rafa/spo_gtfs_sample_20fps_tile__GIST.gif", fps = 20)
