Census Data in R, Part 1: Getting Tables with tidycensus

ACS and decennial data for North Carolina counties and census tracts

Use the tidycensus package to pull ACS and 2020 decennial census data for North Carolina counties and Mecklenburg County census tracts.

About this tutorial

Time: about 25 minutes. Prerequisites: R and RStudio (or Positron) installed, basic familiarity with the tidyverse, and a Census API key. If you haven’t read Census 101 yet, skim it first.

By the end, you will be able to:

store your Census API key so that R can find it
pull ACS estimates for every county in North Carolina
pull ACS estimates for census tracts in one or more counties
find the variable IDs you need
compute percentages and their margins of error
pull 2020 decennial census counts

Setup

tidycensus is an R package written by Kyle Walker. It talks to the Census API for you and returns tidy data frames. Install it along with the tidyverse if you don’t have them yet:

install.packages(c("tidycensus", "tidyverse"))

Then load both packages:

library(tidycensus)
library(tidyverse)

Store your API key

Run the line below once, using the key from your email in place of the placeholder. The install = TRUE argument saves the key to a file named .Renviron in your home folder. R reads that file each time it starts, so you never have to paste the key into a script again.

census_api_key("YOUR_KEY_GOES_HERE", install = TRUE)

Restart R after you run it. You can confirm that R found the key with this line:

Sys.getenv("CENSUS_API_KEY")

Keep the key out of your scripts

After you have run census_api_key() once, delete that line from your script. That way you won’t push your key to GitHub by accident.

Your first request: income by county

get_acs() pulls data from the American Community Survey. Let’s get median household income (variable B19013_001) for every county in North Carolina:

nc_income <- get_acs(
  geography = "county",
  state = "NC",
  variables = "B19013_001",
  year = 2024
)

Getting data from the 2020-2024 5-year ACS

nc_income

# A tibble: 100 × 5
   GEOID NAME                             variable   estimate   moe
   <chr> <chr>                            <chr>         <dbl> <dbl>
 1 37001 Alamance County, North Carolina  B19013_001    65651  1854
 2 37003 Alexander County, North Carolina B19013_001    65354  4962
 3 37005 Alleghany County, North Carolina B19013_001    47172  4894
 4 37007 Anson County, North Carolina     B19013_001    47302  5753
 5 37009 Ashe County, North Carolina      B19013_001    56866  4577
 6 37011 Avery County, North Carolina     B19013_001    61428  7429
 7 37013 Beaufort County, North Carolina  B19013_001    58357  3205
 8 37015 Bertie County, North Carolina    B19013_001    48750  8791
 9 37017 Bladen County, North Carolina    B19013_001    42422  6049
10 37019 Brunswick County, North Carolina B19013_001    77024  2734
# ℹ 90 more rows

You should get a tibble with 100 rows, one per county, and five columns:

Column	What it holds
`GEOID`	The county’s FIPS code, such as `37119` for Mecklenburg
`NAME`	A readable name, such as “Mecklenburg County, North Carolina”
`variable`	The variable ID you requested
`estimate`	The ACS estimate
`moe`	The margin of error, at a 90 percent confidence level

A few defaults are worth knowing about:

get_acs() uses the 5-year ACS unless you set survey = "acs1". With year = 2024, you get the 2020 to 2024 estimates.
Always set year yourself. The default changes as new data comes out. Setting it explicitly keeps your results reproducible.
You can write states and counties as names ("North Carolina"), abbreviations ("NC"), or FIPS codes ("37").

Now use regular dplyr to explore:

nc_income |>
  arrange(desc(estimate)) |>
  slice_head(n = 5)

# A tibble: 5 × 5
  GEOID NAME                             variable   estimate   moe
  <chr> <chr>                            <chr>         <dbl> <dbl>
1 37183 Wake County, North Carolina      B19013_001   105768  1626
2 37179 Union County, North Carolina     B19013_001   102900  1953
3 37037 Chatham County, North Carolina   B19013_001    94317  4141
4 37053 Currituck County, North Carolina B19013_001    93511  8425
5 37029 Camden County, North Carolina    B19013_001    91078 14321

Naming variables and requesting several at once

Variable IDs are hard to remember. You can give them names inside c(). You can also request several variables in one call:

nc_profile <- get_acs(
  geography = "county",
  state = "NC",
  variables = c(
    median_income = "B19013_001",
    median_age = "B01002_001",
    total_pop = "B01003_001"
  ),
  year = 2024
)

Getting data from the 2020-2024 5-year ACS

This result is in long format. Each county appears three times, once per variable. Long data works well with ggplot2 and group_by(). If you’d rather have one row per county, ask for wide output:

nc_profile_wide <- get_acs(
  geography = "county",
  state = "NC",
  variables = c(
    median_income = "B19013_001",
    median_age = "B01002_001",
    total_pop = "B01003_001"
  ),
  output = "wide",
  year = 2024
)

Getting data from the 2020-2024 5-year ACS

names(nc_profile_wide)

[1] "GEOID"          "NAME"           "median_incomeE" "median_incomeM"
[5] "median_ageE"    "median_ageM"    "total_popE"     "total_popM"

In wide output, tidycensus adds an E suffix to estimates and an M suffix to margins of error. You get columns such as median_incomeE and median_incomeM.

Finding variables

The 5-year ACS has more than 25,000 variables. load_variables() downloads the full list so you can search it:

acs_vars <- load_variables(2024, "acs5")

In RStudio, run View(acs_vars) and type in the search box to filter. You can also search with code. The label column describes what each variable measures. The concept column names the table it belongs to.

acs_vars |>
  filter(str_detect(concept, regex("median household income", ignore_case = TRUE))) |>
  select(name, label, concept) |>
  head(12)

# A tibble: 12 × 3
   name        label                                                     concept
   <chr>       <chr>                                                     <chr>  
 1 B19013A_001 Estimate!!Median household income in the past 12 months … Median…
 2 B19013B_001 Estimate!!Median household income in the past 12 months … Median…
 3 B19013C_001 Estimate!!Median household income in the past 12 months … Median…
 4 B19013D_001 Estimate!!Median household income in the past 12 months … Median…
 5 B19013E_001 Estimate!!Median household income in the past 12 months … Median…
 6 B19013F_001 Estimate!!Median household income in the past 12 months … Median…
 7 B19013G_001 Estimate!!Median household income in the past 12 months … Median…
 8 B19013H_001 Estimate!!Median household income in the past 12 months … Median…
 9 B19013I_001 Estimate!!Median household income in the past 12 months … Median…
10 B19013_001  Estimate!!Median household income in the past 12 months … Median…
11 B19019_001  Estimate!!Total:                                          Median…
12 B19019_002  Estimate!!Total:!!1-person households                     Median…

Look closely at the results. You will see B19013_001 along with B19013A_001 through B19013I_001. The lettered versions repeat the table for specific racial and ethnic groups. For example, B19013B covers households headed by someone who is Black alone. Tables with PR in the ID cover Puerto Rico only.

Reading labels

Labels use !! to separate levels of a table. For example, Estimate!!Total:!!Income in the past 12 months below poverty level: means “of the total, the number with income below the poverty level.” The first row of most tables (the variable ending in _001) is the total for that table. You will often use it as a denominator.

Computing percentages

Counts are hard to compare across counties of very different sizes. Mecklenburg has more than a million people. Tyrrell County has about 3,000. Percentages make fairer comparisons.

Table B17001 counts people by poverty status. B17001_002 is the number of people below the poverty level. B17001_001 is the number of people for whom poverty status is determined. The summary_var argument fetches the denominator alongside your variable:

nc_poverty <- get_acs(
  geography = "county",
  state = "NC",
  variables = c(below_poverty = "B17001_002"),
  summary_var = "B17001_001",
  year = 2024
) |>
  mutate(
    pct_poverty = 100 * estimate / summary_est,
    pct_poverty_moe = 100 * moe_prop(estimate, summary_est, moe, summary_moe)
  )

Getting data from the 2020-2024 5-year ACS

nc_poverty |>
  select(NAME, pct_poverty, pct_poverty_moe) |>
  arrange(desc(pct_poverty))

# A tibble: 100 × 3
   NAME                              pct_poverty pct_poverty_moe
   <chr>                                   <dbl>           <dbl>
 1 Robeson County, North Carolina           27.2            1.16
 2 Scotland County, North Carolina          25.8            2.06
 3 Watauga County, North Carolina           25.5            2.09
 4 Richmond County, North Carolina          25.2            3.65
 5 Washington County, North Carolina        24.1            4.71
 6 Hertford County, North Carolina          22.9            2.82
 7 Halifax County, North Carolina           22.8            2.30
 8 Hyde County, North Carolina              22.2            8.20
 9 Bladen County, North Carolina            22.0            3.16
10 Lenoir County, North Carolina            22.0            2.37
# ℹ 90 more rows

The summary_var values arrive in the summary_est and summary_moe columns. moe_prop() is a tidycensus helper that computes the margin of error for a proportion. When you derive a new number from ACS estimates, you should derive a new MOE for it too.

Going smaller: census tracts

To look inside a county, switch geography to "tract". Tract requests need a state. They usually include one or more counties as well.

meck_income <- get_acs(
  geography = "tract",
  state = "NC",
  county = "Mecklenburg",
  variables = c(median_income = "B19013_001"),
  year = 2024
)

Getting data from the 2020-2024 5-year ACS

meck_income

# A tibble: 305 × 5
   GEOID       NAME                                      variable estimate   moe
   <chr>       <chr>                                     <chr>       <dbl> <dbl>
 1 37119000101 Census Tract 1.01; Mecklenburg County; N… median_…   119044 19607
 2 37119000102 Census Tract 1.02; Mecklenburg County; N… median_…   134506 13081
 3 37119000103 Census Tract 1.03; Mecklenburg County; N… median_…   130114 25459
 4 37119000104 Census Tract 1.04; Mecklenburg County; N… median_…   107949 16118
 5 37119000301 Census Tract 3.01; Mecklenburg County; N… median_…    91053 27082
 6 37119000302 Census Tract 3.02; Mecklenburg County; N… median_…    81641 26694
 7 37119000401 Census Tract 4.01; Mecklenburg County; N… median_…    79458 25974
 8 37119000402 Census Tract 4.02; Mecklenburg County; N… median_…    76504 12793
 9 37119000501 Census Tract 5.01; Mecklenburg County; N… median_…    85833 44076
10 37119000502 Census Tract 5.02; Mecklenburg County; N… median_…    93399 12184
# ℹ 295 more rows

That should return 305 tracts. To cover the area around Davidson, pass several counties at once:

lkn_income <- get_acs(
  geography = "tract",
  state = "NC",
  county = c("Mecklenburg", "Iredell", "Cabarrus"),
  variables = c(median_income = "B19013_001"),
  year = 2024
)

Getting data from the 2020-2024 5-year ACS

Watch the margins of error

Tract estimates rest on small samples, so their MOEs are much larger than county MOEs. A quick way to spot shaky estimates is to compare each MOE to its estimate:

meck_income |>
  mutate(moe_ratio = moe / estimate) |>
  arrange(desc(moe_ratio)) |>
  select(NAME, estimate, moe, moe_ratio) |>
  head(10)

# A tibble: 10 × 4
   NAME                                                 estimate   moe moe_ratio
   <chr>                                                   <dbl> <dbl>     <dbl>
 1 Census Tract 56.22; Mecklenburg County; North Carol…    14821 13439     0.907
 2 Census Tract 49; Mecklenburg County; North Carolina     40000 33460     0.836
 3 Census Tract 58.52; Mecklenburg County; North Carol…    70972 44217     0.623
 4 Census Tract 56.09; Mecklenburg County; North Carol…    24700 15011     0.608
 5 Census Tract 9; Mecklenburg County; North Carolina     102946 60559     0.588
 6 Census Tract 38.02; Mecklenburg County; North Carol…    57500 32766     0.570
 7 Census Tract 56.26; Mecklenburg County; North Carol…    87908 48056     0.547
 8 Census Tract 62.18; Mecklenburg County; North Carol…   171250 92762     0.542
 9 Census Tract 30.08; Mecklenburg County; North Carol…   110250 56616     0.514
10 Census Tract 5.01; Mecklenburg County; North Caroli…    85833 44076     0.514

A tract whose MOE is 30 or 40 percent of its estimate deserves a caution note in your presentation. You may also notice some NA estimates. In Mecklenburg, they belong to tracts 9801, 9802, and 9803. Tracts numbered in the 9800s are special land-use tracts that cover places like airports and large parks. Few or no people live in them, so the Bureau can’t compute an income.

Top-coded values

Sort meck_income by estimate from highest to lowest and look at the top tract. Its estimate is exactly 250001, and its MOE is NA. That isn’t a real income. The ACS top-codes medians, so 250001 means “$250,000 or more.” Keep this in mind if you average incomes or rank the wealthiest tracts.

The 2020 decennial census

get_decennial() works much like get_acs(). There are two differences. You need to name a summary file with sumfile, and the results come back in a column called value because decennial data has no margin of error.

Total population by county from the Demographic and Housing Characteristics file:

nc_pop_2020 <- get_decennial(
  geography = "county",
  state = "NC",
  variables = c(total_pop = "P1_001N"),
  year = 2020,
  sumfile = "dhc"
)

Getting data from the 2020 decennial Census

Using the Demographic and Housing Characteristics File

Note: 2020 decennial Census data use differential privacy, a technique that
introduces errors into data to preserve respondent confidentiality.
ℹ Small counts should be interpreted with caution.
ℹ See https://www.census.gov/library/fact-sheets/2021/protecting-the-confidentiality-of-the-2020-census-redistricting-data.html for additional guidance.
This message is displayed once per session.

The Redistricting file (sumfile = "pl") has race and ethnicity counts at every level of geography. Here is the share of each Mecklenburg tract’s population that identifies as Black alone:

meck_black <- get_decennial(
  geography = "tract",
  state = "NC",
  county = "Mecklenburg",
  variables = c(black_alone = "P1_004N"),
  summary_var = "P1_001N",
  year = 2020,
  sumfile = "pl"
) |>
  mutate(pct_black = 100 * value / summary_value)

Getting data from the 2020 decennial Census

Using the PL 94-171 Redistricting Data Summary File

You can browse decennial variables the same way as ACS variables, with load_variables(2020, "pl") or load_variables(2020, "dhc").

Which one should I use?

Use the decennial census for population counts and basic demographics in 2020. Use the ACS for anything about income, poverty, education, employment, or housing costs. Use it too when you need something more recent than 2020.

Saving your data

Save your results as a CSV so your teammates can use them:

write_csv(nc_poverty, "nc_poverty_2024.csv")

When you read the file back, tell readr that GEOID is text. Otherwise it may guess that the column is a number and drop leading zeros:

nc_poverty <- read_csv(
  "nc_poverty_2024.csv",
  col_types = cols(GEOID = col_character())
)

Try it yourself

Find the variable for the median gross rent (hint: search the concept column for “median gross rent”). Pull it for every North Carolina county. Which five counties have the highest rent?
Use summary_var to compute the percentage of households in each Mecklenburg tract that have no vehicle available. Table B08201 is a good place to start. Which tracts have the largest MOEs?
Compare total population in 2020 (decennial) with the 2024 ACS 5-year estimate of total population for each county. Which counties seem to be growing fastest? Why should you be cautious about this comparison?