In mfherman/nycgeo: Spatial Data Files for NYC

knitr::opts_chunk$set(
  collapse = TRUE,
  comment = "#>",
  fig.path = "man/figures/README-",
  out.width = "100%",
  dpi = 150,
  cache = TRUE,
  cache.path = "cache/readme/"
  )

nycgeo

The nycgeo package contains spatial data files for various geographic and administrative boundaries in New York City as well as tools for working with NYC spatial data. Data is in the sf (simple features) format and includes boundaries for boroughs (counties), public use microdata areas (PUMAs), community districts (CDs), neighborhood tabulation areas (NTAs), census tracts, census blocks, city council districts, school districts, and police precincts.

Additionally, selected demographic, social, and economic estimates from the U.S. Census Bureau American Community Survey can be added to the geographic boundaries in nycgeo, allowing for contextualization and easy choropleth mapping. Finally, nycgeo makes it simple to access a subset of spatial data in a particular geographic area, such as all census tracts in Brooklyn and Queens.

Why `nycgeo`?

The spatial files contained in the nycgeo package are available on websites such as the New York City Department of City Planning's Bytes of the Big Apple and the U.S. Census Bureau TIGER/Line® Shapefiles, but this package aims to make accessing the spatial data more convenient. Instead of downloading and converting shapefiles each time you need them, nycgeo provides the files in a consistent format (sf) with added metadata that enable joins with non-spatial data.

Other R packages share some features with nycgeo. In particular, the wonderful tidycensus package can access the Census Bureau's API and download ACS estimates as well as TIGER/Line® Shapefiles (via tigris).

One difference between the boundaries included here and the TIGER/Line® Shapefiles available through tigris is that these boundaries are clipped to the shoreline, allowing for better mapping of New York City. Additionally, nycgeo contains boundaries for geographic areas that are not available from the Census Bureau. This includes neighborhood tabulation areas (NTAs) and community districts (CDs).

Finally, all spatial data included in the package uses the NAD83 / New York Long Island (ftUS) State Plane projected coordinate system (EPSG 2263), which is the standard projection used by New York City government agencies.

Installation

You can install nycgeo from GitHub with:

# install.packages("remotes")
remotes::install_github("mfherman/nycgeo")

To get the most out of nycgeo, you should also install and load the sf package when you use nycgeo. If you haven't attached sf, you will get this friendly reminder when you load nycgeo:

library(nycgeo)

Depending on your operating system and available libraries, sf can be tricky to install the first time. The sf website is a good place to start if you're having trouble. If you're using macOS, this is a good guide to installing the required libraries.

# install.packages("sf")
library(sf)

Examples

Basic Usage

library(nycgeo)
library(sf)
library(tidyverse)

The most basic usage of nycgeo is to get boundaries in the sf format. Use nyc_boundaries() to get your desired geography. To make best use of the package, you should also load the sf package when using nycgeo. For these examples, I'll also load tidyverse as this will allow us to take advantage of pretty tibble printing and will come in handy when we want to manipulate and map the spatial data later.

library(nycgeo)
library(sf)
library(tidyverse)

nyc_boundaries(geography = "tract")

Filter by geography

If you don't need census tracts for the entire city, you can use the filter_by and region arguments of nyc_boundaries() to specify the area you are interested in. For example, the following code returns only census tracts in Brooklyn and Queens.

bk_qn_tracts <- nyc_boundaries(
  geography = "tract",
  filter_by = "borough",
  region = c("brooklyn", "queens")
  )

ggplot(bk_qn_tracts) +
  geom_sf() +
  theme_minimal()

Note, you can select multiple regions by passing a character vector to the region argument, but you can only choose a single geography to filter_by. Additionally, you can only filter by a geography that is larger than or equal to the boundaries you request. For example, it is not possible to filter PUMAs by NTAs because NTAs are smaller than PUMAs.

Adding American Community Survey Data

nycgeo includes selected estimates from the American Community Survey as datasets. You can access these datasets directly or have them appended to the spatial data. To print a tibble of ACS data, simply call the data you want.

nta_acs_data

To add census estimates to an sf object, use add_acs_data = TRUE to an nyc_boundaries()call. For example, here we get all NTAs in Manhattan with ACS data appended. One convenience of having the ACS data joined to the sf object is that you can very simply make a choropleth map. Here we do it with ggplot2, but you could use tmap, leaflet or any other spatial package that works with sf objects.

mn_ntas <- nyc_boundaries(
  geography = "nta",
  filter_by = "borough",
  region = "manhattan",
  add_acs_data = TRUE
  )

ggplot(mn_ntas) +
  geom_sf(aes(fill = pop_ba_above_pct_est)) +
  scale_fill_viridis_c(
    name = "Bachelor's or above",
    labels = scales::percent_format(),
    option = "magma"
    ) +
  theme_void() +
  theme(panel.grid = element_line(color = "transparent")) +
  labs(title = "Which neighborhoods in Manhattan are most educated?")

Joining with other data

One use case of nycgeo() is if you have non-spatial data that relates to census tracts, NTAs, or other geographies and need to join that data with spatial boundaries to plot or otherwise analyze. This non-spatial data may be coded in a variety of ways and might not have names or IDs that match your spatial data. The sf data provided in nycgeo seeks to have a variety of geographic metadata that will match whatever labels your non-spatial data has.

In this example, we have non-spatial data from the NYC Neighborhood Health Atlas at the NTA-level from which we would like to make a choropleth map. To do this, we import the .csv file and then join it to the spatial NTA object matching on NTA IDs. Then, we can map it as in the above example.

nta_health <- read_csv("https://raw.githubusercontent.com/mfherman/nycgeo/master/inst/extdata/nta-health.csv") %>% 
  select(NTA_Code, BlackCarbon)

nyc_boundaries(geography = "nta") %>% 
  left_join(nta_health, by = c("nta_id" = "NTA_Code")) %>% 
  ggplot() +
  geom_sf(aes(fill = BlackCarbon)) +
  scale_fill_viridis_c(name = "Black carbon (absorbance units)", option = "inferno") +
  theme_void() +
  theme(panel.grid = element_line(color = "transparent")) +
  labs(title = "Which neighborhoods have high levels of black carbon pollution?")

Finding which districts a set of points lies within

Point-in-polygon operations are common tasks for spatial analysis. Given a set of points we want to find out which polygon contains each point. A real-world application of this would be counting the number of schools in each community district.

We start with a (non-spatial) data frame of all schools in New York, but with columns for latitude and longitude. Then we use those latitudes and longitudes to convert the data frame to an sf object. From there, we can use the nyc_point_poly() function to find which community district (CD) each point (school) is in and then count by CD to get the total number of schools in each CD.

nyc_schools <- read_csv("https://raw.githubusercontent.com/mfherman/nycgeo/master/inst/extdata/nyc-schools.csv")

schools_sf <- nyc_schools %>% 
  st_as_sf(
    coords = c("longitude", "latitude"),
    crs = 4326,
    stringsAsFactors = FALSE
    )

nyc_point_poly(schools_sf, "cd") %>% 
  st_set_geometry(NULL) %>% 
  count(cd_name, borough_cd_id)