This map shows the percentage of the county population who speak another language other than English in their homes. District 28, in Davidson County, leads all districts with 41.4%.
# ----------------------------------------------------------
# ACS LANGUAGE AT HOME ANALYSIS (2019 vs 2024)
# FINAL FINAL SCRIPT (ZERO WARNINGS GUARANTEED)
# ----------------------------------------------------------
# ---------------------------
# Install packages if needed
# ---------------------------
packages <- c("tidyverse", "tidycensus", "sf", "leaflet", "htmlwidgets", "plotly")
installed <- rownames(installed.packages())
for (p in packages) {
if (!(p %in% installed)) install.packages(p)
}
# ---------------------------
# Load libraries
# ---------------------------
library(tidyverse)
library(tidycensus)
library(sf)
library(leaflet)
library(htmlwidgets)
library(plotly)
# ---------------------------
# Census API key
# ---------------------------
census_api_key("904ea9b9c6b5f8f7d6b1d4846a59ee8cc3b48d1f", install = FALSE)
# ---------------------------
# Define years + variable
# ---------------------------
year_old <- 2019
year_new <- 2024
variable <- "DP02_0115PE"
target_label <- "Percent speaking language other than English at home"
# ---------------------------
# Get 2024 data
# ---------------------------
mydata_new <- get_acs(
geography = "county subdivision",
state = "TN",
variables = variable,
year = year_new,
survey = "acs5",
geometry = TRUE
) %>%
rename(Area = NAME) %>%
transmute(
GEOID,
Area,
geometry,
Estimate_new = estimate,
MoE_new = moe
)
# ---------------------------
# Get 2019 data
# ---------------------------
mydata_old <- get_acs(
geography = "county subdivision",
state = "TN",
variables = variable,
year = year_old,
survey = "acs5",
geometry = FALSE
) %>%
rename(Area = NAME) %>%
transmute(
GEOID,
Area,
Estimate_old = estimate,
MoE_old = moe
)
# ---------------------------
# Join datasets
# ---------------------------
mydata_joined <- mydata_new %>%
left_join(mydata_old, by = c("GEOID", "Area"))
# ---------------------------
# Filter counties
# ---------------------------
search_terms <- c(
"Davidson",
"Rutherford",
"Williamson",
"Cheatham",
"Robertson",
"Sumner",
"Wilson"
)
pattern <- paste(search_terms, collapse = "|")
filtereddata <- mydata_joined %>%
filter(str_detect(Area, regex(pattern, ignore_case = TRUE)))
# ---------------------------
# Compute change (NO moe_sum)
# ---------------------------
filtereddata <- filtereddata %>%
mutate(
Diff = Estimate_new - Estimate_old,
# SAFE MOE formula (no warnings ever)
Diff_MOE = sqrt((MoE_new^2) + (MoE_old^2)),
# Safe significance check
Sig_90 = ifelse(
is.na(MoE_new) | is.na(MoE_old),
FALSE,
abs(Diff) > Diff_MOE
),
SigDirection = case_when(
Sig_90 & Diff > 0 ~ "Significant increase",
Sig_90 & Diff < 0 ~ "Significant decrease",
TRUE ~ "No statistically significant change"
)
)
# ---------------------------
# Prepare map data
# ---------------------------
mapdata <- filtereddata %>%
st_as_sf() %>%
st_transform(4326)
# ---------------------------
# Create graph
# ---------------------------
plotdf <- filtereddata %>%
st_drop_geometry() %>%
mutate(
y_ordered = reorder(Area, Diff),
sig_class = case_when(
SigDirection == "Significant increase" ~ "Increase (significant)",
SigDirection == "Significant decrease" ~ "Decrease (significant)",
TRUE ~ "Not statistically significant"
)
)
mygraph <- plot_ly(
data = plotdf,
x = ~Diff,
y = ~as.character(y_ordered),
type = "scatter",
mode = "markers",
color = ~sig_class
)
mygraph
# ---------------------------
# Create map
# ---------------------------
mapdata$popup <- paste0(
"<strong>", mapdata$Area, "</strong><br/>",
"2019: ", round(mapdata$Estimate_old, 2), "%<br/>",
"2024: ", round(mapdata$Estimate_new, 2), "%<br/>",
"<b>Change:</b> ", round(mapdata$Diff, 2), "%"
)
max_abs <- max(abs(mapdata$Diff), na.rm = TRUE)
if (!is.finite(max_abs)) max_abs <- 0
bins <- pretty(c(-max_abs, max_abs), n = 6)
pal <- colorBin(
palette = "RdBu",
domain = mapdata$Diff,
bins = bins,
reverse = TRUE
)
DivisionMap <- leaflet(mapdata) %>%
addProviderTiles(providers$CartoDB.Positron) %>%
addPolygons(
fillColor = ~pal(Diff),
fillOpacity = 0.6,
color = "black",
weight = 1,
popup = ~popup
)
DivisionMap
# ---------------------------
# Save outputs
# ---------------------------
saveWidget(as_widget(mygraph), "ACS_Language_Graph.html", selfcontained = TRUE)
saveWidget(DivisionMap, "ACS_Language_Map.html", selfcontained = TRUE)
```{# ———————————————————-} # GRAPH SCRIPT (CORRECT — ONE ROW PER SUBDIVISION) # ———————————————————-
library(tidyverse) library(sf) library(ggplot2)
plotdf <- filtereddata %>% st_drop_geometry() %>%
# 🔥 THIS IS THE FIX: use GEOID (true unique ID) group_by(GEOID, Area) %>% summarise( Diff = first(Diff), SigDirection = first(SigDirection), .groups = “drop” ) %>%
# sort like example arrange(desc(Diff)) %>%
# keep visible number (like example) slice(1:30) %>%
mutate( sig_class = ifelse( SigDirection == “Significant increase”, “Increase (significant)”, “Not statistically significant” ) )
plotdf\(Area <- factor(plotdf\)Area, levels = rev(plotdf$Area))
mygraph <- ggplot(plotdf, aes(x = Diff, y = Area)) + geom_vline(xintercept = 0, color = “black”, linewidth = 0.8) + geom_point(aes(color = sig_class), size = 3) + scale_color_manual(values = c( “Increase (significant)” = “#b2182b”, “Not statistically significant” = “#bdbdbd” )) + labs( title = “Estimated Change in Pct. speaking non-English at home (2019 to 2024)”, subtitle = “County subdivisions. Color indicates statistical significance (90%).”, x = “Estimated change”, y = ““, color =”” ) + theme_minimal() + theme( axis.text.y = element_text(size = 10), panel.grid.major.y = element_blank() )
mygraph ```