-
Notifications
You must be signed in to change notification settings - Fork 1
/
3.1_analysis-travel_time_r5.R
326 lines (258 loc) · 13.7 KB
/
3.1_analysis-travel_time_r5.R
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
###################################################################################################
### The purpose of this script is to calculate a travel time matrix for each different mode ###
### combination. r5r is used for the calculations ###
###################################################################################################
# increase the memory available to Java. Needs to be done at the beginning of the script
options(java.parameters = "-Xmx20G") # 3 gegabytes
library(r5r)
library(sf)
library(tidyverse)
# ------------------------------------- Define Variables ------------------------------------- #
# city
#city <- "San Francisco"
#city <- "Mexico City"
city <- "Minneapolis"
#city <- "Cairo"
# --- define paths
# path to folder with gtfs, road network and elevation data
data_path <- paste0("../data_raw/", city, "/GTFS/")
# path where pbf files are stored
pbf_path <- paste0("../data_raw/", city, "/PBFs/")
# are we running the analysis for congested (yes) or freeflow (no) speeds
#congested <- "yes"
congested <- "no"
# ------------------------------------- Load in the data ------------------------------------- #
# city geometry
city_geom <- st_read(paste0("../data/", city, "/level_i_ii/zones_w_micromobility_providers.geojson"))
# prep data for inputting to r5r
# Function takes a base layer, gets the geometry centroid, and renames the id column that we pass to it into a standard name
prep_base_layer = function(layer, id_col){
id_col = sym(id_col)
layer = layer %>% select(!! id_col) %>%
st_centroid() %>% # we route from zone centroids
rename(id = !! id_col)
}
# apply the function
city_geom_r5 <- prep_base_layer(layer = city_geom, id_col = "cell_id")
# ------------------------------------- PBF data
# We have different PBFs for each city (freeflow / congested). We need to make sure we using the correct one.
# PBFs for any city are stored in <City>/PBFs.
# We copy the correct one into the <City>/GTFS directory, and delete PBFs / .mapdb / .dat files from previous runs
pbf_files <- tribble(
~city, ~freeflow_pbf_file, ~congested_pbf_file,
# San Francisco
"San Francisco", "SF_original3.osm.pbf" , "SF_real_speed3.osm.pbf", # docked_Bay.Wheels
#"San Francisco", "dockless_Spin.San.Francisco", TRUE,
# Minneapolis
"Minneapolis", "MN_original.osm.pbf", NA, # docked_Lyft.Scooters.Minneapolis
# "Minneapolis", "dockless_Spin.Minneapolis", TRUE,
# Mexico City
"Mexico City", "planet_mx.osm.pbf", "mx_real1.osm.pbf", # docked_ECOBICI
# Cairo
"Cairo", "cairo.osm.pbf", "cairo_real1.osm.pbf" # docked_ECOBICI
)
# --- identify pbf file names from table
freeflow_pbf <- pbf_files$freeflow_pbf_file[pbf_files$city == city]
congested_pbf <- pbf_files$congested_pbf_file[pbf_files$city == city]
# --- delete files in GTFS directory
# remove existing pbf file
# system(paste0("rm ", data_path, "*.pbf"))
# gsub to replace spaces with '\' https://stackoverflow.com/questions/38909016/r-how-to-replace-space-in-string-with-a-single-backslash-and-space
system(paste0("rm ", gsub(" ", "\\\ ", data_path, fixed = TRUE), "*.pbf"))
# remove other files created by r5 during graph building
#system(paste0("rm ", data_path, "*.p"))
system(paste0("rm ", gsub(" ", "\\\ ", data_path, fixed = TRUE), "*.p"))
system(paste0("rm ", gsub(" ", "\\\ ", data_path, fixed = TRUE), "*.dat"))
system(paste0("rm ", gsub(" ", "\\\ ", data_path, fixed = TRUE), "*.map"))
system(paste0("rm ", gsub(" ", "\\\ ", data_path, fixed = TRUE), "*.mapdb"))
# --- copy pbf file we want into GTFS directory (directory where graph is built)
# check which file to move (depends on whether we want congestion or not)
if(congested == "yes"){
# move the augmented pbf to the gtfs directory
system(paste0("cp ", gsub(" ", "\\\ ", pbf_path, fixed = TRUE), congested_pbf, " ",
gsub(" ", "\\\ ", data_path, fixed = TRUE)))
} else {
# mve the unedited pbf into the pbf directory
system(paste0("cp ", gsub(" ", "\\\ ", pbf_path, fixed = TRUE), freeflow_pbf, " ",
gsub(" ", "\\\ ", data_path, fixed = TRUE)))
}
# ------------------------------------- Setup r5 ------------------------------------- #
# # need to edit pbf to avoid a turn-restrictions related bug in r5
# # https://github.com/ipeaGIT/r5r/issues/141#issuecomment-767648579
# r5r_core <- setup_r5(data_path = data_path,
# verbose = TRUE,
# use_elevation = TRUE, # we need to download a tiff file and then turn this to TRUE
# overwrite = TRUE) # turn to true once we have elevation
# ------------------------------------- Define combinations ------------------------------------- #
# --- tibble with all mode combinations
combinations <- tribble(
~combination, ~mode, ~egress_mode, ~max_walk,
# public transport
1, c("WALK", "TRANSIT"), "WALK", 500,
# public transport + micromobility as a first mile option only
2.1, c("WALK", "TRANSIT", "BICYCLE"), "WALK", 500,
# public transport + micromobility as a last mile option only
2.2, c("WALK", "TRANSIT"), "BICYCLE", 500,
# public transport + micromobility as a first and last mile option only
2.3, c("WALK", "TRANSIT", "BICYCLE"), "BICYCLE", 500,
# cycling / micromobility only
3, "BICYCLE", "WALK", 500,
# car
4, c("WALK", "CAR"), "WALK", 500,
# walk only
5, "WALK", "WALK", 5000
)
# ------------------------------------- Define Routing Parameters ------------------------------------- #
max_walk_dist <- 1000 # meters
max_trip_duration <- 75 # minutes
# This needs to be within the start_date and end_date of the calendar.txt (0_edit_gtfs_calendar.R can help with this)
departure_datetime <- as.POSIXct("24-01-2022 07:30:00",
format = "%d-%m-%Y %H:%M:%S")
time_window <- 60 # in minutes (adding this value to departure datetime, gives you the end of the time window)
percentiles <- 75
# ------------------------------------- Define Routing Parameters ------------------------------------- #
# --- calculate a travel time matrix
# ttm_transit <- travel_time_matrix(r5r_core = r5r_core,
# origins = city_geom_r5 ,
# destinations = city_geom_r5,
# time_window = time_window,
# percentiles = percentiles,
# mode = combinations$mode[1][[1]],
# mode_egress = combinations$egress_mode[1][[1]],
# departure_datetime = departure_datetime,
# max_walk_dist = max_walk_dist,
# max_trip_duration = max_trip_duration,
# max_lts = 2,
# # to suppress r5 output. Change to true if debugging
# verbose = FALSE,
# # slow down function by ~20%
# progress = TRUE)
# FUNCTION: Calculate travel time for all combinations and return the result as 1 dataframe
tt_matrix = function(data_path,
combinations,
zone_layer,
elevation,
# bike speed for normal and electric vehicles
bike_speed_e,
bike_speed_c,
max_lts){
# stop any running r5 instances
# r5r::stop_r5()
# setup r5
print("Setting up r5...")
r5r_core <- setup_r5(data_path = data_path,
verbose = TRUE,
use_elevation = elevation, # we need to download a tiff file and then turn this to TRUE
overwrite = TRUE) # turn to true once we have elevation
# empty list to store results for each combination
results <- vector(mode = "list", length = nrow(combinations))
print("Calculating travel times...")
# calculate travel time matrix for each combination
for(i in 1:nrow(combinations)){
#status updates
print(paste0("CALCULATING TRAVEL TIME FOR combination: ", combinations$combination[i], " ....."))
# calculate a travel time matrix
ttm <- r5r::travel_time_matrix(r5r_core = r5r_core,
origins = zone_layer ,
destinations = zone_layer,
time_window = time_window,
percentiles = percentiles,
mode = combinations$mode[i][[1]],
mode_egress = combinations$egress_mode[i][[1]],
departure_datetime = departure_datetime,
max_walk_dist = combinations$max_walk[i][[1]], #max_walk_dist,
max_trip_duration = max_trip_duration,
# if elevation is TRUE, we are looking at classic bikes, otherwise electric
bike_speed = ifelse(elevation == TRUE, bike_speed_c, bike_speed_e),
max_lts = max_lts,
# number of threads
#n_threads = 2,
# to suppress r5 output. Change to true if debugging
verbose = FALSE,
# slow down function by ~20%
progress = TRUE)
# add column with combination name / number
# if elevation = TRUE, we are looking at classic bikes
if (elevation == TRUE){
ttm$combination <- stringr::str_c(combinations$combination[i], "_classic")} else{
# if elevation = FALSE, we are looking at electric bikes
ttm$combination <- stringr::str_c(combinations$combination[i], "_electric")
}
# add ttm to results list
results[[i]] <- ttm
#status updates
print(paste0("COMPLETED combination: ", combinations$combination[i], " ....."))
}
# stop r5
r5r::stop_r5()
# combine list into 1 dataframe
results <- bind_rows(results)
# pivot wider to get 1 row per OD pair
results <- results %>% pivot_wider(names_from = combination,
values_from = travel_time)
return(results)
}
# apply the function
# classic
tt_matrix_city_c <- tt_matrix(zone_layer = city_geom_r5,
combinations = combinations,
data_path = data_path,
elevation = TRUE,
# bike speed for normal and electric vehicles
bike_speed_e = 22.5,
bike_speed_c = 16.6,
max_lts = 2)
# temporary as I don't have the time to run the below function now:
#write_csv(tt_matrix_city_c, paste0("../data/", city, "/level_i_ii/freeflow_classic_temp.csv"))
#tt_matrix_city_c <- read_csv(paste0("../data/", city, "/level_i_ii/freeflow_classic_temp.csv"))
#write_csv(tt_matrix_city_c, paste0("../data/", city, "/level_i_ii/congested_classic_temp.csv"))
#tt_matrix_city_c <- read_csv(paste0("../data/", city, "/level_i_ii/congested_classic_temp.csv"))
# temporary - END
# electric
tt_matrix_city_e <- tt_matrix(zone_layer = city_geom_r5,
combinations = combinations,
data_path = data_path,
elevation = FALSE,
# bike speed for normal and electric vehicles
bike_speed_e = 22.5,
bike_speed_c = 16.6,
max_lts = 2)
tt_matrix_city <- tt_matrix_city_c %>%
full_join(tt_matrix_city_e, by = c("fromId", "toId"))
# save
#write_csv(tt_matrix_city, paste0("../data/", city, "/level_i_ii/travel_time_matrix_freeflow.csv"))
#write_csv(tt_matrix_city, paste0("../data/", city, "/level_i_ii/travel_time_matrix.csv"))
# --- test free flow vs congested speeds
# travel time matrix with freeflow car speeds
tt_matrix_city_free <- read_csv(paste0("../data/", city, "/level_i_ii/travel_time_matrix_freeflow.csv"))
#
# travel time matrix using edited pbf (more realistic travel speeds)
tt_matrix_city <- read_csv(paste0("../data/", city, "/level_i_ii/travel_time_matrix.csv"))
#
tt_matrix_city <- tt_matrix_city %>%
mutate(across(fromId:toId, as.numeric))
#
# congested joined with free
tt_matrix <- tt_matrix_city %>% left_join(tt_matrix_city_free, by = c("fromId", "toId"))
# #
# percentage change in speeds
tt_matrix <- tt_matrix %>% mutate(perc_change_1 = ((`1_classic.x` - `1_classic.y`) / `1_classic.y`) * 100,
perc_change_2.1 = ((`2.1_classic.x` - `2.1_classic.y`) / `2.1_classic.y`) * 100,
perc_change_2.2 = ((`2.2_classic.x` - `2.2_classic.y`) / `2.2_classic.y`) * 100,
perc_change_2.3 = ((`2.3_classic.x` - `2.3_classic.y`) / `2.3_classic.y`) * 100,
perc_change_3 = ((`3_classic.x` - `3_classic.y`) / `3_classic.y`) * 100,
perc_change_4 = ((`4_classic.x` - `4_classic.y`) / `4_classic.y`) * 100,
perc_change_5 = ((`5_classic.x` - `5_classic.y`) / `5_classic.y`) * 100
)
# #
# plot
tt_matrix %>% select(fromId, toId, contains("perc_change")) %>%
pivot_longer(cols = perc_change_1:perc_change_4) %>%
ggplot(aes(x=value, fill=name)) +
geom_histogram(alpha=0.7, position = 'identity', binwidth = 10) +
facet_wrap(~ name) +
#scale_fill_manual(values=c("#69b3a2", "#404080")) +
theme_minimal() +
labs(fill="") +
# axis limits
scale_x_continuous(limits = c(NA, 200))