The combination of spatial open data, standardised spatial web services, and recent innovations in R’s spatial packages ecosystem make it an exiting domain to play with. In this walkthrough I will demonstrate:

querying OpenStreetMap for geographic features;

processing them with the R sf and dplyr packages;

and packages; interactively displaying them on historical maps using Leaflet and WMS-webservices.

The geographic features we are interested in are the 19th century “gangen” (Dutch for “passages”) in Leuven, Belgium. Around 1880 Leuven counted about 140 of these gangen, housing 13% of the population of Leuven at that time. They were called passages, because a small passage led you from the public street to terrains behind the ‘regular’ houses. Those terrains were subdivided into small housing units rented by the proletariat and laboring classes in Leuven.

Lets see what spatial open data and webservices can learn us about these gangen today.

Fetch spatial objects from OpenStreetMap Load the required libraries: library(osmdata) # fetching OpenStreetmap-data library(leaflet) # display interactive maps library(stringr) # optional, tidyverse string-operations library(dplyr) # data-manipulation library(sf) # spatial features library library(forcats) # optional, tidyverse factor-manipulation The gangen were generally called ‘-gang’, e.g. Zandgang, Mathildegang, Blauwputgang, etc., and walking through Leuven you can still find streetnames ending in -gang, denoting the location of historical gangen. To get the spatial location and layout of every street in Leuven ending in -gang, we use OpenStreetMap (OSM). osmdata is an R package for accessing the data underlying (OSM). It allows us to download all spatial features in a bounding box around Leuven from OSM, and filter on type ‘highway’ to get only the streets: q.leuven <- opq(bbox = 'Leuven, Belgium') %>% # set bounding box via name add_osm_feature(key = 'highway') %>% # get streets add_osm_feature(key = 'name') %>% # include the name osmdata_sf() # return as a sf-object q.leuven ## Object of class 'osmdata' with: ## $bbox : 50.8242096,4.640295,50.9440707,4.7705305 ## $overpass_call : The call submitted to the overpass API ## $meta : metadata including timestamp and version numbers ## $osm_points : 'sf' Simple Features Collection with 27656 points ## $osm_lines : 'sf' Simple Features Collection with 5725 linestrings ## $osm_polygons : 'sf' Simple Features Collection with 109 polygons ## $osm_multilines : NULL ## $osm_multipolygons : 'sf' Simple Features Collection with 4 multipolygons The returned osmdata-object contains about 35.000 spatial features of various types (points, lines, polygons). We are presently only interested in the spatial line representations of the gangen. So we select the osm_lines part of the osmdata-object with the base R $ -operator, select() from all the variables returned by OSM only the name-variable, and filter() down the list of spatial lines to only those where the name contains “gang” somewhere. gangen = q.leuven$osm_lines %>% # use only the lines in the osmdata-object select(name) %>% # keep only the name-variabele # keep only the lines with 'gang' in the name: filter(str_detect(name, 'gang') == TRUE) # fix plot issue: https://stackoverflow.com/a/50846998/125085 names(gangen$geometry) <- NULL gangen <- gangen %>% # Drop unused levels in name-variable w/t forcats::fct_drop(): mutate(name = fct_drop(name)) %>% arrange(name) # sort alphabetically on the name gangen ## Simple feature collection with 58 features and 1 field ## geometry type: LINESTRING ## dimension: XY ## bbox: xmin: 4.686551 ymin: 50.86957 xmax: 4.71337 ymax: 50.88845 ## epsg (SRID): 4326 ## proj4string: +proj=longlat +datum=WGS84 +no_defs ## First 10 features: ## name geometry ## 1 Beursgang LINESTRING (4.692325 50.880... ## 2 Beursgang LINESTRING (4.692272 50.880... ## 3 Beursgang LINESTRING (4.69243 50.8803... ## 4 Blauwputgang LINESTRING (4.711557 50.881... ## 5 Blauwputgang LINESTRING (4.711378 50.882... ## 6 Blauwputgang LINESTRING (4.711617 50.881... ## 7 Bloemenberggang LINESTRING (4.688926 50.882... ## 8 Boghe of Brandgang LINESTRING (4.696834 50.876... ## 9 Busleidengang LINESTRING (4.699747 50.881... ## 10 Busleidengang LINESTRING (4.699971 50.881... We now have a sf -collection with streets having ‘gang’ in their name. The objectalso doubles as a familiar R-dataframe with two columns/variables: the OSM name of the spatial line and geometry , which holds the spatial definition of the line. Using the leaflet-package, which is a wrapper for the JS library in R, we can display this sf -object on an interactive map: leaflet(gangen) %>% addTiles() %>% addPolylines( # use linestring name to show as a label on mouse-over: label = ~name, # optional: highlight linestring in red on mouse-over: highlight = highlightOptions(color = "red")) The remaining 26 gangen are all in their expected place. There is however an issue: when a gang is not a simple street but a complex of streets, it is represented by multiple spatial strings. That is why there are e.g. three records for Beursgang in the sf-object, and you can also see this clearly if you mouse-over the Matildegang, the gang the most south-east in Leuven (different lines should be highlighted in red). This is to be expected. The simple features specification implemented in the sf -package and other spatial libaries is build on basic spatial ‘buildingblocks’, such as points (POINT), lines (LINESTRING), and polygons (POLYGON). And osmdata returns you spatial objects in proper buildingblock-form. So if you want to use multiple simple lines to spatially represent a complex of streets, you need to combine lines into a composit spatial set of lines (MULTILINESTRING). Luckly, the sf -packages nicely integrates into the tidyverse-ecosystem, which allows us to use common dplyr -verbs to do specific spatial operations. We group the records/lines in the sf-objects on the name of the street with group_by() and aggregate the individual records/lines with the same name using summarise() . This gives us a nice sf-dataframe with the expected 26 records of the remaining gangen in Leuven. gangen.merged <- gangen %>% group_by(name) %>% summarise() gangen.merged$name ## [1] Beursgang Blauwputgang Bloemenberggang ## [4] Boghe of Brandgang Busleidengang Cansgang ## [7] Daneelsgang Elzasgang Hallengang ## [10] Hanengang Horengang Korbeeklogang ## [13] Mathildegang Oppendorpgang Paardengang ## [16] Peterseliegang Pioengang Puttegang ## [19] Rapengang Straatjesgang Theresianengang ## [22] Valkerijgang Vestinggang Wolvengang ## [25] Wolvenpoortgang Zandgang Zongang ## 27 Levels: Beursgang Blauwputgang Bloemenberggang ... Zongang Displaying this aggregated sf-object on a map, shows on mouse-over that the lines are properly joined into the expected gangen-complexes. leaflet(gangen.merged) %>% addTiles() %>% addPolylines( label = ~name, highlight = highlightOptions(color = "red"))

Historical maps using WMS-tileservers Next to contemporary map-layers, it would be interesting to display the (remaining) gangen on historical maps. Conveniently, the regional Agentschap Informatie Vlaanderen makes tile-sets of historical maps available through publicly available WMS-webservices. For a 19th century phenomenon as the gangen, we will use the historical Ferraris map, and the Vandermaelen map. The Ferraris map was created between 1770 and 1778, and is the first systematic, large scale mapping in present-day Belgium and Western Europe. The map reflects the detailed state of the Southern Netherlands just before the start of the Industrial Revolution and the end of the Ancien Régime, when the landscape was still similar to that of the Middle Ages. Directly in Leaflet, we use addWMSTiles() to request the Ferraris map tiles from the WMS-webservice: leaflet() %>% setView(4.69839084, 50.880089, zoom = 14) %>% # center on Leuven addWMSTiles( baseUrl = "https://geoservices.informatievlaanderen.be/raadpleegdiensten/histcart/wms", layers = list("ferraris"), # request Ferraris-tiles from available WMS tilesets options = WMSTileOptions(format = "image/png"), attribution = "Source: GIS Geoservices Informatie Vlaanderen") %>% addControl("Leuven on the Ferraris map (1777)", position = "topright") We can do the same for the Vandermaelen map. This “Carte topographique de la Belgique” of Belgium was published by Belgian cartographer Philippe Vandermaelen between 1846 and 1854. With the relatively late industrialisation of Leuven, the map provides an nice background to contextualise the location of the gangen in the second half of the century. leaflet() %>% setView(4.69839084, 50.880089, zoom = 14) %>% addWMSTiles( "https://geoservices.informatievlaanderen.be/raadpleegdiensten/histcart/wms", layers = list("vandermaelen"), options = WMSTileOptions(format = "image/png"), attribution = "Source: GIS Geoservices Informatie Vlaanderen") %>% addControl("Leuven on the Vandermaelen map (1846-1854)", position = "topright") Very visible on the Vandermaelen map is the grid-pattern of (intended) streets in the eastern part of Leuven, laid out in the 1838 “plan-Laenen” of city architect François-Henri Laenen.