class: title-slide, center, middle background-image: url(images/ouaerial.jpeg) background-size: cover # .heat[.fancy[Census and USGS Data]] ## .heat[.fancy[Ani Ruhil]] --- name: agenda ## .fancy[ .heat[ Agenda ]] Brief overview of useful Census and other data sources you may need - using the Census API with `tidycensus` and `censusapi` - data from the World Bank - data from the USGS --- class: inverse, middle, center # .heat[ .fancy[ `tidycensus` ]] --- Kyle Walker has authored `tidycensus` to ease working with Census data and it is built to work well with the `tidyverse` and mapping libraries. * get a Census API key [from here](http://api.census.gov/data/key_signup.html) * install the `tidycensus` package ```r library(tidycensus) library(tidyverse) census_api_key("INSERT YOUR API KEY HERE") ``` Two functions of primary value -- (1) `get_decennial()` for decennial census data from 1990, 2000, and 2010 (2) `get_acs()` for data from the American Community Survey (ACS) series Typical starting point will be finding the table you want and an easy way to do this is by pulling a list of all tables available for a sopecific census product. Below you see the code for grabbing a list of the Summary File 1 tables from the 2010 decennial census, followed by the 1-year 2016 ACS and the 5-year ACS (2012-2016). ```r load_variables(year = 2010, dataset = "sf1", cache = TRUE) -> my.tables.d load_variables(year = 2018, dataset = "acs1", cache = TRUE) -> my.tables.a1 load_variables(year = 2018, dataset = "acs5", cache = TRUE) -> my.tables.a5 ``` --- Each of these will show you what variables are available, and make it easier for you to choose the ones you want to work with. Assume we want total population, variable `P0010001`, from the decennial census ```r get_decennial( year = 2010, geography = "state", variables = c(totpopn = "P001001") ) -> state.popn.d get_decennial( year = 2010, geography = "county", variables = c(totpopn = "P001001"), state = "OH" ) -> county.popn.d get_decennial( geography = "tract", variables = c(totpopn = "P001001"), state = "OH", county = "Athens" ) -> tract.popn.d ``` and then from the ACS. Note that the default for each is 2010 and the most recent ACS (2012-2016). ```r get_acs( year = 2018, geography = "state", variables = c(totpopn = "B01003_001") ) -> state.popn.acs get_acs( year = 2018, geography = "county", variables = c(totpopn = "B01003_001"), state = "OH" ) -> county.popn.acs get_acs( year = 2018, geography = "tract", variables = c(totpopn = "B01003_001"), state = "OH", county = "Athens" ) -> tract.popn.acs ``` --- You can also get an `entire table` instead of having to list variables, such as shown below, with the default decennial census for the package's current version (2010). ```r get_decennial( year = 2010, geography = "county", table = "P012", summary_var = "P001001" ) -> county.table.d get_decennial( year = 2010, geography = "state", table = "P012", summary_var = "P001001" ) -> state.table.d ``` Population data for `all Census tracts in the country`? ```r library(tidycensus) library(purrr) unique(fips_codes$state)[1:51] -> us map_df( us, function(x) { get_acs(geography = "tract", variables = "B01003_001", state = x, geometry = TRUE) } ) -> totalpop ``` --- class: inverse, middle, center # .salt[.fancy[ `censusapi` ]] --- This package will allow you to grab a vast array of Census products, and very niftily as well I might add. ```r library(censusapi) listCensusApis() -> apis ``` Application Programming Interfaces (APIs) could have slightly varying parameters so it pays to check the API documentation [available here](https://www.census.gov/data/developers/data-sets.html). It is easy to find variable names via the built-in function `listCensusMetadata`. As an example, the bureau's small area health insurance estimates are shown below, as well as the small area income and poverty estimates. ```r listCensusMetadata( name = "timeseries/healthins/sahie", type = "variables" ) -> sahie_vars listCensusMetadata( name = "timeseries/poverty/saipe", type = "variables" ) -> saipe_vars ``` Want to see what geographies are available? Switch `type =` to `geography`: ```r listCensusMetadata( name = "timeseries/healthins/sahie", type = "geography" ) -> sahie_geos listCensusMetadata( name = "timeseries/poverty/saipe", type = "geography" ) -> saipe_geos ``` --- Grab the most recent county-level data but note that the latest SAHIE are for 2015 while the latest SAIPE are for 2016. Every variable is downloaded as a `chr` so you will need to flip what should be numeric into numeric. ```r getCensus( name = "timeseries/healthins/sahie", vars = c("NAME", "IPRCAT", "IPR_DESC", "PCTUI_PT"), region = "county:*", regionin = "state:39", time = 2015, key = CENSUS_KEY ) -> sahie_counties head(sahie_counties, n = 12L) getCensus( name = "timeseries/poverty/saipe", vars = c("NAME", "SAEPOVRTALL_PT", "SAEMHI_PT"), region = "county:*", regionin = "state:39", time = 2016, key = CENSUS_KEY) -> saipe_counties head(saipe_counties, n = 12L) saipe_counties %>% mutate( prate = as.numeric(saipe_counties$SAEPOVRTALL_PT), mdhinc = as.numeric(saipe_counties$SAEMHI_PT) ) -> saipe_counties ``` --- If you want data for a lot of geographies, the package will let you do it in seconds. For example, if you want tract-level data for `all tracts`, you can get it as shown below: ```r fips tracts <- NULL for (f in fips) { stateget <- paste("state:", f, sep = "") temp <- getCensus(name = "sf3", vintage = 1990, vars = c("P0070001", "P0070002", "P114A001"), region = "tract:*", regionin = stateget, key = CENSUS_KEY) tracts <- rbind(tracts, temp) } head(tracts) ``` Notice that you had to specify the `region`, and since tracts are nested within states, you had to add `regionin = stateget`. --- class: inverse, middle, center # .fancy[.large[ World Bank data with `data360r` and `wbstats` ]] --- Two packages are available to access data gathered by the World Bank, [`data360r`](https://github.com/mrpsonglao/data360r) and [`wbstats`](https://github.com/GIST-ORNL/wbstats). Install both packages, and then we can start with `data360r`. We can start by seeing what indicators are available, some basic country-level information, and what data-sets are available. ```r library(data360r) #get all indicator metadata in Govdata360 get_metadata360(site = "gov", metadata_type = "indicators") -> df_indicators #get all country metadata in TCdata360 get_metadata360(metadata_type = 'countries') -> df_countries #get all dataset metadata in TCdata360 get_metadata360(metadata_type = 'datasets') -> df_datasets ``` --- Once you have identified a particular table, note it's ID, and then pull it. Say I want indicator 90 = `Can a married woman register a business in the same way as a married man?` I can do the same thing for more than one indicator by specifying the IDs. ```r get_data360(indicator_id = c(90)) -> df_ind90 get_data360(indicator_id = c(28130, 28131)) -> df_indtwo ``` If I only want all data for a specific country or just the indicators we pull in `df_ind1` for a specific country you could do: ```r get_data360(country_iso3 = "IND") -> df_allone search_360("woman business", search_type = "indicator", limit_results = 5) -> df_ind1 get_data360(indicator_id = df_ind1$id, country_iso3 = "IND") -> df_allindtwo ``` --- Now an example with two measures -- legal age of marriage for boys and girls. Note that the package allows you to specify the `long` format (preferred) than the default `wide` format you see earlier results being returned in. Note also the use of `timeframes = c()` that allows you to specify the specific time period you want the indicators for. ```r search_360("marriage", search_type = "indicator") -> marriage.indicators DT::datatable(marriage.indicators) ```
--- ## .fancy[ Legal age of marriage? ] .pull-left[ ```r get_data360(indicator_id = c(204, 205), timeframes = c(2016), output_type = 'long' ) -> df_marriage ggplot(df_marriage, aes(x = Observation, group = Indicator, fill = Indicator)) + geom_bar() + theme(legend.position = "none") + facet_wrap(~ Indicator, ncol = 1) + labs( x = "Legal Age of Marriage", y = "Frequency", title = "Legal Age of Marriage", subtitle = " for Boys vs. Girls (2016)", caption = "Source: World Bank Data") ``` ] .pull-right[ <img src="Module06_files/figure-html/wbo4b2-1.png" width="100%" style="display: block; margin: auto;" /> ] --- The `wbstats` package does pretty much the same thing. Let us see the core functionality by loading the library and then seeing what is available in terms of indicators, topics, and so on. We can then set the most current list of information in `wb_cachelist` to be used via `new_cache`. Doing so speeds up the operations and ensures that you are getting the most uptodate data. ```r library(wbstats) str(wb_cachelist, max.level = 1) ``` ``` ## List of 7 ## $ countries :'data.frame': 304 obs. of 18 variables: ## $ indicators :'data.frame': 16978 obs. of 7 variables: ## $ sources :'data.frame': 43 obs. of 8 variables: ## $ datacatalog:'data.frame': 238 obs. of 29 variables: ## $ topics :'data.frame': 21 obs. of 3 variables: ## $ income :'data.frame': 7 obs. of 3 variables: ## $ lending :'data.frame': 4 obs. of 3 variables: ``` ```r new_cache <- wbcache() ``` --- What indicators are available? ```r wbsearch(pattern = "corruption") -> corruption_vars DT::datatable(corruption_vars) ```
--- If I want information from a particular source, say Bloomberg, ```r wbsearch( pattern = "Bloomberg", fields = "sourceOrg" ) -> blmbrg_vars DT::datatable(blmbrg_vars) ```
--- Searching for indicators tied to multiple subjects is easy as well: ```r wbsearch( pattern = "poverty | unemployment | employment" ) -> povemply_vars DT::datatable(povemply_vars) ```
--- Once we identify what we want, downloading the data is easy as well, needing us to specify just the indicator(s) and then the start and end dates, and then specific country codes if you want data for specific countries. Below I am pulling total population. .pull-left[ ```r wb( indicator = "SP.POP.TOTL", startdate = 1960, enddate = 2016 ) -> pop_data1 wb( country = c("ABW","AFG", "SSF", "ECA", "IND", "CHN"), indicator = "SP.POP.TOTL", startdate = 1960, enddate = 2016 ) -> pop_data2 wb( country = c("ABW","AFG", "SSF", "ECA", "IND", "CHN"), indicator = c("SP.POP.TOTL", "NY.GDP.MKTP.CD"), startdate = 1960, enddate = 2016 ) -> pop_data3 DT::datatable(pop_data3) ``` ] .pull-right[
] --- ```r wb( country = c("ABW","AFG", "SSF", "ECA", "IND", "CHN"), indicator = c("SP.POP.TOTL", "NY.GDP.MKTP.CD"), startdate = 1960, enddate = 2016, return_wide = TRUE ) -> pop_data4 DT::datatable(pop_data4) ```
--- By default `wb()` will return the data in long format but not necessarily in a tidy format. If you want the data returned on call in a wide format, specify `return_wide = TRUE` and you will have tidy data. If you will be working with dates, whether for plotting purposes or otherwise, then activate the `POSIXct = TRUE` switch. Otherwise you will have to do this manually. ```r wb( country = c("IND", "CHN"), indicator = c("SP.POP.TOTL"), startdate = 1960, enddate = 2016, return_wide = TRUE, POSIXct = TRUE ) -> pop_data5 DT::datatable(pop_data5) ```
--- .pull-left[ ```r library(scales) ggplot(pop_data5, aes(x = date_ct, y = SP.POP.TOTL)) + geom_line(aes(color = country)) + scale_y_continuous(labels = comma) + scale_x_date(date_breaks = "10 years") + theme(legend.position = "bottom") + labs(x = "Date", y = "Total Population") + theme_minimal() ``` ] .pull-right[ <img src="Module06_files/figure-html/wbstats06b2-1.png" width="100%" style="display: block; margin: auto;" /> ] --- class: inverse, middle, center # .fancy[.salt[ USGS data with `dataRetrieval` ]] --- The `dataRetrieval` package gives you easy access to water data gathered and warehoused by the USGS, USDA, EPA, and other entities. The package has an excellent [tutorial available here](https://owi.usgs.gov/R/dataRetrieval.html#1) so I will not go into too many details and nuances here. Start by installing the package and then loading it. ```r library(dataRetrieval) ``` You will need to know the site(s) you are interested in, as well as the parameters and statistics of interest. The package comes with a built-in data-set that shows you the parameters available, and the complete list of statistics [is available here](https://help.waterdata.usgs.gov/code/stat_cd_nm_query?stat_nm_cd=%25&fmt=html). Sites can be located [from this inventory](https://waterdata.usgs.gov/nwis/inventory). ```r parameterCdFile -> parameterCdFile names(parameterCdFile) ``` ``` ## [1] "parameter_cd" "parameter_group_nm" "parameter_nm" "casrn" ## [5] "srsname" "parameter_units" ``` --- If you are curious about a specific parameter, you can see what all is available for it. I'll look for anything related to the keyword `storm`, and also what parameter units are available. ```r parameterCdFile[ grep("storm", parameterCdFile$parameter_nm, ignore.case = TRUE), ] -> stormq unique(stormq$parameter_units) ``` ``` ## [1] "nu" "hours" "minutes" "Mgal" "ft3/s" "mgd" "in" ``` Let us do a quick grab of some data for the Hocking River at Athens, Ohio. .pull-left[ ```r "03159500" -> siteNo "00065" -> pCode "2014-10-01" -> start.date "2018-02-26" -> end.date readNWISuv( siteNumbers = siteNo, parameterCd = pCode, startDate = start.date, endDate = end.date ) -> hocking ``` ] .pull-right[ ``` ## Observations: 54,883 ## Variables: 6 ## $ agency_cd <chr> "USGS", "USGS", "USGS", "USGS", "USGS", "USGS", "USGS", "USGS"… ## $ site_no <chr> "03159500", "03159500", "03159500", "03159500", "03159500", "0… ## $ dateTime <dttm> 2014-10-01 04:00:00, 2014-10-01 04:30:00, 2014-10-01 05:00:00… ## $ X_00065_00000 <dbl> 2.91, 2.91, 2.91, 2.91, 2.91, 2.91, 2.91, 2.91, 2.91, 2.91, 2.… ## $ X_00065_00000_cd <chr> "A", "A", "A", "A", "A", "A", "A", "A", "A", "A", "A", "A", "A… ## $ tz_cd <chr> "UTC", "UTC", "UTC", "UTC", "UTC", "UTC", "UTC", "UTC", "UTC",… ``` ] --- Since the column names are based on parameter codes and hence cryptic, you can clean them up, and also see other attributes embedded in the data-set. ```r names(hocking) ``` ``` ## [1] "agency_cd" "site_no" "dateTime" "X_00065_00000" ## [5] "X_00065_00000_cd" "tz_cd" ``` ```r renameNWISColumns(hocking) -> hocking names(hocking) ``` ``` ## [1] "agency_cd" "site_no" "dateTime" "GH_Inst" "GH_Inst_cd" "tz_cd" ``` ```r names(attributes(hocking)) ``` ``` ## [1] "names" "row.names" "class" "url" "siteInfo" ## [6] "variableInfo" "disclaimer" "statisticInfo" "queryTime" ``` --- If I wanted data for multiple sites, I could find the site numbers and then grab the data. ```r c("03158200", "03159246") -> sites "00065" -> pcode "2014-10-01" -> start.date "2018-02-26" -> end.date readNWISuv( siteNumbers = sites, parameterCd = pcode, startDate = start.date, endDate = end.date ) -> hocking2 renameNWISColumns(hocking2) -> hocking2 ``` --- Now a simple `time-series plot of gage height` for both sites. Note that although I asked for data going far back, not all sites have data for all time periods; it helps to check the site inventory first. .pull-left[ ```r attr(hocking2, "variableInfo") -> parameterInfo ifelse( hocking2$site_no == "03158200", "Monday Creek at Doanville", "Sunday Creek Below Millfield" ) -> hocking2$station as.Date( as.character(hocking2$dateTime), format = "%Y-%m-%d" ) -> hocking2$mydates ggplot( data = hocking2, aes(x = mydates, y = GH_Inst, color = station) ) + geom_line() + labs( x = "", y = parameterInfo$variableDescription) + scale_x_date(date_breaks = "24 weeks") + theme_minimal() + theme(legend.position = "bottom") ``` ] .pull-right[ <img src="Module06_files/figure-html/usgs07b-1.png" width="100%" style="display: block; margin: auto;" /> ] --- class: right, middle <img class="circle" src="https://github.com/aniruhil.png" width="175px"/> # Find me at... [<svg style="height:0.8em;top:.04em;position:relative;" viewBox="0 0 512 512"><path d="M459.37 151.716c.325 4.548.325 9.097.325 13.645 0 138.72-105.583 298.558-298.558 298.558-59.452 0-114.68-17.219-161.137-47.106 8.447.974 16.568 1.299 25.34 1.299 49.055 0 94.213-16.568 130.274-44.832-46.132-.975-84.792-31.188-98.112-72.772 6.498.974 12.995 1.624 19.818 1.624 9.421 0 18.843-1.3 27.614-3.573-48.081-9.747-84.143-51.98-84.143-102.985v-1.299c13.969 7.797 30.214 12.67 47.431 13.319-28.264-18.843-46.781-51.005-46.781-87.391 0-19.492 5.197-37.36 14.294-52.954 51.655 63.675 129.3 105.258 216.365 109.807-1.624-7.797-2.599-15.918-2.599-24.04 0-57.828 46.782-104.934 104.934-104.934 30.213 0 57.502 12.67 76.67 33.137 23.715-4.548 46.456-13.32 66.599-25.34-7.798 24.366-24.366 44.833-46.132 57.827 21.117-2.273 41.584-8.122 60.426-16.243-14.292 20.791-32.161 39.308-52.628 54.253z"/></svg> @aruhil](http://twitter.com/aruhil) [<svg style="height:0.8em;top:.04em;position:relative;" viewBox="0 0 512 512"><path d="M326.612 185.391c59.747 59.809 58.927 155.698.36 214.59-.11.12-.24.25-.36.37l-67.2 67.2c-59.27 59.27-155.699 59.262-214.96 0-59.27-59.26-59.27-155.7 0-214.96l37.106-37.106c9.84-9.84 26.786-3.3 27.294 10.606.648 17.722 3.826 35.527 9.69 52.721 1.986 5.822.567 12.262-3.783 16.612l-13.087 13.087c-28.026 28.026-28.905 73.66-1.155 101.96 28.024 28.579 74.086 28.749 102.325.51l67.2-67.19c28.191-28.191 28.073-73.757 0-101.83-3.701-3.694-7.429-6.564-10.341-8.569a16.037 16.037 0 0 1-6.947-12.606c-.396-10.567 3.348-21.456 11.698-29.806l21.054-21.055c5.521-5.521 14.182-6.199 20.584-1.731a152.482 152.482 0 0 1 20.522 17.197zM467.547 44.449c-59.261-59.262-155.69-59.27-214.96 0l-67.2 67.2c-.12.12-.25.25-.36.37-58.566 58.892-59.387 154.781.36 214.59a152.454 152.454 0 0 0 20.521 17.196c6.402 4.468 15.064 3.789 20.584-1.731l21.054-21.055c8.35-8.35 12.094-19.239 11.698-29.806a16.037 16.037 0 0 0-6.947-12.606c-2.912-2.005-6.64-4.875-10.341-8.569-28.073-28.073-28.191-73.639 0-101.83l67.2-67.19c28.239-28.239 74.3-28.069 102.325.51 27.75 28.3 26.872 73.934-1.155 101.96l-13.087 13.087c-4.35 4.35-5.769 10.79-3.783 16.612 5.864 17.194 9.042 34.999 9.69 52.721.509 13.906 17.454 20.446 27.294 10.606l37.106-37.106c59.271-59.259 59.271-155.699.001-214.959z"/></svg> aniruhil.org](https://aniruhil.org) [<svg style="height:0.8em;top:.04em;position:relative;" viewBox="0 0 512 512"><path d="M476 3.2L12.5 270.6c-18.1 10.4-15.8 35.6 2.2 43.2L121 358.4l287.3-253.2c5.5-4.9 13.3 2.6 8.6 8.3L176 407v80.5c0 23.6 28.5 32.9 42.5 15.8L282 426l124.6 52.2c14.2 6 30.4-2.9 33-18.2l72-432C515 7.8 493.3-6.8 476 3.2z"/></svg> ruhil@ohio.edu](mailto:ruhil@ohio.edu)