Below is an interactive dot plot and map showing how Davidson County and bordering counties compare on median household income.
Below is the code used to produce the graphic and map.
############################################################
# STEP 1: Load Packages
############################################################
if (!require("tidycensus")) install.packages("tidycensus")
if (!require("tidyverse")) install.packages("tidyverse")
if (!require("kableExtra")) install.packages("kableExtra")
if (!require("sf")) install.packages("sf")
if (!require("leaflet")) install.packages("leaflet")
if (!require("htmlwidgets")) install.packages("htmlwidgets")
if (!require("plotly")) install.packages("plotly")
library(tidycensus)
library(tidyverse)
library(kableExtra)
library(sf)
library(leaflet)
library(htmlwidgets)
library(plotly)
############################################################
# STEP 2: Verify Census API Key
############################################################
if (!nzchar(Sys.getenv("CENSUS_API_KEY"))) {
stop(
paste(
"\nNo Census API key was found.",
"\n\nTo obtain a free Census API key, visit:",
"\nhttps://api.census.gov/data/key_signup.html",
"\n\nAfter obtaining your key, run the",
"\none-time setup code shown above.",
"\n"
)
)
}
############################################################
# STEP 3: Configure Your Analysis
############################################################
# ACS year
acs_year <- 2024
# State abbreviation
state_abbrev <- "TN"
# ACS variable
#
# Examples:
#
# DP05_0001 = Total population
# DP03_0062 = Median Household Income
# DP03_0128P = Poverty Rate
# DP02_0068P = Bachelor's Degree or Higher Rate
# DP05_0018 = Median Age
# DP04_0046P = Homeownership Rate
# DP04_0134 = Median Gross Rent
# DP04_0142P = Rent Burden Rate
# DP05_0090P = Hispanic or Latino Population Percentage
# DP02_0094P = Foreign-Born Percentage
# DP03_0009P = Unemployment Rate
# DP03_0025 = Average Commute Time
my_variable <- "DP03_0062"
# Geography level
#
# Recommended ACS geography options:
#
# "county"
# "county subdivision"
# "tract"
geography_level <- "county"
# County filter
#
# Use NULL to include all counties in the state.
#
# Examples:
#
# selected_counties <- NULL
# selected_counties <- "Rutherford"
# selected_counties <- c(
# "Davidson",
# "Rutherford",
# "Williamson",
# "Cheatham",
# "Robertson",
# "Sumner",
# "Wilson"
# )
selected_counties <- c(
"Davidson",
"Rutherford",
"Williamson",
"Cheatham",
"Robertson",
"Sumner",
"Wilson"
)
############################################################
# STEP 4: Explore ACS Codebooks (Optional)
############################################################
ProfileTables <- load_variables(
acs_year,
"acs5/profile"
)
# DetailedTables <- load_variables(
# acs_year,
# "acs5"
# )
# SubjectTables <- load_variables(
# acs_year,
# "acs5/subject"
# )
############################################################
# STEP 5: Download ACS Data
############################################################
if (is.null(selected_counties)) {
mydata <- get_acs(
geography = geography_level,
state = state_abbrev,
year = acs_year,
survey = "acs5",
variables = my_variable,
geometry = TRUE
)
} else {
mydata <- get_acs(
geography = geography_level,
state = state_abbrev,
county = selected_counties,
year = acs_year,
survey = "acs5",
variables = my_variable,
geometry = TRUE
)
}
############################################################
# STEP 6: Keep Needed Variables and Geometry
############################################################
mydata <- mydata |>
select(
Area = NAME,
Variable = variable,
Estimate = estimate,
Margin_of_Error = moe,
geometry
)
############################################################
# STEP 7: Transform Coordinates for Leaflet
############################################################
mydata <- st_transform(
mydata,
crs = 4326
)
############################################################
# STEP 8: Download County Boundaries (Optional Overlay)
############################################################
if (is.null(selected_counties)) {
county_boundaries <- get_acs(
geography = "county",
state = state_abbrev,
year = acs_year,
survey = "acs5",
variables = my_variable,
geometry = TRUE
)
} else {
county_boundaries <- get_acs(
geography = "county",
state = state_abbrev,
county = selected_counties,
year = acs_year,
survey = "acs5",
variables = my_variable,
geometry = TRUE
)
}
county_boundaries <- st_transform(
county_boundaries,
4326
)
############################################################
# STEP 9: Sort Results from Highest to Lowest
############################################################
mydata <- mydata |>
arrange(desc(Estimate))
############################################################
# STEP 10: Display Results as a Table
############################################################
ResultsTable <- mydata |>
st_drop_geometry() |>
kbl(
col.names = c(
"Area",
"Variable",
"Estimate",
"Margin of Error"
),
format.args = list(big.mark = ",")
) |>
kable_styling(
bootstrap_options = c("striped", "hover"),
full_width = FALSE
)
############################################################
# STEP 11: Create Pop-Up Content
############################################################
mydata$popup <- paste0(
mydata$Area,
"\n\nEstimate: ",
format(mydata$Estimate, big.mark = ","),
"\nPlus/Minus: ",
format(mydata$Margin_of_Error, big.mark = ",")
)
############################################################
# STEP 12: Build Color Palette
############################################################
qs <- unique(
quantile(
mydata$Estimate,
probs = seq(0, 1, length.out = 6),
na.rm = TRUE
)
)
pal <- colorBin(
palette = "Blues",
domain = mydata$Estimate,
bins = qs,
pretty = FALSE
)
############################################################
# STEP 13: Create Choropleth Map
############################################################
AreaMap <- leaflet(mydata) |>
addProviderTiles(
providers$Esri.WorldTopoMap
) |>
addPolygons(
fillColor = ~pal(Estimate),
fillOpacity = 0.5,
color = "black",
weight = 1,
popup = ~popup
)
# Optional county boundary overlay for
# sub-county geographies
if (geography_level != "county") {
AreaMap <- AreaMap |>
addPolygons(
data = county_boundaries,
fill = FALSE,
color = "black",
weight = 2,
opacity = 1,
group = "County Boundaries"
)
}
AreaMap <- AreaMap |>
addLegend(
pal = pal,
values = ~Estimate,
title = "Estimate",
labFormat = labelFormat(
big.mark = ","
)
)
############################################################
# STEP 14: Create Plotly Dot Plot
############################################################
plotdf <- mydata |>
st_drop_geometry() |>
mutate(
point_color = pal(Estimate),
y_ordered = reorder(
Area,
Estimate
),
hover_text = paste0(
Area,
"\nEstimate: ",
format(Estimate, big.mark = ","),
"\nPlus/Minus: ",
format(Margin_of_Error, big.mark = ",")
)
)
DotPlot <- plot_ly(
data = plotdf,
x = ~Estimate,
y = ~as.character(y_ordered),
type = "scatter",
mode = "markers",
showlegend = FALSE,
marker = list(
color = ~point_color,
size = 8,
line = list(
color = "rgba(120,120,120,0.9)",
width = 0.5
)
),
error_x = list(
type = "data",
array = ~Margin_of_Error,
arrayminus = ~Margin_of_Error,
color = "rgba(0,0,0,0.65)",
thickness = 1
),
text = ~hover_text,
hovertemplate = "%{text}"
) |>
layout(
title = list(
text = paste0(
"ACS Estimates by Geography",
"\nError bars show ACS margins of error."
)
),
xaxis = list(
title = "Estimate",
tickformat = ",.0f",
automargin = TRUE
),
yaxis = list(
title = "",
automargin = TRUE,
categoryorder = "array",
categoryarray = levels(plotdf$y_ordered)
),
margin = list(
l = 200,
r = 20,
b = 60,
t = 60,
pad = 2
)
)
############################################################
# STEP 15: Export Dot Plot (Optional)
############################################################
saveWidget(
widget = as_widget(DotPlot),
file = "ACSGraph.html",
selfcontained = TRUE
)
############################################################
# STEP 16: Export Map (Optional)
############################################################
saveWidget(
widget = AreaMap,
file = "ACSMap.html",
selfcontained = TRUE
)
############################################################
# STEP 17: Show output
############################################################
DotPlot # The dot plot with error bars
AreaMap # The map