Author: Paulo van Breugel
Updated on: 24-03-19
R has advanced spatial tools, including tools to create and visualize spatial data. One of these tools is tmap, a R-library for drawing thematic maps. The syntax is based on A Layered Grammar of Graphics by Hadley Wickham and resembles the syntax of ggplot2. For those used to working with GIS tools such as QGIS, the way of working will feel very natural. You create a map by simply combining different layers one by one. The difference is that in most GIS tools you'll use the GUI (GRASS GIS being an exception in that you can use the GUI or command line), while in R you do this on the command line.
The tmap library is largely focussed on vector maps, but it can also handle raster layers, using the excellent raster library in the backend. If you only need to quickly view a raster layer, this library is all you need. However, tmap makes it easier to combine raster layers with vector maps and to fine-tune the layout. In the example below, I'll create a figure with maps for each month with the average rainfall.
Let's first load the libraries we need to create the plot.
# Libraries library(sp) library(sf) library(raster) library(tmap) library(tmaptools)
Next step is to load a vector layer with the national boundaries. We'll use the 'World' layer from the tmap package.
The default projection of this layer is Eckhart IV (World) (EPSG: 54012). Although an excellent choice to get a good representation of the relative size of the continents, we'll use the Robison projection, which provides a compromise between a good representation of area and of the form.
PROJ <- "+proj=robin +lon_0=0 +x_0=0 +y_0=0 +ellps=WGS84 +datum=WGS84 +units=m +no_defs" World <- st_transform(World, crs = PROJ)
Precipitation raster data
precdata <- getData("worldclim",var="prec",res=10)
We need to reproject this data as well, using the *set_projection* function from tmap. You can optionally aggregate the raster values to obtain lower resolution raster layers for quicker printing.
# Aggregate (ommit for final map) precdata <- aggregate(precdata, fact=2, fun = mean) # Reproject the data precdata <- set_projection(precdata, "robin")
Plotting the data
We can now plot the data. Note that to plot all layers in a stack, you simply omit the layer ('col' parameter) in the *tm_raster* function. So, if you want to plot the rainfall for January only, set 'col = prec1'. I'll leave it to you to find out more about the individual commands using the excellent R help files.
tm1 <- tm_shape(precdata) + tm_raster(palette = get_brewer_pal("YlGnBu", n = 6), breaks = c(0, 25, 50, 100, 250, 500, 1500), alpha = 0.6, title = "Rainfall (mm)") + tm_shape(World) + tm_borders("grey20", lwd = 0.5) + tm_style("white", earth.boundary = c(-180, -87, 180, 87), earth.boundary.color = "white", legend.show = TRUE) + tm_format("World", inner.margins = 0.02, frame = FALSE) + tm_layout(frame=FALSE, legend.outside = TRUE, legend.outside.position = "right", legend.outside.size = 0.1, legend.title.size = 0.8, legend.text.size = 0.6, legend.format = list(scientific = TRUE, format = "f")) + tm_facets(ncol=3) tmap_save(tm = tm1, filename = "monthly rainfall3.png", width = 18, height = 11, units="cm", dpi = 150)
And if you have followed the code above, you should get the same map as below. Tmap offers many ways to tweak the outcome. For example, you may want to change the size of the legend and fonts depending on where / how you want to publish the figure. Or you can place the legend below, above or inside the figures. It is all possible. Of course, with great flexibility comes more complexity, but I am sure you'll figure it out with the help of the manual pages or the many examples and Q&A online.
If you have questions
If you have questions or comments about the text, let me know. You can use this contact form. Please make sure to include the page title ("plotting multiple maps in R") or page name ("plotting multiple maps in R").