Choices for Visual Mapping in R

R is a powerful and free program with many packages capable of producing attractive choropleths. This introduction provides an overview of some of the most popular packages for creating choropleth maps. Rather than being an exhaustive look at each, it will provide an impression of each with a view to helping beginners get a sense of the code required for each package.

The example data used each case is data from the 2004 US Presidential election. This data is included in the GWmodel package as a shapefile, also obtainable at http://www.personal.psu.edu/users/a/c/acr181/election.html. I chose this package and it’s data as I will be using it for future demonstrations of geographically-weighted statistics.

In all examples, the shapefile is saved locally in my workspace in a folder called ‘Projectmaps’ with the file itself named US3.

If you wish to recreate these maps with your own data, you should as a minimum:

Very often the data you want to map will be separate from the geometry you need to map it on. In such cases you will probably want to merge your separate files, e.g with a spatial merge.

R’s Base Graphics

Despite the numerous packages for creating choropleths, the basic graphics contained in R itself have the ability to do this efficiently. Because of the simple argumentation with basic plotting, maps can be produced with little code. The longest part is defining your own colour palette based on your data and passing this to to the cols argument.

library(RColorBrewer)
library(rgdal)
US3 <- readOGR(dsn="Projectmaps", layer = "US3")
mycolours <- brewer.pal(8, "Blues")
mybreaks <- c(0, 20, 40, 50, 60, 80, 100)
cut(US3$Bush_pct, mybreaks)
mycolourscheme <- mycolours[findInterval(US3$Bush_pct, vec = mybreaks)]
plot(US3, col = mycolourscheme, main = "Percentage Vote Share for Bush - 2004", cex = 5)
legend(-119, 31.5, legend = levels(cut(US3$Bush_pct, mybreaks)), fill = mycolourscheme, cex = 0.8, title = "% vote for Bush")

spplot

Spplot is an extension of the basic plot but specifically for spatial objects. I wasn’t a fan of the first output above with R’s basic graphics, but I found that removing the polygon borders improved the aesthetic. Where your spatial areas are larger - e.g US states or British counties, polygon borders are generally a good thing, giving important context to your data. With spplot you need to pass your colour palette to col.regionsin order to make use of RColorBrewer.

library(sp)
library(RColorBrewer)
mycolours <- brewer.pal(8, "YlOrRd")
spplot(US3,"Bush_pct", par.settings = list(axis.line = list(col ="transparent")), main = "Percentage vote for Bush, 2004", cuts = 5, col ="transparent", col.regions = mycolours)

## GISTools

The GISTools package is even easier to use than R’s basic graphics. Shading in particular is made easier via the auto.shading command in combination with the RColorBrewer package. The placing of the legend is usually done by specification of Long lat coordinates (or easting nothings). However, as pointed out to me by Chris Brunsdon(the author of this package), choro.legend borrows from R’s base legend function, so you can instead use the argument px= to force placement of the legend using a word e.g ‘top’, ‘bottom’, or ‘topleft’. This is demonstrated below. Alternatively I could have specified the same coordinates used in plot number 1.

library(GISTools)
myshading = auto.shading(US3$Bush_pct, n=7,
                         cols=brewer.pal(7, "Blues"))
choropleth(US3, US3$Bush_pct,shading=myshading, main = "Percentage Vote Share for Bush - 2004")
choro.legend(px='bottomleft', sh=myshading,fmt="%4.1f",cex=0.8)

TMap

Tmap is very underrated. It has similar syntax to the popular ggplot2 but will also produce a reasonable map with only a couple of lines of code. A default colour scheme will be given where this is not specified and passed to tm_polygons and a legend will also be created by default. For smaller maps tm_text("AREA LABEL HERE") can easily overlay county or state names, and these labels can be told to vary according to the size of the polygon with the size argument in the tm_text later. Tmap has a number of styles you can play around with, including ‘classic’ and ‘natural’. We could include a title with the tm_layout layer, but you may feel that the legend title suffices.

Although I went straight for the tm_shape command, an even simpler way you can to use qmap if you don’t have many additional commands and arguments to specify.

library(tmap)
tm_shape(US3) + tm_polygons(col='Bush_pct', title = "Percentage vote share for Bush - 2004", palette = "Spectral") + tm_style_classic() + tm_scale_bar(position = c("right", "top")) 

Ggplot2

Ggplot2 is undoubedly the most popular package for producing choropleths in R, and very powerful. It is used for all sorts of plotting like the base graphics package, and unlike tmap is not designed solely for mapping. Aside from the base graphics method, ggplot2 involves the most steps to produce a basic choropleth. In particular it needs to be transformed to a data frame using fortify. Once converted to a dataframe with the fortify command, we give each county an id based on row numbers so we can then join the bulk of our data back to the dataframe with merge.

library(ggplot2)
US3@data$id <- rownames(US3@data)
newUS <- fortify(US3, region = "id")
newdf <- merge(newUS, US3@data, by = "id")
myplot <- ggplot() +
  geom_polygon(data = newdf, aes(fill = Bush_pct, 
                                    x = long, 
                                    y = lat, 
                                    group = group)) + ggtitle("Vote Percentage for Bush, 2004") + theme(plot.title = element_text(hjust =0.5))
myplot

The first effort with ggplot (not rendered here) was pretty horrible, but by adding coord_map() we can attempt to fix the distortion of our plot. If you’re not a fan of ggplots distinctive automatic grey background we can do something about that too with the theme_bw or theme_nothing functions. The final product is below. Take note of the excellent viridis package for providing additional colour scaling.

library(ggplot2)
library(scales)
library(ggmap)
library(viridis)
US3@data$id <- rownames(US3@data)
newUS <- fortify(US3, region = "id")
newdf <- merge(newUS, US3@data, by = "id")
Myplot <- ggplot() +

  geom_polygon(data = newdf, aes(fill = Bush_pct, 
                                    x = long, 
                                    y = lat, 
                                    group = group)) +
  theme_nothing(legend = TRUE) + coord_map() + ggtitle("Vote Percentage for Bush, 2004") + theme(plot.title = element_text(hjust =0.5))

NicerPlot <- Myplot + scale_fill_viridis(option = "magma", direction = -1)
NicerPlot

Leaflet

Finally, we can have a look at the possibilities of interactive mapping with leaflet. This package has grown significantly in popularity in recent years and has fast become common currency amongst companies wishing to dynamically visualize its data. It is an excellent option to consider where the patterns in your data are large and complex and where you have constituent polygons of varying sizes.

library(rgdal)
library(RColorBrewer)
library(leaflet)
US3 <- readOGR(dsn="Projectmaps", layer = "US3")
## OGR data source with driver: ESRI Shapefile 
## Source: "Projectmaps", layer: "US3"
## with 3111 features
## It has 71 fields
bins <- c(0, 20, 30, 40, 50, 60, 70, 80, 90, Inf)
pal <- colorBin("YlOrRd", domain = US3$Bush_pct, bins = bins)

labels <- sprintf(
  "<strong>%s</strong><br/>%g Bush",
  US3$NAME, US3$Bush_pct
) %>% lapply(htmltools::HTML)

leaflet(US3) %>%
  setView(-94, 39, 5) %>%
  addProviderTiles("MapBox", options = providerTileOptions(
    id = "mapbox.light",
    accessToken = Sys.getenv('MAPBOX_ACCESS_TOKEN'))) %>%
addPolygons(
  fillColor = ~pal(Bush_pct),
  weight = 1,
  opacity = 1,
  color = "white",
  dashArray = "3",
  fillOpacity = 0.7,
  highlight = highlightOptions(
    weight = 5,
    color = "#666",
    dashArray = "",
    fillOpacity = 0.7,
    bringToFront = TRUE),
  label = labels,
  labelOptions = labelOptions(
    style = list("font-weight" = "normal", padding = "3px 8px"),
    textsize = "15px",
    direction = "auto")) %>%
addLegend(pal = pal, values = ~Bush_pct, opacity = 0.7, title = NULL,
                position = "bottomright")

Final Thoughts

This has been a simple introduction to approaches for creating choropleth maps in R with code which users can easily adapt to their own projects. Don’t worry if the code doesn’t make sense, try and recreate the example or apply the existing code to your own data. Some final thoughts: