-
-
Save ryanburge/190eb8bfb5665554b072ea5d8b88cb95 to your computer and use it in GitHub Desktop.
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
library(rio) | |
library(janitor) | |
hp <- read_csv("E://data/hps.csv") %>% clean_names() | |
gg <- hp %>% | |
mutate(msa = est_msa) %>% | |
mutate(msa = frcode(msa == 35620 ~ "NYC", | |
msa == 31080 ~ "Charlotte", | |
msa == 16980 ~ "Chicago", | |
msa == 19100 ~ "Dallas", | |
msa == 26420 ~ "Houston", | |
msa == 47900 ~ "DC", | |
msa == 33100 ~ "Miami", | |
msa == 37980 ~ "Philadelphia", | |
msa == 12060 ~ "Atlanta", | |
msa == 38060 ~ "Phoenix", | |
msa == 14460 ~ "Boston", | |
msa == 41860 ~ "San Francisco", | |
msa == 40140 ~ "Riverside", | |
msa == 19820 ~ "Detroit", | |
msa == 42660 ~ "Seattle")) %>% | |
mutate(att = support3) %>% | |
mutate(att = frcode(att == 1 ~ "Less than Once a Year", | |
att == 2 ~ "1-3/Yr.", | |
att == 3 ~ "4-11/Yr.", | |
att == 4 ~ "12+/Yr.")) %>% | |
group_by(msa) %>% | |
ct(att, wt = pweight, show_na = FALSE) | |
lvl <- gg %>% | |
filter(att == "Less than Once a Year") %>% | |
select(msa, sort = pct) | |
gg <- left_join(gg, lvl) | |
gg$msa <- fct_reorder(gg$msa, gg$sort, .desc = TRUE) | |
gg %>% | |
filter(msa != "NA") %>% | |
mutate(lab = round(pct, 2)) %>% | |
ggplot(., aes(x = 1, y = pct, fill = fct_rev(att))) + | |
geom_col(color = "black") + | |
coord_flip() + | |
facet_wrap(~ msa, ncol =1, strip.position = "left") + | |
theme_rb() + | |
scale_fill_manual(values = c(moma.colors("ustwo", 4))) + | |
theme(legend.position = "bottom") + | |
scale_y_continuous(labels = percent) + | |
theme(strip.text.y.left = element_text(angle=0)) + | |
guides(fill = guide_legend(reverse=T, nrow = 1)) + | |
theme(axis.title.y=element_blank(), axis.text.y=element_blank(), axis.ticks.y=element_blank()) + | |
theme(panel.grid.minor.y=element_blank(), panel.grid.major.y=element_blank()) + | |
geom_text(aes(label = ifelse(pct >.03, paste0(lab*100, '%'), '')), position = position_stack(vjust = 0.5), size = 4, family = "font", color = "black") + | |
geom_text(aes(label = ifelse(att == "Less than Once a Year", paste0(lab*100, '%'), '')), position = position_stack(vjust = 0.5), size = 4, family = "font", color = "white") + | |
# geom_text(aes(label = ifelse(att == "12+/Yr.", paste0(lab*100, '%'), '')), position = position_stack(vjust = 0.5), size = 4, family = "font", color = "white") + | |
theme(plot.title = element_text(size = 16)) + | |
theme(strip.text.y.left = element_text(angle = 0, hjust = 1)) + | |
labs(x = "", y = "", title = "Religious Attendance in the 15 Largest Metro Areas", caption = "@ryanburge\nData: Household Pulse Survey, 2024") | |
save("att_msa_pulse.png", wd = 9, ht = 6) | |
gg <- hp %>% | |
mutate(msa = est_msa) %>% | |
mutate(msa = frcode(msa == 35620 ~ "NYC", | |
msa == 31080 ~ "Charlotte", | |
msa == 16980 ~ "Chicago", | |
msa == 19100 ~ "Dallas", | |
msa == 26420 ~ "Houston", | |
msa == 47900 ~ "DC", | |
msa == 33100 ~ "Miami", | |
msa == 37980 ~ "Philadelphia", | |
msa == 12060 ~ "Atlanta", | |
msa == 38060 ~ "Phoenix", | |
msa == 14460 ~ "Boston", | |
msa == 41860 ~ "San Francisco", | |
msa == 40140 ~ "Riverside", | |
msa == 19820 ~ "Detroit", | |
msa == 42660 ~ "Seattle")) %>% | |
mutate(educ = eeduc) %>% | |
mutate(educ = frcode(educ == 1 | educ == 2 | educ == 3 ~ "HS\nor\nLess", | |
educ == 4 | educ == 5 ~ "Some\nColl.", | |
educ == 6 ~ "4 Yr.", | |
educ == 7 ~ "Grad\nSchool")) %>% | |
mutate(att = support3) %>% | |
mutate(att = case_when(att == 1 ~ 1, | |
att <= 4 ~ 0)) %>% | |
group_by(msa, educ) %>% | |
mean_ci(att, wt = pweight, ci = .84) %>% filter(msa != "NA") | |
gg %>% | |
mutate(lab = round(mean, 2)) %>% | |
ggplot(., aes(x = educ, y = mean, fill = msa, group = msa)) + | |
geom_col(color = "black") + | |
facet_wrap(~ msa) + | |
theme_rb() + | |
error_bar() + | |
y_pct() + | |
scale_fill_manual(values = c(moma.colors("Warhol", 15))) + | |
lab_bar(top = FALSE, pos = .065, sz = 4, type = lab) + | |
geom_text(aes(y = .065, label = ifelse(msa == "NYC" | msa == "DC" | msa == "Boston", paste0(lab*100, '%'), "")), position = position_dodge(width = .9), size = 4, family = "font", color = "white") + | |
geom_text(aes(y = .065, label = ifelse(msa == "San Francisco" | msa == "Riverside" | msa == "Seattle", paste0(lab*100, '%'), "")), position = position_dodge(width = .9), size = 4, family = "font", color = "white") + | |
labs(x = "Highest Level of Education", y = "", title = "Share Never Attending Religious Services", caption = "@ryanburge\nData: Household Pulse Survey, 2024") | |
save("msa_ed_never_att.png", wd = 6, ht = 8) | |
gg <- hp %>% | |
filter(income > 0) %>% | |
mutate(att = support3) %>% | |
mutate(att = case_when(att == 1 ~ 1, | |
att <= 4 ~ 0)) %>% | |
mutate(msa = est_msa) %>% | |
mutate(msa = frcode(msa == 35620 ~ "NYC", | |
msa == 31080 ~ "Charlotte", | |
msa == 16980 ~ "Chicago", | |
msa == 19100 ~ "Dallas", | |
msa == 26420 ~ "Houston", | |
msa == 47900 ~ "DC", | |
msa == 33100 ~ "Miami", | |
msa == 37980 ~ "Philadelphia", | |
msa == 12060 ~ "Atlanta", | |
msa == 38060 ~ "Phoenix", | |
msa == 14460 ~ "Boston", | |
msa == 41860 ~ "San Francisco", | |
msa == 40140 ~ "Riverside", | |
msa == 19820 ~ "Detroit", | |
msa == 42660 ~ "Seattle")) %>% | |
group_by(msa, income) %>% | |
mean_ci(att, wt = pweight, ci = .84) %>% filter(msa != "NA") | |
gg %>% | |
mutate(lab = round(mean, 2)) %>% | |
ggplot(., aes(x = income, y = mean, color = msa, group = msa)) + | |
error_bar() + | |
geom_point(stroke = .75, shape = 21, fill = "white") + | |
geom_smooth(se = FALSE, method = lm, aes(color = msa)) + | |
facet_wrap(~ msa) + | |
theme_rb() + | |
y_pct() + | |
scale_x_continuous(breaks = c(2,4,6,8), labels = c("$25K", "$50K", "$100K", "$200K")) + | |
scale_color_manual(values = c(moma.colors("Warhol", 15))) + | |
theme(axis.text.x = element_text(size = 8)) + | |
# lab_bar(top = FALSE, pos = .065, sz = 4, type = lab) + | |
# geom_text(aes(y = .065, label = ifelse(msa == "NYC" | msa == "DC" | msa == "Boston", paste0(lab*100, '%'), "")), position = position_dodge(width = .9), size = 4, family = "font", color = "white") + | |
# geom_text(aes(y = .065, label = ifelse(msa == "San Francisco" | msa == "Riverside" | msa == "Seattle", paste0(lab*100, '%'), "")), position = position_dodge(width = .9), size = 4, family = "font", color = "white") + | |
labs(x = "Household Income", y = "", title = "Share Never Attending Religious Services", caption = "@ryanburge\nData: Household Pulse Survey, 2024") | |
save("msa_income_never_att.png", wd = 6, ht = 8) | |
gg <- hp %>% | |
mutate(hisp = rhispanic) %>% | |
mutate(hisp = case_when(hisp == 2 ~ 1, | |
hisp == 1 ~ 0)) %>% | |
mutate(race = rrace) %>% | |
mutate(race = frcode(race == 1 & hisp == 0 ~ "White", | |
race == 2 ~ "Black", | |
hisp == 1 ~ "Hisp.", | |
race == 3 ~ "Asian")) %>% | |
mutate(att = support3) %>% | |
mutate(att = case_when(att == 1 ~ 1, | |
att <= 4 ~ 0)) %>% | |
mutate(msa = est_msa) %>% | |
mutate(msa = frcode(msa == 35620 ~ "NYC", | |
msa == 31080 ~ "Charlotte", | |
msa == 16980 ~ "Chicago", | |
msa == 19100 ~ "Dallas", | |
msa == 26420 ~ "Houston", | |
msa == 47900 ~ "DC", | |
msa == 33100 ~ "Miami", | |
msa == 37980 ~ "Philadelphia", | |
msa == 12060 ~ "Atlanta", | |
msa == 38060 ~ "Phoenix", | |
msa == 14460 ~ "Boston", | |
msa == 41860 ~ "San Francisco", | |
msa == 40140 ~ "Riverside", | |
msa == 19820 ~ "Detroit", | |
msa == 42660 ~ "Seattle")) %>% | |
group_by(msa, race) %>% | |
mean_ci(att, wt = pweight, ci = .84) %>% | |
filter(msa != "NA") %>% | |
filter(race != "NA") | |
gg %>% | |
mutate(lab = round(mean, 2)) %>% | |
ggplot(., aes(x = race, y = mean, fill = msa, group = msa)) + | |
geom_col(color = "black") + | |
facet_wrap(~ msa) + | |
theme_rb() + | |
error_bar() + | |
y_pct() + | |
scale_fill_manual(values = c(moma.colors("Warhol", 15))) + | |
lab_bar(top = FALSE, pos = .065, sz = 4, type = lab) + | |
theme(axis.text.x = element_text(size = 8)) + | |
geom_text(aes(y = .065, label = ifelse(msa == "NYC" | msa == "DC" | msa == "Boston", paste0(lab*100, '%'), "")), position = position_dodge(width = .9), size = 4, family = "font", color = "white") + | |
geom_text(aes(y = .065, label = ifelse(msa == "San Francisco" | msa == "Riverside" | msa == "Seattle", paste0(lab*100, '%'), "")), position = position_dodge(width = .9), size = 4, family = "font", color = "white") + | |
labs(x = "Race of Respondent", y = "", title = "Share Never Attending Religious Services", caption = "@ryanburge\nData: Household Pulse Survey, 2024") | |
save("msa_race_never_att.png", wd = 6, ht = 8) | |
gg <- hp %>% | |
mutate(birthyr = tbirth_year) %>% | |
mutate(age = 2024 - birthyr) %>% | |
mutate(age2 = frcode(age <= 35 ~ "18-35", | |
age >= 36 & age <= 49 ~ "36-49", | |
age >= 50 & age <= 64 ~ "50-64", | |
age >= 65 ~ "65+")) %>% | |
mutate(att = support3) %>% | |
mutate(att = case_when(att == 1 ~ 1, | |
att <= 4 ~ 0)) %>% | |
mutate(msa = est_msa) %>% | |
mutate(msa = frcode(msa == 35620 ~ "NYC", | |
msa == 31080 ~ "Charlotte", | |
msa == 16980 ~ "Chicago", | |
msa == 19100 ~ "Dallas", | |
msa == 26420 ~ "Houston", | |
msa == 47900 ~ "DC", | |
msa == 33100 ~ "Miami", | |
msa == 37980 ~ "Philadelphia", | |
msa == 12060 ~ "Atlanta", | |
msa == 38060 ~ "Phoenix", | |
msa == 14460 ~ "Boston", | |
msa == 41860 ~ "San Francisco", | |
msa == 40140 ~ "Riverside", | |
msa == 19820 ~ "Detroit", | |
msa == 42660 ~ "Seattle")) %>% | |
group_by(msa, age2) %>% | |
mean_ci(att, wt = pweight, ci = .84) %>% | |
filter(msa != "NA") %>% | |
filter(age2 != "NA") | |
gg %>% | |
mutate(lab = round(mean, 2)) %>% | |
ggplot(., aes(x = age2, y = mean, fill = msa, group = msa)) + | |
geom_col(color = "black") + | |
facet_wrap(~ msa) + | |
theme_rb() + | |
error_bar() + | |
y_pct() + | |
scale_fill_manual(values = c(moma.colors("Warhol", 15))) + | |
lab_bar(top = FALSE, pos = .065, sz = 4, type = lab) + | |
theme(axis.text.x = element_text(size = 8)) + | |
geom_text(aes(y = .065, label = ifelse(msa == "NYC" | msa == "DC" | msa == "Boston", paste0(lab*100, '%'), "")), position = position_dodge(width = .9), size = 4, family = "font", color = "white") + | |
geom_text(aes(y = .065, label = ifelse(msa == "San Francisco" | msa == "Riverside" | msa == "Seattle", paste0(lab*100, '%'), "")), position = position_dodge(width = .9), size = 4, family = "font", color = "white") + | |
labs(x = "Age of Respondent", y = "", title = "Share Never Attending Religious Services", caption = "@ryanburge\nData: Household Pulse Survey, 2024") | |
save("msa_age2_never_att.png", wd = 6, ht = 8) | |
regg <- hp %>% | |
filter(income > 0) %>% | |
mutate(att = support3) %>% | |
mutate(att = case_when(att == 1 ~ 1, | |
att <= 4 ~ 0)) %>% | |
mutate(msa = est_msa) %>% | |
mutate(msa = frcode(msa == 35620 ~ "NYC", | |
msa == 31080 ~ "Charlotte", | |
msa == 16980 ~ "Chicago", | |
msa == 19100 ~ "Dallas", | |
msa == 26420 ~ "Houston", | |
msa == 47900 ~ "DC", | |
msa == 33100 ~ "Miami", | |
msa == 37980 ~ "Philadelphia", | |
msa == 12060 ~ "Atlanta", | |
msa == 38060 ~ "Phoenix", | |
msa == 14460 ~ "Boston", | |
msa == 41860 ~ "San Francisco", | |
msa == 40140 ~ "Riverside", | |
msa == 19820 ~ "Detroit", | |
msa == 42660 ~ "Seattle")) %>% | |
mutate(birthyr = tbirth_year) %>% | |
mutate(age = 2024 - birthyr) %>% | |
mutate(hisp = rhispanic) %>% | |
mutate(hisp = case_when(hisp == 2 ~ 1, | |
hisp == 1 ~ 0)) %>% | |
mutate(race = rrace) %>% | |
mutate(white = case_when(race == 1 & hisp == 0 ~ 1, | |
TRUE ~ 0)) %>% | |
mutate(educ = eeduc) %>% | |
select(att, msa, age, white, income, educ) | |
list_of_msa_datasets <- split(regg, regg$msa) | |
names(list_of_msa_datasets) <- gsub(" ", "_", names(list_of_msa_datasets)) | |
for(msa in names(list_of_msa_datasets)) { | |
assign(msa, list_of_msa_datasets[[msa]]) | |
} | |
library(jtools) | |
coef_names <- c("Education" = "educ", | |
"Age" = "age", | |
"Income" = "income", | |
"White" = "white") | |
library(patchwork) | |
plots_list <- list() | |
# Loop through each MSA dataset | |
for (msa_name in names(list_of_msa_datasets)) { | |
# Fit the model | |
model <- glm(att ~ age + white + income + educ, family = "binomial", data = list_of_msa_datasets[[msa_name]]) | |
# Create the plot with plot_summs | |
plot <- plot_summs(model, scale = TRUE, robust = "HC3", coefs = coef_names) + | |
theme_rb() + | |
labs(x = "", y = "", title = msa_name) # Use MSA name as the title | |
# Add the plot to the list | |
plots_list[[msa_name]] <- plot | |
} | |
# Combine all plots using patchwork | |
combined_plot <- reduce(plots_list, `+`) | |
# To arrange the plots in a more structured manner, you can specify the layout | |
# For example, to arrange in a 4x4 grid (adjust based on the number of MSAs) | |
combined_plot + plot_layout(ncol = 4, nrow = 4) | |
# Display the combined plot | |
print(combined_plot) | |
combined_plot <- combined_plot + | |
plot_annotation(title = "What Factors Increase the Likelihood of Never Attending?", | |
theme = theme(plot.title = element_text(hjust = 0.5))) | |
ggsave("E://graphs23/patched.png", combined_plot, width = 10, height = 10) | |
regg <- hp %>% | |
filter(income > 0) %>% | |
mutate(att = support3) %>% | |
mutate(att = case_when(att == 1 ~ 1, | |
att <= 4 ~ 0)) %>% | |
mutate(msa = est_msa) %>% | |
mutate(msa = frcode(msa == 35620 ~ "NYC", | |
msa == 31080 ~ "Charlotte", | |
msa == 16980 ~ "Chicago", | |
msa == 19100 ~ "Dallas", | |
msa == 26420 ~ "Houston", | |
msa == 47900 ~ "DC", | |
msa == 33100 ~ "Miami", | |
msa == 37980 ~ "Philadelphia", | |
msa == 12060 ~ "Atlanta", | |
msa == 38060 ~ "Phoenix", | |
msa == 14460 ~ "Boston", | |
msa == 41860 ~ "San Francisco", | |
msa == 40140 ~ "Riverside", | |
msa == 19820 ~ "Detroit", | |
msa == 42660 ~ "Seattle")) %>% | |
mutate(birthyr = tbirth_year) %>% | |
mutate(age = 2024 - birthyr) %>% | |
mutate(hisp = rhispanic) %>% | |
mutate(hisp = case_when(hisp == 2 ~ 1, | |
hisp == 1 ~ 0)) %>% | |
mutate(race = rrace) %>% | |
mutate(race = frcode(race == 1 & hisp == 0 ~ "White", | |
race == 2 ~ "Black", | |
hisp == 1 ~ "Hisp.", | |
race == 3 ~ "Asian")) %>% | |
mutate(educ = eeduc) %>% | |
select(att, msa, age, race, income, educ) %>% | |
filter(age <= 75) %>% | |
filter(msa != "NA") | |
int <- glm(att ~ race*age + income + educ, family = 'binomial', data = regg) | |
library(interactions) | |
graph <- interact_plot(int, pred = age, modx = race, interval = TRUE, int.width = .76) | |
graph + | |
theme_rb(legend = TRUE) + | |
y_pct() + | |
labs(x = "Age", y = "Estimate of Share Never Attending", title = "The Interaction of Race and Age on Never Attending", | |
subtitle = "In Major Metro Areas - Controls for Education and Income", | |
caption = "@ryanburge\nData: Household Pulse Survey, 2024") | |
ggsave("E://graphs23/never_attend_pulse_interact.png", bg = "white", width = 7, height = 6) | |
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment