Main merge temp

code/EconomicVar_Tract_2010_ACS.R

@@ -0,0 +1,60 @@
+# Author: Ashlynn Wimer
+# Date: 07/28/23
+# About: This script uses tidycensus to grab tract level economic Decennial Census
+# and ACS data.
+
+
+# Libraries
+library(tidycensus)
+library(tidyr)
+library(dplyr)
+library(purrr)
+
+# Get a list of all relevant states' FIPS codes.
+states_fips <- tigris::fips_codes |> 
+  filter(!state_code %in% c('60', '72', '66', '69', '78', '74')) |>
+  select(state_code) |> unique() 
+
+states_fips <- states_fips$state_code
+
+# Acquire data
+econ2010 <- map_dfr(.x = states_fips, ~ get_acs(geography = 'tract', 
+        variables = c(
+          NMMLFU    = 'B12006_006',  NMMLF  = 'B12006_004',
+          NMFLFU    = 'B12006_011',  NMFLF  = 'B12006_009',
+          MMLFU     = 'B12006_017',  MMLF   = 'B12006_015',
+          MFLFU     = 'B12006_022',  MFLF   = 'B12006_020',
+          SMLFU     = 'B12006_028',  SMLF   = 'B12006_026',
+          SFLFU     = 'B12006_033',  SFLF   = 'B12006_031',
+          WMLFU     = 'B12006_039',  WMLF   = 'B12006_037',
+          WFLFU     = 'B12006_044',  WFLF   = 'B12006_042',
+          DMLFU     = 'B12006_050',  DMLF   = 'B12006_048',
+          DFLFU     = 'B12006_055',  DFLF   = 'B12006_053',
+          pci       = 'B19301_001',  cntPov = 'B06012_002', 
+          cntPovUni = 'B06012_001'), 
+        year = 2010, state=.x)) 
+
+# Clean the data
+econ2010 <- econ2010 |> 
+  select(GEOID, variable, estimate) |>
+  spread("variable", "estimate") |>
+  mutate(unempP = round(100 * (NMMLFU + NMFLFU + MMLFU +  
+                               MFLFU + SMLFU + SFLFU +  
+                               WMLFU + WFLFU + DMLFU +  
+                               DFLFU) / 
+                          (NMMLFU + NMMLF + NMFLFU + NMFLF + 
+                           SMLFU  + SMLF  + SFLFU  + SFLF + 
+                           MMLFU  + MMLF  + MFLFU  + MFLF + 
+                           WMLFU  + WMLF  + WFLFU  + WFLF + 
+                           DMLFU  + DFLF  + DFLFU  + DFLF), 
+                        2),
+         povP   = round(100 * cntPov / cntPovUni, 2)
+         )|>
+  select(GEOID, unempP, povP, pciE = pci)
+
+# Final cleaning
+econ2010 <- econ2010 |> 
+  map_df(~ ifelse(is.nan(.x), NA, .x))
+
+# Save
+write.csv(econ2010, "../data_final/EC03_2010_T_ACS.csv", row.names = FALSE)

code/access_OTP.R

@@ -0,0 +1,170 @@
+# Author: Ashlynn Wimer
+# Date: July 17th, 2023
+# About: This script prepares the nearest distance metrics to Opioid Treatment 
+# Programs at the tract and zip level, merges those metrics into the drive time
+# metrics calculated in a separate script, cleans the files, and then saves them.
+
+# It is essentially just a modification of access_mentalhealth.R, a script 
+# written by Susan Paykin in January 2021. 
+
+# Libraries #
+library(dplyr)
+library(stringr)
+library(sf)
+library(units)
+
+# Read in point dataset
+opioid_treatment_progs <- read.csv("Access Metrics - Health Resources/Tract/Driving/OpioidTreatmentProgram_Geocoded.csv")
+
+# Convert to spatial data
+otp.sf <- st_as_sf(opioid_treatment_progs,
+                   coords = c('Longitude', 'Latitude'),
+                   crs = 4326)
+
+### Nearest Resource Analysis
+
+# Read in location Data
+zips <- read_sf("../data_final/geometryFiles/tl_2018_zcta/zctas2018.shp")
+tracts <- read_sf("../data_final/geometryFiles/tl_2018_tract/tracts2018.shp")
+
+zips <- zips |> st_transform(3857)
+tracts <- tracts |> st_transform(3857)
+otp.sf <- otp.sf |> st_transform(3857)
+
+### Nearest OTP - ZCTAs ----
+
+# Create centroids for ZCTAs
+zipCentroids <- st_centroid(zips)
+
+# Identify the OTP that is nearest to the zip centroid
+nearestOTP_index <- st_nearest_feature(zipCentroids, otp.sf)
+nearestOTP <- otp.sf[nearestOTP_index,]
+
+# Calculate distance
+minDistZipsOTP <- st_distance(zipCentroids, nearestOTP, by_element = TRUE)
+head(minDistZipsOTP)
+
+# Convert to miles
+minDistZipsOTP_mi <- set_units(minDistZipsOTP, "mi")
+
+# Merge data - everything is ordered so cbind works.
+minDistZipsOTP_sf <- cbind(zips, minDistZipsOTP_mi)
+head(minDistZipsOTP_sf)
+
+# Clean up
+minDistZipsOTP_sf <- minDistZipsOTP_sf |> 
+  select(ZCTA = ZCTA5CE10, 
+         minDisOTP = minDistZipsOTP_mi) |>
+  mutate(ZCTA = str_pad(ZCTA, width=5, side='left', pad='0')) |>
+  st_drop_geometry()
+
+### Nearest OTP - Tracts ----
+
+# Create centroids for tracts
+tractCentroids <- st_centroid(tracts)
+
+# Identify the OTP that is nearest to the tract centroid
+nearestOTP_index <- st_nearest_feature(tractCentroids, otp.sf)
+nearestOTP <- otp.sf[nearestOTP_index,]
+
+# Calculate distance
+minDistTractsOTP <- st_distance(tractCentroids, nearestOTP, by_element = TRUE)
+head(minDistTractsOTP)
+
+# Convert units to miles.
+minDistTractsOTP_mi <- set_units(minDistTractsOTP, "mi")
+
+# Merge data - everything is ordered so cbind works.
+minDistTractsOTP_sf <- cbind(tracts, minDistTractsOTP_mi)
+head(minDistTractsOTP_sf)
+
+# Clean up data
+minDistTractsOTP_sf <- minDistTractsOTP_sf |>
+  select(GEOID, STATEFP, COUNTYFP, 
+         TRACTCE, minDisOTP = minDistTractsOTP_mi) |>
+  mutate(GEOID = str_pad(GEOID, width=11, side='left', pad='0')) |>
+  st_drop_geometry()
+
+### Merge travel access metrics -- Tracts #####
+
+# Read in driving times, rename variables and clean GEOID
+tract_drive <- read.csv('Access Metrics - Health Resources/Tract/Driving/opioidtreatmentprogram_drive_tract.csv') |>
+  rename(
+    GEOID = origin,
+    countDrive = count.within.30,
+    timeDrive = minutes
+    ) |> 
+  mutate(GEOID = str_pad(GEOID, width=11, side='left', pad='0'))
+head(tract_drive)
+
+# Merge
+otp_accessT <- left_join(tract_drive, minDistTractsOTP_sf, by='GEOID')
+head(otp_accessT)
+
+otp_accessT <- otp_accessT |> 
+  select(GEOID, minDisOTP, timeDrive, countDrive)
+head(otp_accessT)
+
+# Save file
+write.csv(otp_accessT, "../data_final/Access07_T.csv", row.names=FALSE)
+
+##### Merge travel access metrics - Zip Code #####
+
+# Read in driving metrics, rename variables, fix ZCTAs
+otp_drive <- read.csv('Access Metrics - Health Resources/ZIP Code/Driving/opioidtreatmentprogram_drive_zip.csv') |>
+  rename(ZCTA = origin,
+         countDrive = count.within.30,
+         timeDrive = minutes) |>
+  mutate(ZCTA = str_pad(ZCTA, width=5, side='left', pad='0'))
+
+# Merge
+otp_accessZ <- left_join(otp_drive, minDistZipsOTP_sf, by='ZCTA')
+head(otp_accessZ)
+
+otp_accessZ <- otp_accessZ |>
+  select(ZCTA, minDisOTP, countDrive, timeDrive)
+
+head(otp_accessZ)
+
+# Save file
+write.csv(otp_accessZ, "../data_final/Access07_Z.csv", row.names=FALSE)
+
+##### County Access #####
+
+tract_OTP <- read.csv("../data_final/Access07_T.csv")
+
+tract_OTP$GEOID <- str_pad(tract_OTP$GEOID, width=11, side='left', pad = '0')
+
+tract_OTP$COUNTYFP <- substr(tract_OTP$GEOID, 1, 5)
+
+county_OTP <- tract_OTP |>
+  group_by(COUNTYFP) |>
+  summarise(cntT = n(),
+            cntTimeDrive = sum(timeDrive <= 30),
+            avTimeDrive = round(mean(timeDrive), 2)
+            ) |>
+  mutate(pctTimeDrive = round(cntTimeDrive / cntT, 2))
+
+head(county_OTP)
+
+##### State Access #####
+
+tract_OTP$STATEFP <- substr(tract_OTP$GEOID, 1, 2)
+
+state_OTP <- tract_OTP |>
+  group_by(STATEFP) |>
+  summarise(cntT = n(),
+            cntTimeDrive = sum(timeDrive <= 30, na.rm = TRUE),
+            avTimeDrive = round(mean(timeDrive, na.rm=TRUE), 2)
+            ) |>
+  mutate(pctTimeDrive = round(cntTimeDrive / cntT, 2))
+
+head(state_OTP)
+
+### Save state and county files
+
+# County
+write.csv(county_OTP, "../data_final/Access07_C.csv", row.names=FALSE)
+
+# State
+write.csv(state_OTP, '../data_final/Access07_S.csv', row.names=FALSE)