Cycling uptake functions about cyclingpotential-hack

[Robinlovelace]

Output of reproducible code above:

# from code/reproducible-example.R

# remotes::install_github("itsleeds/pct")

library(dplyr)
#> 
#> Attaching package: 'dplyr'
#> The following objects are masked from 'package:stats':
#> 
#>     filter, lag
#> The following objects are masked from 'package:base':
#> 
#>     intersect, setdiff, setequal, union
library(pct)
library(stplanr)
library(sf)
#> Linking to GEOS 3.8.0, GDAL 3.0.4, PROJ 7.0.0

# Inter-district travel ---------------------------------------------------


u3 = "https://github.com/U-Shift/cyclingpotential-hack/releases/download/1.0/routes_fast.geojson"
route_segments_fast = sf::read_sf(u3)
routes_fast = route_segments_fast %>%
  group_by(DICOFREor11, DICOFREde11) %>%
  summarise(
    Origem = first(DICOFREor11),
    Destino = first(DICOFREde11),
    Bike = mean(Bike),
    All = mean(Total),
    Length_fast_m = sum(distances),
    Hilliness_average = mean(gradient_segment),
    Hilliness_90th_percentile = quantile(gradient_segment, probs = 0.9)
  )
#> `summarise()` regrouping output by 'DICOFREor11' (override with `.groups` argument)

unique(sf::st_geometry_type(routes_fast))
#> [1] MULTILINESTRING
#> 18 Levels: GEOMETRY POINT LINESTRING POLYGON MULTIPOINT ... TRIANGLE
nrow(routes_fast)
#> [1] 332
routes_fast$pcycle_current = routes_fast$Bike / routes_fast$All
plot(routes_fast["pcycle_current"])

m_pct = pct::model_pcycle_pct_2020(
  pcycle = routes_fast$pcycle_current,
  distance = routes_fast$Length_fast_m,
  # gradient = routes_fast$Hilliness_average,
  gradient = routes_fast$Hilliness_average,
  weights = routes_fast$All
)
m_pct
#> 
#> Call:  stats::glm(formula = pcycle ~ distance + sqrt(distance) + I(distance^2) + 
#>     gradient + distance * gradient + sqrt(distance) * gradient, 
#>     family = "quasibinomial", weights = weights)
#> 
#> Coefficients:
#>             (Intercept)                 distance           sqrt(distance)  
#>               5.881e+00                1.628e-03               -2.774e-01  
#>           I(distance^2)                 gradient        distance:gradient  
#>              -1.086e-08               -1.830e+02               -1.996e-02  
#> sqrt(distance):gradient  
#>               4.096e+00  
#> 
#> Degrees of Freedom: 331 Total (i.e. Null);  325 Residual
#> Null Deviance:       1394 
#> Residual Deviance: 1341  AIC: NA

pcycle_pct_govtarget = pct::uptake_pct_govtarget_2020(
  distance = routes_fast$Length_fast_m,
  gradient = routes_fast$Hilliness_average
)
#> Distance assumed in m, switching to km

pcycle_pct_godutch = pct::uptake_pct_godutch_2020(
  distance = routes_fast$Length_fast_m,
  gradient = routes_fast$Hilliness_average
)
#> Distance assumed in m, switching to km

plot(
  routes_fast$Length_fast_m,
  routes_fast$pcycle_current,
  cex = routes_fast$All / mean(routes_fast$All),
  ylim = c(0, 0.5)
  )
points(routes_fast$Length_fast_m, m_pct$fitted.values, col = "red")
points(routes_fast$Length_fast_m, pcycle_pct_godutch, col = "green")
points(routes_fast$Length_fast_m, pcycle_pct_govtarget, col = "grey")

routes_fast$slc_godutch = routes_fast$All * pcycle_pct_godutch
length(unique(routes_fast$geometry))
#> [1] 332

rnet_fast = overline(sf::st_cast(routes_fast, "LINESTRING"), attrib = "slc_godutch")
#> Warning in st_cast.sf(routes_fast, "LINESTRING"): repeating attributes for all
#> sub-geometries for which they may not be constant
rnet_fast$slc_godutch = round(rnet_fast$slc_godutch)
summary(rnet_fast$slc_godutch)
#>    Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
#>     5.0   139.5   287.0   398.5   527.0  3786.0
rnet_99th_percentile = quantile(rnet_fast$slc_godutch, probs = 0.99)
rnet_fast$slc_godutch[rnet_fast$slc_godutch > rnet_99th_percentile] = rnet_99th_percentile
mapview::mapview(rnet_fast, alpha = 0.5, lwd = rnet_fast$slc_godutch / 100)
#> Warning in if ("gl" %in% names(list(...)) & isTRUE(list(...)$gl) &
#> inherits(sf::st_geometry(x), : the condition has length > 1 and only the first
#> element will be used

^{Created on 2020-09-02 by the reprex package (v0.3.0)}

[Robinlovelace]

Here is a reproducible script to show how the PCT uptake function works.

Can we generate different uptake scenarios for Lisbon?

Or other places?

u3 = "https://github.com/U-Shift/cyclingpotential-hack/releases/download/1.0/routes_fast.geojson"
route_segments_fast = sf::read_sf(u3)
routes_fast = route_segments_fast %>%
  group_by(DICOFREor11, DICOFREde11) %>%
  summarise(
    Origem = first(DICOFREor11),
    Destino = first(DICOFREde11),
    Bike = mean(Bike),
    All = mean(Total),
    Length_fast_m = sum(distances),
    Hilliness_average = mean(gradient_segment),
    Hilliness_90th_percentile = quantile(gradient_segment, probs = 0.9)
  )

unique(sf::st_geometry_type(routes_fast))
nrow(routes_fast)
routes_fast$pcycle_current = routes_fast$Bike / routes_fast$All
plot(routes_fast["pcycle_current"])

m_pct = pct::model_pcycle_pct_2020(
  pcycle = routes_fast$pcycle_current,
  distance = routes_fast$Length_fast_m,
  # gradient = routes_fast$Hilliness_average,
  gradient = routes_fast$Hilliness_average,
  weights = routes_fast$All
)
m_pct

pcycle_pct_govtarget = pct::uptake_pct_govtarget_2020(
  distance = routes_fast$Length_fast_m,
  gradient = routes_fast$Hilliness_average
)

pcycle_pct_godutch = pct::uptake_pct_godutch_2020(
  distance = routes_fast$Length_fast_m,
  gradient = routes_fast$Hilliness_average
)

plot(
  routes_fast$Length_fast_m,
  routes_fast$pcycle_current,
  cex = routes_fast$All / mean(routes_fast$All),
  ylim = c(0, 0.5)
  )
points(routes_fast$Length_fast_m, m_pct$fitted.values, col = "red")
points(routes_fast$Length_fast_m, pcycle_pct_godutch, col = "green")
points(routes_fast$Length_fast_m, pcycle_pct_govtarget, col = "grey")

routes_fast$slc_godutch = routes_fast$All * pcycle_pct_godutch
length(unique(routes_fast$geometry))

rnet_fast = overline(sf::st_cast(routes_fast, "LINESTRING"), attrib = "slc_godutch")
rnet_fast$slc_godutch = round(rnet_fast$slc_godutch)
summary(rnet_fast$slc_godutch)
rnet_99th_percentile = quantile(rnet_fast$slc_godutch, probs = 0.99)
rnet_fast$slc_godutch[rnet_fast$slc_godutch > rnet_99th_percentile] = rnet_99th_percentile
mapview::mapview(rnet_fast, alpha = 0.5, lwd = rnet_fast$slc_godutch / 100)