library(tidycensus)
library(tidyverse)
library(janitor)
library(stringr)
library(openxlsx)
library(tigris)
library(sf)
library(ggplot2)
library(scales)
options(tigris_use_cache = TRUE)
# ---------------------------------------------------------------------------------
# parameters and file paths
# ---------------------------------------------------------------------------------
acs_year <- 2024
acs_survey <- "acs5"
pov_min <- 131
pov_max <- 165
data_dir <- "data"
output_dir <- "output"
dir.create(output_dir, showWarnings = FALSE, recursive = TRUE)
# Geocorr files
lower_crosswalk_path <- file.path(data_dir, "geocorr2022_2613906081.csv")
upper_crosswalk_path <- file.path(data_dir, "geocorr2022_2613903963.csv")
# Output files
excel_output_path <- file.path(
output_dir,
"tx_snap_lege_districts.xlsx"
)
lower_map_path <- file.path(
output_dir,
"tx_snap_lower.png"
)
upper_map_path <- file.path(
output_dir,
"tx_snap_upper.png"
)VHK grant data - Children 0-3 on SNAP, 131-165% of Poverty
READ ME
This analysis uses ACS PUMS microdata from tidycensus to estimate Texas households participating in SNAP with children and income between 131% and 165% of poverty, used here as an ACS proxy for the SNAP/FPL policy band. The ACS PUMS variable POVPIP measures income as a percentage of the Census Bureau poverty threshold. This is similar to, but not identical to, the federal poverty guidelines used for administrative program eligibility. Because Texas SNAP eligibility is commonly discussed using FPL thresholds, this analysis sometimes uses “131–165% FPL” as shorthand while interpreting estimates as ACS poverty-threshold-based approximations.
Because ACS PUMS does not directly identify Texas House or Senate districts, estimates are first calculated at the PUMA level and then allocated to state legislative districts using Geocorr 2022 population-weighted PUMA-to-district crosswalks. Results are produced separately for Texas House districts (lower) and Texas Senate districts (upper).
The main outcome, pct_snap_hh_0_3_131_165, represents the percent of SNAP households in each district that have at least one child age 0–3 and household income between 131% and 165% of poverty. Count columns are weighted ACS PUMS estimates rounded to whole numbers. Percent columns are shown as percentages. Because PUMAs and legislative districts do not perfectly align, district-level estimates should be interpreted as approximations.
You can view the 2024 ACS PUMS DATA DICTIONARY here.
Setup
# ---------------------------------------------------------------------------------
# pull Texas ACS PUMS microdata
# ---------------------------------------------------------------------------------
# Variables:
# SERIALNO = household ID
# SPORDER = person number within household
# PUMA = Public Use Microdata Area
# AGEP = age
# POVPIP = income-to-poverty ratio, as percent of Census poverty threshold
# Used here as an ACS proxy for the SNAP/FPL policy band.
# FS = household received Food Stamps/SNAP in past 12 months
# PWGTP = person weight
# WGTP = household weight
vars <- c(
"SERIALNO",
"SPORDER",
"PUMA",
"AGEP",
"POVPIP",
"FS",
"PWGTP",
"WGTP"
)
tx_pums <- get_pums(
variables = vars,
state = "TX",
survey = acs_survey,
year = acs_year,
recode = FALSE
)Cleaning & Analysis
# ---------------------------------------------------------------------------------
# clean PUMS variables
# ---------------------------------------------------------------------------------
tx_clean <- tx_pums %>%
clean_names() %>%
mutate(
puma = str_pad(as.character(puma), width = 5, pad = "0"),
age = as.numeric(agep),
povpip = as.numeric(povpip),
child_0_3 = age >= 0 & age <= 3,
child_u18 = age < 18,
is_snap_hh = as.character(fs) == "1",
income_131_165 = povpip >= pov_min & povpip <= pov_max
)
# ---------------------------------------------------------------------------------
# household-level indicators
# ---------------------------------------------------------------------------------
hh_level <- tx_clean %>%
group_by(serialno, puma) %>%
summarise(
has_0_3 = any(child_0_3, na.rm = TRUE),
has_u18 = any(child_u18, na.rm = TRUE),
is_snap_hh = any(is_snap_hh, na.rm = TRUE),
povpip = first(na.omit(povpip)),
income_131_165 = povpip >= pov_min & povpip <= pov_max,
hh_wt = first(wgtp),
unwt_persons = n(),
unwt_0_3 = sum(child_0_3, na.rm = TRUE),
unwt_u18 = sum(child_u18, na.rm = TRUE),
.groups = "drop"
)
# ---------------------------------------------------------------------------------
# PUMA-level household estimates
# ---------------------------------------------------------------------------------
puma_hh <- hh_level %>%
group_by(puma) %>%
summarise(
total_hh =
sum(hh_wt, na.rm = TRUE),
snap_hh =
sum(hh_wt[is_snap_hh], na.rm = TRUE),
# Households on SNAP with at least one child age 0-3, regardless of poverty level
snap_hh_0_3 =
sum(
hh_wt[
is_snap_hh &
has_0_3
],
na.rm = TRUE
),
# Households on SNAP with income 131-165% of poverty, regardless of child status
snap_hh_131_165 =
sum(
hh_wt[
is_snap_hh &
income_131_165
],
na.rm = TRUE
),
# Main universe: households on SNAP with at least one child age 0-3
# and income 131-165% of poverty
snap_hh_0_3_131_165 =
sum(
hh_wt[
is_snap_hh &
has_0_3 &
income_131_165
],
na.rm = TRUE
),
snap_hh_u18_131_165 =
sum(
hh_wt[
is_snap_hh &
has_u18 &
income_131_165
],
na.rm = TRUE
),
snap_hh_u18 =
sum(
hh_wt[
is_snap_hh &
has_u18
],
na.rm = TRUE
),
unwt_total_hh =
n(),
unwt_snap_hh =
sum(is_snap_hh, na.rm = TRUE),
.groups = "drop"
)
# ---------------------------------------------------------------------------------
# PUMA-level population estimates
# ---------------------------------------------------------------------------------
puma_pop <- tx_clean %>%
group_by(puma) %>%
summarise(
# total pop all ages
pop_total =
sum(pwgtp, na.rm = TRUE),
# pop all ages on SNAP
pop_snap =
sum(pwgtp[is_snap_hh], na.rm = TRUE),
# pop ages 0-3 on SNAP, regardless of poverty level
pop_0_3_snap =
sum(pwgtp[is_snap_hh & child_0_3], na.rm = TRUE),
# Population ages 0-3 living in SNAP households with household income 131-165% of poverty
pop_0_3_snap_131_165 =
sum(pwgtp[is_snap_hh & child_0_3 & income_131_165], na.rm = TRUE),
.groups = "drop"
)
puma_estimates <- puma_hh %>%
left_join(puma_pop, by = "puma")Checks
# ---------------------------------------------------------------------------------
# QA checks for statewide totals
# ---------------------------------------------------------------------------------
qa_statewide <- puma_estimates %>%
summarise(
total_hh = sum(total_hh, na.rm = TRUE),
snap_hh = sum(snap_hh, na.rm = TRUE),
snap_hh_0_3 = sum(snap_hh_0_3, na.rm = TRUE),
snap_hh_131_165 = sum(snap_hh_131_165, na.rm = TRUE),
snap_hh_0_3_131_165 = sum(snap_hh_0_3_131_165, na.rm = TRUE),
pop_total = sum(pop_total, na.rm = TRUE),
pop_snap = sum(pop_snap, na.rm = TRUE),
pop_0_3_snap = sum(pop_0_3_snap, na.rm = TRUE),
pop_0_3_snap_131_165 = sum(pop_0_3_snap_131_165, na.rm = TRUE)
)
qa_statewide %>%
print(width = Inf)# A tibble: 1 × 9
total_hh snap_hh snap_hh_0_3 snap_hh_131_165 snap_hh_0_3_131_165 pop_total
<dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
1 10992823 1257339 274113 127583 29249 30188424
pop_snap pop_0_3_snap pop_0_3_snap_131_165
<dbl> <dbl> <dbl>
1 4461936 393901 41084# these should all be TRUE
qa_nested <- qa_statewide %>%
transmute(
snap_hh_lte_total_hh = snap_hh <= total_hh,
snap_hh_0_3_lte_snap_hh = snap_hh_0_3 <= snap_hh,
snap_hh_131_165_lte_snap_hh = snap_hh_131_165 <= snap_hh,
snap_hh_0_3_131_165_lte_snap_hh_0_3 = snap_hh_0_3_131_165 <= snap_hh_0_3,
snap_hh_0_3_131_165_lte_snap_hh_131_165 = snap_hh_0_3_131_165 <= snap_hh_131_165,
pop_snap_lte_pop_total = pop_snap <= pop_total,
pop_0_3_snap_lte_pop_snap = pop_0_3_snap <= pop_snap,
pop_0_3_snap_131_165_lte_pop_0_3_snap = pop_0_3_snap_131_165 <= pop_0_3_snap
)
qa_nested %>%
print(width = Inf)# A tibble: 1 × 8
snap_hh_lte_total_hh snap_hh_0_3_lte_snap_hh snap_hh_131_165_lte_snap_hh
<lgl> <lgl> <lgl>
1 TRUE TRUE TRUE
snap_hh_0_3_131_165_lte_snap_hh_0_3 snap_hh_0_3_131_165_lte_snap_hh_131_165
<lgl> <lgl>
1 TRUE TRUE
pop_snap_lte_pop_total pop_0_3_snap_lte_pop_snap
<lgl> <lgl>
1 TRUE TRUE
pop_0_3_snap_131_165_lte_pop_0_3_snap
<lgl>
1 TRUE These estimates all look to be accurate:
** Note: Income-to-poverty ratios using Census Bureau Poverty Thresholds are used as a proxy for Federal Poverty Line (FPL).
| Measure | Estimate |
|---|---|
| Total Texas households | 10,992,823 |
| SNAP households | 1,257,339 |
| SNAP households with kids ages 0–3 | 274,113 |
| SNAP households at 131–165% FPL | 127,583 |
| SNAP households with kids 0–3 and 131–165% FPL | 29,249 |
| Total Texas population | 30,188,424 |
| People in SNAP households | 4,461,936 |
| Kids ages 0–3 in SNAP households | 393,901 |
| Kids ages 0–3 in SNAP households at 131–165% FPL | 41,084 |
Geographic Crosswalk
# ---------------------------------------------------------------------------------
# function to read Geocorr crosswalks
# ---------------------------------------------------------------------------------
read_geocorr_crosswalk <- function(path, chamber = c("lower", "upper")) {
chamber <- match.arg(chamber)
raw <- read_csv(path, show_col_types = FALSE) %>%
clean_names() %>%
filter(state != "State code")
if (chamber == "lower") {
out <- raw %>%
transmute(
puma = str_pad(as.character(puma22), width = 5, pad = "0"),
district = str_pad(as.character(sldl24), width = 3, pad = "0"),
district_label = paste0("HD-", district),
chamber = "lower",
geography = "Texas House District",
state = stab,
puma_name = puma22name,
geocorr_pop20 = as.numeric(pop20),
afact = as.numeric(afact)
)
}
if (chamber == "upper") {
out <- raw %>%
transmute(
puma = str_pad(as.character(puma22), width = 5, pad = "0"),
district = str_pad(as.character(sldu24), width = 3, pad = "0"),
district_label = paste0("SD-", district),
chamber = "upper",
geography = "Texas Senate District",
state = stab,
puma_name = puma22name,
geocorr_pop20 = as.numeric(pop20),
afact = as.numeric(afact)
)
}
out %>%
filter(
!is.na(puma),
!is.na(district),
!is.na(afact)
)
}
puma_to_lower <- read_geocorr_crosswalk(
path = lower_crosswalk_path,
chamber = "lower"
)
puma_to_upper <- read_geocorr_crosswalk(
path = upper_crosswalk_path,
chamber = "upper"
)
# ---------------------------------------------------------------------------------
# check crosswalk allocation factors
# ---------------------------------------------------------------------------------
# afact should sum to about 1 within each PUMA
lower_crosswalk_check <- puma_to_lower %>%
group_by(puma) %>%
summarise(
afact_sum = sum(afact, na.rm = TRUE),
n_districts = n(),
.groups = "drop"
) %>%
arrange(desc(abs(afact_sum - 1)))
upper_crosswalk_check <- puma_to_upper %>%
group_by(puma) %>%
summarise(
afact_sum = sum(afact, na.rm = TRUE),
n_districts = n(),
.groups = "drop"
) %>%
arrange(desc(abs(afact_sum - 1)))
print(lower_crosswalk_check)# A tibble: 217 × 3
puma afact_sum n_districts
<chr> <dbl> <int>
1 02509 1.000 6
2 04621 1.000 3
3 00100 1.00 3
4 00400 1.00 3
5 01903 1.000 4
6 02302 1.00 4
7 02501 1.000 3
8 02513 1.00 5
9 02515 1.00 4
10 02600 1.00 6
# ℹ 207 more rowsprint(upper_crosswalk_check)# A tibble: 217 × 3
puma afact_sum n_districts
<chr> <dbl> <int>
1 00600 1.00 3
2 01909 1.000 2
3 02306 1.000 2
4 02312 1.000 3
5 02502 1.00 3
6 02505 1.000 3
7 02600 1.00 4
8 04602 1.00 3
9 04614 1.000 3
10 04622 1.00 3
# ℹ 207 more rows# ---------------------------------------------------------------------------------
# function to allocate PUMA estimates to legislative districts
# ---------------------------------------------------------------------------------
allocate_to_districts <- function(puma_estimates, crosswalk) {
allocated_puma_district <- puma_estimates %>%
left_join(crosswalk, by = "puma") %>%
mutate(
total_hh_alloc =
total_hh * afact,
snap_hh_alloc =
snap_hh * afact,
snap_hh_0_3_alloc =
snap_hh_0_3 * afact,
snap_hh_131_165_alloc =
snap_hh_131_165 * afact,
snap_hh_0_3_131_165_alloc =
snap_hh_0_3_131_165 * afact,
# snap_hh_u18_131_165_alloc =
# snap_hh_u18_131_165 * afact,
#
# snap_hh_u18_alloc =
# snap_hh_u18 * afact,
pop_0_3_snap_alloc =
pop_0_3_snap * afact,
pop_0_3_snap_131_165_alloc =
pop_0_3_snap_131_165 * afact,
pop_snap_alloc =
pop_snap * afact,
pop_total_alloc =
pop_total * afact
)
district_summary <- allocated_puma_district %>%
group_by(
geography,
chamber,
district,
district_label
) %>%
summarise(
pumas =
paste(sort(unique(puma)), collapse = "; "),
puma_names =
paste(sort(unique(puma_name)), collapse = "; "),
n_pumas =
n_distinct(puma),
geocorr_pop20 =
sum(geocorr_pop20, na.rm = TRUE),
total_hh =
sum(total_hh_alloc, na.rm = TRUE),
snap_hh =
sum(snap_hh_alloc, na.rm = TRUE),
snap_hh_0_3 =
sum(snap_hh_0_3_alloc, na.rm = TRUE),
snap_hh_131_165 =
sum(snap_hh_131_165_alloc, na.rm = TRUE),
snap_hh_0_3_131_165 =
sum(snap_hh_0_3_131_165_alloc, na.rm = TRUE),
# snap_hh_u18_131_165 =
# sum(snap_hh_u18_131_165_alloc, na.rm = TRUE),
#
# snap_hh_u18 =
# sum(snap_hh_u18_alloc, na.rm = TRUE),
#
pop_0_3_snap =
sum(pop_0_3_snap_alloc, na.rm = TRUE),
pop_0_3_snap_131_165 =
sum(pop_0_3_snap_131_165_alloc, na.rm = TRUE),
pop_snap =
sum(pop_snap_alloc, na.rm = TRUE),
pop_total =
sum(pop_total_alloc, na.rm = TRUE),
.groups = "drop"
) %>%
mutate(
pct_snap_hh_0_3_131_165 =
snap_hh_0_3_131_165 / snap_hh,
pct_snap_pop_0_3_131_165 =
pop_0_3_snap_131_165 / pop_snap,
# pct_hh_u18_131_165_snap =
# snap_hh_u18_131_165 / total_hh,
#
# pct_hh_u18_snap =
# snap_hh_u18 / total_hh,
pov_level = "131-165% of Census Poverty Threshold",
note = "ACS PUMS household estimates allocated from PUMAs using Geocorr population-weighted allocation factors."
) %>%
mutate(
across(
c(
total_hh,
snap_hh,
snap_hh_0_3,
snap_hh_131_165,
snap_hh_0_3_131_165,
pop_0_3_snap,
pop_0_3_snap_131_165,
pop_snap,
pop_total,
geocorr_pop20
),
~round(.x, 0)
),
across(
starts_with("pct_"),
~round(.x * 100, 2)
)
) %>%
select(
geography,
chamber,
district,
district_label,
pct_snap_hh_0_3_131_165, ## primary measure (%) - households
snap_hh_0_3_131_165, ## primary measure (#) - households
pct_snap_pop_0_3_131_165, ## primary measure (%) - person-level
pop_0_3_snap_131_165, ## primary measure (#) - person-level
pop_0_3_snap,
snap_hh_0_3,
snap_hh_131_165,
pop_snap,
pop_total,
total_hh,
snap_hh,
pumas,
puma_names,
n_pumas,
geocorr_pop20,
pov_level,
note
) %>%
arrange(desc(pct_snap_hh_0_3_131_165))
list(
allocated_puma_district = allocated_puma_district,
district_summary = district_summary
)
}
lower_results <- allocate_to_districts(
puma_estimates = puma_estimates,
crosswalk = puma_to_lower
)
upper_results <- allocate_to_districts(
puma_estimates = puma_estimates,
crosswalk = puma_to_upper
)
lower_df <- lower_results$district_summary
upper_df <- upper_results$district_summary
# ---------------------------------------------------------------------------------
# save final Excel file
# ---------------------------------------------------------------------------------
wb <- createWorkbook()
addWorksheet(wb, "lower")
writeData(wb, "lower", lower_df)
freezePane(wb, "lower", firstRow = TRUE)
setColWidths(wb, "lower", cols = 1:ncol(lower_df), widths = "auto")
addWorksheet(wb, "upper")
writeData(wb, "upper", upper_df)
freezePane(wb, "upper", firstRow = TRUE)
setColWidths(wb, "upper", cols = 1:ncol(upper_df), widths = "auto")
# ---------------------------------------------------------------------------------
# add README sheet with notes and column definitions
# ---------------------------------------------------------------------------------
readme_overview <- tibble::tribble(
~field, ~description,
"Purpose",
"This workbook estimates Texas households participating in SNAP with at least one child age 0-3 and household income between 131% and 165% of Census Bureau Poverty Thresholds (used as proxy for FPL), by Texas legislative district.",
"Data source",
"2024 ACS PUMS microdata pulled with tidycensus, plus MCDC Geocorr 2022 PUMA-to-state-legislative-district crosswalks using 2022 PUMAs and 2024 legislative districts.",
"Geography method",
"ACS PUMS does not directly identify Texas House or Senate districts. Estimates are first calculated by PUMA and then allocated to legislative districts using Geocorr population-weighted allocation factors.",
"Unit of analysis",
"Households are used as the analytic unit. This is a household-level estimate, not an IPUMS family-unit estimate.",
"Interpretation note",
"District estimates are approximations because PUMAs and legislative districts do not perfectly align."
)
readme_columns <- tibble::tribble(
~field, ~description,
"geography",
"Legislative district type: Texas House District or Texas Senate District.",
"chamber",
"Legislative chamber: lower for Texas House, upper for Texas Senate.",
"district",
"Legislative district number, padded to three digits.",
"district_label",
"District label with chamber prefix, such as HD-001 or SD-001.",
"pct_snap_hh_0_3_131_165",
"Among SNAP households in the district, the percent with at least one child age 0-3 and household income between 131% and 165% of poverty. This is a household-level measure. Percent values are stored as percentages, not proportions; for example, 12.5 means 12.5%.",
"snap_hh_0_3_131_165",
"Estimated number of SNAP households with at least one child age 0-3 and household income between 131% and 165% of poverty threshold. This is a household-level estimate using household weights.",
"pct_snap_pop_0_3_131_165",
"Among all people living in SNAP households in the district, the percent who are ages 0-3 and live in households with income between 131% and 165% of poverty threshold. This is a person-level measure. Percent values are stored as percentages, not proportions; for example, 12.5 means 12.5%.",
"pop_0_3_snap_131_165",
"Estimated number of children ages 0-3 living in SNAP households with household income between 131% and 165% of poverty threshold. This is a person-level estimate using child records and person weights.",
"pop_0_3_snap",
"Estimated number of children ages 0-3 living in SNAP households, regardless of regardless of household income-to-poverty ratio.",
"snap_hh_0_3",
"Estimated number of SNAP households with at least one child age 0-3, regardless of household income-to-poverty ratio.",
"snap_hh_131_165",
"Estimated number of SNAP households with household income between 131% and 165% of poverty threshold, regardless of whether children are present.",
"pop_snap",
"Estimated total population of all ages living in SNAP households.",
"pop_total",
"Estimated total population of all ages in the district, calculated from ACS PUMS PUMA-level estimates allocated to districts.",
"total_hh",
"Estimated total number of households in the district, calculated from ACS PUMS PUMA-level estimates allocated to districts.",
"snap_hh",
"Estimated total number of SNAP households in the district, calculated from ACS PUMS PUMA-level estimates allocated to districts.",
"pumas",
"PUMA codes that overlap the legislative district and contribute to the allocated estimate.",
"puma_names",
"Names of PUMAs that overlap the legislative district and contribute to the allocated estimate.",
"n_pumas",
"Number of PUMAs that overlap the legislative district.",
"geocorr_pop20",
"2020 Census population from the Geocorr PUMA-to-district overlap file. This is a crosswalk/reference population used for context and QA; it is not an ACS PUMS estimate.",
"pov_level",
"Income to poverty ratio used in the analysis, measured as a person's or family's total annual income as a percentage of their official federal poverty threshold. Here, 131-165% of poverty",
"note",
"Brief methodology note describing that ACS PUMS household and person estimates were allocated from PUMAs to districts using Geocorr population-weighted allocation factors."
)
addWorksheet(wb, "README")
# Styles
header_style <- createStyle(
textDecoration = "bold",
fontSize = 14,
fgFill = "#D9EAF7",
halign = "left"
)
table_header_style <- createStyle(
textDecoration = "bold",
fgFill = "#EDEDED",
border = "Bottom"
)
wrap_style <- createStyle(
wrapText = TRUE,
valign = "top"
)
# Analysis Overview
writeData(wb, "README", "Analysis Overview", startRow = 1, startCol = 1)
addStyle(wb, "README", header_style, rows = 1, cols = 1:2, gridExpand = TRUE)
writeData(wb, "README", readme_overview, startRow = 2, startCol = 1)
addStyle(wb, "README", table_header_style, rows = 2, cols = 1:2, gridExpand = TRUE)
addStyle(
wb, "README", wrap_style,
rows = 3:(nrow(readme_overview) + 2),
cols = 1:2,
gridExpand = TRUE
)
# Column Definitions
col_def_start <- nrow(readme_overview) + 5
writeData(wb, "README", "Column Definitions", startRow = col_def_start, startCol = 1)
addStyle(wb, "README", header_style, rows = col_def_start, cols = 1:2, gridExpand = TRUE)
writeData(wb, "README", readme_columns, startRow = col_def_start + 1, startCol = 1)
addStyle(
wb, "README", table_header_style,
rows = col_def_start + 1,
cols = 1:2,
gridExpand = TRUE
)
addStyle(
wb, "README", wrap_style,
rows = (col_def_start + 2):(col_def_start + nrow(readme_columns) + 1),
cols = 1:2,
gridExpand = TRUE
)
# formatting
setColWidths(wb, "README", cols = 1, widths = 28)
setColWidths(wb, "README", cols = 2, widths = 120)
freezePane(wb, "README", firstRow = FALSE)
saveWorkbook(wb, excel_output_path, overwrite = TRUE)Mapping
# ---------------------------------------------------------------------------------
# pull district geometries for maps
# ---------------------------------------------------------------------------------
tx_lower_geom <- state_legislative_districts(
state = "TX",
house = "lower",
year = 2024,
cb = TRUE
) %>%
st_transform(3083) %>%
mutate(
district = str_pad(as.character(SLDLST), width = 3, pad = "0")
)
tx_upper_geom <- state_legislative_districts(
state = "TX",
house = "upper",
year = 2024,
cb = TRUE
) %>%
st_transform(3083) %>%
mutate(
district = str_pad(as.character(SLDUST), width = 3, pad = "0")
)
lower_map_data <- tx_lower_geom %>%
left_join(
lower_df %>%
select(
district,
district_label,
pct_snap_hh_0_3_131_165
),
by = "district"
)
upper_map_data <- tx_upper_geom %>%
left_join(
upper_df %>%
select(
district,
district_label,
pct_snap_hh_0_3_131_165
),
by = "district"
)
# ---------------------------------------------------------------------------------
# lower chamber map
# ---------------------------------------------------------------------------------
lower_map <- ggplot(lower_map_data) +
geom_sf(
aes(fill = pct_snap_hh_0_3_131_165),
color = "white",
linewidth = 0.1
) +
scale_fill_viridis_c(
option = "rocket",
direction = -1,
labels = label_percent(scale = 1),
na.value = "grey95"
) +
labs(
title = "Texas House Districts",
subtitle = "Percent of SNAP households with children ages 0-3 and income 131-165% of poverty",
fill = "Percent",
caption = "Source: 2024 ACS PUMS via tidycensus; PUMA estimates allocated to districts using MCDC Geocorr 2022."
) +
theme_void() +
theme(
plot.title = element_text(face = "bold"),
plot.subtitle = element_text(size = 10),
legend.position = "right"
)
ggsave(
filename = lower_map_path,
plot = lower_map,
width = 10,
height = 8,
dpi = 300
)
lower_map
# ---------------------------------------------------------------------------------
# upper chamber map
# ---------------------------------------------------------------------------------
upper_map <- ggplot(upper_map_data) +
geom_sf(
aes(fill = pct_snap_hh_0_3_131_165),
color = "white",
linewidth = 0.15
) +
scale_fill_viridis_c(
option = "rocket",
direction = -1,
labels = label_percent(scale = 1),
na.value = "grey95"
) +
labs(
title = "Texas Senate Districts",
subtitle = "Percent of SNAP households with children ages 0-3 and income 131-165% of poverty",
fill = "Percent",
caption = "Source: 2024 ACS PUMS via tidycensus; PUMA estimates allocated to districts using MCDC Geocorr 2022."
) +
theme_void() +
theme(
plot.title = element_text(face = "bold"),
plot.subtitle = element_text(size = 10),
legend.position = "right"
)
ggsave(
filename = upper_map_path,
plot = upper_map,
width = 10,
height = 8,
dpi = 300
)
upper_map
State-Level Estimates
statewide_pct <- puma_estimates %>%
summarise(
snap_hh_0_3_131_165 = sum(snap_hh_0_3_131_165, na.rm = TRUE),
snap_hh_0_3 = sum(snap_hh_0_3, na.rm = TRUE),
snap_hh_131_165 = sum(snap_hh_131_165, na.rm = TRUE),
snap_hh = sum(snap_hh, na.rm = TRUE),
pop_0_3_snap_131_165 = sum(pop_0_3_snap_131_165, na.rm = TRUE),
pop_0_3_snap = sum(pop_0_3_snap, na.rm = TRUE),
pop_snap = sum(pop_snap, na.rm = TRUE),
pop_total = sum(pop_total, na.rm = TRUE),
pct_snap_hh_0_3_131_165 =
snap_hh_0_3_131_165 / snap_hh * 100,
pct_snap_pop_0_3_131_165 =
pop_0_3_snap_131_165 / pop_snap * 100
) %>%
mutate(
across(
c(
snap_hh_0_3_131_165,
snap_hh_0_3,
snap_hh_131_165,
snap_hh,
pop_0_3_snap_131_165,
pop_0_3_snap,
pop_snap,
pop_total
),
~round(.x, 0)
),
across(
starts_with("pct_"),
~round(.x, 2)
)
)
statewide_pct %>%
print(width = Inf)# A tibble: 1 × 10
snap_hh_0_3_131_165 snap_hh_0_3 snap_hh_131_165 snap_hh pop_0_3_snap_131_165
<dbl> <dbl> <dbl> <dbl> <dbl>
1 29249 274113 127583 1257339 41084
pop_0_3_snap pop_snap pop_total pct_snap_hh_0_3_131_165
<dbl> <dbl> <dbl> <dbl>
1 393901 4461936 30188424 2.33
pct_snap_pop_0_3_131_165
<dbl>
1 0.92Statewide, an estimated 29,249 households receiving SNAPhave at least one child age 0–3 and income between 131% and 165% of poverty, representing 2.33% of all households receiving SNAP.
| Measure | Estimate | Interpretation |
|---|---|---|
| SNAP households with children ages 0–3 and income 131–165% of poverty | 29,249 | Estimated number of SNAP households with at least one child age 0–3 and household income between 131% and 165% of poverty. |
| SNAP households with children ages 0–3 | 274,113 | Estimated number of SNAP households with at least one child age 0–3, regardless of income band. |
| SNAP households with income 131–165% of poverty | 127,583 | Estimated number of SNAP households with household income between 131% and 165% of poverty, regardless of whether children are present. |
| Total SNAP households | 1,257,339 | Estimated total number of households receiving SNAP. |
| Children ages 0–3 in SNAP households with income 131–165% of poverty | 41,084 | Estimated number of children ages 0–3 living in SNAP households with income between 131% and 165% of poverty. |
| Children ages 0–3 in SNAP households | 393,901 | Estimated number of children ages 0–3 living in SNAP households, regardless of income band. |
| Total population living in SNAP households | 4,461,936 | Estimated number of people of all ages living in SNAP households. |
| Total population | 30,188,424 | Estimated total Texas population in the ACS PUMS data. |
| Percent of SNAP households with children ages 0–3 and income 131–165% of poverty | 2.33% | Among all SNAP households, the share with at least one child age 0–3 and income between 131% and 165% of poverty. |
| Percent of SNAP population that is ages 0–3 and in households with income 131–165% of poverty | 0.92% | Among all people living in SNAP households, the share who are ages 0–3 and live in households with income between 131% and 165% of poverty. |
Additional Checks (crosswalk matches and reliability)
# anti-join checks for PUMA to district crosswalks
lower_pums_missing_from_crosswalk <- puma_estimates %>%
distinct(puma) %>%
anti_join(
puma_to_lower %>% distinct(puma),
by = "puma"
)
upper_pums_missing_from_crosswalk <- puma_estimates %>%
distinct(puma) %>%
anti_join(
puma_to_upper %>% distinct(puma),
by = "puma"
)
lower_crosswalk_missing_from_pums <- puma_to_lower %>%
distinct(puma) %>%
anti_join(
puma_estimates %>% distinct(puma),
by = "puma"
)
upper_crosswalk_missing_from_pums <- puma_to_upper %>%
distinct(puma) %>%
anti_join(
puma_estimates %>% distinct(puma),
by = "puma"
)
crosswalk_match_summary <- tibble(
check = c(
"lower: PUMS PUMAs missing from crosswalk",
"upper: PUMS PUMAs missing from crosswalk",
"lower: Crosswalk PUMAs missing from PUMS",
"upper: Crosswalk PUMAs missing from PUMS"
),
n_missing = c(
nrow(lower_pums_missing_from_crosswalk),
nrow(upper_pums_missing_from_crosswalk),
nrow(lower_crosswalk_missing_from_pums),
nrow(upper_crosswalk_missing_from_pums)
)
) %>%
mutate(
passes = n_missing == 0
)
crosswalk_match_summary %>%
print(width = Inf)# A tibble: 4 × 3
check n_missing passes
<chr> <int> <lgl>
1 lower: PUMS PUMAs missing from crosswalk 0 TRUE
2 upper: PUMS PUMAs missing from crosswalk 0 TRUE
3 lower: Crosswalk PUMAs missing from PUMS 0 TRUE
4 upper: Crosswalk PUMAs missing from PUMS 0 TRUE All Texas PUMAs in the ACS PUMS data are represented in both the House and Senate district crosswalks, and all PUMAs in the crosswalks are present in the PUMS data.
# statewide reliability check using ACS PUMS household replicate weights
pums_rep_cache_path <- file.path(
data_dir,
"tx_pums_2024_acs5_statewide_hh_repwts_with_puma.rds"
)
rep_vars <- c(
"SERIALNO",
"PUMA",
"AGEP",
"POVPIP",
"FS",
"WGTP"
)
if (file.exists(pums_rep_cache_path)) {
tx_pums_reps <- readRDS(pums_rep_cache_path)
} else {
tx_pums_reps <- get_pums(
variables = rep_vars,
state = "TX",
survey = acs_survey,
year = acs_year,
recode = FALSE,
rep_weights = "housing"
)
saveRDS(tx_pums_reps, pums_rep_cache_path)
}# check that object exists
exists("tx_pums_reps")[1] TRUE# check dimensions
dim(tx_pums_reps)[1] 1326295 89# check column names
names(tx_pums_reps) [1] "SERIALNO" "SPORDER" "PUMA" "AGEP" "POVPIP" "FS"
[7] "STATE" "WGTP" "PWGTP" "WGTP1" "WGTP2" "WGTP3"
[13] "WGTP4" "WGTP5" "WGTP6" "WGTP7" "WGTP8" "WGTP9"
[19] "WGTP10" "WGTP11" "WGTP12" "WGTP13" "WGTP14" "WGTP15"
[25] "WGTP16" "WGTP17" "WGTP18" "WGTP19" "WGTP20" "WGTP21"
[31] "WGTP22" "WGTP23" "WGTP24" "WGTP25" "WGTP26" "WGTP27"
[37] "WGTP28" "WGTP29" "WGTP30" "WGTP31" "WGTP32" "WGTP33"
[43] "WGTP34" "WGTP35" "WGTP36" "WGTP37" "WGTP38" "WGTP39"
[49] "WGTP40" "WGTP41" "WGTP42" "WGTP43" "WGTP44" "WGTP45"
[55] "WGTP46" "WGTP47" "WGTP48" "WGTP49" "WGTP50" "WGTP51"
[61] "WGTP52" "WGTP53" "WGTP54" "WGTP55" "WGTP56" "WGTP57"
[67] "WGTP58" "WGTP59" "WGTP60" "WGTP61" "WGTP62" "WGTP63"
[73] "WGTP64" "WGTP65" "WGTP66" "WGTP67" "WGTP68" "WGTP69"
[79] "WGTP70" "WGTP71" "WGTP72" "WGTP73" "WGTP74" "WGTP75"
[85] "WGTP76" "WGTP77" "WGTP78" "WGTP79" "WGTP80" # Check that household replicate weights came through
names(tx_pums_reps) %>%
str_subset("^WGTP[0-9]+$") [1] "WGTP1" "WGTP2" "WGTP3" "WGTP4" "WGTP5" "WGTP6" "WGTP7" "WGTP8"
[9] "WGTP9" "WGTP10" "WGTP11" "WGTP12" "WGTP13" "WGTP14" "WGTP15" "WGTP16"
[17] "WGTP17" "WGTP18" "WGTP19" "WGTP20" "WGTP21" "WGTP22" "WGTP23" "WGTP24"
[25] "WGTP25" "WGTP26" "WGTP27" "WGTP28" "WGTP29" "WGTP30" "WGTP31" "WGTP32"
[33] "WGTP33" "WGTP34" "WGTP35" "WGTP36" "WGTP37" "WGTP38" "WGTP39" "WGTP40"
[41] "WGTP41" "WGTP42" "WGTP43" "WGTP44" "WGTP45" "WGTP46" "WGTP47" "WGTP48"
[49] "WGTP49" "WGTP50" "WGTP51" "WGTP52" "WGTP53" "WGTP54" "WGTP55" "WGTP56"
[57] "WGTP57" "WGTP58" "WGTP59" "WGTP60" "WGTP61" "WGTP62" "WGTP63" "WGTP64"
[65] "WGTP65" "WGTP66" "WGTP67" "WGTP68" "WGTP69" "WGTP70" "WGTP71" "WGTP72"
[73] "WGTP73" "WGTP74" "WGTP75" "WGTP76" "WGTP77" "WGTP78" "WGTP79" "WGTP80"file.exists(pums_rep_cache_path)[1] TRUEfile.info(pums_rep_cache_path)$size[1] 89129209# # statewide CV / reliability calculation
library(data.table)
Attaching package: 'data.table'The following objects are masked from 'package:lubridate':
hour, isoweek, isoyear, mday, minute, month, quarter, second, wday,
week, yday, yearThe following objects are masked from 'package:dplyr':
between, first, lastThe following object is masked from 'package:purrr':
transposeThe following object is masked from 'package:base':
%notin%rep_cols <- paste0("WGTP", 1:80)
# make sure replicate data has PUMA
if (!"PUMA" %in% names(tx_pums_reps)) {
stop("tx_pums_reps does not include PUMA. Re-run the replicate-weight pull with PUMA included.")
}
# Target households from the main analysis
# This should reproduce 29,249 estimate exactly
target_hh <- hh_level %>%
mutate(
target_29k =
is_snap_hh &
has_0_3 &
income_131_165
) %>%
select(
serialno,
puma,
target_29k,
hh_wt
)
main_estimate_check <- target_hh %>%
summarise(
estimate = sum(hh_wt[target_29k], na.rm = TRUE)
)
main_estimate_check# A tibble: 1 × 1
estimate
<dbl>
1 29249# pull one household level replicate weight row per SERIALNO/PUMA
dt_reps <- as.data.table(tx_pums_reps)
dt_reps[, puma := str_pad(as.character(PUMA), width = 5, pad = "0")]
dt_reps[, serialno := as.character(SERIALNO)]
# convert replicate weights to numeric
for (w in rep_cols) {
dt_reps[, (w) := as.numeric(get(w))]
}
hh_rep_weights <- dt_reps[
,
lapply(.SD, function(x) x[1]),
by = .(serialno, puma),
.SDcols = rep_cols
]
# join target universe from main analysis to replicate weights
target_hh_dt <- as.data.table(
target_hh %>%
mutate(
serialno = as.character(serialno),
puma = str_pad(as.character(puma), width = 5, pad = "0")
)
)
setkey(target_hh_dt, serialno, puma)
setkey(hh_rep_weights, serialno, puma)
cv_data <- hh_rep_weights[target_hh_dt]
# check for missing replicate weights after join
missing_rep_join <- cv_data[
target_29k == TRUE & is.na(WGTP1),
.N
]
if (missing_rep_join > 0) {
warning("Some target households are missing replicate weights after join.")
}
# full estimate, using original main analysis household weight
estimate <- cv_data[
target_29k == TRUE,
sum(hh_wt, na.rm = TRUE)
]
# replicate estimates
rep_estimates <- sapply(
rep_cols,
function(w) {
cv_data[
target_29k == TRUE,
sum(get(w), na.rm = TRUE)
]
}
)
# ACS PUMS SE formula
se <- sqrt((4 / 80) * sum((rep_estimates - estimate)^2, na.rm = TRUE))
cv <- se / estimate
reliability <- case_when(
is.na(cv) ~ "Not available",
cv < 0.15 ~ "Reliable",
cv < 0.30 ~ "Use caution",
TRUE ~ "Unreliable"
)
statewide_cv_29k <- tibble(
measure = "Statewide SNAP households with child age 0-3 and income 131-165% of poverty",
estimate = round(estimate, 0),
standard_error = round(se, 0),
cv_percent = round(cv * 100, 2),
reliability = reliability
)
statewide_cv_29k# A tibble: 1 × 5
measure estimate standard_error cv_percent reliability
<chr> <dbl> <dbl> <dbl> <chr>
1 Statewide SNAP households with… 29249 1178 4.03 Reliable Statewide estimate: 29,249 households on sNAPn with kids 0-3 and 131-165% above poverty
Standard error: 1,178
Coefficient of variation: 4.03%
Reliability flag: Reliable
# District-level SE / CV for primary household measure only
# Primary measure: snap_hh_0_3_131_165
rep_cols <- paste0("WGTP", 1:80)
# PUMA-level full estimate using the same target universe as the statewide 29,249 estimate
puma_cv_full <- cv_data[
,
.(
snap_hh_0_3_131_165 = sum(
hh_wt[target_29k == TRUE],
na.rm = TRUE
)
),
by = puma
]
# PUMA-level replicate estimates
puma_cv_reps <- rbindlist(
lapply(
seq_along(rep_cols),
function(i) {
w <- rep_cols[i]
cv_data[
,
.(
rep = i,
snap_hh_0_3_131_165_rep = sum(
get(w)[target_29k == TRUE],
na.rm = TRUE
)
),
by = puma
]
}
)
)
# Check: PUMA totals should still add to 29,249
puma_cv_full[
,
round(sum(snap_hh_0_3_131_165, na.rm = TRUE), 0)
][1] 29249# allocate PUMA estimates to districts and calculate CV
calc_district_cv <- function(crosswalk) {
# full district estimate
full_district <- puma_cv_full %>%
as_tibble() %>%
left_join(crosswalk, by = "puma") %>%
group_by(
geography,
chamber,
district,
district_label
) %>%
summarise(
estimate = sum(snap_hh_0_3_131_165 * afact, na.rm = TRUE),
.groups = "drop"
)
# rep district estimates
rep_district <- puma_cv_reps %>%
as_tibble() %>%
left_join(crosswalk, by = "puma") %>%
group_by(
geography,
chamber,
district,
district_label,
rep
) %>%
summarise(
rep_estimate = sum(snap_hh_0_3_131_165_rep * afact, na.rm = TRUE),
.groups = "drop"
)
# calc SE, CV, and reliability
district_cv <- rep_district %>%
left_join(
full_district,
by = c(
"geography",
"chamber",
"district",
"district_label"
)
) %>%
group_by(
geography,
chamber,
district,
district_label
) %>%
summarise(
estimate = first(estimate),
standard_error = sqrt(
(4 / 80) *
sum(
(rep_estimate - estimate)^2,
na.rm = TRUE
)
),
.groups = "drop"
) %>%
mutate(
cv = standard_error / estimate,
reliability = case_when(
is.na(cv) ~ "Not available",
estimate < 1 ~ "Estimate is zero or near zero",
cv < 0.15 ~ "Reliable",
cv < 0.30 ~ "Use caution",
TRUE ~ "Unreliable"
),
estimate = round(estimate, 0),
standard_error = round(standard_error, 0),
cv_percent = round(cv * 100, 2)
) %>%
select(
geography,
chamber,
district,
district_label,
estimate,
standard_error,
cv_percent,
reliability
)
district_cv
}# calculate House and Senate CVs
lower_cv <- calc_district_cv(puma_to_lower)Warning in left_join(., crosswalk, by = "puma"): Detected an unexpected many-to-many relationship between `x` and `y`.
ℹ Row 1 of `x` matches multiple rows in `y`.
ℹ Row 97 of `y` matches multiple rows in `x`.
ℹ If a many-to-many relationship is expected, set `relationship =
"many-to-many"` to silence this warning.upper_cv <- calc_district_cv(puma_to_upper)Warning in left_join(., crosswalk, by = "puma"): Detected an unexpected many-to-many relationship between `x` and `y`.
ℹ Row 1 of `x` matches multiple rows in `y`.
ℹ Row 63 of `y` matches multiple rows in `x`.
ℹ If a many-to-many relationship is expected, set `relationship =
"many-to-many"` to silence this warning.Warning above is ok, just means one PUMA can map to multiple lege districts
#reliability distribution
lower_cv %>%
count(reliability) %>%
print(width = Inf)# A tibble: 2 × 2
reliability n
<chr> <int>
1 Unreliable 103
2 Use caution 47upper_cv %>%
count(reliability) %>%
print(width = Inf)# A tibble: 3 × 2
reliability n
<chr> <int>
1 Reliable 6
2 Unreliable 2
3 Use caution 23# least reliable House estimates
lower_cv %>%
arrange(desc(cv_percent)) %>%
print(n = 25, width = Inf)# A tibble: 150 × 8
geography chamber district district_label estimate standard_error
<chr> <chr> <chr> <chr> <dbl> <dbl>
1 Texas House District lower 073 HD-073 35 38
2 Texas House District lower 066 HD-066 20 21
3 Texas House District lower 061 HD-061 15 12
4 Texas House District lower 020 HD-020 8 7
5 Texas House District lower 136 HD-136 118 90
6 Texas House District lower 052 HD-052 34 26
7 Texas House District lower 049 HD-049 28 21
8 Texas House District lower 044 HD-044 25 18
9 Texas House District lower 122 HD-122 69 44
10 Texas House District lower 070 HD-070 39 24
11 Texas House District lower 115 HD-115 100 62
12 Texas House District lower 016 HD-016 144 88
13 Texas House District lower 089 HD-089 41 25
14 Texas House District lower 051 HD-051 187 106
15 Texas House District lower 047 HD-047 45 25
16 Texas House District lower 014 HD-014 93 52
17 Texas House District lower 127 HD-127 136 74
18 Texas House District lower 048 HD-048 32 18
19 Texas House District lower 006 HD-006 86 46
20 Texas House District lower 057 HD-057 63 33
21 Texas House District lower 010 HD-010 146 76
22 Texas House District lower 101 HD-101 96 50
23 Texas House District lower 015 HD-015 100 52
24 Texas House District lower 086 HD-086 58 30
25 Texas House District lower 098 HD-098 130 67
cv_percent reliability
<dbl> <chr>
1 110. Unreliable
2 102. Unreliable
3 82.3 Unreliable
4 80.9 Unreliable
5 76.6 Unreliable
6 76.2 Unreliable
7 75.0 Unreliable
8 69.0 Unreliable
9 64 Unreliable
10 62.2 Unreliable
11 61.8 Unreliable
12 60.8 Unreliable
13 60.5 Unreliable
14 56.6 Unreliable
15 55.9 Unreliable
16 55.5 Unreliable
17 54.6 Unreliable
18 54.1 Unreliable
19 53.1 Unreliable
20 52.4 Unreliable
21 52.2 Unreliable
22 52.2 Unreliable
23 51.8 Unreliable
24 51.7 Unreliable
25 51.4 Unreliable
# ℹ 125 more rows# least reliable Senate estimates
upper_cv %>%
arrange(desc(cv_percent)) %>%
print(n = 25, width = Inf)# A tibble: 31 × 8
geography chamber district district_label estimate standard_error
<chr> <chr> <chr> <chr> <dbl> <dbl>
1 Texas Senate District upper 008 SD-008 185 74
2 Texas Senate District upper 025 SD-025 258 101
3 Texas Senate District upper 014 SD-014 726 217
4 Texas Senate District upper 007 SD-007 868 228
5 Texas Senate District upper 004 SD-004 731 190
6 Texas Senate District upper 012 SD-012 500 126
7 Texas Senate District upper 011 SD-011 725 173
8 Texas Senate District upper 031 SD-031 1002 223
9 Texas Senate District upper 018 SD-018 754 163
10 Texas Senate District upper 005 SD-005 421 90
11 Texas Senate District upper 030 SD-030 382 81
12 Texas Senate District upper 024 SD-024 1097 233
13 Texas Senate District upper 003 SD-003 1015 204
14 Texas Senate District upper 029 SD-029 835 162
15 Texas Senate District upper 016 SD-016 1197 230
16 Texas Senate District upper 001 SD-001 840 161
17 Texas Senate District upper 002 SD-002 421 78
18 Texas Senate District upper 023 SD-023 1499 272
19 Texas Senate District upper 019 SD-019 1182 212
20 Texas Senate District upper 013 SD-013 1308 233
21 Texas Senate District upper 026 SD-026 1086 191
22 Texas Senate District upper 017 SD-017 1454 256
23 Texas Senate District upper 021 SD-021 1435 245
24 Texas Senate District upper 009 SD-009 997 169
25 Texas Senate District upper 015 SD-015 902 141
cv_percent reliability
<dbl> <chr>
1 40.0 Unreliable
2 39.0 Unreliable
3 29.9 Use caution
4 26.3 Use caution
5 26.0 Use caution
6 25.2 Use caution
7 23.9 Use caution
8 22.3 Use caution
9 21.7 Use caution
10 21.4 Use caution
11 21.2 Use caution
12 21.2 Use caution
13 20.0 Use caution
14 19.4 Use caution
15 19.2 Use caution
16 19.2 Use caution
17 18.6 Use caution
18 18.2 Use caution
19 17.9 Use caution
20 17.8 Use caution
21 17.6 Use caution
22 17.6 Use caution
23 17.1 Use caution
24 16.9 Use caution
25 15.7 Use caution
# ℹ 6 more rows# attach to final district tables
lower_df_with_cv <- lower_df %>%
left_join(
lower_cv %>%
select(
district,
standard_error,
cv_percent,
reliability
),
by = "district"
)
upper_df_with_cv <- upper_df %>%
left_join(
upper_cv %>%
select(
district,
standard_error,
cv_percent,
reliability
),
by = "district"
)# save new Excel copy with CV/reliability columns added
excel_output_path_with_cv <- file.path(
output_dir,
"tx_snap_lege_districts_with_reliability.xlsx"
)
wb_cv <- createWorkbook()
# Lower chamber with reliability columns
addWorksheet(wb_cv, "lower")
writeData(wb_cv, "lower", lower_df_with_cv)
freezePane(wb_cv, "lower", firstRow = TRUE)
setColWidths(wb_cv, "lower", cols = 1:ncol(lower_df_with_cv), widths = "auto")
# Upper chamber with reliability columns
addWorksheet(wb_cv, "upper")
writeData(wb_cv, "upper", upper_df_with_cv)
freezePane(wb_cv, "upper", firstRow = TRUE)
setColWidths(wb_cv, "upper", cols = 1:ncol(upper_df_with_cv), widths = "auto")
# Add statewide reliability summary
addWorksheet(wb_cv, "statewide_reliability")
writeData(wb_cv, "statewide_reliability", statewide_cv_29k)
freezePane(wb_cv, "statewide_reliability", firstRow = TRUE)
setColWidths(
wb_cv,
"statewide_reliability",
cols = 1:ncol(statewide_cv_29k),
widths = "auto"
)
# Add reliability summary tabs
lower_reliability_summary <- lower_cv %>%
count(reliability, name = "n_districts")
upper_reliability_summary <- upper_cv %>%
count(reliability, name = "n_districts")
addWorksheet(wb_cv, "lower_reliability_summary")
writeData(wb_cv, "lower_reliability_summary", lower_reliability_summary)
freezePane(wb_cv, "lower_reliability_summary", firstRow = TRUE)
setColWidths(
wb_cv,
"lower_reliability_summary",
cols = 1:ncol(lower_reliability_summary),
widths = "auto"
)
addWorksheet(wb_cv, "upper_reliability_summary")
writeData(wb_cv, "upper_reliability_summary", upper_reliability_summary)
freezePane(wb_cv, "upper_reliability_summary", firstRow = TRUE)
setColWidths(
wb_cv,
"upper_reliability_summary",
cols = 1:ncol(upper_reliability_summary),
widths = "auto"
)
# README / notes
readme_cv <- tibble::tribble(
~field, ~description,
"Purpose",
"This workbook adds CV-based reliability fields to the legislative district estimates.",
"Primary reliability measure",
"Reliability is calculated for snap_hh_0_3_131_165: estimated SNAP households with at least one child age 0-3 and household income between 131% and 165% of poverty.",
"standard_error",
"Standard error for the primary household count estimate, calculated using ACS PUMS household replicate weights.",
"cv_percent",
"Coefficient of variation for the primary household count estimate, shown as a percent. CV = standard error / estimate.",
"reliability",
"CV-based reliability flag. CV <15% = Reliable; 15%-29.9% = Use caution; 30%+ = Unreliable. These are rule-of-thumb thresholds, not official Census suppression rules.",
"Important caveat",
"District estimates also include uncertainty from allocating PUMA-level estimates to legislative districts using Geocorr allocation factors. That allocation uncertainty is not fully captured by the CV."
)
addWorksheet(wb_cv, "README")
writeData(wb_cv, "README", readme_cv)
freezePane(wb_cv, "README", firstRow = TRUE)
setColWidths(wb_cv, "README", cols = 1, widths = 28)
setColWidths(wb_cv, "README", cols = 2, widths = 120)
# Save workbook
saveWorkbook(wb_cv, excel_output_path_with_cv, overwrite = TRUE)
excel_output_path_with_cv[1] "output/tx_snap_lege_districts_with_reliability.xlsx"