::p_load(tmap, tidyverse, sf) pacman
In-Class Exercise 3: Analytical Mapping
0.1 Overview
…
0.2 Installing and loading packages
The tmap package documentation can be found here.
0.3 Importing data
<- read_rds("data/rds/NGA_wp.rds") NGA_wp
0.4 Visualizing distribution of functional and non-functional water point using Choropleth Maps
<- tm_shape(NGA_wp) +
p1 tm_fill("wp_functional",
n = 10,
#style is used for data classification method
style = "equal",
#colour palettes are always plural
palette = "Blues") +
# line width
tm_borders(lwd = 0.1,
# opacity
alpha = 1) +
# main.title places the title outside the plot
tm_layout(main.title = "Distribution of functional water point",
# legend.outside puts legends in- or outside of your plot
legend.outside = FALSE)
# things to take note: tm_fill() and tm_borders() combined is tm_polygon()
p1
<- tm_shape(NGA_wp) +
p1 tm_fill("wp_functional",
n = 10,
#style is used for data classification method
style = "equal",
#colour palettes are always plural
palette = "Blues") +
# line width
tm_borders(lwd = 0.1,
# opacity
alpha = 1) +
# main.title places the title outside the plot
tm_layout(main.title = "Distribution of functional water points by LGAs",
# legend.outside puts legends in- or outside of your plot
legend.outside = FALSE)
# things to take note: tm_fill() and tm_borders() combined is tm_polygon()
# p1 is a map object
p1
<- tm_shape(NGA_wp) +
p2 tm_fill("total_wp",
n = 10,
#style is used for data classification method
style = "equal",
#colour palettes are always plural
palette = "Blues") +
# line width
tm_borders(lwd = 0.1,
# opacity
alpha = 1) +
# main.title places the title outside the plot
tm_layout(main.title = "Distribution of total water points by LGAs",
# legend.outside puts legends in- or outside of your plot
legend.outside = FALSE)
p2
# show both together on one row
tmap_arrange(p2, p1, nrow=1)
# You can even show both maps as one by using Rates: usually we do this to see if both categories tell the same story
<- NGA_wp %>%
NGA_wp
# use mutate() to calc % (i.e., rate) of functional and non-functional water points
mutate(pct_functional = wp_functional/total_wp) %>%
mutate(pct_nonfunctional = wp_nonfunctional/total_wp)
tm_shape(NGA_wp) +
tm_fill("pct_functional",
n = 10,
style = "equal",
palette = "Blues",
legend.hist = TRUE) +
tm_borders(lwd = 0.1,
alpha = 1) +
tm_layout(main.title = "Rate map of functional water points by LGAs",
legend.outside = TRUE)
0.5 Visualization using a Percentile Map
Tells you which areas are the top 10%. The six specific categories are: 0-1%, 1-10%, 10-50%, 50-90%, 90-99%, and 99-100%.
# Data Preparation
# Exclude NA values
<- NGA_wp %>%
NGA_wp drop_na()
# create customised classification and extract values
<- c(0, .01, .1, .5, .9, .99, 1)
percent <- NGA_wp["pct_functional"] %>%
var
#
st_set_geometry(NULL)
quantile(var[,1], percent)
0% 1% 10% 50% 90% 99% 100%
0.0000000 0.0000000 0.2169811 0.4791667 0.8611111 1.0000000 1.0000000
0.5.1 Writing functions to do the same functionality for specific data sets…
# R function to extract a variable (i.e., wp_nonfunctional as a vector out of an sf data.frame)
<- function(vname, df) {
get.var <- df[vname] %>%
v st_set_geometry(NULL)
<- unname(v[,1])
v
# return vector with values (without a col name)
return(v)
}
# percentile mapping function
<- function(vnam, df, legtitle=NA, mtitle="Percentile Map"){
percentmap <- c(0,.01,.1,.5,.9,.99,1)
percent <- get.var(vnam, df)
var <- quantile(var, percent)
bperc tm_shape(df) +
tm_polygons() +
tm_shape(df) +
tm_fill(vnam,
title=legtitle,
breaks=bperc,
palette="Blues",
labels=c("< 1%", "1% - 10%", "10% - 50%", "50% - 90%", "90% - 99%", "> 99%")) +
tm_borders() +
tm_layout(main.title = mtitle,
title.position = c("right","bottom"))
}
# test the function
percentmap("total_wp", NGA_wp)
0.6 Visualizing using a Box Plot
# boxplot is as augmented quartile map with an additional lower and upper category.
ggplot(data = NGA_wp,
aes(x = "",
y = wp_nonfunctional)) +
geom_boxplot()
# creating a boxbreaks function
# arguments - v: vector with observations, mult: multiplier for IQR (default 1.5)
# returns - bb: vector with 7 breakpoints compute quartile and fences
<- function(v,mult=1.5) {
boxbreaks <- unname(quantile(v))
qv <- qv[4] - qv[2]
iqr <- qv[4] + mult * iqr
upfence <- qv[2] - mult * iqr
lofence # initialize break points vector
<- vector(mode="numeric",length=7)
bb # logic for lower and upper fences
if (lofence < qv[1]) { # no lower outliers
1] <- lofence
bb[2] <- floor(qv[1])
bb[else {
} 2] <- lofence
bb[1] <- qv[1]
bb[
}if (upfence > qv[5]) { # no upper outliers
7] <- upfence
bb[6] <- ceiling(qv[5])
bb[else {
} 6] <- upfence
bb[7] <- qv[5]
bb[
}3:5] <- qv[2:4]
bb[return(bb)
}
# creating the get.var function
# arguments - vname: variable name (in quotes), df: name of sf data.frame
# returns - v: vector with values (without col name)
<- function(vname,df) {
get.var <- df[vname] %>% st_set_geometry(NULL)
v <- unname(v[,1])
v return(v)
}
# test function
<- get.var("wp_nonfunctional", NGA_wp)
var boxbreaks(var)
[1] -56.5 0.0 14.0 34.0 61.0 131.5 278.0
0.6.1 Boxmap function
<- function(vnam, df,
boxmap legtitle=NA,
mtitle="Box Map",
mult=1.5){
<- get.var(vnam,df)
var <- boxbreaks(var)
bb tm_shape(df) +
tm_polygons() +
tm_shape(df) +
tm_fill(vnam,title=legtitle,
breaks=bb,
palette="Set3",
labels = c("lower outlier",
"< 25%",
"25% - 50%",
"50% - 75%",
"> 75%",
"upper outlier")) +
tm_borders() +
tm_layout(main.title = mtitle,
title.position = c("left",
"top"))
}
tmap_mode("plot")
boxmap("wp_nonfunctional", NGA_wp)
0.7 Recode to zero
This code chunk is used to recode LGAs with zero total water points into NA.
<- NGA_wp %>%
NGA_wp mutate(wp_functional = na_if(
< 0)) total_wp, total_wp