Finding structural holes and create maps

Loading libraries

Load amenity/shop/transit facility type to all others betweenness edges information, bind df, summarise and store.

# Load the files of betweenness count
path = "D:/EconNet/Paris/Betweenness Count/individual/"
files = list.files(path=path, pattern=NULL, all.files=FALSE,
           full.names=FALSE)

# bind all the dataframes together
df = data.frame()
for (i in files) {
  df = rbind(df, read.csv(paste0(path,'/',i)))
}

# summarise by osmids
cnt_osm = df %>%
  group_by(osmid_from, osmid_to) %>%
  summarise(n = sum(n))
## `summarise()` has grouped output by 'osmid_from'. You can override using the
## `.groups` argument.
cnt_osm = cnt_osm[order(cnt_osm$n, decreasing = TRUE), ]

# write the total file to csv
write.csv(cnt_osm,"D:/EconNet/Paris/Betweenness Count/total_kernel_edges.csv")

head(cnt_osm)
## # A tibble: 6 x 3
## # Groups:   osmid_from [6]
##   osmid_from   osmid_to       n
##        <dbl>      <dbl>   <int>
## 1 5763285653  693550585 1700301
## 2 5763283950 5763283952 1696237
## 3 5763283952 5763285653 1696237
## 4  693550585   25554699 1685632
## 5   25554699   15644057 1685631
## 6   15644057 6369010761 1682661

Create one- and bi-directional keys, for bi-directional summarise both way betweenness edge count as one by the new key.

# Determine min and max osm per edge, for determining betweenness count from both directions
from = unlist(as.list(cnt_osm[,'osmid_from']))
to = unlist(as.list(cnt_osm[,'osmid_to']))
min = c()
max = c()
onedir = c()
bidir = c()
for (i in 1:nrow(cnt_osm)) {
  min = append(min, min(from[i], to[i]))
  max = append(max, max(from[i], to[i]))
  onedir = append(onedir, paste0(from[i],'-',to[i]))
  bidir = append(bidir, paste0(min[i],'-',max[i]))
}
cnt_osm$min_osm = min
cnt_osm$max_osm = max
cnt_osm$key_onedir = onedir
cnt_osm$key_bidir = bidir


# Sum by the key of min-max osmid-edges
cnt2 = cnt_osm %>%
  group_by(key_bidir) %>%
  summarise(n = sum(n))
cnt2 = cnt2[order(cnt2$n, decreasing = TRUE), ]


head(cnt2)
## # A tibble: 6 x 2
##   key_bidir                   n
##   <chr>                   <int>
## 1 5763283950-5763283952 1707257
## 2 5763283952-5763285653 1707257
## 3 693550585-5763285653  1700301
## 4 25554699-693550585    1685632
## 5 15644057-25554699     1685631
## 6 15644057-6369010761   1682661

Note that the echo = FALSE parameter was added to the code chunk to prevent printing of the R code that generated the plot.

Load the cleaned edges from before and create both one- and bidirectional keys from the edge direction, get distinct information.

# Read the cleaned edges
ce = read.csv("D:/EconNet/Paris/R_cleaned/edgesC_clean.csv")
ce = ce %>% select(-X)

# Add the counts to the edge-table
for (i in 1:nrow(ce)) {
  ce[i,'min_osm'] = min(ce[i,'u'], ce[i,'v'])
  ce[i,'max_osm'] = max(ce[i,'u'], ce[i,'v'])
  ce[i,'key_bidir'] = paste0(ce[i,'min_osm'],'-',ce[i,'max_osm'])
  ce[i,'key_onedir'] = paste0(ce[i,'u'],'-',ce[i,'v'])
}

shp_ce_sf = st_read("D:/EconNet/Paris/Cycling_City/BikeRoutes.shp")
## Reading layer `BikeRoutes' from data source 
##   `D:\EconNet\Paris\Cycling_City\BikeRoutes.shp' using driver `ESRI Shapefile'
## Simple feature collection with 48039 features and 23 fields
## Geometry type: LINESTRING
## Dimension:     XY
## Bounding box:  xmin: 2.225318 ymin: 48.81586 xmax: 2.469711 ymax: 48.90192
## Geodetic CRS:  WGS 84
# Add the counts to the edge-table
for (i in 1:nrow(shp_ce_sf)) {
  shp_ce_sf[i,'min_osm'] = min(shp_ce_sf$u[i], shp_ce_sf$v[i])
  shp_ce_sf[i,'max_osm'] = max(shp_ce_sf$u[i], shp_ce_sf$v[i])
  shp_ce_sf[i,'key_bidir'] = paste0(shp_ce_sf$min_osm[i],'-',
                                    shp_ce_sf$max_osm[i],'-')
  shp_ce_sf[i,'key_onedir'] = paste0(shp_ce_sf$u[i],'-',shp_ce_sf$v[i])
}

ce = ce %>% distinct(key_onedir, .keep_all = TRUE)
shp_ce_sf = shp_ce_sf %>% distinct(key_onedir, .keep_all = TRUE)
ce = merge(shp_ce_sf[c('key_onedir','geometry')], ce, by ='key_onedir')

ce1 = ce %>% distinct(key_bidir, .keep_all = TRUE)

head(ce)
## Simple feature collection with 6 features and 24 fields
## Geometry type: LINESTRING
## Dimension:     XY
## Bounding box:  xmin: 2.361818 ymin: 48.87632 xmax: 2.36405 ymax: 48.87679
## Geodetic CRS:  WGS 84
##                key_onedir           u           v      osmid
## 1  10000480189-5520660517 10000480189  5520660517 1096400097
## 2  10000480189-5900331100 10000480189  5900331100   24396069
## 3  10000480189-7660847829 10000480189  7660847829  112618158
## 4   10000480189-920531606 10000480189   920531606   78443912
## 5 10000522753-10000522756 10000522753 10000522756 1091942513
## 6 10000522756-10000522753 10000522756 10000522753 1091942513
##                    name    highway maxspeed oneway length lanes clsp     x1_m
## 1 Place Raoul Follereau pedestrian       11      0 38.586     1   10 453331.0
## 2 Place Raoul Follereau pedestrian       11      0 44.561     1   10 453331.0
## 3 Place Raoul Follereau pedestrian       30      0 29.370     1   30 453331.0
## 4          Passage L/10 pedestrian       30      0 30.933     1   30 453331.0
## 5                  <NA>    service       30      0 23.566     1   30 453206.7
## 6                  <NA>    service       30      0 23.566     1   30 453219.1
##      y1_m     x2_m    y2_m        hwys Fow Flane  Fsp Ftype perc_length
## 1 5413929 453295.9 5413944  pedestrian   0     0 0.00  0.10    42.44460
## 2 5413929 453370.1 5413908  pedestrian   0     0 0.00  0.10    49.01710
## 3 5413929 453345.8 5413955  pedestrian   0     0 0.04  0.10    33.59928
## 4 5413929 453315.5 5413903  pedestrian   0     0 0.04  0.10    35.38735
## 5 5413938 453219.1 5413918 residential   0     0 0.04  0.05    25.73407
## 6 5413918 453206.7 5413938 residential   0     0 0.04  0.05    25.73407
##       min_osm     max_osm               key_bidir
## 1  5520660517 10000480189  5520660517-10000480189
## 2  5900331100 10000480189  5900331100-10000480189
## 3  7660847829 10000480189  7660847829-10000480189
## 4   920531606 10000480189   920531606-10000480189
## 5 10000522753 10000522756 10000522753-10000522756
## 6 10000522753 10000522756 10000522753-10000522756
##                         geometry
## 1 LINESTRING (2.363034 48.876...
## 2 LINESTRING (2.36405 48.8763...
## 3 LINESTRING (2.363712 48.876...
## 4 LINESTRING (2.363306 48.876...
## 5 LINESTRING (2.361819 48.876...
## 6 LINESTRING (2.361819 48.876...

Merge on one- or bidirectional keys edge information, geometries and edge counts. Create an importance factor based on street penalties and form with edge betweenness the structural holes factor.

First the bi-directional which we see as main, then the one-directional, which we seem as minor

## Simple feature collection with 6 features and 19 fields
## Geometry type: LINESTRING
## Dimension:     XY
## Bounding box:  xmin: 2.34814 ymin: 48.85752 xmax: 2.351181 ymax: 48.85854
## Geodetic CRS:  WGS 84
##                u          v  min_osm.x  max_osm.x       n     osmid
## 21339 5763285653  693550585  693550585 5763285653 1700301 627023260
## 19168 5763283950 5763283952 5763283950 5763283952 1696237 815695339
## 19170 5763283952 5763285653 5763283952 5763285653 1696237 815695339
## 9445   693550585   25554699   25554699  693550585 1685632 545612376
## 2932    25554699   15644057   15644057   25554699 1685631 545612376
## 2934    15644057 6369010761   15644057 6369010761 1682661 860156627
##                   name        hwys Ftype clsp  Fsp oneway   Fow length
## 21339 Rue Saint-Martin residential  0.05   30 0.04      1 0.025  3.908
## 19168    Rue de Rivoli    cycleway  0.00   10 0.00      0 0.000 81.743
## 19170    Rue de Rivoli    cycleway  0.00   10 0.00      0 0.000 44.249
## 9445     Rue de Rivoli    cycleway  0.00   10 0.00      1 0.025 45.883
## 2932     Rue de Rivoli    cycleway  0.00   10 0.00      1 0.025 71.377
## 2934     Rue de Rivoli    cycleway  0.00   10 0.00      1 0.025  6.124
##       perc_length lanes Flane                       geometry importance_factor
## 21339    4.374224     1  0.00 LINESTRING (2.349631 48.858...           0.11930
## 19168   81.743000     1  0.00 LINESTRING (2.351181 48.857...           0.00000
## 19170   44.249000     1  0.00 LINESTRING (2.350171 48.857...           0.00000
## 9445    47.030075     1  0.00 LINESTRING (2.349653 48.858...           0.02500
## 2932    73.161425     1  0.00 LINESTRING (2.34909 48.8582...           0.02500
## 2934     6.590955     2  0.05 LINESTRING (2.348216 48.858...           0.07625
##       improve_importance_factor
## 21339                 202845.91
## 19168                      0.00
## 19170                      0.00
## 9445                   42140.80
## 2932                   42140.77
## 2934                  128302.90

Lets pre-process the cycleways for the graph and show a base map

## Warning: Transforming 'ext' to Web Mercator (EPSG: 3857), since 'ext' has a
## different CRS. The CRS of the returned basemap will be Web Mercator, which is
## the default CRS used by the supported tile services.
## Loading basemap 'streets' from map service 'osm'...
## # A tibble: 1 x 7
##   format width height colorspace matte filesize density
##   <chr>  <int>  <int> <chr>      <lgl>    <int> <chr>  
## 1 PNG      712    381 sRGB       FALSE   508473 72x72

## Focus on maps with bi-directional count, because cities often re-design the whole street, thus both ways of the cycle-network on that street. Results will be shown in graphs. First the edge betweenness count graph.

The maps are interactive. You can zoom and click on streets for more information. The files are saved to an HTML widget with the saveWidget function.

# Order to light on the focus area
st_bind = ce1_merge[order(ce1_merge$n),]

# Set a color palette. For the edge count this is a heatmap-like
pal <- colorNumeric(
  palette = c('#fff782','Orange','Red'),
  domain = st_bind$n)

# The bike routes are set as green
palb <- colorFactor(palette = "Green",
                    domain = NULL)

# Leaflets allow for smooth integration with sf.           
sf_lines = leaflet() %>%
  addTiles() %>%
  addPolylines(data = st_bind, color = ~pal(st_bind$n),
               # Log is used to balance the weight element
               weight = log(st_bind$n),
               # popups give information when a road is selected
               popup = paste0(st_bind$name,
                              "<br> <b> Type: </b>",st_bind$hwys,' (',st_bind$clsp,' km/h)',
                              "<br> <b> Times on route: </b>",st_bind$n,
                              "<br> <b> Road penalty: </b>",st_bind$importance_factor,
                              "<br> <b> Struc. hole score: </b>", round(st_bind$improve_importance_factor))) %>%
  addPolylines(data = st_bind_bike, color = 'green') %>%
  # Add legends to the graph
  addLegend(pal = pal, values = st_bind$n, position = "topright", title = 'Times on Route count') %>%
  addLegend(pal = palb, values = 'dedicated cycleway', position = "topright")

# Save the graph as an interactive leaflet widget
saveWidget(sf_lines, file = "D:/EconNet/Paris/Leaflet maps/Betweenness_edges.html")

sf_lines

Also for the second graph, which indicates the penalties resulting from its street design, and where the roads are situated in relation to the dedicated cycleways

# Order to light on the focus area
st_bind = st_bind[order(st_bind$importance_factor),]

# Do the same for the importance factor
pal <- colorNumeric(
  palette = c('#b0d6e3','#ff41c3'),
  domain = st_bind$importance_factor)

# Leaflets allow for smooth integration with sf.           
sf_lines = leaflet() %>%
  addTiles() %>%
  addPolylines(data = st_bind, 
               color = ~pal(st_bind$importance_factor),
               # *20 is used to balance the weight element
               weight = st_bind$importance_factor*20,
               # popups give information when a road is selected
               popup = paste0(st_bind$name,
                              "<br> <b> Type: </b>",st_bind$hwys,' (',st_bind$clsp,' km/h)',
                              "<br> <b> Times on route: </b>",st_bind$n,
                              "<br> <b> Road penalty: </b>",st_bind$importance_factor,
                              "<br> <b> Struc. hole score: </b>", round(st_bind$improve_importance_factor))) %>%
  addPolylines(data = st_bind_bike, color = 'green') %>%
  # Add legends to the graph
  addLegend(pal = pal, values = st_bind$importance_factor, position = "topright",title = 'Road penalty score') %>%
  addLegend(pal = palb, values = 'dedicated cycleway', position = "topright")

# Save the graph as an interactive leaflet widget
saveWidget(sf_lines, file = "D:/EconNet/Paris/Leaflet maps/Importance_factor.html")

sf_lines

The third graph brings the edge betweenness count and the importance factor togheter in a improve importance factor, which roads need to be imroved first as a result of both these metrics. These we consider the structural holes of the network.

# Order to light on the focus area
st_bind = st_bind[order(st_bind$improve_importance_factor),]

# As for calculating the betweenness edges or structural holes
pal <- colorNumeric(
  palette = "Reds",
  domain = st_bind$improve_importance_factor)

# Leaflets allow for smooth integration with sf
sf_lines = leaflet() %>%
  addTiles() %>%
  addPolylines(data = st_bind, 
               color = ~pal(st_bind$improve_importance_factor),
               # log is used to balance the weight element
               weight = log(st_bind$improve_importance_factor),
               # popups give information when a road is selected
               popup = paste0(st_bind$name,
               "<br> <b> Type: </b>",st_bind$hwys,' (',st_bind$clsp,' km/h)',
               "<br> <b> Times on route: </b>",st_bind$n,
               "<br> <b> Road penalty: </b>",st_bind$importance_factor,
               "<br> <b> Struc. hole score: </b>", round(st_bind$improve_importance_factor))) %>%
  addPolylines(data = st_bind_bike, color = 'green') %>%
  # Add legends to the graph
  addLegend(pal = pal, values = st_bind$improve_importance_factor, position = "topright",title = 'Structural hole score') %>%
  addLegend(pal = palb, values = 'dedicated cycleway', position = "topright")

# Save the graph as an interactive leaflet widget
saveWidget(sf_lines, file = "D:/EconNet/Paris/Leaflet maps/Structural_holes.html")

sf_lines