A user on Stackoverflow recently asked about how to convert a data set containing crime statistics (e.g. auto accidents, shootings, etc) into rates per million population by U.S. state.
To calculate crime rates one needs to merge a source of U.S. state level population data with the data frame containing the event by state data. Fortunately, United States taxpayers have paid to have this information collected by funding of the U.S. Census Bureau.
To illustrate this approach one can use the 2019 U.S. population estimates from the U.S. Census Bureau that are included in a data set of state level population estimates from 2010 – 2019.
The numerators in our rates will be derived from the June 20, 2020 COVID-19 statistics courtesy of our friends at the Johns Hopkins University Center for Systems Science and Engineering.
First, we download and load the US population data, which looks like this.
popData <- "https://www2.census.gov/programs-surveys/popest/tables/2010-2019/state/totals/nst-est2019-01.xlsx" download.file(popData, "./data/nst-est2019-01.xlsx",mode="wb")
Next we load the
readxl package, read the population data, assign column names and keep the July 2019 estimate. We use the
range= argument to skip the first few rows of data that are summary rows for the total U.S. and aggregated regions of states, as well as the last row which contains data for Puerto Rico.
library(readxl) colNames <- c("state","census2010","est_base",paste0("est_",2010:2019)) usPopEst <- read_excel("./data/nst-est2019-01.xlsx", range="A10:M60", col_names = colNames)[c(1,13)] usPopEst$state <- gsub("\\.","",usPopEst$state)
Now, we read the COVID-19 data as of June 20th, 2020 from my forked version of the JHU CCSE COVID-19 Github repository.
# read covid-19 data for June 20th covidFile <- "https://raw.githubusercontent.com/lgreski/COVID-19/master/csse_covid_19_data/csse_covid_19_daily_reports/06-20-2020.csv" covidData <- read.csv(covidFile, header=TRUE)
Since the COVID-19 data is at a county / zip code level, we need to aggregate it to the state level. Once aggregated, we merge the confirmed cases and deaths columns with the population estimates and calculate incidence rates per million population. The
dplyr package provides the R functions we need to group, aggregate, and calculate the desired rates.
# aggregate stats to state level library(dplyr) covidData %>% filter(Country_Region == "US") %>% rename(state = Province_State) %>% group_by(state) %>% summarise(Confirmed = sum(Confirmed), Deaths = sum(Deaths)) %>% inner_join(.,usPopEst) %>% mutate(confirmed_per_million = Confirmed / (est_2019 / 1000000), death_per_million = Deaths / (est_2019 / 1000000)) -> summedStates head(summedStates)
…and the output:
> head(summedStates) # A tibble: 6 x 6 state Confirmed Deaths est_2019 confirmed_per_million death_per_million
1 Alabama 29549 838 4903185 6026. 171. 2 Alaska 741 12 731545 1013. 16.4 3 Arizona 50127 1346 7278717 6887. 185. 4 Arkansas 15142 224 3017804 5018. 74.2 5 California 175213 5494 39512223 4434. 139. 6 Colorado 30333 1647 5758736 5267. 286. >