Get and save to DB the Fairfax tax administration real estate database.

get_db_conn <- function(db_name = "sdad",
               db_host = "postgis1",
               db_port = "5432",
               db_user = Sys.getenv(x = "db_userid"), # requires you to setup environmental vars (above)
               db_pass = Sys.getenv(x = "db_pwd")) {
                   RPostgreSQL::dbConnect(
                  drv = RPostgreSQL::PostgreSQL(),
                  dbname = "sdad",
                  host = "10.250.124.195",
                  port = 5432,
                  user = db_user,
                  password = db_pass)
}

dc_dbWriteTable <-
  function(
    db_conn,
    schema_name,
    table_name,
    table_data,
    table_owner = "data_commons"
  ) {
    # check for geometry/geography columns
    tf <- sapply(table_data, {function(x) inherits(x, 'sfc')})
    # if TRUE, use sf
    if (TRUE %in% tf) {
      sf_write_result <- sf::st_write(obj = table_data, dsn = db_conn, layer = c(schema_name, table_name), row.names = FALSE)
      print(sf_write_result)
      # if FALSE, use DBI
    } else {
      write_result <- DBI::dbWriteTable(conn = db_conn, name = c(schema_name, table_name), value = table_data, row.names = FALSE)
      print(write_result)
    }
    # change table owner
    chgown_result <- DBI::dbSendQuery(conn = db_conn, statement = paste0("ALTER TABLE ", schema_name, ".", table_name, " OWNER TO ", table_owner))
    print(chgown_result)
  }
library(dplyr)
library(sf)

# upload the data
fairfax_housing_units <- read.csv('~/Github/vdh/fairfax_data/Tax_Administration_s_Real_Estate_Fairfax.csv') %>% dplyr::select(PARID,LIVUNIT)
fairfax_parcels_blocks <- readRDS("~/Github/vdh/fairfax_data/fairfax_parcels_blocks_wgs84.RDS") %>% dplyr::select(PIN,GEOID20,HOUSING20,POP20)

# merge the two data
fairfax_parcels_blocks <- merge(fairfax_housing_units, fairfax_parcels_blocks, by.x='PARID', by.y='PIN')

# Remove all NA units
fairfax_parcels_blocks <- fairfax_parcels_blocks[!is.na(fairfax_parcels_blocks$LIVUNIT),]
fairfax_parcels_blocks <- fairfax_parcels_blocks[!is.na(fairfax_parcels_blocks$GEOID20),]

con <- get_db_conn()
dc_dbWriteTable(con, "dc_working", "fairfax_parcels_blocks", fairfax_parcels_blocks)
DBI::dbDisconnect(con)

Create a demographic multiplier for each parcel

fairfax_parcels_blocks_sf <- fairfax_parcels_blocks

# set as data.table to perform easy group by aggregation
library(data.table)
setDT(fairfax_parcels_blocks_sf)

# get count of parcels per block group: group by block group (first 12 integers of Full_Block) and get count of units per block group (sum(Total_Units))
# remove rows where nchar(substr)!=12
faifax_bg_prcel_cnt <- fairfax_parcels_blocks_sf[, .(cnt = sum(LIVUNIT)), substr(GEOID20, 1, 12)][nchar(substr)==12]

# update column names
colnames(faifax_bg_prcel_cnt) <- c("bg_geoid", "prcl_cnt")

# set va_arl_block_parcels_sf back to sf for geo functions
fairfax_parcels_blocks_sf <- sf::st_as_sf(fairfax_parcels_blocks_sf)

# create block group geoid (to make merge easier)
fairfax_parcels_blocks_sf$bg_geoid <- substr(fairfax_parcels_blocks_sf$GEOID20, 1, 12)

# merge on bg_geoid
fairfax_parcels_blocks_cnts_sf <- merge(fairfax_parcels_blocks_sf, faifax_bg_prcel_cnt, by = "bg_geoid")

# create parcel-level demographic multiplier by dividing Total_Units per parcel by total count of parcels in the block group (prcl_cnt) 
fairfax_parcels_blocks_cnts_sf$mult <- fairfax_parcels_blocks_cnts_sf$LIVUNIT/fairfax_parcels_blocks_cnts_sf$prcl_cnt

Get ACS demographic data

library(tidyr)

# create named acs variable list
named_acs_var_list <- c(
        total_pop = "B01001_001",
        wht_alone = "B02001_002",
        afr_amer_alone = "B02001_003",
        amr_ind_alone = "B02001_004",
        asian_alone = "B02001_005",
        hispanic = "B03002_012",
        male = "B01001_002",
        female = "B01001_026",
        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")


# ACS Data: fairfax county
library(tidycensus)
census_api_key(Sys.getenv('census_api_key'))

acs_data <- data.table::setDT(
  tidycensus::get_acs(
    year = 2019,
    state = "51",
    county = "059",
    geography = "block group",
    output = "wide",
    variables = named_acs_var_list
  )
)


# compute age group
acs_data <- acs_data %>%
        mutate(total_pop=total_popE, wht_alone = wht_aloneE, afr_amer_alone = afr_amer_aloneE, amr_ind_alone = amr_ind_aloneE,
               asian_alone = asian_aloneE, hispanic = hispanicE, male = maleE, female = femaleE,
               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,
               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,
               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) %>%
        dplyr::select(GEOID, NAME,total_pop,wht_alone,afr_amer_alone,amr_ind_alone,asian_alone,hispanic,male,female,pop_under_20,pop_20_64,pop_65_plus) %>%
        gather(measure, value, -c(GEOID, NAME)) %>%
        rename(geoid = GEOID, region_name = NAME) %>%
        select(geoid,region_name,measure,value)

# reshape ACS data

# rename columns
colnames(acs_data) <- c("bg_geoid", "name", "variable", "estimate")

maps parcelle with country block

# merge with ACS data and generate parcel demographic estimates by multiplying ACS estimate by parcel multipliers
fairfax_parcels_blocks_dmgs_sf <- merge(fairfax_parcels_blocks_cnts_sf, acs_data, by='bg_geoid', allow.cartesian=TRUE)
fairfax_parcels_blocks_dmgs_sf$prcl_estimate <- fairfax_parcels_blocks_dmgs_sf$mult * fairfax_parcels_blocks_dmgs_sf$estimate

# save to DB
con <- get_db_conn()
dc_dbWriteTable(con, "dc_working", "fairfax_parcel_demographics", fairfax_parcels_blocks_dmgs_sf)
DBI::dbDisconnect(con)

Create "wide" table of demographic counts per parcel

# switch to data.table
fairfax_parcels_blocks_cnts_dmgs_dt <- data.table::as.data.table(fairfax_parcels_blocks_dmgs_sf)

# drop geometry column - huge because so many repeats and not needed here
fairfax_parcels_blocks_cnts_dmgs_dt$geometry <- NULL

# filter to needed columns
fairfax_parcels_blocks_cnts_dmgs_dt <-fairfax_parcels_blocks_cnts_dmgs_dt[, .(parid = PARID, geoid = GEOID20, measure = variable, value = prcl_estimate, mult=mult)]

# Cast long file to wide
fairfax_parcels_blocks_cnts_dmgs_dt_wide <- data.table::dcast(fairfax_parcels_blocks_cnts_dmgs_dt, parid + geoid + mult ~ measure, value.var = "value", fun.aggregate = sum)

# add back parcel geo
fairfax_parcel_geo <- fairfax_parcels_blocks_dmgs_sf[, c("PARID")]
fairfax_parcel_geo_sf <- sf::st_as_sf(setDF(fairfax_parcel_geo))
fairfax_parcel_geo_sf_unq <- unique(fairfax_parcel_geo_sf)

colnames(fairfax_parcel_geo_sf_unq) <- c("parid", "geometry")
fairfax_parcels_blocks_cnts_dmgs_dt_wide_geo <- merge(fairfax_parcels_blocks_cnts_dmgs_dt_wide, fairfax_parcel_geo_sf_unq, by = "parid")

# add race percentage columns
fairfax_parcels_blocks_cnts_dmgs_dt_wide_geo[, wht_alone_pct := round(100*(wht_alone/total_pop), 2)]
fairfax_parcels_blocks_cnts_dmgs_dt_wide_geo[, afr_amer_alone_pct := round(100*(afr_amer_alone/total_pop), 2)]

# back to sf before writing to DB
fairfax_parcels_blocks_cnts_dmgs_dt_wide_geo <- sf::st_as_sf(fairfax_parcels_blocks_cnts_dmgs_dt_wide_geo)

# save to DB
con <- get_db_conn()
dc_dbWriteTable(con, "dc_working", "fairfax_parcels_blocks_cnts_dmgs_dt_wide_geo", fairfax_parcels_blocks_cnts_dmgs_dt_wide_geo)
DBI::dbDisconnect(con)


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