For information on MOE calculations: https://www.census.gov/content/dam/Census/library/publications/2020/acs/acs_general_handbook_2020_ch08.pdf

load packages

library(sf)
library(tidyverse)
library(tmap)
library(tmaptools)
library(tigris)
library(tidycensus)
library(rmapshaper)
library(matrixStats)
library(SpatialAcc)
library(reticulate)
library(dplyr)
library(tidygeocoder)
library(readxl)
library(units)
library(DBI)

read in ACS variables

# gets race, age, etc data
acs_vars <- function(year, geography)
{
        if (year > 2015)
        {
                dcmdva.ct.demographics <- get_acs(
                        geography = geography,
                        year = year,
                        variables = c(pop_tot = "B01001_001", pop_race_tot = "B02001_001", pop_white = "B02001_002",
                                      pop_black = "B02001_003", pop_native = "B02001_004",
                                      pop_asian = "B02001_005", pop_pacific_islander = "B02001_006",
                                      pop_other = "B02001_007", pop_two_or_more_races = "B02001_008",
                                      pop_eth_tot = "B03003_001", pop_hispanic_or_latino = "B03003_003",
                                      male_under_5 = "B01001_003", male_5_9 = "B01001_004",
                                      male_10_14 = "B01001_005", male_15_17 = "B01001_006",
                                      male_18_19 = "B01001_007", male_20 = "B01001_008",
                                      male_21 = "B01001_009", male_22_24 = "B01001_010",
                                      male_25_29 = "B01001_011", male_30_34 = "B01001_012",
                                      male_35_39 = "B01001_013", male_40_44 = "B01001_014",
                                      male_45_49 = "B01001_015", male_50_54 = "B01001_016",
                                      male_55_59 = "B01001_017", male_60_61 = "B01001_018",
                                      male_62_64 = "B01001_019", male_65_66 = "B01001_020",
                                      male_67_69 = "B01001_021", male_70_74 = "B01001_022",
                                      male_75_79 = "B01001_023", male_80_84 = "B01001_024",
                                      male_over_85 = "B01001_025",
                                      female_under_5 = "B01001_027", female_5_9 = "B01001_028",
                                      female_10_14 = "B01001_029", female_15_17 = "B01001_030",
                                      female_18_19 = "B01001_031", female_20 = "B01001_032",
                                      female_21 = "B01001_033", female_22_24 = "B01001_034",
                                      female_25_29 = "B01001_035", female_30_34 = "B01001_036",
                                      female_35_39 = "B01001_037", female_40_44 = "B01001_038",
                                      female_45_49 = "B01001_039", female_50_54 = "B01001_040",
                                      female_55_59 = "B01001_041", female_60_61 = "B01001_042",
                                      female_62_64 = "B01001_043", female_65_66 = "B01001_044",
                                      female_67_69 = "B01001_045", female_70_74 = "B01001_046",
                                      female_75_79 = "B01001_047", female_80_84 = "B01001_048",
                                      female_over_85 = "B01001_049",
                                      vet_denom = "B21001_001", num_vet = "B21001_002", total_hh = "C16002_001",
                                      limited_english_1 = "C16002_004", limited_english_2 = "C16002_007",
                                      limited_english_3 = "C16002_010", limited_english_4 = "C16002_013"),
                        state = c("VA", "MD", "DC"),
                        survey = "acs5",
                        output = "wide",
                        geometry = TRUE) # need for density calculation
        }
        else
        {
                dcmdva.ct.demographics <- get_acs(
                geography = geography,
                year = year,
                variables = c(pop_tot = "B01001_001", pop_race_tot = "B02001_001", pop_white = "B02001_002",
                              pop_black = "B02001_003", pop_native = "B02001_004",
                              pop_asian = "B02001_005", pop_pacific_islander = "B02001_006",
                              pop_other = "B02001_007", pop_two_or_more_races = "B02001_008",
                              pop_eth_tot = "B03001_001", pop_hispanic_or_latino = "B03001_003",
                              male_under_5 = "B01001_003", male_5_9 = "B01001_004",
                              male_10_14 = "B01001_005", male_15_17 = "B01001_006",
                              male_18_19 = "B01001_007", male_20 = "B01001_008",
                              male_21 = "B01001_009", male_22_24 = "B01001_010",
                              male_25_29 = "B01001_011", male_30_34 = "B01001_012",
                              male_35_39 = "B01001_013", male_40_44 = "B01001_014",
                              male_45_49 = "B01001_015", male_50_54 = "B01001_016",
                              male_55_59 = "B01001_017", male_60_61 = "B01001_018",
                              male_62_64 = "B01001_019", male_65_66 = "B01001_020",
                              male_67_69 = "B01001_021", male_70_74 = "B01001_022",
                              male_75_79 = "B01001_023", male_80_84 = "B01001_024",
                              male_over_85 = "B01001_025",
                              female_under_5 = "B01001_027", female_5_9 = "B01001_028",
                              female_10_14 = "B01001_029", female_15_17 = "B01001_030",
                              female_18_19 = "B01001_031", female_20 = "B01001_032",
                              female_21 = "B01001_033", female_22_24 = "B01001_034",
                              female_25_29 = "B01001_035", female_30_34 = "B01001_036",
                              female_35_39 = "B01001_037", female_40_44 = "B01001_038",
                              female_45_49 = "B01001_039", female_50_54 = "B01001_040",
                              female_55_59 = "B01001_041", female_60_61 = "B01001_042",
                              female_62_64 = "B01001_043", female_65_66 = "B01001_044",
                              female_67_69 = "B01001_045", female_70_74 = "B01001_046",
                              female_75_79 = "B01001_047", female_80_84 = "B01001_048",
                              female_over_85 = "B01001_049",
                              vet_denom = "B21001_001", num_vet = "B21001_002"),
                state = c("VA", "MD", "DC"),
                survey = "acs5",
                output = "wide",
                geometry = TRUE) # need for density calculation
        }
        return(dcmdva.ct.demographics)
}
# counties
dcmdva.ct.demographics.2019 <- acs_vars(2019, "county") %>% mutate(year = 2019)
dcmdva.ct.demographics.2018 <- acs_vars(2018, "county") %>% mutate(year = 2018)
dcmdva.ct.demographics.2017 <- acs_vars(2017, "county") %>% mutate(year = 2017)
dcmdva.ct.demographics.2016 <- acs_vars(2016, "county") %>% mutate(year = 2016)
dcmdva.ct.demographics.2015 <- acs_vars(2015, "county") %>% mutate(year = 2015)
dcmdva.ct.demographics.2014 <- acs_vars(2014, "county") %>% mutate(year = 2014)
dcmdva.ct.demographics.2013 <- acs_vars(2013, "county") %>% mutate(year = 2013)
dcmdva.ct.demographics.2012 <- acs_vars(2012, "county") %>% mutate(year = 2012)
dcmdva.ct.demographics.2011 <- acs_vars(2011, "county") %>% mutate(year = 2011)
dcmdva.ct.demographics.2010 <- acs_vars(2010, "county") %>% mutate(year = 2010)
dcmdva.ct.demographics.2009 <- acs_vars(2009, "county") %>% mutate(year = 2009)
# tracts
dcmdva.tr.demographics.2019 <- acs_vars(2019, "tract") %>% mutate(year = 2019)
dcmdva.tr.demographics.2018 <- acs_vars(2018, "tract") %>% mutate(year = 2018)
dcmdva.tr.demographics.2017 <- acs_vars(2017, "tract") %>% mutate(year = 2017)
dcmdva.tr.demographics.2016 <- acs_vars(2016, "tract") %>% mutate(year = 2016)
dcmdva.tr.demographics.2015 <- acs_vars(2015, "tract") %>% mutate(year = 2015)
dcmdva.tr.demographics.2014 <- acs_vars(2014, "tract") %>% mutate(year = 2014)
dcmdva.tr.demographics.2013 <- acs_vars(2013, "tract") %>% mutate(year = 2013)
dcmdva.tr.demographics.2012 <- acs_vars(2012, "tract") %>% mutate(year = 2012)
dcmdva.tr.demographics.2011 <- acs_vars(2011, "tract") %>% mutate(year = 2011)
dcmdva.tr.demographics.2010 <- acs_vars(2010, "tract") %>% mutate(year = 2010)
dcmdva.tr.demographics.2009 <- acs_vars(2009, "tract") %>% mutate(year = 2009)
# block groups
dcmdva.bg.demographics.2019 <- acs_vars(2019, "block group") %>% mutate(year = 2019)
dcmdva.bg.demographics.2018 <- acs_vars(2018, "block group") %>% mutate(year = 2018)
dcmdva.bg.demographics.2017 <- acs_vars(2017, "block group") %>% mutate(year = 2017)
dcmdva.bg.demographics.2016 <- acs_vars(2016, "block group") %>% mutate(year = 2016)
dcmdva.bg.demographics.2015 <- acs_vars(2015, "block group") %>% mutate(year = 2015)
dcmdva.bg.demographics.2014 <- acs_vars(2014, "block group") %>% mutate(year = 2014)
dcmdva.bg.demographics.2013 <- acs_vars(2013, "block group") %>% mutate(year = 2013)
### No Block Group measurements from 2012-2009
# dcmdva.bg.demographics.2012 <- acs_vars(2012, "block group") %>% mutate(year = 2012) # no bg data for 2012 - 2009
# dcmdva.bg.demographics.2011 <- acs_vars(2011, "block group") %>% mutate(year = 2011)
# dcmdva.bg.demographics.2010 <- acs_vars(2010, "block group") %>% mutate(year = 2010)
# dcmdva.bg.demographics.2009 <- acs_vars(2009, "block group") %>% mutate(year = 2009)

County level - might be some NAs that need handling

# split into a and b because of the different columns across the years (more variables b/w 2019-2016)
dcmdva.ct.demo.a <- rbind(dcmdva.ct.demographics.2019, dcmdva.ct.demographics.2018,
                          dcmdva.ct.demographics.2017, dcmdva.ct.demographics.2016)
dcmdva.ct.demo.a$density <- dcmdva.ct.demo.a$pop_totE / (dcmdva.ct.demo.a %>% st_area() %>% set_units(mi^2))
dcmdva.ct.demo.a.2 <- st_drop_geometry(dcmdva.ct.demo.a) %>%
        mutate(pop_total = pop_totE, pop_hispanic_or_latino = pop_hispanic_or_latinoE,
               perc_hispanic_or_latino = pop_hispanic_or_latinoE/pop_eth_totE, pop_white = pop_whiteE,
               perc_white = pop_whiteE/pop_race_totE, pop_black = pop_blackE, perc_black = pop_blackE/pop_race_totE,
               pop_native = pop_nativeE, perc_native = pop_nativeE/pop_race_totE,
               pop_AAPI = pop_asianE + pop_pacific_islanderE, perc_AAPI = (pop_asianE + pop_pacific_islanderE)/pop_race_totE,
               pop_other = pop_otherE, perc_other = pop_otherE/pop_race_totE, pop_two_or_more_races = pop_two_or_more_racesE,
               perc_two_or_more_races = pop_two_or_more_racesE/pop_race_totE,
               pop_under_20 = male_under_5E + male_5_9E + male_10_14E + male_15_17E + male_18_19E + 
                       female_under_5E + female_5_9E + female_10_14E + female_15_17E + female_18_19E,
               perc_under_20 = (male_under_5E + male_5_9E + male_10_14E + male_15_17E + male_18_19E + 
                       female_under_5E + female_5_9E + female_10_14E + female_15_17E + female_18_19E)/pop_totE,
               pop_20_64 = male_20E + male_21E + male_22_24E + male_25_29E + male_30_34E + male_35_39E +
                       male_40_44E +  male_45_49E + male_50_54E + male_55_59E + male_60_61E + male_62_64E + 
                       female_20E + female_21E + female_22_24E + female_25_29E + female_30_34E + female_35_39E +
                       female_40_44E +  female_45_49E + female_50_54E + female_55_59E + female_60_61E + female_62_64E,
               perc_20_64 = (male_20E + male_21E + male_22_24E + male_25_29E + male_30_34E + male_35_39E +
                       male_40_44E +  male_45_49E + male_50_54E + male_55_59E + male_60_61E + male_62_64E + 
                       female_20E + female_21E + female_22_24E + female_25_29E + female_30_34E + female_35_39E +
                       female_40_44E +  female_45_49E + female_50_54E + female_55_59E + female_60_61E + female_62_64E)/pop_totE,
               pop_65_plus = male_65_66E + male_67_69E + male_70_74E + male_75_79E + male_80_84E + male_over_85E +
                       female_65_66E + female_67_69E + female_70_74E + female_75_79E + female_80_84E + male_over_85E,
               perc_65_plus = (male_65_66E + male_67_69E + male_70_74E + male_75_79E + male_80_84E + male_over_85E +
                       female_65_66E + female_67_69E + female_70_74E + female_75_79E + female_80_84E + male_over_85E)/pop_totE,
               pop_veteran = num_vetE, perc_veteran = num_vetE/vet_denomE,
               hh_limited_english = limited_english_1E + limited_english_2E + limited_english_3E + limited_english_4E,
               perc_hh_limited_english = (limited_english_1E + limited_english_2E + limited_english_3E + limited_english_4E)/total_hhE,
               pop_white_MOE = pop_whiteM, pop_black_MOE = pop_blackM, pop_native_MOE = pop_nativeM,
               pop_AAPI_MOE = sqrt(pop_asianM^2 + pop_pacific_islanderM^2), pop_other_MOE = pop_otherM, pop_two_or_more_races_MOE = pop_two_or_more_racesM,
               pop_hispanic_or_latino_MOE = pop_hispanic_or_latinoM,
               pop_under_20_MOE = sqrt(sum(male_under_5M^2, male_5_9M^2, male_10_14M^2, male_15_17M^2, male_18_19M^2, female_under_5M^2,
                                            female_5_9M^2, female_10_14M^2, female_15_17M^2,  female_18_19M^2, na.rm=T)),
               pop_20_64_MOE = sqrt(sum(male_20M^2, male_21M^2, male_22_24M^2, male_25_29M^2, male_30_34M^2,  male_35_39M^2, male_40_44M^2,
                                            male_45_49M^2, male_50_54M^2, male_55_59M^2, male_60_61M^2, male_62_64M^2 + female_20M^2 + 
                                            female_21M^2, female_22_24M^2, female_25_29M^2, female_30_34M^2, female_35_39M^2, female_40_44M^2,
                                            female_45_49M^2, female_50_54M^2, female_55_59M^2, female_60_61M^2, female_62_64M^2, na.rm=T)),
               pop_65_plus_MOE = sqrt(sum(male_65_66M^2, male_67_69M^2, male_70_74M^2, male_75_79M^2,
                                          male_80_84M^2, male_over_85M^2, female_65_66M^2, female_67_69M^2,
                                          female_70_74M^2, female_75_79M^2, female_80_84M^2, male_over_85M^2, na.rm=T)), num_vet_MOE = num_vetM,
               hh_limited_english_MOE = sqrt(sum(limited_english_1M^2, limited_english_2M^2,
                                                 limited_english_3M^2, limited_english_4M^2, na.rm=T))) %>%
        mutate(perc_under_20_MOE = pop_under_20_MOE / pop_totE, perc_20_64_MOE = pop_20_64_MOE / pop_totE,  perc_65_plus_MOE = pop_65_plus_MOE / pop_totE,
               perc_white_MOE = pop_white_MOE / pop_race_totE, perc_black_MOE = pop_black_MOE / pop_race_totE,
               perc_native_MOE = pop_native_MOE / pop_race_totE, perc_AAPI_MOE = pop_AAPI_MOE / pop_race_totE,
               perc_hispanic_or_latino_MOE = pop_hispanic_or_latino_MOE / pop_eth_totE,
               perc_other_MOE = pop_other_MOE / pop_race_totE, perc_two_or_more_races_MOE = pop_two_or_more_races_MOE / pop_race_totE,
               perc_veteran_MOE = sqrt(num_vet_MOE^2 - (perc_veteran * vet_denomM^2)) / vet_denomE,
               perc_hh_limited_english_MOE = sqrt(hh_limited_english_MOE^2 - (perc_hh_limited_english * total_hhM^2)) / total_hhE) %>%
        select(GEOID, NAME, year, pop_total, pop_hispanic_or_latino, perc_hispanic_or_latino, pop_white, perc_white, pop_black,
               perc_black, pop_native, perc_native, pop_AAPI, perc_AAPI, pop_other, perc_other, pop_two_or_more_races, perc_two_or_more_races,
               pop_under_20, perc_under_20, pop_20_64, perc_20_64, pop_65_plus, perc_65_plus, density, pop_veteran, perc_veteran,
               hh_limited_english, perc_hh_limited_english, pop_white_MOE, perc_white_MOE, pop_black_MOE, perc_black_MOE, pop_native_MOE,
               perc_native_MOE, pop_AAPI_MOE, perc_AAPI_MOE, pop_other_MOE, perc_other_MOE, pop_two_or_more_races_MOE,
               perc_two_or_more_races_MOE, pop_hispanic_or_latino_MOE, perc_hispanic_or_latino_MOE,
               pop_under_20_MOE, perc_under_20_MOE, pop_20_64_MOE, perc_20_64_MOE, pop_65_plus_MOE, perc_65_plus_MOE,
               num_vet_MOE, perc_veteran_MOE, hh_limited_english_MOE, perc_hh_limited_english_MOE) %>%
        gather(measure, value, -c(GEOID, NAME, year)) %>%
        rename(geoid = GEOID, region_name = NAME) %>%
        mutate(region_type = "county",
               measure_type = case_when(grepl('MOE$', measure) ~ "MOE",
                                        grepl('^perc', measure) ~ "percent",
                                        measure == "density" ~ "density",
                                        TRUE ~ "count"),
               measure_units = ifelse(measure == "density", "1/mi^2", NA)) %>%
        relocate("geoid", "region_type", "region_name", "year", "measure", "value", "measure_type", "measure_units") %>%
        drop_na(value)
dcmdva.ct.demo.b <- rbind(dcmdva.ct.demographics.2015, dcmdva.ct.demographics.2014, dcmdva.ct.demographics.2013, dcmdva.ct.demographics.2012,
                          dcmdva.ct.demographics.2011, dcmdva.ct.demographics.2010, dcmdva.ct.demographics.2009)
dcmdva.ct.demo.b$density <- dcmdva.ct.demo.b$pop_totE / (dcmdva.ct.demo.b %>% st_area() %>% set_units(mi^2))
dcmdva.ct.demo.b.2 <- st_drop_geometry(dcmdva.ct.demo.b) %>%
        mutate(pop_total = pop_totE, pop_hispanic_or_latino = pop_hispanic_or_latinoE,
               perc_hispanic_or_latino = pop_hispanic_or_latinoE/pop_eth_totE, pop_white = pop_whiteE,
               perc_white = pop_whiteE/pop_race_totE, pop_black = pop_blackE, perc_black = pop_blackE/pop_race_totE,
               pop_native = pop_nativeE, perc_native = pop_nativeE/pop_race_totE,
               pop_AAPI = pop_asianE + pop_pacific_islanderE, perc_AAPI = (pop_asianE + pop_pacific_islanderE)/pop_race_totE,
               pop_other = pop_otherE, perc_other = pop_otherE/pop_race_totE, pop_two_or_more_races = pop_two_or_more_racesE,
               perc_two_or_more_races = pop_two_or_more_racesE/pop_race_totE,
               pop_under_20 = male_under_5E + male_5_9E + male_10_14E + male_15_17E + male_18_19E + 
                       female_under_5E + female_5_9E + female_10_14E + female_15_17E + female_18_19E,
               perc_under_20 = (male_under_5E + male_5_9E + male_10_14E + male_15_17E + male_18_19E + 
                       female_under_5E + female_5_9E + female_10_14E + female_15_17E + female_18_19E)/pop_totE,
               pop_20_64 = male_20E + male_21E + male_22_24E + male_25_29E + male_30_34E + male_35_39E +
                       male_40_44E +  male_45_49E + male_50_54E + male_55_59E + male_60_61E + male_62_64E + 
                       female_20E + female_21E + female_22_24E + female_25_29E + female_30_34E + female_35_39E +
                       female_40_44E +  female_45_49E + female_50_54E + female_55_59E + female_60_61E + female_62_64E,
               perc_20_64 = (male_20E + male_21E + male_22_24E + male_25_29E + male_30_34E + male_35_39E +
                       male_40_44E +  male_45_49E + male_50_54E + male_55_59E + male_60_61E + male_62_64E + 
                       female_20E + female_21E + female_22_24E + female_25_29E + female_30_34E + female_35_39E +
                       female_40_44E +  female_45_49E + female_50_54E + female_55_59E + female_60_61E + female_62_64E)/pop_totE,
               pop_65_plus = male_65_66E + male_67_69E + male_70_74E + male_75_79E + male_80_84E + male_over_85E +
                       female_65_66E + female_67_69E + female_70_74E + female_75_79E + female_80_84E + male_over_85E,
               perc_65_plus = (male_65_66E + male_67_69E + male_70_74E + male_75_79E + male_80_84E + male_over_85E +
                       female_65_66E + female_67_69E + female_70_74E + female_75_79E + female_80_84E + male_over_85E)/pop_totE,
               pop_veteran = num_vetE, perc_veteran = num_vetE/vet_denomE,
               pop_white_MOE = pop_whiteM, pop_black_MOE = pop_blackM, pop_native_MOE = pop_nativeM,
               pop_AAPI_MOE = sqrt(pop_asianM^2 + pop_pacific_islanderM^2), pop_other_MOE = pop_otherM, pop_two_or_more_races_MOE = pop_two_or_more_racesM,
               pop_hispanic_or_latino_MOE = pop_hispanic_or_latinoM,
               pop_under_20_MOE = sqrt(male_under_5M^2 + male_5_9M^2 + male_10_14M^2 + male_15_17M^2 +
                                               male_18_19M^2 + female_under_5M^2 + female_5_9M^2 + female_10_14M^2 +
                                               female_15_17M^2 +  female_18_19M^2),
               pop_20_64_MOE = sqrt(male_20M^2 + male_21M^2 + male_22_24M^2 + male_25_29M^2 + male_30_34M^2 + male_35_39M^2 + male_40_44M^2 +
                                            male_45_49M^2 + male_50_54M^2 + male_55_59M^2 + male_60_61M^2 + male_62_64M^2 + female_20M^2 + 
                                            female_21M^2 + female_22_24M^2 + female_25_29M^2 + female_30_34M^2 + female_35_39M^2 + female_40_44M^2 +
                                            female_45_49M^2 + female_50_54M^2 + female_55_59M^2 + female_60_61M^2 + female_62_64M^2),
               pop_65_plus_MOE = sqrt(male_65_66M^2 + male_67_69M^2 + male_70_74M^2 + male_75_79M^2 + male_80_84M^2 + male_over_85M^2 + female_65_66M^2 +
                                          female_67_69M^2 + female_70_74M^2 + female_75_79M^2 + female_80_84M^2 + male_over_85M^2),
               num_vet_MOE = num_vetM) %>%
        mutate(perc_under_20_MOE = pop_under_20_MOE / pop_totE, perc_20_64_MOE = pop_20_64_MOE / pop_totE,  perc_65_plus_MOE = pop_65_plus_MOE / pop_totE,
               perc_white_MOE = pop_white_MOE / pop_race_totE, perc_black_MOE = pop_black_MOE / pop_race_totE,
               perc_native_MOE = pop_native_MOE / pop_race_totE, perc_AAPI_MOE = pop_AAPI_MOE / pop_race_totE,
               perc_other_MOE = pop_other_MOE / pop_race_totE, perc_two_or_more_races_MOE = pop_two_or_more_races_MOE / pop_race_totE,
               perc_veteran_MOE = sqrt(num_vet_MOE^2 - (perc_veteran * vet_denomM^2)) / vet_denomE) %>%
        select(GEOID, NAME, year, pop_total, pop_hispanic_or_latino, perc_hispanic_or_latino, pop_white, perc_white, pop_black,
               perc_black, pop_native, perc_native, pop_AAPI, perc_AAPI, pop_other, perc_other, pop_two_or_more_races, perc_two_or_more_races,
               pop_under_20, perc_under_20, pop_20_64, perc_20_64, pop_65_plus, perc_65_plus, density, pop_veteran, perc_veteran,
               pop_white_MOE, perc_white_MOE, pop_black_MOE, perc_black_MOE, pop_native_MOE, perc_native_MOE,
               pop_AAPI_MOE, perc_AAPI_MOE, pop_other_MOE, perc_other_MOE, pop_two_or_more_races_MOE, perc_two_or_more_races_MOE,
               pop_hispanic_or_latino_MOE, pop_under_20_MOE, perc_under_20_MOE, pop_20_64_MOE, perc_20_64_MOE, pop_65_plus_MOE,
               perc_65_plus_MOE, num_vet_MOE, perc_veteran_MOE) %>%
        gather(measure, value, -c(GEOID, NAME, year)) %>%
        rename(geoid = GEOID, region_name = NAME) %>%
        mutate(region_type = "county",
               measure_type = case_when(grepl('MOE$', measure) ~ "MOE",
                                        grepl('^perc', measure) ~ "percent",
                                        measure == "density" ~ "density",
                                        TRUE ~ "count"),
               measure_units = ifelse(measure == "density", "1/mi^2", NA)) %>%
        relocate("geoid", "region_type", "region_name", "year", "measure", "value", "measure_type", "measure_units") %>%
        drop_na(value)

dcmdva.ct.demo.final <- rbind(dcmdva.ct.demo.a.2, dcmdva.ct.demo.b.2)
ncr.ct.demo.final <- dcmdva.ct.demo.final[dcmdva.ct.demo.final$geoid %in% c("51013", "51059", "51107", "51510", "51600", "51153", "51683", "51685", "51610", "11001", "24031", "24033", "24017", "24021"),]

Tract level

dcmdva.tr.demo.a <- rbind(dcmdva.tr.demographics.2019, dcmdva.tr.demographics.2018,
                          dcmdva.tr.demographics.2017, dcmdva.tr.demographics.2016)
dcmdva.tr.demo.a$density <- dcmdva.tr.demo.a$pop_totE / (dcmdva.tr.demo.a %>% st_area() %>% set_units(mi^2))
dcmdva.tr.demo.a$density_MOE <- dcmdva.tr.demo.a$pop_totM / (dcmdva.tr.demo.a %>% st_area() %>% set_units(mi^2))
dcmdva.tr.demo.a.2 <- st_drop_geometry(dcmdva.tr.demo.a) %>%
        mutate(pop_total = pop_totE, pop_hispanic_or_latino = pop_hispanic_or_latinoE,
               perc_hispanic_or_latino = pop_hispanic_or_latinoE/pop_eth_totE, pop_white = pop_whiteE,
               perc_white = pop_whiteE/pop_race_totE, pop_black = pop_blackE, perc_black = pop_blackE/pop_race_totE,
               pop_native = pop_nativeE, perc_native = pop_nativeE/pop_race_totE,
               pop_AAPI = pop_asianE + pop_pacific_islanderE, perc_AAPI = (pop_asianE + pop_pacific_islanderE)/pop_race_totE,
               pop_other = pop_otherE, perc_other = pop_otherE/pop_race_totE, pop_two_or_more_races = pop_two_or_more_racesE,
               perc_two_or_more_races = pop_two_or_more_racesE/pop_race_totE,
               pop_under_20 = male_under_5E + male_5_9E + male_10_14E + male_15_17E + male_18_19E + 
                       female_under_5E + female_5_9E + female_10_14E + female_15_17E + female_18_19E,
               perc_under_20 = (male_under_5E + male_5_9E + male_10_14E + male_15_17E + male_18_19E + 
                       female_under_5E + female_5_9E + female_10_14E + female_15_17E + female_18_19E)/pop_totE,
               pop_20_64 = male_20E + male_21E + male_22_24E + male_25_29E + male_30_34E + male_35_39E +
                       male_40_44E +  male_45_49E + male_50_54E + male_55_59E + male_60_61E + male_62_64E + 
                       female_20E + female_21E + female_22_24E + female_25_29E + female_30_34E + female_35_39E +
                       female_40_44E +  female_45_49E + female_50_54E + female_55_59E + female_60_61E + female_62_64E,
               perc_20_64 = (male_20E + male_21E + male_22_24E + male_25_29E + male_30_34E + male_35_39E +
                       male_40_44E +  male_45_49E + male_50_54E + male_55_59E + male_60_61E + male_62_64E + 
                       female_20E + female_21E + female_22_24E + female_25_29E + female_30_34E + female_35_39E +
                       female_40_44E +  female_45_49E + female_50_54E + female_55_59E + female_60_61E + female_62_64E)/pop_totE,
               pop_65_plus = male_65_66E + male_67_69E + male_70_74E + male_75_79E + male_80_84E + male_over_85E +
                       female_65_66E + female_67_69E + female_70_74E + female_75_79E + female_80_84E + male_over_85E,
               perc_65_plus = (male_65_66E + male_67_69E + male_70_74E + male_75_79E + male_80_84E + male_over_85E +
                       female_65_66E + female_67_69E + female_70_74E + female_75_79E + female_80_84E + male_over_85E)/pop_totE,
               pop_veteran = num_vetE, perc_veteran = num_vetE/vet_denomE,
               hh_limited_english = limited_english_1E + limited_english_2E + limited_english_3E + limited_english_4E,
               perc_hh_limited_english = (limited_english_1E + limited_english_2E + limited_english_3E + limited_english_4E)/total_hhE,
               pop_white_MOE = pop_whiteM, pop_black_MOE = pop_blackM, pop_native_MOE = pop_nativeM,
               pop_AAPI_MOE = sqrt(pop_asianM^2 + pop_pacific_islanderM^2), pop_other_MOE = pop_otherM, pop_two_or_more_races_MOE = pop_two_or_more_racesM,
               pop_hispanic_or_latino_MOE = pop_hispanic_or_latinoM,
               pop_under_20_MOE = sqrt(sum(male_under_5M^2, male_5_9M^2, male_10_14M^2, male_15_17M^2, male_18_19M^2, female_under_5M^2,
                                            female_5_9M^2, female_10_14M^2, female_15_17M^2,  female_18_19M^2, na.rm=T)),
               pop_20_64_MOE = sqrt(sum(male_20M^2, male_21M^2, male_22_24M^2, male_25_29M^2, male_30_34M^2,  male_35_39M^2, male_40_44M^2,
                                            male_45_49M^2, male_50_54M^2, male_55_59M^2, male_60_61M^2, male_62_64M^2 + female_20M^2 + 
                                            female_21M^2, female_22_24M^2, female_25_29M^2, female_30_34M^2, female_35_39M^2, female_40_44M^2,
                                            female_45_49M^2, female_50_54M^2, female_55_59M^2, female_60_61M^2, female_62_64M^2, na.rm=T)),
               pop_65_plus_MOE = sqrt(sum(male_65_66M^2, male_67_69M^2, male_70_74M^2, male_75_79M^2,
                                          male_80_84M^2, male_over_85M^2, female_65_66M^2, female_67_69M^2,
                                          female_70_74M^2, female_75_79M^2, female_80_84M^2, male_over_85M^2, na.rm=T)), num_vet_MOE = num_vetM,
               hh_limited_english_MOE = sqrt(sum(limited_english_1M^2, limited_english_2M^2,
                                                 limited_english_3M^2, limited_english_4M^2, na.rm=T)), pop_total_MOE = pop_totM) %>%
        mutate(perc_under_20_MOE = sqrt(pop_under_20_MOE^2 - (perc_under_20 * pop_totM)^2) / pop_totE,
               perc_20_64_MOE = sqrt(pop_20_64_MOE^2 - (perc_20_64 * pop_totM)^2) / pop_totE,
               perc_65_plus_MOE = sqrt(pop_65_plus_MOE^2 - (perc_65_plus * pop_totM)^2) / pop_totE,
               perc_white_MOE = sqrt(pop_white_MOE^2 - (perc_white * pop_race_totM)^2) / pop_race_totE,
               perc_black_MOE = sqrt(pop_black_MOE^2 - (perc_black * pop_race_totM)^2) / pop_race_totE,
               perc_native_MOE = sqrt(pop_native_MOE^2 - (perc_native * pop_race_totM)^2) / pop_race_totE,
               perc_AAPI_MOE = sqrt(pop_AAPI_MOE^2 - (perc_AAPI * pop_race_totM)^2) / pop_race_totE,
               perc_hispanic_or_latino_MOE = sqrt(pop_hispanic_or_latino_MOE^2 - (perc_hispanic_or_latino * pop_eth_totM)^2) / pop_eth_totE,
               perc_other_MOE = sqrt(pop_other_MOE^2 - (perc_other * pop_race_totM)^2) / pop_race_totE,
               perc_two_or_more_races_MOE = sqrt(pop_two_or_more_races_MOE^2 - (perc_two_or_more_races * pop_race_totM)^2) / pop_race_totE,
               perc_veteran_MOE = sqrt(num_vet_MOE^2 - (perc_veteran * vet_denomM^2)) / vet_denomE,
               perc_hh_limited_english_MOE = sqrt(hh_limited_english_MOE^2 - (perc_hh_limited_english * total_hhM^2)) / total_hhE) %>%
        select(GEOID, NAME, year, pop_total, pop_hispanic_or_latino, perc_hispanic_or_latino, pop_white, perc_white, pop_black,
               perc_black, pop_native, perc_native, pop_AAPI, perc_AAPI, pop_other, perc_other, pop_two_or_more_races, perc_two_or_more_races,
               pop_under_20, perc_under_20, pop_20_64, perc_20_64, pop_65_plus, perc_65_plus, density, pop_veteran, perc_veteran,
               hh_limited_english, perc_hh_limited_english, pop_white_MOE, perc_white_MOE, pop_black_MOE, perc_black_MOE, pop_native_MOE,
               perc_native_MOE, pop_AAPI_MOE, perc_AAPI_MOE, pop_other_MOE, perc_other_MOE, pop_two_or_more_races_MOE,
               perc_two_or_more_races_MOE, pop_hispanic_or_latino_MOE, perc_hispanic_or_latino_MOE,
               pop_under_20_MOE, perc_under_20_MOE, pop_20_64_MOE, perc_20_64_MOE, pop_65_plus_MOE, perc_65_plus_MOE,
               num_vet_MOE, perc_veteran_MOE, hh_limited_english_MOE, perc_hh_limited_english_MOE, pop_total_MOE, density_MOE) %>%
        gather(measure, value, -c(GEOID, NAME, year)) %>%
        rename(geoid = GEOID, region_name = NAME) %>%
        mutate(region_type = "tract",
               measure_type = case_when(grepl('MOE$', measure) ~ "MOE",
                                        grepl('^perc', measure) ~ "percent",
                                        measure == "density" ~ "density",
                                        TRUE ~ "count"),
               measure_units = ifelse(measure == "density", "1/mi^2", NA)) %>%
        relocate("geoid", "region_type", "region_name", "year", "measure", "value", "measure_type", "measure_units") %>% drop_na(value)
dcmdva.tr.demo.b <- rbind(dcmdva.tr.demographics.2015, dcmdva.tr.demographics.2014, dcmdva.tr.demographics.2013, dcmdva.tr.demographics.2012,
                          dcmdva.tr.demographics.2011, dcmdva.tr.demographics.2010, dcmdva.tr.demographics.2009)
dcmdva.tr.demo.b$density <- dcmdva.tr.demo.b$pop_totE / (dcmdva.tr.demo.b %>% st_area() %>% set_units(mi^2))
dcmdva.tr.demo.b$density_MOE <- dcmdva.tr.demo.b$pop_totM / (dcmdva.tr.demo.b %>% st_area() %>% set_units(mi^2))
dcmdva.tr.demo.b.2 <- st_drop_geometry(dcmdva.tr.demo.b) %>%
        mutate(pop_total = pop_totE, pop_hispanic_or_latino = pop_hispanic_or_latinoE,
               perc_hispanic_or_latino = pop_hispanic_or_latinoE/pop_eth_totE, pop_white = pop_whiteE,
               perc_white = pop_whiteE/pop_race_totE, pop_black = pop_blackE, perc_black = pop_blackE/pop_race_totE,
               pop_native = pop_nativeE, perc_native = pop_nativeE/pop_race_totE,
               pop_AAPI = pop_asianE + pop_pacific_islanderE, perc_AAPI = (pop_asianE + pop_pacific_islanderE)/pop_race_totE,
               pop_other = pop_otherE, perc_other = pop_otherE/pop_race_totE, pop_two_or_more_races = pop_two_or_more_racesE,
               perc_two_or_more_races = pop_two_or_more_racesE/pop_race_totE,
               pop_under_20 = male_under_5E + male_5_9E + male_10_14E + male_15_17E + male_18_19E + 
                       female_under_5E + female_5_9E + female_10_14E + female_15_17E + female_18_19E,
               perc_under_20 = (male_under_5E + male_5_9E + male_10_14E + male_15_17E + male_18_19E + 
                       female_under_5E + female_5_9E + female_10_14E + female_15_17E + female_18_19E)/pop_totE,
               pop_20_64 = male_20E + male_21E + male_22_24E + male_25_29E + male_30_34E + male_35_39E +
                       male_40_44E +  male_45_49E + male_50_54E + male_55_59E + male_60_61E + male_62_64E + 
                       female_20E + female_21E + female_22_24E + female_25_29E + female_30_34E + female_35_39E +
                       female_40_44E +  female_45_49E + female_50_54E + female_55_59E + female_60_61E + female_62_64E,
               perc_20_64 = (male_20E + male_21E + male_22_24E + male_25_29E + male_30_34E + male_35_39E +
                       male_40_44E +  male_45_49E + male_50_54E + male_55_59E + male_60_61E + male_62_64E + 
                       female_20E + female_21E + female_22_24E + female_25_29E + female_30_34E + female_35_39E +
                       female_40_44E +  female_45_49E + female_50_54E + female_55_59E + female_60_61E + female_62_64E)/pop_totE,
               pop_65_plus = male_65_66E + male_67_69E + male_70_74E + male_75_79E + male_80_84E + male_over_85E +
                       female_65_66E + female_67_69E + female_70_74E + female_75_79E + female_80_84E + male_over_85E,
               perc_65_plus = (male_65_66E + male_67_69E + male_70_74E + male_75_79E + male_80_84E + male_over_85E +
                       female_65_66E + female_67_69E + female_70_74E + female_75_79E + female_80_84E + male_over_85E)/pop_totE,
               pop_veteran = num_vetE, perc_veteran = num_vetE/vet_denomE,
               pop_white_MOE = pop_whiteM, pop_black_MOE = pop_blackM, pop_native_MOE = pop_nativeM,
               pop_AAPI_MOE = sqrt(pop_asianM^2 + pop_pacific_islanderM^2), pop_other_MOE = pop_otherM, pop_two_or_more_races_MOE = pop_two_or_more_racesM,
               pop_hispanic_or_latino_MOE = pop_hispanic_or_latinoM,
               pop_under_20_MOE = sqrt(sum(male_under_5M^2, male_5_9M^2, male_10_14M^2, male_15_17M^2, male_18_19M^2, female_under_5M^2,
                                            female_5_9M^2, female_10_14M^2, female_15_17M^2,  female_18_19M^2, na.rm=T)),
               pop_20_64_MOE = sqrt(sum(male_20M^2, male_21M^2, male_22_24M^2, male_25_29M^2, male_30_34M^2,  male_35_39M^2, male_40_44M^2,
                                            male_45_49M^2, male_50_54M^2, male_55_59M^2, male_60_61M^2, male_62_64M^2 + female_20M^2 + 
                                            female_21M^2, female_22_24M^2, female_25_29M^2, female_30_34M^2, female_35_39M^2, female_40_44M^2,
                                            female_45_49M^2, female_50_54M^2, female_55_59M^2, female_60_61M^2, female_62_64M^2, na.rm=T)),
               pop_65_plus_MOE = sqrt(sum(male_65_66M^2, male_67_69M^2, male_70_74M^2, male_75_79M^2,
                                          male_80_84M^2, male_over_85M^2, female_65_66M^2, female_67_69M^2,
                                          female_70_74M^2, female_75_79M^2, female_80_84M^2, male_over_85M^2, na.rm=T)),
               num_vet_MOE = num_vetM, pop_total_MOE = pop_totM) %>%
        mutate(perc_under_20_MOE = sqrt(pop_under_20_MOE^2 - (perc_under_20 * pop_totM)^2) / pop_totE,
               perc_20_64_MOE = sqrt(pop_20_64_MOE^2 - (perc_20_64 * pop_totM)^2) / pop_totE,
               perc_65_plus_MOE = sqrt(pop_65_plus_MOE^2 - (perc_65_plus * pop_totM)^2) / pop_totE,
               perc_white_MOE = sqrt(pop_white_MOE^2 - (perc_white * pop_race_totM)^2) / pop_race_totE,
               perc_black_MOE = sqrt(pop_black_MOE^2 - (perc_black * pop_race_totM)^2) / pop_race_totE,
               perc_native_MOE = sqrt(pop_native_MOE^2 - (perc_native * pop_race_totM)^2) / pop_race_totE,
               perc_AAPI_MOE = sqrt(pop_AAPI_MOE^2 - (perc_AAPI * pop_race_totM)^2) / pop_race_totE,
               perc_hispanic_or_latino_MOE = sqrt(pop_hispanic_or_latino_MOE^2 - (perc_hispanic_or_latino * pop_eth_totM)^2) / pop_eth_totE,
               perc_other_MOE = sqrt(pop_other_MOE^2 - (perc_other * pop_race_totM)^2) / pop_race_totE,
               perc_two_or_more_races_MOE = sqrt(pop_two_or_more_races_MOE^2 - (perc_two_or_more_races * pop_race_totM)^2) / pop_race_totE,
               perc_veteran_MOE = sqrt(num_vet_MOE^2 - (perc_veteran * vet_denomM^2)) / vet_denomE) %>%
        select(GEOID, NAME, year, pop_total, pop_hispanic_or_latino, perc_hispanic_or_latino, pop_white, perc_white, pop_black,
               perc_black, pop_native, perc_native, pop_AAPI, perc_AAPI, pop_other, perc_other, pop_two_or_more_races, perc_two_or_more_races,
               pop_under_20, perc_under_20, pop_20_64, perc_20_64, pop_65_plus, perc_65_plus, density, pop_veteran, perc_veteran,
               pop_white_MOE, perc_white_MOE, pop_black_MOE, perc_black_MOE, pop_native_MOE,
               perc_native_MOE, pop_AAPI_MOE, perc_AAPI_MOE, pop_other_MOE, perc_other_MOE, pop_two_or_more_races_MOE,
               perc_two_or_more_races_MOE, pop_hispanic_or_latino_MOE, perc_hispanic_or_latino_MOE,
               pop_under_20_MOE, perc_under_20_MOE, pop_20_64_MOE, perc_20_64_MOE, pop_65_plus_MOE, perc_65_plus_MOE,
               num_vet_MOE, perc_veteran_MOE, pop_total_MOE, density_MOE) %>%
        gather(measure, value, -c(GEOID, NAME, year)) %>%
        rename(geoid = GEOID, region_name = NAME) %>%
        mutate(region_type = "tract",
               measure_type = case_when(grepl('MOE$', measure) ~ "MOE",
                                        grepl('^perc', measure) ~ "percent",
                                        measure == "density" ~ "density",
                                        TRUE ~ "count"),
               measure_units = ifelse(measure == "density", "1/mi^2", NA)) %>%
        relocate("geoid", "region_type", "region_name", "year", "measure", "value", "measure_type", "measure_units")
dcmdva.tr.demo.final <- rbind(dcmdva.tr.demo.a.2, dcmdva.tr.demo.b.2)
ncr.tr.demo.final <- dcmdva.tr.demo.final[substr(dcmdva.tr.demo.final$geoid, 1, 5) %in% c("51013", "51059", "51107", "51510", "51600", "51153", "51683", "51685", "51610", "11001", "24031", "24033", "24017", "24021"),]

Block Group level

dcmdva.bg.demo.a <- rbind(dcmdva.bg.demographics.2019, dcmdva.bg.demographics.2018,
                          dcmdva.bg.demographics.2017, dcmdva.bg.demographics.2016)
dcmdva.bg.demo.a$density <- dcmdva.bg.demo.a$pop_totE / (dcmdva.bg.demo.a %>% st_area() %>% set_units(mi^2))
dcmdva.bg.demo.a$density_MOE <- dcmdva.bg.demo.a$pop_totM / (dcmdva.bg.demo.a %>% st_area() %>% set_units(mi^2))
dcmdva.bg.demo.a.2 <- st_drop_geometry(dcmdva.bg.demo.a) %>%
        mutate(pop_total = pop_totE, pop_hispanic_or_latino = pop_hispanic_or_latinoE,
               perc_hispanic_or_latino = pop_hispanic_or_latinoE/pop_eth_totE, pop_white = pop_whiteE,
               perc_white = pop_whiteE/pop_race_totE, pop_black = pop_blackE, perc_black = pop_blackE/pop_race_totE,
               pop_native = pop_nativeE, perc_native = pop_nativeE/pop_race_totE,
               pop_AAPI = pop_asianE + pop_pacific_islanderE, perc_AAPI = (pop_asianE + pop_pacific_islanderE)/pop_race_totE,
               pop_other = pop_otherE, perc_other = pop_otherE/pop_race_totE, pop_two_or_more_races = pop_two_or_more_racesE,
               perc_two_or_more_races = pop_two_or_more_racesE/pop_race_totE,
               pop_under_20 = male_under_5E + male_5_9E + male_10_14E + male_15_17E + male_18_19E + 
                       female_under_5E + female_5_9E + female_10_14E + female_15_17E + female_18_19E,
               perc_under_20 = (male_under_5E + male_5_9E + male_10_14E + male_15_17E + male_18_19E + 
                       female_under_5E + female_5_9E + female_10_14E + female_15_17E + female_18_19E)/pop_totE,
               pop_20_64 = male_20E + male_21E + male_22_24E + male_25_29E + male_30_34E + male_35_39E +
                       male_40_44E +  male_45_49E + male_50_54E + male_55_59E + male_60_61E + male_62_64E + 
                       female_20E + female_21E + female_22_24E + female_25_29E + female_30_34E + female_35_39E +
                       female_40_44E +  female_45_49E + female_50_54E + female_55_59E + female_60_61E + female_62_64E,
               perc_20_64 = (male_20E + male_21E + male_22_24E + male_25_29E + male_30_34E + male_35_39E +
                       male_40_44E +  male_45_49E + male_50_54E + male_55_59E + male_60_61E + male_62_64E + 
                       female_20E + female_21E + female_22_24E + female_25_29E + female_30_34E + female_35_39E +
                       female_40_44E +  female_45_49E + female_50_54E + female_55_59E + female_60_61E + female_62_64E)/pop_totE,
               pop_65_plus = male_65_66E + male_67_69E + male_70_74E + male_75_79E + male_80_84E + male_over_85E +
                       female_65_66E + female_67_69E + female_70_74E + female_75_79E + female_80_84E + male_over_85E,
               perc_65_plus = (male_65_66E + male_67_69E + male_70_74E + male_75_79E + male_80_84E + male_over_85E +
                       female_65_66E + female_67_69E + female_70_74E + female_75_79E + female_80_84E + male_over_85E)/pop_totE,
               pop_veteran = num_vetE, perc_veteran = num_vetE/vet_denomE,
               hh_limited_english = limited_english_1E + limited_english_2E + limited_english_3E + limited_english_4E,
               perc_hh_limited_english = (limited_english_1E + limited_english_2E + limited_english_3E + limited_english_4E)/total_hhE,
               pop_white_MOE = pop_whiteM, pop_black_MOE = pop_blackM, pop_native_MOE = pop_nativeM,
               pop_AAPI_MOE = sqrt(pop_asianM^2 + pop_pacific_islanderM^2), pop_other_MOE = pop_otherM, pop_two_or_more_races_MOE = pop_two_or_more_racesM,
               pop_hispanic_or_latino_MOE = pop_hispanic_or_latinoM,
               pop_under_20_MOE = sqrt(sum(male_under_5M^2, male_5_9M^2, male_10_14M^2, male_15_17M^2, male_18_19M^2, female_under_5M^2,
                                            female_5_9M^2, female_10_14M^2, female_15_17M^2,  female_18_19M^2, na.rm=T)),
               pop_20_64_MOE = sqrt(sum(male_20M^2, male_21M^2, male_22_24M^2, male_25_29M^2, male_30_34M^2,  male_35_39M^2, male_40_44M^2,
                                            male_45_49M^2, male_50_54M^2, male_55_59M^2, male_60_61M^2, male_62_64M^2 + female_20M^2 + 
                                            female_21M^2, female_22_24M^2, female_25_29M^2, female_30_34M^2, female_35_39M^2, female_40_44M^2,
                                            female_45_49M^2, female_50_54M^2, female_55_59M^2, female_60_61M^2, female_62_64M^2, na.rm=T)),
               pop_65_plus_MOE = sqrt(sum(male_65_66M^2, male_67_69M^2, male_70_74M^2, male_75_79M^2,
                                          male_80_84M^2, male_over_85M^2, female_65_66M^2, female_67_69M^2,
                                          female_70_74M^2, female_75_79M^2, female_80_84M^2, male_over_85M^2, na.rm=T)), num_vet_MOE = num_vetM,
               hh_limited_english_MOE = sqrt(sum(limited_english_1M^2, limited_english_2M^2,
                                                 limited_english_3M^2, limited_english_4M^2, na.rm=T)), pop_total_MOE = pop_totM) %>%
        mutate(perc_under_20_MOE = sqrt(pop_under_20_MOE^2 - (perc_under_20 * pop_totM)^2) / pop_totE,
               perc_20_64_MOE = sqrt(pop_20_64_MOE^2 - (perc_20_64 * pop_totM)^2) / pop_totE,
               perc_65_plus_MOE = sqrt(pop_65_plus_MOE^2 - (perc_65_plus * pop_totM)^2) / pop_totE,
               perc_white_MOE = sqrt(pop_white_MOE^2 - (perc_white * pop_race_totM)^2) / pop_race_totE,
               perc_black_MOE = sqrt(pop_black_MOE^2 - (perc_black * pop_race_totM)^2) / pop_race_totE,
               perc_native_MOE = sqrt(pop_native_MOE^2 - (perc_native * pop_race_totM)^2) / pop_race_totE,
               perc_AAPI_MOE = sqrt(pop_AAPI_MOE^2 - (perc_AAPI * pop_race_totM)^2) / pop_race_totE,
               perc_hispanic_or_latino_MOE = sqrt(pop_hispanic_or_latino_MOE^2 - (perc_hispanic_or_latino * pop_eth_totM)^2) / pop_eth_totE,
               perc_other_MOE = sqrt(pop_other_MOE^2 - (perc_other * pop_race_totM)^2) / pop_race_totE,
               perc_two_or_more_races_MOE = sqrt(pop_two_or_more_races_MOE^2 - (perc_two_or_more_races * pop_race_totM)^2) / pop_race_totE,
               perc_veteran_MOE = sqrt(num_vet_MOE^2 - (perc_veteran * vet_denomM^2)) / vet_denomE,
               perc_hh_limited_english_MOE = sqrt(hh_limited_english_MOE^2 - (perc_hh_limited_english * total_hhM^2)) / total_hhE) %>%
        select(GEOID, NAME, year, pop_total, pop_hispanic_or_latino, perc_hispanic_or_latino, pop_white, perc_white, pop_black,
               perc_black, pop_native, perc_native, pop_AAPI, perc_AAPI, pop_other, perc_other, pop_two_or_more_races, perc_two_or_more_races,
               pop_under_20, perc_under_20, pop_20_64, perc_20_64, pop_65_plus, perc_65_plus, density, pop_veteran, perc_veteran,
               hh_limited_english, perc_hh_limited_english, pop_white_MOE, perc_white_MOE, pop_black_MOE, perc_black_MOE, pop_native_MOE,
               perc_native_MOE, pop_AAPI_MOE, perc_AAPI_MOE, pop_other_MOE, perc_other_MOE, pop_two_or_more_races_MOE,
               perc_two_or_more_races_MOE, pop_hispanic_or_latino_MOE, perc_hispanic_or_latino_MOE,
               pop_under_20_MOE, perc_under_20_MOE, pop_20_64_MOE, perc_20_64_MOE, pop_65_plus_MOE, perc_65_plus_MOE,
               num_vet_MOE, perc_veteran_MOE, hh_limited_english_MOE, perc_hh_limited_english_MOE, pop_total_MOE, density_MOE) %>%
        gather(measure, value, -c(GEOID, NAME, year)) %>%
        rename(geoid = GEOID, region_name = NAME) %>%
        mutate(region_type = "block group",
               measure_type = case_when(grepl('MOE$', measure) ~ "MOE",
                                        grepl('^perc', measure) ~ "percent",
                                        measure == "density" ~ "density",
                                        TRUE ~ "count"),
               measure_units = ifelse(measure == "density", "1/mi^2", NA)) %>%
        relocate("geoid", "region_type", "region_name", "year", "measure", "value", "measure_type", "measure_units")
dcmdva.bg.demo.b <- rbind(dcmdva.bg.demographics.2015, dcmdva.bg.demographics.2014, dcmdva.bg.demographics.2013)
dcmdva.bg.demo.b$density <- dcmdva.bg.demo.b$pop_totE / (dcmdva.bg.demo.b %>% st_area() %>% set_units(mi^2))
dcmdva.bg.demo.b$density_MOE <- dcmdva.bg.demo.b$pop_totM / (dcmdva.bg.demo.b %>% st_area() %>% set_units(mi^2))
dcmdva.bg.demo.b.2 <- st_drop_geometry(dcmdva.bg.demo.b) %>%
        mutate(pop_total = pop_totE, pop_hispanic_or_latino = pop_hispanic_or_latinoE,
               perc_hispanic_or_latino = pop_hispanic_or_latinoE/pop_eth_totE, pop_white = pop_whiteE,
               perc_white = pop_whiteE/pop_race_totE, pop_black = pop_blackE, perc_black = pop_blackE/pop_race_totE,
               pop_native = pop_nativeE, perc_native = pop_nativeE/pop_race_totE,
               pop_AAPI = pop_asianE + pop_pacific_islanderE, perc_AAPI = (pop_asianE + pop_pacific_islanderE)/pop_race_totE,
               pop_other = pop_otherE, perc_other = pop_otherE/pop_race_totE, pop_two_or_more_races = pop_two_or_more_racesE,
               perc_two_or_more_races = pop_two_or_more_racesE/pop_race_totE,
               pop_under_20 = male_under_5E + male_5_9E + male_10_14E + male_15_17E + male_18_19E + 
                       female_under_5E + female_5_9E + female_10_14E + female_15_17E + female_18_19E,
               perc_under_20 = (male_under_5E + male_5_9E + male_10_14E + male_15_17E + male_18_19E + 
                       female_under_5E + female_5_9E + female_10_14E + female_15_17E + female_18_19E)/pop_totE,
               pop_20_64 = male_20E + male_21E + male_22_24E + male_25_29E + male_30_34E + male_35_39E +
                       male_40_44E +  male_45_49E + male_50_54E + male_55_59E + male_60_61E + male_62_64E + 
                       female_20E + female_21E + female_22_24E + female_25_29E + female_30_34E + female_35_39E +
                       female_40_44E +  female_45_49E + female_50_54E + female_55_59E + female_60_61E + female_62_64E,
               perc_20_64 = (male_20E + male_21E + male_22_24E + male_25_29E + male_30_34E + male_35_39E +
                       male_40_44E +  male_45_49E + male_50_54E + male_55_59E + male_60_61E + male_62_64E + 
                       female_20E + female_21E + female_22_24E + female_25_29E + female_30_34E + female_35_39E +
                       female_40_44E +  female_45_49E + female_50_54E + female_55_59E + female_60_61E + female_62_64E)/pop_totE,
               pop_65_plus = male_65_66E + male_67_69E + male_70_74E + male_75_79E + male_80_84E + male_over_85E +
                       female_65_66E + female_67_69E + female_70_74E + female_75_79E + female_80_84E + male_over_85E,
               perc_65_plus = (male_65_66E + male_67_69E + male_70_74E + male_75_79E + male_80_84E + male_over_85E +
                       female_65_66E + female_67_69E + female_70_74E + female_75_79E + female_80_84E + male_over_85E)/pop_totE,
               pop_veteran = num_vetE, perc_veteran = num_vetE/vet_denomE,
               pop_white_MOE = pop_whiteM, pop_black_MOE = pop_blackM, pop_native_MOE = pop_nativeM,
               pop_AAPI_MOE = sqrt(pop_asianM^2 + pop_pacific_islanderM^2), pop_other_MOE = pop_otherM, pop_two_or_more_races_MOE = pop_two_or_more_racesM,
               pop_hispanic_or_latino_MOE = pop_hispanic_or_latinoM,
               pop_under_20_MOE = sqrt(sum(male_under_5M^2, male_5_9M^2, male_10_14M^2, male_15_17M^2, male_18_19M^2, female_under_5M^2,
                                            female_5_9M^2, female_10_14M^2, female_15_17M^2,  female_18_19M^2, na.rm=T)),
               pop_20_64_MOE = sqrt(sum(male_20M^2, male_21M^2, male_22_24M^2, male_25_29M^2, male_30_34M^2,  male_35_39M^2, male_40_44M^2,
                                            male_45_49M^2, male_50_54M^2, male_55_59M^2, male_60_61M^2, male_62_64M^2 + female_20M^2 + 
                                            female_21M^2, female_22_24M^2, female_25_29M^2, female_30_34M^2, female_35_39M^2, female_40_44M^2,
                                            female_45_49M^2, female_50_54M^2, female_55_59M^2, female_60_61M^2, female_62_64M^2, na.rm=T)),
               pop_65_plus_MOE = sqrt(sum(male_65_66M^2, male_67_69M^2, male_70_74M^2, male_75_79M^2,
                                          male_80_84M^2, male_over_85M^2, female_65_66M^2, female_67_69M^2,
                                          female_70_74M^2, female_75_79M^2, female_80_84M^2, male_over_85M^2, na.rm=T)),
               num_vet_MOE = num_vetM, pop_total_MOE = pop_totM) %>%
        mutate(perc_under_20_MOE = sqrt(pop_under_20_MOE^2 - (perc_under_20 * pop_totM)^2) / pop_totE,
               perc_20_64_MOE = sqrt(pop_20_64_MOE^2 - (perc_20_64 * pop_totM)^2) / pop_totE,
               perc_65_plus_MOE = sqrt(pop_65_plus_MOE^2 - (perc_65_plus * pop_totM)^2) / pop_totE,
               perc_white_MOE = sqrt(pop_white_MOE^2 - (perc_white * pop_race_totM)^2) / pop_race_totE,
               perc_black_MOE = sqrt(pop_black_MOE^2 - (perc_black * pop_race_totM)^2) / pop_race_totE,
               perc_native_MOE = sqrt(pop_native_MOE^2 - (perc_native * pop_race_totM)^2) / pop_race_totE,
               perc_AAPI_MOE = sqrt(pop_AAPI_MOE^2 - (perc_AAPI * pop_race_totM)^2) / pop_race_totE,
               perc_hispanic_or_latino_MOE = sqrt(pop_hispanic_or_latino_MOE^2 - (perc_hispanic_or_latino * pop_eth_totM)^2) / pop_eth_totE,
               perc_other_MOE = sqrt(pop_other_MOE^2 - (perc_other * pop_race_totM)^2) / pop_race_totE,
               perc_two_or_more_races_MOE = sqrt(pop_two_or_more_races_MOE^2 - (perc_two_or_more_races * pop_race_totM)^2) / pop_race_totE,
               perc_veteran_MOE = sqrt(num_vet_MOE^2 - (perc_veteran * vet_denomM^2)) / vet_denomE) %>%
        select(GEOID, NAME, year, pop_total, pop_hispanic_or_latino, perc_hispanic_or_latino, pop_white, perc_white, pop_black,
               perc_black, pop_native, perc_native, pop_AAPI, perc_AAPI, pop_other, perc_other, pop_two_or_more_races, perc_two_or_more_races,
               pop_under_20, perc_under_20, pop_20_64, perc_20_64, pop_65_plus, perc_65_plus, density, pop_veteran, perc_veteran,
               pop_white_MOE, perc_white_MOE, pop_black_MOE, perc_black_MOE, pop_native_MOE,
               perc_native_MOE, pop_AAPI_MOE, perc_AAPI_MOE, pop_other_MOE, perc_other_MOE, pop_two_or_more_races_MOE,
               perc_two_or_more_races_MOE, pop_hispanic_or_latino_MOE, perc_hispanic_or_latino_MOE,
               pop_under_20_MOE, perc_under_20_MOE, pop_20_64_MOE, perc_20_64_MOE, pop_65_plus_MOE, perc_65_plus_MOE,
               num_vet_MOE, perc_veteran_MOE, pop_total_MOE, density_MOE) %>%
        gather(measure, value, -c(GEOID, NAME, year)) %>%
        rename(geoid = GEOID, region_name = NAME) %>%
        mutate(region_type = "block group",
               measure_type = case_when(grepl('MOE$', measure) ~ "MOE",
                                        grepl('^perc', measure) ~ "percent",
                                        measure == "density" ~ "density",
                                        TRUE ~ "count"),
               measure_units = ifelse(measure == "density", "1/mi^2", NA)) %>%
        relocate("geoid", "region_type", "region_name", "year", "measure", "value", "measure_type", "measure_units")
dcmdva.bg.demo.final <- rbind(dcmdva.bg.demo.a.2, dcmdva.bg.demo.b.2)
ncr.bg.demo.final <- dcmdva.bg.demo.final[substr(dcmdva.bg.demo.final$geoid, 1, 5) %in% c("51013", "51059", "51107", "51510", "51600", "51153", "51683", "51685", "51610", "11001", "24031", "24033", "24017", "24021"),]

subset for NCR and send county, tract, and block group data to database

# convert to 0-100 scale
dcmdva.ct.demo.final[grepl("perc", dcmdva.ct.demo.final$measure),]$value <- 100 * dcmdva.ct.demo.final[grepl("perc", dcmdva.ct.demo.final$measure),]$value
ncr.ct.demo.final[grepl("perc", ncr.ct.demo.final$measure),]$value <- 100 * ncr.ct.demo.final[grepl("perc", ncr.ct.demo.final$measure),]$value
dcmdva.tr.demo.final[grepl("perc", dcmdva.tr.demo.final$measure),]$value <- 100 * dcmdva.tr.demo.final[grepl("perc", dcmdva.tr.demo.final$measure),]$value
ncr.tr.demo.final[grepl("perc", ncr.tr.demo.final$measure),]$value <- 100 * ncr.tr.demo.final[grepl("perc", ncr.tr.demo.final$measure),]$value
dcmdva.bg.demo.final[grepl("perc", dcmdva.bg.demo.final$measure),]$value <- 100 * dcmdva.bg.demo.final[grepl("perc", dcmdva.bg.demo.final$measure),]$value
ncr.bg.demo.final[grepl("perc", ncr.bg.demo.final$measure),]$value <- 100 * ncr.bg.demo.final[grepl("perc", ncr.bg.demo.final$measure),]$value
# send to database
source("~/git/vdh/src/helper_functions.R")
con <- get_db_conn()
dc_dbWriteTable(con, "dc_common", "dcmdva_ct_acs_2009_2019_demographics_update", dcmdva.ct.demo.final)
dc_dbWriteTable(con, "dc_common", "ncr_ct_acs_2009_2019_demographics_update", ncr.ct.demo.final)
dc_dbWriteTable(con, "dc_common", "dcmdva_tr_acs_2009_2019_demographics_update", dcmdva.tr.demo.final)
dc_dbWriteTable(con, "dc_common", "ncr_tr_acs_2009_2019_demographics_update", ncr.tr.demo.final)
dc_dbWriteTable(con, "dc_common", "dcmdva_bg_acs_2013_2019_demographics_update", dcmdva.bg.demo.final)
dc_dbWriteTable(con, "dc_common", "ncr_bg_acs_2013_2019_demographics_update", ncr.bg.demo.final)
dbDisconnect(con)

get demographic data at the health district level, first from ACS

acs_vars.hd <- function(year, geography)
{
        if (year > 2015)
        {
                dcmdva.ct.demographics <- get_acs(
                        geography = geography,
                        year = year,
                        variables = c(pop_tot = "B01001_001", pop_race_tot = "B02001_001", pop_white = "B02001_002",
                                      pop_black = "B02001_003", pop_native = "B02001_004",
                                      pop_asian = "B02001_005", pop_pacific_islander = "B02001_006",
                                      pop_other = "B02001_007", pop_two_or_more_races = "B02001_008",
                                      pop_eth_tot = "B03003_001", pop_hispanic_or_latino = "B03003_003",
                                      male_under_5 = "B01001_003", male_5_9 = "B01001_004",
                                      male_10_14 = "B01001_005", male_15_17 = "B01001_006",
                                      male_18_19 = "B01001_007", male_20 = "B01001_008",
                                      male_21 = "B01001_009", male_22_24 = "B01001_010",
                                      male_25_29 = "B01001_011", male_30_34 = "B01001_012",
                                      male_35_39 = "B01001_013", male_40_44 = "B01001_014",
                                      male_45_49 = "B01001_015", male_50_54 = "B01001_016",
                                      male_55_59 = "B01001_017", male_60_61 = "B01001_018",
                                      male_62_64 = "B01001_019", male_65_66 = "B01001_020",
                                      male_67_69 = "B01001_021", male_70_74 = "B01001_022",
                                      male_75_79 = "B01001_023", male_80_84 = "B01001_024",
                                      male_over_85 = "B01001_025",
                                      female_under_5 = "B01001_027", female_5_9 = "B01001_028",
                                      female_10_14 = "B01001_029", female_15_17 = "B01001_030",
                                      female_18_19 = "B01001_031", female_20 = "B01001_032",
                                      female_21 = "B01001_033", female_22_24 = "B01001_034",
                                      female_25_29 = "B01001_035", female_30_34 = "B01001_036",
                                      female_35_39 = "B01001_037", female_40_44 = "B01001_038",
                                      female_45_49 = "B01001_039", female_50_54 = "B01001_040",
                                      female_55_59 = "B01001_041", female_60_61 = "B01001_042",
                                      female_62_64 = "B01001_043", female_65_66 = "B01001_044",
                                      female_67_69 = "B01001_045", female_70_74 = "B01001_046",
                                      female_75_79 = "B01001_047", female_80_84 = "B01001_048",
                                      female_over_85 = "B01001_049",
                                      vet_denom = "B21001_001", num_vet = "B21001_002", total_hh = "C16002_001",
                                      limited_english_1 = "C16002_004", limited_english_2 = "C16002_007",
                                      limited_english_3 = "C16002_010", limited_english_4 = "C16002_013"),
                        state = "VA",
                        survey = "acs5",
                        output = "wide",
                        geometry = FALSE)
        }
        else
        {
                dcmdva.ct.demographics <- get_acs(
                geography = geography,
                year = year,
                variables = c(pop_tot = "B01001_001", pop_race_tot = "B02001_001", pop_white = "B02001_002",
                              pop_black = "B02001_003", pop_native = "B02001_004",
                              pop_asian = "B02001_005", pop_pacific_islander = "B02001_006",
                              pop_other = "B02001_007", pop_two_or_more_races = "B02001_008",
                              pop_eth_tot = "B03001_001", pop_hispanic_or_latino = "B03001_003",
                              male_under_5 = "B01001_003", male_5_9 = "B01001_004",
                              male_10_14 = "B01001_005", male_15_17 = "B01001_006",
                              male_18_19 = "B01001_007", male_20 = "B01001_008",
                              male_21 = "B01001_009", male_22_24 = "B01001_010",
                              male_25_29 = "B01001_011", male_30_34 = "B01001_012",
                              male_35_39 = "B01001_013", male_40_44 = "B01001_014",
                              male_45_49 = "B01001_015", male_50_54 = "B01001_016",
                              male_55_59 = "B01001_017", male_60_61 = "B01001_018",
                              male_62_64 = "B01001_019", male_65_66 = "B01001_020",
                              male_67_69 = "B01001_021", male_70_74 = "B01001_022",
                              male_75_79 = "B01001_023", male_80_84 = "B01001_024",
                              male_over_85 = "B01001_025",
                              female_under_5 = "B01001_027", female_5_9 = "B01001_028",
                              female_10_14 = "B01001_029", female_15_17 = "B01001_030",
                              female_18_19 = "B01001_031", female_20 = "B01001_032",
                              female_21 = "B01001_033", female_22_24 = "B01001_034",
                              female_25_29 = "B01001_035", female_30_34 = "B01001_036",
                              female_35_39 = "B01001_037", female_40_44 = "B01001_038",
                              female_45_49 = "B01001_039", female_50_54 = "B01001_040",
                              female_55_59 = "B01001_041", female_60_61 = "B01001_042",
                              female_62_64 = "B01001_043", female_65_66 = "B01001_044",
                              female_67_69 = "B01001_045", female_70_74 = "B01001_046",
                              female_75_79 = "B01001_047", female_80_84 = "B01001_048",
                              female_over_85 = "B01001_049",
                              vet_denom = "B21001_001", num_vet = "B21001_002"),
                state = "VA",
                survey = "acs5",
                output = "wide",
                geometry = FALSE)
        }
        return(dcmdva.ct.demographics)
}
va.ct.demographics.2019 <- acs_vars.hd(2019, "county") %>% mutate(year = 2019, GEOID = ifelse(GEOID == "51515", "51019", GEOID))
va.ct.demographics.2018 <- acs_vars.hd(2018, "county") %>% mutate(year = 2018, GEOID = ifelse(GEOID == "51515", "51019", GEOID))
va.ct.demographics.2017 <- acs_vars.hd(2017, "county") %>% mutate(year = 2017, GEOID = ifelse(GEOID == "51515", "51019", GEOID))
va.ct.demographics.2016 <- acs_vars.hd(2016, "county") %>% mutate(year = 2016, GEOID = ifelse(GEOID == "51515", "51019", GEOID))
va.ct.demographics.2015 <- acs_vars.hd(2015, "county") %>% mutate(year = 2015, GEOID = ifelse(GEOID == "51515", "51019", GEOID))
va.ct.demographics.2014 <- acs_vars.hd(2014, "county") %>% mutate(year = 2014, GEOID = ifelse(GEOID == "51515", "51019", GEOID))
va.ct.demographics.2013 <- acs_vars.hd(2013, "county") %>% mutate(year = 2013, GEOID = ifelse(GEOID == "51515", "51019", GEOID))
va.ct.demographics.2012 <- acs_vars.hd(2012, "county") %>% mutate(year = 2012, GEOID = ifelse(GEOID == "51515", "51019", GEOID))
va.ct.demographics.2011 <- acs_vars.hd(2011, "county") %>% mutate(year = 2011, GEOID = ifelse(GEOID == "51515", "51019", GEOID))
va.ct.demographics.2010 <- acs_vars.hd(2010, "county") %>% mutate(year = 2010, GEOID = ifelse(GEOID == "51515", "51019", GEOID))
va.ct.demographics.2009 <- acs_vars.hd(2009, "county") %>% mutate(year = 2009, GEOID = ifelse(GEOID == "51515", "51019", GEOID))
health_district <- read.csv("/project/biocomplexity/sdad/projects_data/vdh/va_county_to_hd.csv")
health_district$county_id <- as.character(health_district$county_id)
con <- get_db_conn()
health_district_geoids <- st_read(con, query = "SELECT * FROM dc_geographies.va_hd_vdh_2021_health_district_geo_names")
DBI::dbDisconnect(con)
health_district_2 <- left_join(health_district, health_district_geoids, by = c("health_district" = "region_name"))

get virginia health district data, formatting

va.ct <- get_acs(geography = "county",
                 year = 2019,
                 variables = c(pop_tot = "B01001_001"),
                 state = "VA",
                 survey = "acs5",
                 output = "wide",
                 geometry = T)
health_district_3 <- merge(va.ct, health_district_2, by.x = "GEOID", by.y = "county_id") %>% mutate(area = st_area(geometry)) %>% group_by(geoid, health_district) %>% summarize(area = sum(area)/(1609.34)^2)
va.ct.demo.a <- rbind(va.ct.demographics.2019, va.ct.demographics.2018,
                      va.ct.demographics.2017, va.ct.demographics.2016)
va.hd.demo.a <- merge(va.ct.demo.a, health_district_2[, c("county_id", "health_district", "geoid")], by.x = "GEOID", by.y = "county_id") %>%
        group_by(geoid, health_district, year) %>%
        summarize(pop_total = sum(pop_totE), pop_white = sum(pop_whiteE), pop_black = sum(pop_blackE), pop_native = sum(pop_nativeE),
                  pop_AAPI = sum(pop_asianE) + sum(pop_pacific_islanderE), pop_other = sum(pop_otherE),
                  pop_race_tot = sum(pop_race_totE), pop_eth_tot = sum(pop_eth_totE),
                  pop_two_or_more_races = sum(pop_two_or_more_racesE), pop_hispanic_or_latino = sum(pop_hispanic_or_latinoE),
                  pop_under_20 = sum(male_under_5E + male_5_9E + male_10_14E + male_15_17E + male_18_19E + 
                          female_under_5E + female_5_9E + female_10_14E + female_15_17E + female_18_19E),
                  pop_20_64 = sum(male_20E + male_21E + male_22_24E + male_25_29E + male_30_34E + male_35_39E +
                       male_40_44E +  male_45_49E + male_50_54E + male_55_59E + male_60_61E + male_62_64E + 
                       female_20E + female_21E + female_22_24E + female_25_29E + female_30_34E + female_35_39E +
                       female_40_44E +  female_45_49E + female_50_54E + female_55_59E + female_60_61E + female_62_64E),
                  pop_65_plus = sum(male_65_66E + male_67_69E + male_70_74E + male_75_79E + male_80_84E + male_over_85E +
                       female_65_66E + female_67_69E + female_70_74E + female_75_79E + female_80_84E + male_over_85E),
                  pop_veteran = sum(num_vetE), vet_denom = sum(vet_denomE),
                  hh_limited_english = sum(limited_english_1E + limited_english_2E + limited_english_3E + limited_english_4E),
                  total_hh = sum(total_hhE),
                  pop_white_MOE = sqrt(sum(pop_whiteM^2)), pop_black_MOE = sqrt(sum(pop_blackM^2)), pop_other_MOE = sqrt(sum(pop_otherM^2)),
                  pop_native_MOE = sqrt(sum(pop_nativeM^2)), pop_AAPI_MOE = sqrt(sum(pop_asianM^2 + pop_pacific_islanderM^2)),
                  pop_eth_tot_MOE = sqrt(sum(pop_eth_totM^2, na.rm=T)), pop_two_or_more_races_MOE = sqrt(sum(pop_two_or_more_racesM^2)),
                  pop_hispanic_or_latino_MOE = sqrt(sum(pop_hispanic_or_latinoM^2)),
                  pop_under_20_MOE = sqrt(sum(male_under_5M^2, male_5_9M^2, male_10_14M^2, male_15_17M^2, male_18_19M^2, female_under_5M^2,
                                            female_5_9M^2, female_10_14M^2, female_15_17M^2,  female_18_19M^2, na.rm=T)),
                  pop_20_64_MOE = sqrt(sum(male_20M^2, male_21M^2, male_22_24M^2, male_25_29M^2, male_30_34M^2,  male_35_39M^2, male_40_44M^2,
                                            male_45_49M^2, male_50_54M^2, male_55_59M^2, male_60_61M^2, male_62_64M^2 + female_20M^2 + 
                                            female_21M^2, female_22_24M^2, female_25_29M^2, female_30_34M^2, female_35_39M^2, female_40_44M^2,
                                            female_45_49M^2, female_50_54M^2, female_55_59M^2, female_60_61M^2, female_62_64M^2, na.rm=T)),
                  pop_65_plus_MOE = sqrt(sum(male_65_66M^2, male_67_69M^2, male_70_74M^2, male_75_79M^2,
                                          male_80_84M^2, male_over_85M^2, female_65_66M^2, female_67_69M^2,
                                          female_70_74M^2, female_75_79M^2, female_80_84M^2, male_over_85M^2, na.rm=T)),
                  num_vet_MOE = sqrt(sum(num_vetM^2)), vet_denom_MOE = sqrt(sum(vet_denomM^2)),
                  hh_limited_english_MOE = sqrt(sum(limited_english_1M^2, limited_english_2M^2,
                                                 limited_english_3M^2, limited_english_4M^2, na.rm=T)),
                  total_hh_MOE = sqrt(sum(total_hhM^2))) %>%
        merge(health_district_3, by = c("geoid", "health_district")) %>%
        mutate(perc_hispanic_or_latino = 100*pop_hispanic_or_latino/pop_eth_tot, perc_white = 100*pop_white/pop_race_tot,
               perc_black = 100*pop_black/pop_race_tot, perc_native = 100*pop_native/pop_race_tot, perc_AAPI = 100*pop_AAPI/pop_race_tot,
               perc_other = 100*pop_other/pop_race_tot, perc_two_or_more_races = 100*pop_two_or_more_races/pop_race_tot,
               perc_under_20 = 100*pop_under_20/pop_total, perc_20_64 = 100*pop_20_64/pop_total, perc_65_plus = 100*pop_65_plus/pop_total,
               perc_veteran = 100*pop_veteran/vet_denom, perc_hh_limited_english = 100*hh_limited_english/total_hh,
               density = pop_total / area,
               perc_hispanic_or_latino_MOE = 100*pop_hispanic_or_latino_MOE/pop_eth_tot, perc_white_MOE = 100*pop_white_MOE/pop_race_tot,
               perc_black_MOE = 100*pop_black_MOE/pop_race_tot, perc_native_MOE = 100*pop_native_MOE/pop_race_tot,
               perc_AAPI_MOE = 100*pop_AAPI_MOE/pop_race_tot, perc_other_MOE = 100*pop_other_MOE/pop_race_tot,
               perc_two_or_more_races_MOE = 100*pop_two_or_more_races_MOE/pop_race_tot,
               perc_under_20_MOE = 100*pop_under_20_MOE/pop_total, perc_20_64_MOE = 100*pop_20_64_MOE/pop_total,
               perc_65_plus_MOE = 100*pop_65_plus_MOE/pop_total,
               perc_hh_limited_english_MOE = 100*sqrt(hh_limited_english_MOE^2 - (perc_hh_limited_english * total_hh_MOE)^2) / total_hh,
               perc_veteran_MOE = 100*sqrt(num_vet_MOE^2 - (perc_veteran * vet_denom_MOE)^2) / vet_denom) %>%
        select(geoid, health_district, year, pop_total, pop_hispanic_or_latino, perc_hispanic_or_latino, pop_white, perc_white, pop_black,
               perc_black, pop_native, perc_native, pop_AAPI, perc_AAPI, pop_other, perc_other, pop_two_or_more_races, perc_two_or_more_races,
               pop_under_20, perc_under_20, pop_20_64, perc_20_64, pop_65_plus, perc_65_plus, density, pop_veteran, perc_veteran,
               hh_limited_english, perc_hh_limited_english, pop_white_MOE, perc_white_MOE, pop_black_MOE, perc_black_MOE, pop_native_MOE,
               perc_native_MOE, pop_AAPI_MOE, perc_AAPI_MOE, pop_other_MOE, perc_other_MOE, pop_two_or_more_races_MOE,
               perc_two_or_more_races_MOE, pop_hispanic_or_latino_MOE, perc_hispanic_or_latino_MOE,
               pop_under_20_MOE, perc_under_20_MOE, pop_20_64_MOE, perc_20_64_MOE, pop_65_plus_MOE, perc_65_plus_MOE,
               num_vet_MOE, perc_veteran_MOE, hh_limited_english_MOE, perc_hh_limited_english_MOE) %>%
        gather(measure, value, -c(geoid, health_district, year)) %>%
        rename(region_name = health_district) %>%
        mutate(region_type = "health district",
               measure_type = case_when(grepl('MOE$', measure) ~ "MOE",
                                        grepl('^perc', measure) ~ "percent",
                                        measure == "density" ~ "density",
                                        TRUE ~ "count"),
               measure_units = ifelse(measure == "density", "1/mi^2", NA)) %>%
        relocate("geoid", "region_type", "region_name", "year", "measure", "value", "measure_type", "measure_units") %>%
        drop_na(value)
va.ct.demo.b <- rbind(va.ct.demographics.2015, va.ct.demographics.2014, va.ct.demographics.2013, va.ct.demographics.2012,
                          va.ct.demographics.2011, va.ct.demographics.2010, va.ct.demographics.2009)
va.hd.demo.b <- merge(va.ct.demo.b, health_district_2[, c("county_id", "health_district", "geoid")], by.x = "GEOID", by.y = "county_id") %>%
        group_by(geoid, health_district, year) %>%
        summarize(pop_total = sum(pop_totE), pop_white = sum(pop_whiteE), pop_black = sum(pop_blackE), pop_native = sum(pop_nativeE),
                  pop_AAPI = sum(pop_asianE) + sum(pop_pacific_islanderE), pop_other = sum(pop_otherE),
                  pop_race_tot = sum(pop_race_totE), pop_eth_tot = sum(pop_eth_totE),
                  pop_two_or_more_races = sum(pop_two_or_more_racesE), pop_hispanic_or_latino = sum(pop_hispanic_or_latinoE),
                  pop_under_20 = sum(male_under_5E + male_5_9E + male_10_14E + male_15_17E + male_18_19E + 
                          female_under_5E + female_5_9E + female_10_14E + female_15_17E + female_18_19E),
                  pop_20_64 = sum(male_20E + male_21E + male_22_24E + male_25_29E + male_30_34E + male_35_39E +
                       male_40_44E +  male_45_49E + male_50_54E + male_55_59E + male_60_61E + male_62_64E + 
                       female_20E + female_21E + female_22_24E + female_25_29E + female_30_34E + female_35_39E +
                       female_40_44E +  female_45_49E + female_50_54E + female_55_59E + female_60_61E + female_62_64E),
                  pop_65_plus = sum(male_65_66E + male_67_69E + male_70_74E + male_75_79E + male_80_84E + male_over_85E +
                       female_65_66E + female_67_69E + female_70_74E + female_75_79E + female_80_84E + male_over_85E),
                  pop_veteran = sum(num_vetE), vet_denom = sum(vet_denomE),
                  pop_white_MOE = sqrt(sum(pop_whiteM^2)), pop_black_MOE = sqrt(sum(pop_blackM^2)), pop_other_MOE = sqrt(sum(pop_otherM^2)),
                  pop_native_MOE = sqrt(sum(pop_nativeM^2)), pop_AAPI_MOE = sqrt(sum(pop_asianM^2 + pop_pacific_islanderM^2)),
                  pop_eth_tot_MOE = sqrt(sum(pop_eth_totM^2, na.rm=T)), pop_two_or_more_races_MOE = sqrt(sum(pop_two_or_more_racesM^2)),
                  pop_hispanic_or_latino_MOE = sqrt(sum(pop_hispanic_or_latinoM^2)),
                  pop_under_20_MOE = sqrt(sum(male_under_5M^2, male_5_9M^2, male_10_14M^2, male_15_17M^2, male_18_19M^2, female_under_5M^2,
                                            female_5_9M^2, female_10_14M^2, female_15_17M^2,  female_18_19M^2, na.rm=T)),
                  pop_20_64_MOE = sqrt(sum(male_20M^2, male_21M^2, male_22_24M^2, male_25_29M^2, male_30_34M^2,  male_35_39M^2, male_40_44M^2,
                                            male_45_49M^2, male_50_54M^2, male_55_59M^2, male_60_61M^2, male_62_64M^2 + female_20M^2 + 
                                            female_21M^2, female_22_24M^2, female_25_29M^2, female_30_34M^2, female_35_39M^2, female_40_44M^2,
                                            female_45_49M^2, female_50_54M^2, female_55_59M^2, female_60_61M^2, female_62_64M^2, na.rm=T)),
                  pop_65_plus_MOE = sqrt(sum(male_65_66M^2, male_67_69M^2, male_70_74M^2, male_75_79M^2,
                                          male_80_84M^2, male_over_85M^2, female_65_66M^2, female_67_69M^2,
                                          female_70_74M^2, female_75_79M^2, female_80_84M^2, male_over_85M^2, na.rm=T)),
                  num_vet_MOE = sqrt(sum(num_vetM^2)), vet_denom_MOE = sqrt(sum(vet_denomM^2))) %>%
        merge(health_district_3, by = c("geoid", "health_district")) %>%
        mutate(perc_hispanic_or_latino = 100*pop_hispanic_or_latino/pop_eth_tot, perc_white = 100*pop_white/pop_race_tot,
               perc_black = 100*pop_black/pop_race_tot, perc_native = 100*pop_native/pop_race_tot, perc_AAPI = 100*pop_AAPI/pop_race_tot,
               perc_other = 100*pop_other/pop_race_tot, perc_two_or_more_races = 100*pop_two_or_more_races/pop_race_tot,
               perc_under_20 = 100*pop_under_20/pop_total, perc_20_64 = 100*pop_20_64/pop_total, perc_65_plus = 100*pop_65_plus/pop_total,
               perc_veteran = 100*pop_veteran/vet_denom, density = pop_total / area,
               perc_hispanic_or_latino_MOE = 100*pop_hispanic_or_latino_MOE/pop_eth_tot, perc_white_MOE = 100*pop_white_MOE/pop_race_tot,
               perc_black_MOE = 100*pop_black_MOE/pop_race_tot, perc_native_MOE = 100*pop_native_MOE/pop_race_tot,
               perc_AAPI_MOE = 100*pop_AAPI_MOE/pop_race_tot, perc_other_MOE = 100*pop_other_MOE/pop_race_tot,
               perc_two_or_more_races_MOE = 100*pop_two_or_more_races_MOE/pop_race_tot,
               perc_under_20_MOE = 100*pop_under_20_MOE/pop_total, perc_20_64_MOE = 100*pop_20_64_MOE/pop_total,
               perc_65_plus_MOE = 100*pop_65_plus_MOE/pop_total,
               perc_veteran_MOE = 100*sqrt(num_vet_MOE^2 - (perc_veteran * vet_denom_MOE)^2) / vet_denom) %>%
        select(geoid, health_district, year, pop_total, pop_hispanic_or_latino, perc_hispanic_or_latino, pop_white, perc_white, pop_black,
               perc_black, pop_native, perc_native, pop_AAPI, perc_AAPI, pop_other, perc_other, pop_two_or_more_races, perc_two_or_more_races,
               pop_under_20, perc_under_20, pop_20_64, perc_20_64, pop_65_plus, perc_65_plus, density, pop_veteran, perc_veteran,
               pop_white_MOE, perc_white_MOE, pop_black_MOE, perc_black_MOE, pop_native_MOE,
               perc_native_MOE, pop_AAPI_MOE, perc_AAPI_MOE, pop_other_MOE, perc_other_MOE, pop_two_or_more_races_MOE,
               perc_two_or_more_races_MOE, pop_hispanic_or_latino_MOE, perc_hispanic_or_latino_MOE,
               pop_under_20_MOE, perc_under_20_MOE, pop_20_64_MOE, perc_20_64_MOE, pop_65_plus_MOE, perc_65_plus_MOE,
               num_vet_MOE, perc_veteran_MOE) %>%
        gather(measure, value, -c(geoid, health_district, year)) %>%
        rename(region_name = health_district) %>%
        mutate(region_type = "health district",
               measure_type = case_when(grepl('MOE$', measure) ~ "MOE",
                                        grepl('^perc', measure) ~ "percent",
                                        measure == "density" ~ "density",
                                        TRUE ~ "count"),
               measure_units = ifelse(measure == "density", "1/mi^2", NA)) %>%
        relocate("geoid", "region_type", "region_name", "year", "measure", "value", "measure_type", "measure_units") %>%
        drop_na(value)
va.hd.demo <- rbind(va.hd.demo.a, va.hd.demo.b)

send to database

# send va hd data to db
con <- get_db_conn()
dc_dbWriteTable(con, "dc_common", "va_hd_acs_2009_2019_demographics_update", va.hd.demo)
DBI::dbDisconnect(con)


uva-bi-sdad/dc.synth.dmg documentation built on June 6, 2022, 8:09 p.m.