Section 20 Maps with leaflet
For spatial data corresponding to geographical locations, it can be helpful to visualise the data on a map. There are various mapping tools available in R. Here, we will discuss the R package leaflet (Cheng, Karambelkar, and Xie 2019), which is an R implementation of the JavaScript library of the same name.
Maps produced using leaflet are best suited for web pages, as you can zoom and scroll a map as you would do with, say, Google maps. These maps can also be used inside shiny apps.
You will need to be online to use leaflet.
20.1 A basic map
We’ll first load the leaflet package.
library(leaflet)To get a basic map, we will normally need to specify the longitude and latitude of the centre of our map, as well the initial zoom level (which may take some trial and error). For example, to obtain a map of Sheffield:
leaflet() %>%
addTiles() %>%
setView(lng = -1.473798,
lat = 53.38252,
zoom = 13)The addTiles() specifies what is drawn on each map ‘tile’, with the default being streets and points of interest provided by OpenStreetMap. For example, to change to an aerial photo (as in Google Earth), we can instead use the addProviderTiles() command, and specify a different ‘tile provider’:
leaflet() %>%
setView(lng = -1.473798, lat = 53.38252, zoom = 13) %>%
addProviderTiles(providers$Esri.WorldImagery)20.2 Indicating points on a map
If we have vectors of coordinates (latitude and longitude), we can add points with the addMarkers() function (there are variants such as addCircles() and addPopups().) We can include the argument popup so that some text is displayed if the mouse is clicked on the marker. Here’s an example where I’ve marked where I work, and my favourite pub.
latitude <- c(53.380909, 53.371566)
longitude <- c(-1.486197, -1.577604)
placeNames <- c("Hicks Building", "The Sportsman")
leaflet() %>%
setView(lng = -1.525, lat = 53.38, zoom =12) %>%
addTiles() %>%
addMarkers(lat = latitude,
lng = longitude,
popup = placeNames)Exercise 20.1 Use leaflet to make a map that displays a region and two points of interest to you, for example, where you live, and your nearest railway station.
20.3 Plotting region boundaries
We can plot region boundaries on maps, and use different colours to make the regions more distinctive, or to indicate different values of some variable in different regions (a “choropleth map”). This can be a little tricky, in that the data sets that indicate region boundaries can be quite complicated.
Here’s an example of plotting ‘county and unitary’ authorities in England and Wales. To start with, we’ll use a base map that uses plainer-looking tiles than OpenStreetMap:
leaflet() %>%
setView(lng = -1.47, lat = 52.3, zoom = 6) %>%
addProviderTiles(providers$Esri.WorldGrayCanvas)20.3.1 shapefiles
Boundary data can be stored in different formats. We will search for a format called a shapefile. If we try a Google search for “england county boundaries shapefile”, we can find this page at data.gov.uk, and a zip file containing the shapefile can be downloaded from here.
20.3.2 Importing and simplifying shapefiles
Shapefiles can be imported using the rgdal package (Bivand, Keitt, and Rowlingson 2020) and the function readOGR(). The argument dsn is the path to the folder containing the shapefile, and the argument layer is the name of the shapefile in that folder (without the .shp extension). We’ll suppose here that the path to the folder is "~/Desktop/Counties": we do
boundaries <- rgdal::readOGR(dsn = "~/Desktop/Counties",
layer = "Counties_and_Unitary_Authorities_(December_2016)_Boundaries")## OGR data source with driver: ESRI Shapefile
## Source: "/Users/st1jeo/Desktop/Counties", layer: "Counties_and_Unitary_Authorities_(December_2016)_Boundaries"
## with 174 features
## It has 10 fieldsThis creates an object of class SpatialPolygonsDataFrame: it’s more complicated than an ordinary dataframe, and we won’t worry too much about its structure. But if we try
str(boundaries@data)## 'data.frame': 174 obs. of 10 variables:
## $ objectid : int 1 2 3 4 5 6 7 8 9 10 ...
## $ ctyua16cd : chr "E06000001" "E06000002" "E06000003" "E06000004" ...
## $ ctyua16nm : chr "Hartlepool" "Middlesbrough" "Redcar and Cleveland" "Stockton-on-Tees" ...
## $ ctyua16nmw: chr NA NA NA NA ...
## $ bng_e : int 447157 451141 464359 444937 428029 354246 362744 369490 332763 511894 ...
## $ bng_n : int 531476 516887 519597 518183 515649 382146 388456 422806 436633 431716 ...
## $ long : num -1.27 -1.21 -1.01 -1.31 -1.57 ...
## $ lat : num 54.7 54.5 54.6 54.6 54.5 ...
## $ st_areasha: num 9.36e+07 5.39e+07 2.45e+08 2.05e+08 1.97e+08 ...
## $ st_lengths: num 71707 43841 97993 119582 107206 ...we can see we have 174 regions (which sounds right!)
The shape file is a large object; it may help to simplify it, so that our code will run a little more quickly. We can do this with the rmapshaper package (Teucher and Russell 2020).
simplifiedBoundaries <-rmapshaper::ms_simplify(boundaries)(it’s worth checking that we get the same-looking map, whether we use boundaries or simplifiedBoundaries)
20.3.3 Plotting the boundaries
We can now add the boundaries to our map, using the command addPolygons():
leaflet() %>%
setView(lng = -1.47, lat = 52.3, zoom = 6) %>%
addProviderTiles(providers$Esri.WorldGrayCanvas) %>%
addPolygons(data = simplifiedBoundaries,
color = "blue",
fillOpacity = 0,
weight = 1,
popup = ~ctyua16nm)For the arguments to addPolygons, we include
fillOpacity: a number between 0 and 1. Set to 0 for no colouring within each region;weight: the thickness of the boundary lines;popup = ~ctyua16nm: if the mouse is clicked inside a region, display the value in thectyua16nmcolumn in thesimplifiedBoundaries@datadata frame.
20.3.4 Using different colours for regions
We can use a different colour for each region, depending on some other data value. Here, we will just use blue for English counties and red for Welsh ones. If we inspect simplifiedBoundaries@data$ctyua16nmw, we can see that the Welsh counties are region numbers 153 to 174; we’ll create c a colour vector manually:
country <- rep("blue", 174)
country[153:174] <- "red"We’ll redraw the map, increasing fillOpacity slightly above 0:
leaflet() %>%
setView(lng = -1.47, lat = 52.3, zoom = 6) %>%
addProviderTiles(providers$Esri.WorldGrayCanvas) %>%
addPolygons(data = simplifiedBoundaries,
fillOpacity = 0.05,
weight = 1,
popup = ~ctyua16nm,
fillColor = country,
color = country)20.4 Further reading
- The leaflet R package is well documented: a manual is available here.
- Although not necessary for using the R package, if you want to understand more about leaflet itself, documentation and tutorials are available at The JavaScript leaflet homepage.
20.5 Acknowledgement
The shapefile was provided by the Office for National Statistics and obtained from
https://data.gov.uk/dataset/11302ddc-65bc-4a8f-96a9-af5c456e442c/counties-and-unitary-authorities-december-2016-full-clipped-boundaries-in-england-and-wales, accessed 13/09/21. Contains public sector information licensed under the Open Government Licence v3.0.