Section 4.3 Bivariate maps
Usually in your thematic maps you are mapping one variable. In the distortion examples above, we incorporate a second variable in order to make the map more legible. However, occasionally, you might want to use your map to illustrate the relationship between two variables. In this case, you may want to create a Bivariate Choropleth Map. Bivariate choropleths follow the same concept as the univariate ones (which display only one variable), except they show two variables at once. This is achieved through the creative use of colour.
We first came across the idea of bivariate choropleth maps from the blog of [Stevens (2015)]. We recommend having a read as he discusses in great detail how to develop the colour schemes necessary for these maps. The idea is that there are two main colours, each one representing gradual change in one of our two variables of interest. Then, these two colours are combined to create an overlapping colour scheme.
The process itself is not too complicated. We start with binning our two variables into classes. If we bin variable 1 into \(n\) classes, and then again variable 2 into \(n\) classes, then when we compare them across one another to create our bivariate map, we will end up with \(n^2\) classes. For example, if we bin both variables into three classes (low, medium, high), then when displaying them together, we will have nine classes to visualise. The blog by [Stevens (2015)] illustrates with nice visuals, so we recommend having a look. However, the practical example is in QGIS, and here we are working in R. There have been adaptations into R (e.g., see [Grossenbacher (2019)]) which we can borrow from here too.
Creating the map boils down to five key steps. First, we take our two variables of interest, and create bins. Second, we create a new variable, which we will use for the shading. Third we create our colour scheme. Fourth, we join this to our
sf object. And, we map! Let’s go through these steps now. Using the example of voting to leave the EU.
It may be an interesting question to look into not only how the voting outcome was distributed (
Pct_Leave) but how this varies with voter turnout (Pct_Turnout). We might be interested in areas with high turnout and high percentage voting to leave, as these may be areas where people felt passionately. On the other hand, other areas may have had low turnout, which may have influenced the result, and in those places, how did people who did turn up vote? These are the kinds of questions we can answer with a bivariate map.
Subsection 4.3.1 Step 1: Bin our variables of interest
We are interested in two key variables:
Pct_Leave - the percentage of people who voted to leave the EU, and Pct_Turnout - the percentage of the electorate who actually voted. These are both numeric continuous variables. In order to create our bivariate choropleth map, we have to bin these values into \(n\) discrete categories. Here let’s go with \(n = 3\text{.}\)
We can use the
cut() and the quantile() functions in base R to class our variable into three quantiles. The quantile() function identifies the sample quantiles in a continuous variable. We need to include the parameters x=: the numeric vector whose sample quantiles are wanted (in this case, our variables Pct_Leave and Pct_Turnout), and probs=: numeric vector of probabilities with values in [0,1]. We can use the sequence generator function seq() to generate these from 0 to 1, by the increment of \(\frac{1}{3}\) for 3 groups. So first, let’s create the breaks for the Pct_Leave variable.
leave_breaks <- quantile(eu_sf$Pct_Leave,
probs = seq(0,1, by = 1/3), # probabilities
na.rm=TRUE, # remove NAs
names=TRUE, # keep names attribute
include.lowest=TRUE) # include value equal to
# lowest ‘breaks’ value
And then the same again but for the turnout variable
Pct_Turnout.
turnout_breaks <- quantile(eu_sf$Pct_Turnout,
probs = seq(0,1, by = 1/3),
na.rm=TRUE,
names=TRUE,
include.lowest=TRUE)
We can have a look at the output:
leave_breaks
## 0% 33.33333% 66.66667% 100% ## 21.38 51.72 59.12 75.56
The results are the cutoff values which we want to use to “cut” our numeric variable. We do this with the
cut() function, , where we specify again what to cut (the variables Pct_Leave and Pct_Turnout), and the breaks at which to cut (the objects we created above, leave_breaks and turnout_breaks):
eu_sf <- eu_sf %>%
mutate(leave_quantiles = cut(Pct_Leave,
breaks = leave_breaks),
turnout_quantiles = cut(Pct_Turnout,
breaks = turnout_breaks))
We have two resulting variables,
leave_quantiles and turnout_quantiles which classify each one of our observations into one of these quartiles for both the variables. In the next step, we use these to create a new variable.
Subsection 4.3.2 Step 2: New variable
In this step we create a new variable called
group, which tells us which quartile the specific Local Authority falls into. By applying the as.numeric() function we translate the ranges of the quartile into their label (i.e., 1st, 2nd, or 3rd quartile). We do this for both variables, and paste them together using the paste() function, and the separator “-”:
eu_sf <- eu_sf %>%
mutate(group = paste(
as.numeric(turnout_quantiles), "-",
as.numeric(leave_quantiles))
)
We now have a new column, called group, which tells us for each Local Authority, which quartile it falls into for each variable. For example, a value of “1 - 1” means the Local Authority belongs to the first quartile in both variables. This area would be considered to have low percentage voting leave, and also low turnout. On the other hand, “3 - 1” means that there was high turnout, but a low percentage voted to leave. We use this variable to assign the appropriate colour for our colour scheme for each of the \(3^2\) (9) combinations.
Subsection 4.3.3 Step 3: Create colour scheme
Picking a colour scheme which reflects both gradual change in each individual variable, and the combined change in both is not an easy task! Luckily, [Stevens (2015)] has created some scale recommendations for us to choose from. We copy two of them below, but you can see the blog for another 2 options.
In this code below, we specify for each of the nine values of the variable created in the previous step (“1 - 1”, “1 - 2”, “1 - 3”, “2 - 1”, ... “3 - 3”) an associated colour using the relevant hex code.
library(tibble)
bivariate_color_scale_1 <- tibble(
"3 - 3" = "#574249", # high - high
"2 - 3" = "#985356",
"1 - 3" = "#c85a5a", # low - high
"3 - 2" = "#627f8c",
"2 - 2" = "#ad9ea5", # medium - medium
"1 - 2" = "#e4acac",
"3 - 1" = "#64acbe", # high - low
"2 - 1" = "#b0d5df",
"1 - 1" = "#e8e8e8" # low - low
) %>%
gather("group", "fill_col")
bivariate_color_scale_2 <- tibble(
"3 - 3" = "#3b4994", # high - high
"2 - 3" = "#5698b9",
"1 - 3" = "#5ac8c8", # low - high
"3 - 2" = "#8c62aa",
"2 - 2" = "#a5add3", # medium - medium
"1 - 2" = "#ace4e4",
"3 - 1" = "#be64ac", # high - low
"2 - 1" = "#dfb0d6",
"1 - 1" = "#e8e8e8" # low - low
) %>%
gather("group", "fill_col")
You can have a look at both of the colour scales we just created, and see which one you like. Feel free to use either in your future bivariate mapping adventures. Or construct your own, perhaps following the two additional ones provided by [Stevens (2015)], or some entirely new ones you may have constructed.
Subsection 4.3.4 Step 4: Join colour scheme
Now that we have a colour scheme, we can join to the spatial object, using
left_join(). The common element is the group variable, so with this approach, we join the relevant colour to each value of group in that column. Let’s join bivariate_color_scale_2, the second of the two we created above.
eu_sf <- left_join(eu_sf, bivariate_color_scale_2, by = "group")
We now have an additional column in our
eu_sf dataframe, called fill_col which we can use to shade each local authority according to the composite variable depicting the percentage who voted to leave and the percentage who turned up to vote.
Subsection 4.3.5 Step 5: Create legend
The legend is a little tricky, as we need to separate out the values into separate Leave and Turnout columns. We can achieve this with the
separate() function:
# separate the groups
bivariate_color_scale <- bivariate_color_scale_2 %>%
separate(group, into = c("Pct_Turnout", "Pct_Leave"), sep = " - ")
Then to create the legend, we actually build a
ggplot() object. This genius bit of code is borrowed from Grossenbacher (2019) implementation of bivariate choropleth maps.
legend <- ggplot() +
geom_tile( data = bivariate_color_scale,
aes(x = Pct_Turnout, y = Pct_Leave, fill = fill_col)) +
scale_fill_identity() +
labs(x = "Higher % turnout ->",
y = "Higher % voted leave ->") +
theme(axis.title = element_text(size = 6)) + # makes text small for
# adding legend to map
coord_fixed() # forces a specified ratio to create quadratic tiles
This code returns the legend as a chart itself. Have a look at what it looks like:
legend

Subsection 4.3.6 Step 6: Map
Now finally we put it all on the map. For the choropleth map, we use our variable
fill_col which contains the matched colour to the group that each observation belongs to. We pass this in the familiar geom_sf() geometry and use the fill= parameter to colour the Local Authorities according to their turnout / voted leave combination. We also have to add the scale_fill_identity() function, as the values in the fill_col variable are actually the hex codes for the colour which we use to shade the Local Authorities.
map <- ggplot(eu_sf) +
geom_sf(aes( fill = fill_col)) +
scale_fill_identity() +
theme_void()
Finally, to display the legend and the map together, we can use the
ggdraw() and draw_plot() functions from the cowplot package.
library(cowplot)
ggdraw() +
draw_plot(map, 0, 0, 1, 1) +
draw_plot(legend, 0.05, 0.075, 0.2, 0.2)

This map may now be able to provide insight into spatial patterns in turnout and voting to leave. For example, in the South you can see lots of pink, representing areas of high turnout and low percentage voting to leave the EU. You can also spot the dark blue areas; these are Local Authorities which saw high voter turnout and a high proportion voting to leave.
Overall, these maps can help visualise two variables on one map, and motivate discussion about relationships between variables in different places.
