#30DayMapChallenge 2021 - Day 09 - Monochrome

bw_quarters.R

# Split NZ into halves or quarters of population
# Using SA2 populations (not as precice as meshblocks, approx enough)
# For #30DayMapChallenge 2021 - Day 09 - Monochrome
#  -- David Friggens, November 2021

# Setup -------------------------------------------------------------------

library(readr)
library(dplyr)
library(sf)
library(rmapshaper)
library(ggplot2)
library(patchwork)

PLOT_FONT <- "Poppins"

COLOUR_MONO1 <- "#333333"
COLOUR_MONO2 <- "#dddddd"
COLOUR_LINE <- COLOUR_MONO1
COLOUR_BG <- "#aaaaaa"

pop_sa2 <-
  read_tsv("population/population_sa2.tsv",
           col_types = "cccccciin") %>% 
  filter(YEAR == 2021) %>% 
  select(area_code = AREA,
         population = `Value  Flags`)

geo_sa2 <- 
  st_read("statsnz/statistical-area-2-2021-centroid-inside.gpkg") %>% 
  st_drop_geometry() %>% 
  # polygons already in NZTM 2193 so will likely only be using e/n, not lat/lon
  select(area_code = SA22021_V1_00,
         easting = EASTING,
         northing = NORTHING,
         latitude = LATITUDE,
         longitude = LONGITUDE)

geo_nz <- 
  st_read("statsnz/territorial-authority-2021-clipped-generalised.gpkg") %>% 
  filter(! TA2021_V1_00 %in% c("067", "999")) %>%  # ignore Chatham Islands thanks
  count() %>% 
  ms_simplify(keep = 0.02)


nz_bbox <- geo_nz %>% st_bbox()
nz_west <- nz_bbox[1]
nz_south <- nz_bbox[2]
nz_east <- nz_bbox[3]
nz_north <- nz_bbox[4]


# Half N-S ----------------------------------------------------------------

midway_ns <- 
  geo_sa2 %>% 
  inner_join(pop_sa2, by = "area_code") %>% 
  summarise(mid_north = round(sum(northing * population) / sum(population), 0)) %>% 
  pull(mid_north)

# geo_north <- 
#   geo_nz %>% 
#   ms_clip(bbox = c(nz_west, midway_ns, nz_east, nz_north))
# geo_south <- 
#   geo_nz %>% 
#   ms_clip(bbox = c(nz_west, nz_south, nz_east, midway_ns))
# 
# ggplot() +
#   geom_sf(data = geo_north,
#           fill = "blue") +
#   geom_sf(data = geo_south,
#           fill = "red")

# Half W-E ----------------------------------------------------------------

midway_we <- 
  geo_sa2 %>% 
  inner_join(pop_sa2, by = "area_code") %>% 
  summarise(mid_east = round(sum(easting * population) / sum(population), 0)) %>% 
  pull(mid_east)

# geo_west <- 
#   geo_nz %>% 
#   ms_clip(bbox = c(nz_west, nz_south, midway_we, nz_north))
# geo_east <- 
#   geo_nz %>% 
#   ms_clip(bbox = c(midway_we, nz_south, nz_east, nz_north))
# 
# ggplot() +
#   geom_sf(data = geo_west,
#           fill = "blue") +
#   geom_sf(data = geo_east,
#           fill = "red")


# Quarter N-S -------------------------------------------------------------

midway_nn <- 
  geo_sa2 %>% 
  filter(northing >= midway_ns) %>% 
  inner_join(pop_sa2, by = "area_code") %>% 
  summarise(mid_north = round(sum(northing * population) / sum(population), 0)) %>% 
  pull(mid_north)
midway_ss <- 
  geo_sa2 %>% 
  filter(northing <= midway_ns) %>% 
  inner_join(pop_sa2, by = "area_code") %>% 
  summarise(mid_north = round(sum(northing * population) / sum(population), 0)) %>% 
  pull(mid_north)

geo_north_north <- 
  geo_nz %>% 
  ms_clip(bbox = c(nz_west, midway_nn, nz_east, nz_north))
geo_mid_north <- 
  geo_nz %>% 
  ms_clip(bbox = c(nz_west, midway_ns, nz_east, midway_nn))
geo_mid_south <- 
  geo_nz %>% 
  ms_clip(bbox = c(nz_west, midway_ss, nz_east, midway_ns))
geo_south_south <- 
  geo_nz %>% 
  ms_clip(bbox = c(nz_west, nz_south, nz_east, midway_ss))

gg_ns4 <- 
  ggplot() +
  geom_sf(data = geo_north_north,
          size = 0.5, colour = COLOUR_LINE,
          fill = COLOUR_MONO1) +
  geom_sf(data = geo_mid_north,
          size = 0.5, colour = COLOUR_LINE,
          fill = COLOUR_MONO2) +
  geom_sf(data = geo_mid_south,
          size = 0.5, colour = COLOUR_LINE,
          fill = COLOUR_MONO1) +
  geom_sf(data = geo_south_south,
          size = 0.5, colour = COLOUR_LINE,
          fill = COLOUR_MONO2) +
  coord_sf(datum = NULL) +
  theme_void() 


# Quarter W-E -------------------------------------------------------------

midway_ee <- 
  geo_sa2 %>% 
  filter(easting >= midway_we) %>% 
  inner_join(pop_sa2, by = "area_code") %>% 
  summarise(mid_east = round(sum(easting * population) / sum(population), 0)) %>% 
  pull(mid_east)
midway_ww <- 
  geo_sa2 %>% 
  filter(easting <= midway_we) %>% 
  inner_join(pop_sa2, by = "area_code") %>% 
  summarise(mid_east = round(sum(easting * population) / sum(population), 0)) %>% 
  pull(mid_east)

geo_east_east <- 
  geo_nz %>% 
  ms_clip(bbox = c(midway_ee, nz_south, nz_east, nz_north))
geo_mid_east <- 
  geo_nz %>% 
  ms_clip(bbox = c(midway_we, nz_south, midway_ee, nz_north))
geo_mid_west <- 
  geo_nz %>% 
  ms_clip(bbox = c(midway_ww, nz_south, midway_we, nz_north))
geo_west_west <- 
  geo_nz %>% 
  ms_clip(bbox = c(nz_west, nz_south, midway_ww, nz_north))

gg_we4 <- 
  ggplot() +
  geom_sf(data = geo_east_east,
          size = 0.5, colour = COLOUR_LINE,
          fill = COLOUR_MONO1) +
  geom_sf(data = geo_mid_east,
          size = 0.5, colour = COLOUR_LINE,
          fill = COLOUR_MONO2) +
  geom_sf(data = geo_mid_west,
          size = 0.5, colour = COLOUR_LINE,
          fill = COLOUR_MONO1) +
  geom_sf(data = geo_west_west,
          size = 0.5, colour = COLOUR_LINE,
          fill = COLOUR_MONO2) +
  coord_sf(datum = NULL) +
  theme_void()

# Join --------------------------------------------------------------------

gg_ns4 + gg_we4 +
  plot_annotation(title = "New Zealand Population in Quarters",
                  subtitle = "#30DayMapChallenge 2021 - Day 09 - Monochrome",
                  caption = "David Friggens @dakvid, Data: Stats NZ",
                  theme = theme(plot.title = element_text(size = 72, family = PLOT_FONT, colour = COLOUR_LINE),
                                plot.subtitle = element_text(size = 28, family = PLOT_FONT, colour = COLOUR_LINE),
                                plot.caption = element_text(size = 24, family = PLOT_FONT, colour = COLOUR_LINE))) &
  theme(plot.background = element_rect(fill = COLOUR_BG, colour = COLOUR_BG),
        panel.background = element_rect(fill = COLOUR_BG, colour = COLOUR_BG))

ggsave("Day_09/Day_09_NZ_population_quarters.png",
       width = 24, height = 20)