load packages
Getting census API key
## Your original .Renviron will be backed up and stored in your R HOME directory if needed.
## Your API key has been stored in your .Renviron and can be accessed by Sys.getenv("CENSUS_API_KEY").
## To use now, restart R or run `readRenviron("~/.Renviron")`
## [1] "aeb63576bd50b4e3e7ce9beabbac2d8ccc2bf6ce"
downloading population data, getting 2020 ACS data for population density
nyc_pop <- get_acs(geography = "tract", variables = "B01003_001", state = "NY", county = c("Bronx", "Kings", "New York", "Queens", "Richmond"), geometry = TRUE, year = 2020)
## Getting data from the 2016-2020 5-year ACS
## Downloading feature geometry from the Census website. To cache shapefiles for use in future sessions, set `options(tigris_use_cache = TRUE)`.
## | | | 0% | | | 1% | |= | 1% | |= | 2% | |== | 2% | |== | 3% | |=== | 4% | |=== | 5% | |==== | 5% | |==== | 6% | |===== | 7% | |===== | 8% | |====== | 8% | |====== | 9% | |======= | 9% | |======= | 10% | |======== | 11% | |======== | 12% | |========= | 12% | |========= | 13% | |========== | 14% | |=========== | 15% | |=========== | 16% | |============ | 17% | |============= | 18% | |============= | 19% | |============== | 20% | |=============== | 21% | |=============== | 22% | |================ | 22% | |================ | 23% | |================= | 24% | |================= | 25% | |================== | 25% | |================== | 26% | |=================== | 27% | |==================== | 28% | |==================== | 29% | |===================== | 30% | |====================== | 32% | |======================= | 33% | |======================== | 34% | |======================== | 35% | |========================= | 35% | |========================= | 36% | |========================== | 37% | |========================== | 38% | |=========================== | 38% | |=========================== | 39% | |============================ | 40% | |============================ | 41% | |============================= | 41% | |============================= | 42% | |============================== | 43% | |=============================== | 44% | |=============================== | 45% | |================================ | 46% | |================================= | 47% | |================================= | 48% | |================================== | 49% | |=================================== | 50% | |==================================== | 51% | |===================================== | 52% | |===================================== | 53% | |====================================== | 54% | |======================================= | 55% | |======================================= | 56% | |======================================== | 57% | |========================================= | 58% | |========================================= | 59% | |========================================== | 59% | |========================================== | 60% | |=========================================== | 61% | |=========================================== | 62% | |============================================ | 62% | |============================================ | 63% | |============================================= | 64% | |============================================= | 65% | |============================================== | 65% | |============================================== | 66% | |=============================================== | 67% | |================================================ | 68% | |================================================ | 69% | |================================================= | 70% | |================================================== | 71% | |================================================== | 72% | |=================================================== | 72% | |=================================================== | 73% | |==================================================== | 74% | |==================================================== | 75% | |===================================================== | 75% | |===================================================== | 76% | |====================================================== | 77% | |====================================================== | 78% | |======================================================= | 78% | |======================================================= | 79% | |======================================================== | 80% | |========================================================= | 81% | |========================================================= | 82% | |========================================================== | 83% | |=========================================================== | 84% | |=========================================================== | 85% | |============================================================ | 86% | |============================================================= | 87% | |============================================================= | 88% | |============================================================== | 88% | |============================================================== | 89% | |=============================================================== | 90% | |=============================================================== | 91% | |================================================================ | 91% | |================================================================ | 92% | |================================================================= | 93% | |================================================================= | 94% | |================================================================== | 94% | |================================================================== | 95% | |=================================================================== | 96% | |==================================================================== | 96% | |==================================================================== | 97% | |===================================================================== | 98% | |===================================================================== | 99% | |======================================================================| 99% | |======================================================================| 100%
loading my download NYC building footprint data
build_footprint <- st_read("/Users/elvagao/Desktop/Rstudio_GTECH/Project_GTECH38520/Building_Footprints/geo_export_1aaeed34-8dab-47fe-a8f7-60738ad1c437.shp")
## Reading layer `geo_export_1aaeed34-8dab-47fe-a8f7-60738ad1c437' from data source `/Users/elvagao/Desktop/Rstudio_GTECH/Project_GTECH38520/Building_Footprints/geo_export_1aaeed34-8dab-47fe-a8f7-60738ad1c437.shp'
## using driver `ESRI Shapefile'
## Simple feature collection with 1082678 features and 16 fields
## Geometry type: POLYGON
## Dimension: XY
## Bounding box: xmin: -74.2555 ymin: 40.49843 xmax: -73.70006 ymax: 40.91505
## Geodetic CRS: WGS 84
tract geometry
nyc_tracts <- tracts(state = "NY", county = c("New York", "Kings", "Queens", "Bronx", "Richmond"), cb = TRUE, year = 2020) %>% st_transform(2263)
simplify geometry and crop to NYC tracts
build_footprint <- st_transform(build_footprint, st_crs(nyc_pop))
fixing the geometries before intersecting (fixed after knowing the error)
build_footprint <- st_make_valid(build_footprint)
try again intersecting again
build_footprint <- st_intersection(build_footprint, st_union(nyc_pop))
## Warning: attribute variables are assumed to be spatially constant throughout
## all geometries
*warning meaning that I am intersecting many building polygons with a single, merged geometry. If a single building is split into many pieces during intersection, the attribute values would stay across those pieces.
count buildings per tract
building_count <- st_join(build_footprint, nyc_pop) %>% group_by(GEOID) %>% summarise(num_buildings = n())
joining back to population data
nyc_merged <- nyc_pop %>% left_join(st_drop_geometry(building_count), by = "GEOID") %>% mutate(num_buildings = replace_na(num_buildings, 0))
suitable score
nyc_merged <- nyc_merged %>% mutate(pop_density_scaled = scale(estimate), building_density_scaled = scale(num_buildings), suitability = -1 * (pop_density_scaled + building_density_scaled))
Define a suitability score
nyc_merged <- nyc_merged %>% mutate(suitability = scale(num_buildings) - scale(estimate))
nyc_merged$suitability[is.infinite(nyc_merged$suitability)] <- 0
NA
nyc_merged <- nyc_merged %>% mutate(suitability = replace_na(suitability, 0))
summary(nyc_merged$suitability)
## V1
## Min. :-7.4073
## 1st Qu.:-0.7976
## Median : 0.1802
## Mean : 0.0000
## 3rd Qu.: 0.8308
## Max. : 5.7789
sum(is.na(nyc_merged$suitability))
## [1] 0
sum(is.nan(nyc_merged$suitability))
## [1] 0
nyc_merged$suitability[is.na(nyc_merged$suitability) | is.nan(nyc_merged$suitability)] <- 0
suitability map, they update he ‘RdYlGn’ to ‘brewer.rd_yl_gn’ version and putting 0 for midpoint for clarity
tmap_mode("view")
## ℹ tmap mode set to "view".
tm_shape(nyc_merged) + tm_polygons(fill = "suitability", fill.scale = tm_scale_intervals(style = "quantile", values = "brewer.rd_yl_gn", midpoint = 0), fill.legend = tm_legend(title = "Suitability")) + tm_borders() + tm_tiles("OpenStreetMap") + tm_title("NYC Safety Zones")
## Registered S3 method overwritten by 'jsonify':
## method from
## print.json jsonlite