library(tidycensus)
library(tidyverse)
library(sf)
library(ggplot2)
v15_Profile <- load_variables(2017 , "acs5/profile", cache = TRUE) #demographic profile tables

View(v15_Profile)
## Warning in system2("/usr/bin/otool", c("-L", shQuote(DSO)), stdout = TRUE):
## running command ''/usr/bin/otool' -L '/Library/Frameworks/R.framework/Resources/
## modules/R_de.so'' had status 1
#Search for variables by keywords in the label 
v15_Profile[grep(x = v15_Profile$label, "MEDHOUSINC"), c("name", "label")]
## # A tibble: 0 x 2
## # … with 2 variables: name <chr>, label <chr>
v15_Profile[grep(x = v15_Profile$label, "Built 2000 to 2009"), c("name", "label")]
## # A tibble: 2 x 2
##   name       label                                                              
##   <chr>      <chr>                                                              
## 1 DP04_0019  Estimate!!YEAR STRUCTURE BUILT!!Total housing units!!Built 2000 to…
## 2 DP04_0019P Percent Estimate!!YEAR STRUCTURE BUILT!!Total housing units!!Built…

Extract from ACS summary file data

Here I get data profile variables from 2017 for Bexar County, TX Census Tracts

The data profile tables are very useful because they contain lots of pre-calculated variables.

Here is a query where we extract several variables from tehe 2017 ACS for Bexar County, Texas. We can also get the spatial data by requesting geometry=TRUE. Useing output="wide" will put each variable in a column of the data set, with each row being a census tract.

sa_acs<-get_acs(geography = "tract",
                state="TX",
                county = c("Bexar"),
                year = 2017,
                variables=c( "DP05_0001E", 
                            "DP03_0062E") ,
                geometry = T, output = "wide")
## Getting data from the 2013-2017 5-year ACS
## Downloading feature geometry from the Census website.  To cache shapefiles for use in future sessions, set `options(tigris_use_cache = TRUE)`.
## Using the ACS Data Profile
#create a county FIPS code - 5 digit
sa_acs$county<-substr(sa_acs$GEOID, 1, 5)

#rename variables and filter missing cases
sa_acs2<-sa_acs%>%
  mutate(totpop= DP05_0001E, hmed=DP03_0062E) %>%
#  st_transform(crs = 102740)%>%
  na.omit()

Write data out to shapefile

#change the directory
sf::st_write(sa_acs2,dsn="~/Desktop/SPRING_2020/gis_class/acs_data", driver="ESRI Shapefile", delete_layer=T, update=T)

Some basic mapping of variables

Here I generate a quantile break for % black in census tracts and compare it to a Jenks break

library(classInt)
library(patchwork)
library(dplyr)

hmed_map<-sa_acs2 %>%
  mutate(cpov=cut(hmed,breaks = quantile(hmed, na.rm=T, p=seq(0,1,length.out = 6)),include.lowest = T),
           jpov = cut(hmed,breaks=data.frame(classIntervals(var=sa_acs2$hmed, n=5, style="pretty")[2])[,1], include.lowest = T))

library(ggsn)





p2<-ggplot(hmed_map, aes(fill = jpov, color = jpov)) + 
  geom_sf() + 
  ggtitle("Medium Household Income", 
          subtitle = "Bexar County Texas, 2017 - Pretty Breaks")+
  scale_fill_brewer(palette = "Reds") + 
  scale_color_brewer(palette = "Reds")+
    theme(axis.text.x = element_blank(), axis.text.y = element_blank())+
  north(hmed_map)+
  scalebar(hmed_map, location="bottomleft", dist=5, transform = T,dist_unit = "km",  model="WGS84", st.size =2)
p2

ggsave(filename="~/Desktop/SPRING_2020/gis_class/acs_data/lap1map1.png")