Downloading and pre-processing census data

Open your browser and navigate to National Historical Geographic Information System (NHGIS) at https://www.nhgis.org/. NHGIS provides normalized Census data across survey years and is generally easier to navigate than the US Census Bureau’s website. Click to search for data, and then set the geographic level to Block Group and the year to the most recent ACS 5-year period (e.g. 2020-2024). Then set the topic to Income - Household and Family Income and Race, Ethnicity and Origins - Hispanic Origin. Make sure the Topic selector is set to OR and not AND.

Navigate to the GIS Files tab, and select the Block Group level data for New York State. Proceed with the wizard, selecting CSV as the output format. After a brief period, both the tabular and spatial data will be available for download- download both files and save them to your data directory.

After unzipping the downloaded files, you’ll notice a range of files with different extensions in each folder. Let’s import both datesets using pandas and geopandas to see what they look like.

import pandas as pd
import geopandas as gpd

census_data = pd.read_csv("./Data/nhgis0080_csv/nhgis0080_ds272_20245_blck_grp.csv")
census_geo = gpd.read_file(
    "./Data/nhgis0080_shape/nhgis0080_shapefile_tl2024_360_blck_grp_2024/NY_blck_grp_2024.shp"
)

We can inspect each dataset now to see what they look like using the head() method.

census_data.head()

census_geo.head()

We can see that the tabular data contains 116 columns, including several identifier and metadata columns, as well as the columns containing the data we selected for download. The spatial data contains fewer columns, but importantly contains a geometry column that contains the spatial information for each block group as well as a matching GISJOIN column that we can use to combine the datasets together. Let’s do that now.

census_data_merged = census_geo.merge(census_data, on="GISJOIN", how="inner").set_index(
    "GISJOIN"
)

census_data_merged

Now we have a merged block group level dataset for New York State that contains data on race/ethnicity and income. However, the column names are not very intuitive. Let’s rename the columns to make them easier to work with using the metadata file that was included in the download. The metadata file contains a mapping of the column names to more descriptive names, so we can use that to rename our columns. If we open the .txt metadata file, we can see the alphanumeric column names in our dataset map to descriptive column names:

The metadata makes clear another important fact about the Census: the data we use are estimates, and each estimate has a margin of error associated with it. The more fine-grained the geographic level, the larger the margin of error. For our purposes, we will use the estimates as-is, but it’s important to be aware of the uncertainty in these estimates when working with Census data.

We’ll use the metadata file to create a mapping of the original column names to more descriptive names, and then use that mapping to rename our columns. This will make it easier to work with the data in the next steps when we start creating maps.

column_xwalk = {
    "AUPFE001": "Total Pop (Race/Ethnicity)",
    "AUPFE002": "Not Hispanic or Latino",
    "AUPFE003": "Not Hispanic or Latino: White alone",
    "AUPFE004": "Not Hispanic or Latino: Black or African American alone",
    "AUPFE005": "Not Hispanic or Latino: American Indian and Alaska Native alone",
    "AUPFE006": "Not Hispanic or Latino: Asian alone",
    "AUPFE007": "Not Hispanic or Latino: Native Hawaiian and Other Pacific Islander alone",
    "AUPFE008": "Not Hispanic or Latino: Some other race alone",
    "AUPFE009": "Not Hispanic or Latino: Two or more races",
    "AUPFE010": "Not Hispanic or Latino: Two or more races: Two races including Some other race",
    "AUPFE011": "Not Hispanic or Latino: Two or more races: Two races excluding Some other race, and three or more races",
    "AUPFE012": "Hispanic or Latino",
    "AUPFE013": "Hispanic or Latino: White alone",
    "AUPFE014": "Hispanic or Latino: Black or African American alone",
    "AUPFE015": "Hispanic or Latino: American Indian and Alaska Native alone",
    "AUPFE016": "Hispanic or Latino: Asian alone",
    "AUPFE017": "Hispanic or Latino: Native Hawaiian and Other Pacific Islander alone",
    "AUPFE018": "Hispanic or Latino: Some other race alone",
    "AUPFE019": "Hispanic or Latino: Two or more races",
    "AUPFE020": "Hispanic or Latino: Two or more races: Two races including Some other race",
    "AUPFE021": "Hispanic or Latino: Two or more races: Two races excluding Some other race, and three or more races",
    "AURTE001": "Total Pop (Income)",
    "AURTE002": "Less than $10,000",
    "AURTE003": "$10,000 to $14,999",
    "AURTE004": "$15,000 to $19,999",
    "AURTE005": "$20,000 to $24,999",
    "AURTE006": "$25,000 to $29,999",
    "AURTE007": "$30,000 to $34,999",
    "AURTE008": "$35,000 to $39,999",
    "AURTE009": "$40,000 to $44,999",
    "AURTE010": "$45,000 to $49,999",
    "AURTE011": "$50,000 to $59,999",
    "AURTE012": "$60,000 to $74,999",
    "AURTE013": "$75,000 to $99,999",
    "AURTE014": "$100,000 to $124,999",
    "AURTE015": "$125,000 to $149,999",
    "AURTE016": "$150,000 to $199,999",
    "AURTE017": "$200,000 or more",
}

Now we’ll use this crosswalk to rename the columns in our merged dataset using the rename method in pandas.

census_data_merged.rename(columns=column_xwalk, inplace=True)

At this point, we’ve renamed the variables of interest to us, but there are still over 100 columns in our dataset, many of which we won’t need. Let’s reduce the columns to only include those we’ll need for our analysis. This includes the columns we renamed, as well as the GISJOIN column and the geometry column, as well as a few other identifiers that will come in handy for further subsetting to NYC.

For the purposes of our analysis, we don’t need the granularity of the various Hispanic or Latino origin categories, so we’ll just keep the total number of Hispanic or Latino residents in each block group. We can drop the other columns to simplify our dataset.

race_ethnicity_keep_columns = list(column_xwalk.values())[:12]
income_keep_columns = list(column_xwalk.values())[21:]

cols_to_keep = (
    race_ethnicity_keep_columns
    + income_keep_columns
    + [
        "geometry",
        "COUNTY",
        "TRACTCE",
        "BLKGRPCE",
        "ALAND",  # land area (square meters)
    ]
)

census_data_merged_cleaned = census_data_merged[cols_to_keep].copy()

census_data_merged_cleaned.head()

Now let’s use the COUNTY column to filter the dataset to only include rows for New York City.

census_nyc = census_data_merged_cleaned[
    census_data_merged_cleaned["COUNTY"].isin(
        [
            "Kings County",
            "Queens County",
            "New York County",
            "Bronx County",
            "Richmond County",
        ]
    )
].copy()

We can print the shape of the resulting GeoDataFrame and plot it to confirm that we have the expected number of block groups and that they look correct: many fewer columns and rows, as well as a map of NYC’s block groups.

print(census_nyc.shape)
census_nyc.plot().set_axis_off()

(6593, 34)

Finally before saving, let’s change the coordinate reference system to EPSG:2263, which is the standard coordinate reference system for New York State Plane. This will make it easier to work with other NYC Open Data products later on.

census_nyc_sp = census_nyc.to_crs(epsg=2263)

census_nyc_sp.plot().set_axis_off()

census_nyc_sp.to_file("./Data/census_nyc_sp.fgb")