Skip to content

Instantly share code, notes, and snippets.

@aoles

aoles/coordinates_bayern-cycling-1000.Rda

Last active June 15, 2020 14:49
Show Gist options
  • Save aoles/ad2b1df588fea320a432e2ad52c5909e to your computer and use it in GitHub Desktop.
Save aoles/ad2b1df588fea320a432e2ad52c5909e to your computer and use it in GitHub Desktop.
Download ZIP
Raw
cyclingWeightings.Rmd
---
title: "Cycling weightings compared"
output: html_document
---
```{r setup, include=FALSE}
knitr::opts_chunk$set(echo = TRUE)
```
## Preprocessing of data
```{r}
library(openrouteservice)
library(leaflet)
load("recommendended-bike-1000.Rda")
n = 1000
lapply(res, function(routes) {
sapply(routes, function(x) x$features[[1]]$properties$summary$distance/1000)
}) %>% data.frame -> distance
lapply(res, function(routes) {
sapply(routes, function(x) x$features[[1]]$properties$summary$duration)
}) %>% data.frame -> duration
lapply(res, function(routes) {
sapply(routes, function(x) x$features[[1]]$properties$ascent)
}) %>% data.frame -> ascent
lapply(res, function(routes) {
lapply(routes, function(x) {
x = x$features[[1]]$properties$extras$suitability$summary
values = sapply(x, `[[`, "value")
amount = sapply(x, `[[`, "amount")
v = rep(0, 8)
v[values-2] <- amount
v
})
}) %>% lapply(function(routes) {
sapply(routes, function(v) {
sum(0:7 * v/100)
})
}) %>% data.frame -> suitability
```
## Route characteristics
```{r}
range(distance$fastest)
```
```{r}
hist(distance$fastest, xlim = c(0, 500), breaks = 50)
```
## Suitability
```{r}
compare <- function(what, ref, f) {
sapply(lapply(what[setdiff(names(what), ref)], f, what[ref]), sum)
}
```
Number of routes that have a better suitability index compared to GraphHopper's `fastest`...
```{r}
compare(suitability, "fastest", `>`)
```
.. and a worse suitability.
```{r}
compare(suitability, "fastest", `<`)
```
Even though current `recommended` is a better choice than `fastest`, all the other alternatives seem to perform better with the proposed `new` being the best.
But how do actually compare the other weightings to `new`?
```{r}
compare(suitability, "new", `<`)
compare(suitability, "new", `>`)
```
Routes which become worse with the new weighting.
```{r}
(worse_suitability <- which(suitability$new < suitability$fastest | suitability$new < suitability$recommended))
```
```{r, out.width = '100%', out.height = '780px', echo = FALSE}
leaflet() %>%
addTiles(group = "OpenStreetMap") %>%
addTiles("https://dev.{s}.tile.openstreetmap.fr/cyclosm/{z}/{x}/{y}.png", group ="CyclOSM" ) %>%
addTiles("https://{s}.tile.thunderforest.com/cycle/{z}/{x}/{y}.png?apikey=13efc496ac0b486ea05691c820824f5f", group ="OpenCycleMap") %>%
addGeoJSON(res$fastest[worse_suitability], fill=FALSE, color = c("#000"), group = "current fastest") %>%
addGeoJSON(res$recommended[worse_suitability], fill=FALSE, color = c("#f00"), group = "current recommended") %>%
addGeoJSON(res$new[worse_suitability], fill=FALSE, color = c("#0f0"), group = "new") %>%
addLayersControl(
baseGroups = c("OpenStreetMap", "CyclOSM", "OpenCycleMap"),
overlayGroups = c("current fastest", "current recommended", "new"),
position = "bottomleft"
) %>%
fitBBox(c(8.974646, 47.265430, 13.849470, 50.567790))
```
## Ascent
Total ascent values generated across routes.
```{r}
sapply(ascent, sum)
```
`new` is second best even though it produces the longest routes.
```{r}
sapply(distance, sum)
```
Number of routes that have a lower ascent (better) compared to GraphHopper's `fastest` ...
```{r}
compare(ascent, "fastest", `<`)
```
... and higher ascent (worse).
```{r}
compare(ascent, "fastest", `>`)
```
Again, the `new` weighting seems to be the best alternative to `fastest`.
## Summary
Percentage of `recommended` routes which are worse than `fastest` in terms of both suitability and ascent.
```{r}
sum(ascent$recommended > ascent$fastest & suitability$recommended < suitability$fastest) / n * 100
```
Percentage of `new` routes which are worse than `fastest` in terms of both suitability and ascent.
```{r}
sum(ascent$new > ascent$fastest & suitability$new < suitability$fastest) / n * 100
```
Percentage of `new` routes which are worse than `recommended` in terms of both suitability and ascent.
```{r}
sum(ascent$new > ascent$recommended & suitability$new < suitability$recommended) / n * 100
```
Raw
generate_routes.R
library("openrouteservice")
set.seed(0L)
options(openrouteservice.url = "http://localhost:8082/ors")
profile = "cycling-regular"
n = 1000
rep = 1
### generate n random coordinates and convert them to start/endpoints
coordinates_file = "coordinates_bayern-cycling-1000.Rda"
if (file.exists(coordinates_file)) {
load(coordinates_file)
} else {
library("sf")
library("geojson")
x = readLines("~/Data/polygons/bayern.geojson")
bbox = geo_bbox(x)
poly = st_read(paste(x, collapse = ""))
pts = vector(mode = "list", n)
i = 0;
while (i < n) {
coords = c(runif(1, bbox[1], bbox[3]), runif(1, bbox[2], bbox[4]))
if (st_within(st_point(coords), poly, sparse=FALSE))
if ( !inherits(try(ors_directions(list(coords, coords), profile=profile), silent=TRUE), "try-error") )
pts[[(i = i+1)]] <- coords
}
# duration matrix with diagonal excluded
res <- ors_matrix(pts, profile = profile)
durations <- res$durations
durations[seq(1, by = n+1, length.out = n)] <- NA
max <- max(res$durations)
d <- vector(mode = "numeric", n)
coordinates <- vector(mode = "list", n)
for (i in 1:n) {
r <- runif(1, max = max)
w <- which.min(abs(durations - r))
d[i] <- durations[w]
index <- arrayInd(w, dim(durations))
a = index[1,1]
b = index[1,2]
durations[a, b] <- NA
durations[b, a] <- NA
coordinates[[i]] <- list(pts[[a]], pts[[b]])
}
coordinates <- coordinates[order(d)]
save(coordinates, file = coordinates_file, compress = "xz")
}
profile_args = list(profile = profile, format = 'geojson', instructions = FALSE, elevation = TRUE, extra_info = "suitability")
profile_args = list(
fastest = c(profile_args, preference = "fastest"),
recommended = c(profile_args, preference = "recommended"),
reduced = c(profile_args, preference = "recommendedr"),
scaled = c(profile_args, preference = "recommendeds"),
new = c(profile_args, preference = "recommendednew")
)
res = lapply(names(profile_args), function(name) {
lapply(1:n, function(id) {
cat(" \r", name, id)
label = sprintf("[BENCHMARK] %s; %d", name, id)
cl = as.call(c(ors_directions, coordinates[id], profile_args[[name]], id = label))
res = try(eval(cl), silent = TRUE)
if ( inherits(res, "try-error") ) {
attr(res, "condition")$message
} else {
cl$id = NULL
query_times = replicate(rep-1, attr(eval(cl), "query_time"))
attr(res, "query_time") = c(attr(res, "query_time"), query_times)
res
}
})
})
names(res) <- names(profile_args)
save(res, file = "recommendended-bike-1000.Rda", compress = "xz")
Raw
recommendended-bike-1000.Rda
This file has been truncated, but you can view the full file.
View raw
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment