This tutorial covers how to create interactive maps in R, focusing on two types of visualizations:
There are many ways to create maps using R, but we will focus on leaflet, a popular open-source JavaScript library that is used by many professional organizations to create interactive maps.
In this tutorial we will require the following packages:
library(leaflet) # The map-making package
library(geojsonio) # A package for geographic and spatial data, requires the latest version of dplyr
library(dplyr) # Used for data manipulation and merging
library(htmltools) # Used for constructing map labels using HTMLIn this tutorial we will use two datasets:
data(quakes)
head(quakes)## lat long depth mag stations
## 1 -20.42 181.62 562 4.8 41
## 2 -20.62 181.03 650 4.2 15
## 3 -26.00 184.10 42 5.4 43
## 4 -17.97 181.66 626 4.1 19
## 5 -20.42 181.96 649 4.0 11
## 6 -19.68 184.31 195 4.0 12happiness <- read.csv('https://remiller1450.github.io/data/HappyPlanet.csv', stringsAsFactors = FALSE)
head(happiness)## Country Region Happiness LifeExpectancy Footprint HLY HPI HPIRank
## 1 Albania 7 5.5 76.2 2.2 41.7 47.91 54
## 2 Algeria 3 5.6 71.7 1.7 40.1 51.23 40
## 3 Angola 4 4.3 41.7 0.9 17.8 26.78 130
## 4 Argentina 1 7.1 74.8 2.5 53.4 58.95 15
## 5 Armenia 7 5.0 71.7 1.4 36.1 48.28 48
## 6 Australia 2 7.9 80.9 7.8 63.7 36.64 102
## GDPperCapita HDI Population
## 1 5316 0.801 3.15
## 2 7062 0.733 32.85
## 3 2335 0.446 16.10
## 4 14280 0.869 38.75
## 5 4945 0.775 3.02
## 6 31794 0.962 20.40Generally speaking, there are three steps to creating a map using leaflet:
leaflet functionaddTiles, addMarkers, and addPolygonsLayers are usually added by “piping” via the %>% operator. Put simply, this means passing the resulting object on the left of the %>% symbol as the first argument to the function on the right of the %>%, thereby chaining together a series of functions that each operate on the results of the prior step.
## Step 1 call leaflet
map <- leaflet() %>%
addTiles() ## Step 2 add a layer using %>%
## Step 3 print the map
mapLeaflet constructs a base map using map tiles, an approach popularized by google maps that is now used by nearly every interactive map online.
Calling addTiles() with no arguments will generate a base map using OpenStreetMap tiles. For most purposes OpenStreetMap tiles are all you’ll need; however, the leaflet package also comes with dozens of different tile options created by third party providers and stored in an object named providers. Provider tile sets can be used to create a base map via the addProviderTiles function:
head(providers) ## Providers contains a list of different tiles## $OpenStreetMap
## [1] "OpenStreetMap"
##
## $OpenStreetMap.Mapnik
## [1] "OpenStreetMap.Mapnik"
##
## $OpenStreetMap.BlackAndWhite
## [1] "OpenStreetMap.BlackAndWhite"
##
## $OpenStreetMap.DE
## [1] "OpenStreetMap.DE"
##
## $OpenStreetMap.CH
## [1] "OpenStreetMap.CH"
##
## $OpenStreetMap.France
## [1] "OpenStreetMap.France"map <- leaflet() %>%
addProviderTiles(providers$Esri.WorldTopoMap)
mapThe simplest addition to base map is a marker, or identifiable point. Markers can be added to the base map using addMarkers function. At minimum a marker requires a longitude and latitude, but it is useful to also include labels, which display text when you hover over the marker, or popups, which display a dialog box when you click on the marker.
fiji_reduced <- quakes[1:10,]
map <- leaflet(data = fiji_reduced) %>%
addTiles() %>%
addMarkers(lng = ~long, lat = ~lat,
label = ~paste(depth, "meters"),
popup = ~paste(mag, "Richter Magnitude", "<br/>", depth, "Depth"))
mapWhat is displayed by a popup or label is an interpretted HTML strings, meaning commands like “<br/>” can be used to separate the text onto a new line, and other HTML commands can be used to modify the text’s appearance.
Additionally, in the example above, you should notice the \(\sim\) character that is used in front of the inputs to lng, lat, label, and popup. This tells addMarkers to look for those variable within the data object specified in the call to leaflet. So, in this example, leaflet looks for the variables “lat” and “lat” in fiji_reduced. Later in this tutorial we’ll see examples that don’t do this.
Question #1
Construct a base map using “Esri.OceanBasemap” tiles. Then, filter the quakes dataset to include only observations reported by more than 40 stations and plot the location of these earthquakes on your map with hoverable labels that display the number of stations reporting (ie: text like “47 Stations Reporting”).
Adding a large number of markers to a map is computationally demanding, and it also makes the map difficult to read. To address these and other concerns, leaflet includes a suite of built in functions includingmarkerClusterOptions, which aggregates nearby markers into a group that will split when you click on it to zoom in.
map <- leaflet(data = quakes) %>%
addTiles() %>%
addMarkers(lng = ~long, lat = ~lat,
label = ~paste(depth, "meters"),
clusterOptions = markerClusterOptions())
mapThe appearance of markers can be changed in two ways. The most straightfoward is to use addCircleMarkers, which plots circular points at the marker locations. These markers have attributes like color, opacity, and fill that can be used to add information to the map. The addCircleMarkers function uses syntax similar to addMarkers, a few of its options are demonstrated in the example below:
quakes$col <- ifelse(quakes$mag > 5, "red", "blue")
map <- leaflet(data = quakes) %>%
addTiles() %>%
addCircleMarkers(lng = ~long, lat = ~lat,
label = ~paste(depth, "meters"),
color = ~col,
fillOpacity = 0.25,
stroke = FALSE)
mapHere, we color earthquakes with magnitudes over 5.0 red, while coloring less intense earthquakes blue. We use the argument fillOpacity = 0.25 to make the markers more transparent, allowing viewers to see relative concentrations, and we add the argument stroke = FALSE to remove the solid outline that by default surrounds circular markers.
Custom markers can be added to a map as icons. Icons are created from an image using its web url or local file path. The example below sets up two icons, a ferry icon and a skull with crossbones icon, using images from an online repository.
oceanIcons <- iconList(
ferry = makeIcon(
"http://globetrotterlife.org/blog/wp-content/uploads/leaflet-maps-marker-icons/ferry-18.png"),
danger = makeIcon(
"http://globetrotterlife.org/blog/wp-content/uploads/leaflet-maps-marker-icons/danger-24.png"))
iconid <- ifelse(quakes$mag < 5, "ferry", "danger")
leaflet(data = quakes[1:4,]) %>% addTiles() %>%
addMarkers(~long, ~lat, icon = oceanIcons[iconid])To use an icon, it must first be setup using makeIcon, and then stored in a special list using iconList before being passed to addMarkers.
In this example you should notice that we do not use a \(\sim\) in front of oceanIcons. That is because oceanIcons exists outside of the quakes data frame.
Question #2
Construct a base map using the default map tiles. Then, using your filtered dataset (which contains earthquakes reported by more than 40 stations), construct a map displaying these earthquakes using clustered circle markers which are colored yellow if the earthquake was reported by more than 70 stations and black otherwise. You do not need to add any labels or popups to this map.
A choropleth is a type of map that displays geographic units like countries, states, census areas, or watersheds with each unit colored based upon a particular variable. The first step in creating a choropleth is to define the boundaries of each geographic unit using spatial polygons.
In this example, we load the spatial polygons (an sp object) of country boundaries from a web url using the geojson_read function.
shapeurl <- "https://raw.githubusercontent.com/johan/world.geo.json/master/countries.geo.json"
WorldCountry <- geojson_read(shapeurl, what = "sp")
Map <- leaflet(WorldCountry) %>% addTiles() %>% addPolygons()
MapLearning how to create your own spatial polygons is beyond the scope of this tutorial, but for most applications you can find a shapefile that someone else has created which contains the spatial polygon boundaries you need. Unfortunately I’m not aware of a single repository of different geo.json files, so finding the boundaries you need can sometimes require a bit of searching.
The addPolygons function will add the area boundaries for each geographic unit defined in sp file. To modify the color, or other attributes, of these polygons we can supply additional arguments:
Map <- leaflet(WorldCountry) %>% addTiles() %>%
addPolygons(
fillColor = "red",
weight = 2,
opacity = 1,
color = "white",
fillOpacity = 0.7)
MapYou should recognize that the arguments within addPolygons expects you to provide a vector of length 1 (like we did in the example above) or a vector that is the same length as the number of spatial polygons (what is needed to make a choropleth). In the example below, we merge data from the Happy Planet with the countries in our sp file, and then use the variable LifeExpectancy to create a choropleth:
## Read in the Happy Planet data
happiness <- read.csv('https://remiller1450.github.io/data/HappyPlanet.csv', stringsAsFactors = FALSE)
## Join the Happy Planet data to the countries in WorldCountry
CountryHappy <- left_join(data.frame(Name = WorldCountry$name), happiness, by = c("Name" ="Country"))
## Create a color palette (using the "magma" color scheme)
pal <- colorBin("magma", domain = CountryHappy$LifeExpectancy)
## Use the palette and LifeExpectancy to color the spatial polygons
Map <- leaflet(WorldCountry) %>% addTiles() %>%
addPolygons(
fillColor = pal(CountryHappy$LifeExpectancy),
weight = 2,
opacity = 1,
color = "white",
fillOpacity = 0.7)
MapIn this example you should notice a few things:
left_join because every spatial polygon (country) needs to have a fillColor.WorldCountry had a matching entry in the Happy Planet data, but addPolygons will use grey to color spatial polygons with an NA value for their fill color.colorBin.To add interactivity to our map we can use the highlight argument to indicate when the mouse is hovering over a particular spatial polygon.
Map <- leaflet(WorldCountry) %>% addTiles() %>%
addPolygons(
fillColor = pal(CountryHappy$LifeExpectancy),
weight = 2,
opacity = 1,
color = "white",
fillOpacity = 0.7,
highlight = highlightOptions(
weight = 3,
color = "grey",
fillOpacity = 0.7,
bringToFront = TRUE))
MapIn addition to adding highlights, we can add labels and popups to our spatial polygons.
myLabels <- paste("<strong>", CountryHappy$Name, "</strong>", "<br/>",
"Life Expectancy:", CountryHappy$LifeExpectancy)
myPopups <- paste("Happy Planet Rank", CountryHappy$HPIRank)
Map <- leaflet(WorldCountry) %>% addTiles() %>%
addPolygons(
fillColor = pal(CountryHappy$LifeExpectancy),
weight = 2,
opacity = 1,
color = "white",
fillOpacity = 0.7,
highlight = highlightOptions(
weight = 3,
color = "grey",
fillOpacity = 0.7,
bringToFront = TRUE),
label = lapply(myLabels, HTML),
popup = myPopups)
MapIn this example you should notice:
addPolygons to read the HTML commands properly we need to apply (using lapply) the function HTML from the htmltools package to each element of the vector we supply to the label argument.addPolygons that we’re using HTML.As a final step, we might consider adding a legend to our choropleth. The addLegend function will add a legend as a new layer:
myLabels <- paste("<strong>", CountryHappy$Name, "</strong>", "<br/>",
"Life Expectancy:", CountryHappy$LifeExpectancy)
myPopups <- paste("Happy Planet Rank", CountryHappy$HPIRank)
Map <- leaflet(WorldCountry) %>% addTiles() %>%
addPolygons(
fillColor = pal(CountryHappy$LifeExpectancy),
weight = 2,
opacity = 1,
color = "white",
fillOpacity = 0.7,
highlight = highlightOptions(weight = 3,
color = "grey",
fillOpacity = 0.7,
bringToFront = TRUE),
label = lapply(myLabels, HTML),
popup = myPopups) %>%
addLegend(pal = pal, values = CountryHappy$LifeExpectancy,
title = "Life Expectancy", position = "bottomright")
MapQuestion #3:
The us_cities data set in the geojsonio package contains state capitals and all cities in the United States with populations larger than 40,000. The file located at “https://raw.githubusercontent.com/PublicaMundi/MappingAPI/master/data/geojson/us-states.json” contains spatial polygon boundaries for all 50 states (and the District of Columbia and Puerto Rico). The USArrests data set in the datasets package contains violant crime rates by state.
For this question you should create a map with the following features:
USArrests data set.UrbanPop below the name on a new line.radius = 2) for each city in the us_cities data set where state capitals are colored “red” and all other cities are colored yellowHints:
bringToFront = FALSE in your highlight to avoid masking the city labels when hovering over a stateshapeurl <- "https://raw.githubusercontent.com/PublicaMundi/MappingAPI/master/data/geojson/us-states.json"
UnitedStates <- geojson_read(shapeurl, what = "sp")
data("USArrests")
data("us_cities")