Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

hamishgibbs add mapping #11

Open
wants to merge 6 commits into
base: master
Choose a base branch
from
Open
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
Binary file added mapping/.DS_Store
Binary file not shown.
5 changes: 5 additions & 0 deletions mapping/add_basemap.R
Original file line number Diff line number Diff line change
@@ -0,0 +1,5 @@
add_basemap <- function(scale = 'medium'){
basemap_data <- rnaturalearth::ne_countries(scale = scale, returnclass = 'sf')

return(geom_sf(data = basemap_data, size = 0.15, fill = '#EBEBEB'))
}
21 changes: 21 additions & 0 deletions mapping/define_fill_continuous_sequential.R
Original file line number Diff line number Diff line change
@@ -0,0 +1,21 @@
define_fill_continuous_sequential <- function(data, attribute, palette){
#Define the breaks, limits, and labels of a continuous fill scale on an entire dataset
#data - data.frame, dataset to define scale attributes
#attribute - character, column to be symbolised by the colour scale
#palette - character, palette from colorspace::hcl_palettes
#Useful for sharing a predefined colour scale between plots.

attr_max <- data %>% pull(!!dplyr::enquo(attribute)) %>% max()
attr_min <- data %>% pull(!!dplyr::enquo(attribute)) %>% min()
scale_breaks <- seq(attr_min, attr_max, length.out = 5)
scale_labels <- as.character(round(scale_breaks, 2))
scale_range <- c(attr_min, attr_max)

scale_colours <- colorspace::sequential_hcl(5, palette = palette)

scale <- ggplot2::scale_fill_gradientn(breaks = as.character(scale_breaks), labels = scale_labels,
limits = scale_range, na.value = "#AAAAAA", colors = scale_colours)

return(scale)

}
9 changes: 9 additions & 0 deletions mapping/get_country_data.R
Original file line number Diff line number Diff line change
@@ -0,0 +1,9 @@
#wrapper for raster::getData to get admin level 0 data by iso3c code

get_country_data <- function(iso){

country <- raster::getData(name = 'GADM', country = iso, level = 0, path = '~/Documents/GADM') %>% sf::st_as_sf()

return(country)

}
105 changes: 105 additions & 0 deletions mapping/map_regions.R
Original file line number Diff line number Diff line change
@@ -0,0 +1,105 @@
#NOTE - rmapshaper::ms_simplify(snap_interval = 0.05) is a manually defined paramter
#if sliver polygons appear or islands are dropped - this will need to be defined
#dynamically or updated manually.

#currently - throwing an error for regional plots when adding legend
#working on: get legend from a plot - apply it to others

#get_country_data stores files locally - need to specify the correct location to store geodata
#in the environment generating the reports.

suppressPackageStartupMessages({
require(dplyr)
require(sf)
require(ggplot2)
require(rgeos)
require(rmapshaper)
require(colorspace)
library(docstring)
require(ggpubr)
})

source('mapping/validate_geometry.R', local = T)
source('mapping/get_country_data.R', local = T)
source('mapping/scale_fill_attribute.R', local = T)
source('mapping/add_basemap.R', local = T)
source('mapping/scale_axes_input_bounds.R', local = T)

.args <- if (interactive()) c(
"mapping/regionref.rds"
) else commandArgs(trailingOnly = TRUE)

#if maps will be produced on a continent level - mapping will start from a reference file like afrkey
#boundary data will be GADM

afrkey <- readRDS(.args[1])

countries <- lapply(afrkey %>% pull(iso), get_country_data)

countries <- lapply(countries, ms_simplify)

countries <- do.call(rbind, countries) %>%
rename(iso = GID_0,
name_0 = NAME_0) %>%
as_tibble() %>%
st_as_sf() %>%
validate_geometry()

countries <- countries %>%
left_join(afrkey, by = c('iso' = 'iso')) %>%
#join attributes here
mutate(dummy = runif(54, 0, 1))

regions <- countries %>%
group_by(region) %>%
group_split(keep = TRUE)

region_polys <- countries %>%
group_by(region) %>%
summarise() %>%
ms_simplify(snap_interval = 0.05) %>%
#join attributes here
mutate(dummy = runif(5, 0, 1))

create_country_plot <- function(map_data, attribute, colour_scale){

p <- ggplot() +
add_basemap() +
geom_sf(data = map_data, aes(fill = map_data %>% pull(!!dplyr::enquo(attribute))), size = 0.2, colour = 'black') +
colour_scale +
scale_x_input_bounds(map_data) +
scale_y_input_bounds(map_data) +
theme_void() +
theme(legend.position = 'none')


return(p)

}

create_region_plot <- function(region_polys, attribute, colour_scale){

p <- ggplot() +
add_basemap() +
geom_sf(data = region_polys, aes(fill = region_polys %>% pull(!!dplyr::enquo(attribute))), size = 0.2, colour = 'black') +
colour_scale +
scale_x_input_bounds(region_polys) +
scale_y_input_bounds(region_polys) +
theme_void() +
theme(legend.position = 'none')

return(p)

}

#define colour scale for whatever attribute by overall dataset (continent) - not region
countries_fill_scale <- scale_fill_attribute(countries, 'dummy')
regional_fill_scale <- scale_fill_attribute(region_polys, 'dummy')

regional_plots <- lapply(regions, create_country_plot, attribute = 'dummy', colour_scale = countries_fill_scale)

region_plot <- create_region_plot(region_polys = region_polys, attribute = 'dummy', regional_fill_scale)

region_panel <- cowplot::plot_grid(plotlist = regional_plots)


5 changes: 0 additions & 5 deletions mapping/region_setup.R
Original file line number Diff line number Diff line change
Expand Up @@ -14,9 +14,4 @@ afrregions <- regionref[region %like% "Africa"][-excl][,.SD,.SDcols=c("alpha-3",
afrregions[`intermediate-region` == "", `intermediate-region` := `sub-region` ]
afrkey <- setkey(afrregions[,.(region = `intermediate-region`, iso=`alpha-3`, code=`country-code`)], iso)

pop.dt <- readRDS("../covidm/data/wpp2019_pop2020.rds")[country_code %in% afrkey$code]
tots.dt <- pop.dt[afrkey, on=.(country_code = code)][, .(tot = sum(f+m)*1000), by=iso]

afrkey[tots.dt, on=.(iso=iso), tot := tot ]

saveRDS(afrkey, tail(.args, 1))
Binary file added mapping/regionref.rds
Binary file not shown.
11 changes: 11 additions & 0 deletions mapping/scale_axes_input_bounds.R
Original file line number Diff line number Diff line change
@@ -0,0 +1,11 @@
scale_x_input_bounds <- function(data){
bounds <- st_bbox(data)
axes_scale <- ggplot2::xlim(bounds$xmin, bounds$xmax)
return(axes_scale)
}

scale_y_input_bounds <- function(data){
bounds <- st_bbox(data)
axes_scale <- ggplot2::ylim(bounds$ymin, bounds$ymax)
return(axes_scale)
}
22 changes: 22 additions & 0 deletions mapping/scale_fill_attribute.R
Original file line number Diff line number Diff line change
@@ -0,0 +1,22 @@
source('mapping/define_fill_continuous_sequential.R', local = T)

scale_fill_attribute <- function(data, attribute){
#Handle picking predefined colour scales for a given attribute.
#data - data.frame, dataset to define scale attributes
#attribute - character, column to be symbolised by the colour scale
#To add a new palette mapping, add a named item to palette_reference.

palette_reference <- list()
palette_reference[['dummy']] = 'Mint'

if(!attribute %in% names(palette_reference)){
stop(paste0('Unknown attribute scale: ', attribute, '. See scale_fill_attribute to add a new attribute scale.'))
}

pal <- palette_reference[[attribute]]

s <- define_fill_continuous_sequential(data = data, attribute = attribute, palette = pal)

return(s)

}
21 changes: 21 additions & 0 deletions mapping/validate_geometry.R
Original file line number Diff line number Diff line change
@@ -0,0 +1,21 @@
#validate the geometry of an sf object and test for invalid geometries
validate_geometry <- function(data){

data <- data %>%
dplyr::mutate(valid_geom = st_is_valid(geometry, reason = T)) %>%
dplyr::mutate(geometry = st_make_valid(geometry)) %>%
dplyr::mutate(valid_geom = st_is_valid(geometry, reason = T))

testthat::test_that('all la geometry is valid.', {

vg <- data %>% pull(valid_geom) %>% unique()

testthat::expect_equal(length(vg), 1)

testthat::expect_equal(vg, 'Valid Geometry')

})

return(data)

}