3D-Grafik der geschlechtsspezifischen Lohnunterschiede in Österreich (https://twitter.com/matschnetzer/status/1603058949191778304)

gpghex.R

librarian::shelf(tidyverse,pdftools,sf,msthemes,rayshader,wesanderson)

# Download PDF from Statistik Austria: https://www.statistik.at/fileadmin/publications/Lohnsteuer-2021.pdf
rawpdf <- pdf_text("Lohnsteuer-2021.pdf")

# Manual editing to obtain the required data from the two-column table
tmppdf <- rawpdf[67:68] |> strsplit('\n') |> unlist()
df <- tmppdf[c(11:55,67:86)] |> as.data.frame() |> setNames("string") |> 
  separate(string, into = c("gem","male","female", "gem2", "male2", "female2"), sep = "\\s\\s+")

# Manual editing of district names to correspond with map
findf <- bind_rows(df |> select(gem,male,female),
                   df |> select(gem = gem2, male = male2, female = female2)) |> 
  mutate(across(c(male,female), ~as.numeric(str_remove(., " "))),
         gpg = (1-female/male)*100,
         gem = str_replace_all(gem, c(" \\("="\\(",
                                      "Wr."="Wiener",
                               "Klagenfurt\\(Stadt\\)"="Klagenfurt Stadt",
                               "Villach\\(Stadt\\)"="Villach Stadt",
                               "St. Veit an der Glan"="Sankt Veit an der Glan",
                               "Krems/Donau\\(Stadt\\)"="Krems an der Donau\\(Stadt\\)",
                               "St. Pölten\\(Stadt\\)"="Sankt Pölten\\(Stadt\\)",
                               "Waidhofen/Ybbs\\(Stadt\\)"="Waidhofen an der Ybbs\\(Stadt\\)",
                               "St. Pölten\\(Land\\)"="Sankt Pölten\\(Land\\)",
                               "Waidhofen/Thaya"="Waidhofen an der Thaya",
                               "Linz\\(Stadt\\)"="Stadt Linz",
                               "Steyr\\(Stadt\\)"="Stadt Steyr",
                               "Wels\\(Stadt\\)"="Stadt Wels",
                               "Braunau am Inn"="Braunau",
                               "Kirchdorf/Krems"="Kirchdorf",
                               "Ried im Innkreis"="Ried",
                               "St. Johann/Pongau"="Sankt Johann im Pongau",
                               "Innsbruck\\(Stadt\\)"="Innsbruck-Stadt",
                               "Innsbruck\\(Land\\)"="Innsbruck-Land")))

# Load map and merge Vienna's districts
map <- st_read("https://raw.githubusercontent.com/ginseng666/GeoJSON-TopoJSON-Austria/master/2021/simplified-99.9/bezirke_999_geo.json") |> 
  mutate(iso = as.numeric(iso),
         iso = ifelse(iso %in% 901:923, 900, iso),
         name = ifelse(iso == 900, "Wien", name)) |> 
  group_by(iso,name) |>  
  summarise(across(geometry, ~ sf::st_union(.)), .groups = "keep") |> 
  summarise(across(geometry, ~ sf::st_combine(.)))

# Merge map and data
shape <- map |> left_join(findf, by = c("name"="gem")) |> st_transform(3857) |> 
  mutate(index_target = 1:n())
target <- st_geometry(shape)

# Create grid with hexagons
grid <- st_make_grid(target,
                     cellsize = 15 * 1000,
                     crs = st_crs(shape),
                     what = "polygons",
                     square = FALSE)
grid <- st_sf(index = 1:length(lengths(grid)), grid)
cent_grid <- st_centroid(grid)
cent_merge <- st_join(cent_grid, shape, left = FALSE)
grid_new <- inner_join(grid, st_drop_geometry(cent_merge))

# Create district borders
bezirke <- grid_new |>  group_by(iso) |> 
  summarise(across(grid, ~ sf::st_union(.)), .groups = "keep") |> 
  summarise(across(grid, ~ sf::st_combine(.)))

# Color palette
pal <- wes_palette("Zissou1", 100, type = "continuous")

# 2D plot
bezplotflat <- ggplot() +
  geom_sf(data = grid_new, size = 0.1, aes(fill = gpg)) +
  geom_sf(data = bezirke, size = 0.3, color = "black", fill = NA) +
  scale_fill_gradientn(colours = pal, name = "", labels = scales::percent_format(scale = 1), 
                       guide = guide_colorbar(barheight = 10, barwidth = 0.3)) +
  coord_sf(expand = FALSE) +
  labs(title = "Geschlechtsspezifische Lohnunterschiede",
       subtitle = "Unterschied des Medianlohns von ganzjährig Vollzeit beschäftigten Frauen und Männern, 2021",
       caption = "Daten: Lohnsteuerstatistik, Statistik Austria. Grafik: @matschnetzer") +
  theme_ms(grid = F, alttf = T) +
  theme(legend.position = "right",
        plot.subtitle = element_text(margin = margin(b=5, t=5, unit="pt"), size = 9, family = "Roboto"),
        plot.caption = element_text(margin = margin(t=5, unit="pt"), size = 7, family = "Roboto"),
        axis.text = element_blank())

ggsave(bezplotflat, file = "gpghexflat.png", width = 7.5, height = 4, dpi = 320)

# 3D animated plot
bezplot3d <- ggplot() +
  geom_sf(data = grid_new, linewidth = 0.3, aes(fill = gpg), color="transparent") +
  scale_fill_gradientn(colours = pal, name = "", guide = guide_colorbar(barheight = 10, barwidth = 0.3)) +
  coord_sf(expand = FALSE) +
  theme_ms(grid = F, alttf = T) +
  theme(legend.position = "none",
        axis.text = element_blank())

plot_gg(bezplot3d, width = 7, height = 4, scale = 300, multicore = TRUE,
        windowsize = c(1200, 700), zoom = 0.45, phi = 50)
render_snapshot()
render_movie(filename = "gpghex.mp4")