############################################################# STEP 0: Save your API key
# The script will not work unless you first paste your
# API key into the API key setup code, shown above, and
# run the code. The code needs to run only once per computer.
############################################################
############################################################# STEP 1: Load Packages
############################################################
if (!require("tidycensus")) install.packages("tidycensus")
## Loading required package: tidycensus
if (!require("tidyverse")) install.packages("tidyverse")
## Loading required package: tidyverse
## ── Attaching core tidyverse packages ──────────────────────── tidyverse 2.0.0 ──
## ✔ dplyr 1.2.1 ✔ readr 2.2.0
## ✔ forcats 1.0.1 ✔ stringr 1.6.0
## ✔ ggplot2 4.0.3 ✔ tibble 3.3.1
## ✔ lubridate 1.9.5 ✔ tidyr 1.3.2
## ✔ purrr 1.2.2
## ── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
## ✖ dplyr::filter() masks stats::filter()
## ✖ dplyr::lag() masks stats::lag()
## ℹ Use the conflicted package (<http://conflicted.r-lib.org/>) to force all conflicts to become errors
if (!require("kableExtra")) install.packages("kableExtra")
## Loading required package: kableExtra
##
## Attaching package: 'kableExtra'
##
## The following object is masked from 'package:dplyr':
##
## group_rows
if (!require("sf")) install.packages("sf")
## Loading required package: sf
## Linking to GEOS 3.14.1, GDAL 3.12.1, PROJ 9.7.1; sf_use_s2() is TRUE
if (!require("leaflet")) install.packages("leaflet")
## Loading required package: leaflet
if (!require("htmlwidgets")) install.packages("htmlwidgets")
## Loading required package: htmlwidgets
if (!require("plotly")) install.packages("plotly")
## Loading required package: plotly
##
## Attaching package: 'plotly'
##
## The following object is masked from 'package:ggplot2':
##
## last_plot
##
## The following object is masked from 'package:stats':
##
## filter
##
## The following object is masked from 'package:graphics':
##
## layout
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\nThen run the one-time setup code at the",
"\ntop of this script to store your key.",
"\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
)
}
## Getting data from the 2020-2024 5-year ACS
## Downloading feature geometry from the Census website. To cache shapefiles for use in future sessions, set `options(tigris_use_cache = TRUE)`.
## Using the ACS Data Profile
############################################################
# 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
)
}
## Getting data from the 2020-2024 5-year ACS
## Downloading feature geometry from the Census website. To cache shapefiles for use in future sessions, set `options(tigris_use_cache = TRUE)`.
## Using the ACS Data Profile
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
############################################################
AreaMap # Name of the map
DotPlot # Name of the dot plot with error bars