geodemographic.R

library(tidycensus)
library(tigris)
library(tidyverse)
library(sf)
options(tigris_use_cache = TRUE)
set.seed(123456)

austin_counties <- counties(year = 2021) %>%
  filter(CBSAFP == "12420") %>%
  pull(COUNTYFP)

austin_inputs <- get_acs(
  geography = "tract",
  variables = c(median_income = "B19013_001",
                median_age = "B01002_001",
                total_population = "B01003_001"),
  state = "TX",
  county = austin_counties,
  output = "wide",
  geometry = TRUE,
  keep_geo_vars = TRUE
) %>%
  mutate(pop_density = total_populationE / (ALAND / 2589988.11 )) %>%
  na.omit()

austin_kmeans <- austin_inputs %>%
  st_drop_geometry() %>%
  select(median_incomeE, median_ageE, pop_density) %>%
  scale() %>%
  kmeans(centers = 6)

austin_clusters <- austin_inputs %>%
  mutate(cluster = as.character(austin_kmeans$cluster))

library(mapview)

mapview(austin_clusters, zcol = "cluster",
        col.regions = RColorBrewer::brewer.pal(6, "Set1"))

ggplot(austin_clusters, aes(fill = cluster)) + 
  geom_sf(size = 0.1) + 
  scale_fill_brewer(palette = "Set1") + 
  theme_void() + 
  labs(fill = "Cluster ")


ggplot(austin_clusters, 
       aes(x = median_incomeE, y = pop_density, color = cluster)) + 
  geom_point() + 
  scale_color_brewer(palette = "Set1") + 
  scale_y_log10() + 
  scale_x_continuous(labels = scales::dollar_format()) + 
  theme_minimal(base_size = 16) + 
  labs(color = "Cluster",
       x = "Median household income",
       y = "Population density (logged)")


library(ggiraph)
library(patchwork)

austin_map <- ggplot(austin_clusters, aes(fill = cluster, data_id = GEOID)) + 
  geom_sf_interactive(size = 0.1) + 
  scale_fill_brewer(palette = "Set1") + 
  theme_void() + 
  labs(fill = "Cluster ")


austin_plot <- ggplot(austin_clusters, 
       aes(x = median_incomeE, y = pop_density, color = cluster, data_id = GEOID)) + 
  geom_point_interactive() + 
  scale_color_brewer(palette = "Set1") + 
  scale_y_log10() + 
  scale_x_continuous(labels = scales::dollar_format()) + 
  theme_minimal(base_size = 12) + 
  labs(color = "Cluster",
       x = "Median household income",
       y = "Population density (logged)")


girafe(ggobj = austin_map + austin_plot, width_svg = 10, height_svg = 5.5)  %>%
  girafe_options(opts_hover(css = "fill:cyan;"), 
                 opts_zoom(min = 1, max = 8))