Report 7

Introduction 

I selected Hamilton County in Tennessee because it is the aera where I grew up. For my analysis, I examined educational attainment, specifically the percentage of residents who hold a bachelor’s degree or higher. I was interested in this variable because education is often associated with income and economic opportunities. I expected areas with higher levels of educational attainment to also have higher household incomes. 

County-Level Analysis 

At the county level, Hamilton County had a higher educational attainment rate than its neighbors. Approximately 37.5% of Hamilton County residents hold a bachelor’s degree or higher, while the surrounding counties have considerably lower rates: Bradley County (24.4%), Marion County (17.1%), Sequatchie County (20.6%), Bledsoe County (12.4%), Rhea County (18.6%), and Meigs County (13.5%). 

This map helped show how Hamilton County compared with nearby counties, but it did not reveal the variation that exists within the county itself. While this data were useful for regional comparisons, they masked important differences between communities within the city of Chattanooga. 

Within-County Analysis 

The sub-county map revealed significant variation across Hamilton County. Some areas had much higher percentages of residents with a bachelor’s degree, while others had considerably lower percentages. One result that stood out to me was Ooltewah (District 9 and where I live), where approximately 38% of residents held a bachelor’s degree or higher. Given Ooltewah’s reputation as a relatively affluent area of Hamilton County, I expected the percentage to be much higher. 

The highest levels of educational attainment were found in Districts 6 and 2, which include affluent communities such as Signal Mountain and Lookout Mountain. District 6 had the highest percentage at 54.4%. In contrast, the lowest educational attainment was found in District 4, which includes much of downtown Chattanooga, at just 20.5%. 

A geographic pattern became apparent when examining the map. In general, educational attainment tended to increase as the distance from downtown Chattanooga increased. Areas farther from downtown Chattanooga, such as Signal Mountain and Lookout Mountain, had higher bachelor’s degree attainment rates, while areas closer to downtown generally had lower educational attainment levels.

To better understand Ooltewah’s results, I also compared educational attainment data with median household income. District 9 had the second-highest median household income in the county at approximately $92,000, trailing only District 7, which had a median household income of $101,390. What surprised me was that Ooltewah’s income level was actually more than Signal Mountain’s ($90,653) and Lookout Mountain’s (73,561), although Ooltewah’s bachelor’s degree attainment level was lower than those two. This challenged my assumption that wealthier areas would consequently have a higher percentage of college graduates. 

One possible explanation is that District 9 includes not only Ooltewah but also the Harrison and Highway 58 areas. Being a Chattanooga native, I can say those areas are generally less affluent than Ooltewah. As a result, combining these communities into a single district may lower the overall percentage of residents with bachelor’s degrees and help explain why District 9 differs from some of the other high-income areas in the county. 

Reflection  

Before beginning the analysis, I expected affluent areas like Ooltewah to have some of the highest percentages of residents with bachelor’s degrees in Hamilton County and to closely resemble communities such as Signal Mountain and Lookout Mountain. While the data partially supported that expectation, the results were not as straightforward as I anticipated. 

The combination of the bachelor degree attainment rate and a high median household income showed that educational attainment alone may not fully explain economic outcomes. Although education and income are often related, the relationship was not as strong as I expected in this case, especially in my area of Ooltewah. 

The sub-county map provided information that was hidden in the county-level map. Looking only at Hamilton County as a whole would have suggested a single overall value for educational attainment, while the sub-county analysis revealed significant differences between communities. 

If I were a journalist, I would pursue a story looking at the relationship between educational attainment and income in Hamilton County. Specifically, I would explore why some areas maintain relatively high household incomes despite not having the highest levels of education. Such a story could provide insight into local industries, career paths, and economic opportunities available to residents. 

 

Code

############################################################
# STEP 1: Load Packages
############################################################

if (!require("tidycensus")) install.packages("tidycensus")
if (!require("tidyverse")) install.packages("tidyverse")
if (!require("kableExtra")) install.packages("kableExtra")
if (!require("sf")) install.packages("sf")
if (!require("leaflet")) install.packages("leaflet")
if (!require("htmlwidgets")) install.packages("htmlwidgets")
if (!require("plotly")) install.packages("plotly")

library(tidycensus)
library(tidyverse)
library(kableExtra)
library(sf)
library(leaflet)
library(htmlwidgets)
library(plotly)

############################################################
# STEP 2: Verify Census API Key
############################################################

if (!nzchar(Sys.getenv("CENSUS_API_KEY"))) {
  
  stop(
    paste(
      "\nNo Census API key was found.",
      "\n\nTo obtain a free Census API key, visit:",
      "\nhttps://api.census.gov/data/key_signup.html",
      "\n\nThen run the one-time setup code at the",
      "\ntop of this script to store your key.",
      "\n"
    )
  )
  
}

############################################################
# STEP 3: Explore 2024 Codebooks and Choose an ACS Variable
############################################################

ProfileTables <- load_variables(2024, "acs5/profile")

# DetailedTables <- load_variables(2024, "acs5")
# SubjectTables  <- load_variables(2024, "acs5/subject")

# Replace the variable code below with the ACS variable
# you want to download.

my_variable <- "DP02_0068P"

# Example variables:
#
# DP05_0001  = Total population
# DP03_0062  = Median Household Income
# DP03_0128P = Poverty Rate
# DP02_0068P = Bachelor's Degree or Higher Rate
# DP05_0018  = Median Age
# DP04_0046P = Homeownership Rate
# DP04_0134  = Median Gross Rent
# DP04_0142P = Pct. renters paying 35%+ HH income in rent (GRAPI)
# DP05_0090P = Hispanic or Latino Population Percentage
# DP02_0094P = Foreign-Born Percentage
# DP03_0009P = Unemployment Rate
# DP03_0025  = Average Commute Time

############################################################
# STEP 4: Download ACS Data
############################################################

# Common ACS geography options:
#
# "state"
# "county"
# "county subdivision"
# "tract"
# "block group"
# "zcta"
#
# Most of this script works for any ACS geography.
# Some geographies may require additional filtering
# or display adjustments.

mydata <- get_acs(
  geography = "county subdivision",
  state = "TN",
  year = 2024,
  survey = "acs5",
  variables = my_variable,
  geometry = TRUE
)

############################################################
# STEP 5: Keep Needed Variables and Geometry
############################################################

mydata <- mydata |>
  select(
    Area = NAME,
    Variable = variable,
    Estimate = estimate,
    Margin_of_Error = moe,
    geometry
  )

############################################################
# STEP 5B: Transform Coordinates for Leaflet
############################################################

mydata <- st_transform(
  mydata,
  crs = 4326
)

############################################################
# STEP 6: Optional Filtering
############################################################

# OPTION A: Filter for specific counties.
#
# Uncomment to keep only county subdivisions whose
# names contain the specified county names. To search
# for just a single county, use just the county name
# without the | character.

mydata <- mydata |>
  filter(
    str_detect(
      Area,
      "Hamilton"
    )
  )

############################################################
# STEP 7: Sort Results from Highest to Lowest
############################################################

mydata <- mydata |>
  arrange(desc(Estimate))

############################################################
# STEP 8: Display Results as a Table
############################################################

ResultsTable <- mydata |>
  st_drop_geometry() |>
  kbl(
    col.names = c(
      "Area",
      "Variable",
      "Estimate",
      "Margin of Error"
    ),
    format.args = list(big.mark = ",")
  ) |>
  kable_styling(
    bootstrap_options = c("striped", "hover"),
    full_width = FALSE
  )

ResultsTable

############################################################
# STEP 9: Create Pop-Up Content
############################################################

mydata$popup <- paste0(
  "", mydata$Area, "
",
"Estimate: ",
format(mydata$Estimate, big.mark = ","),
"
",
"Plus/Minus: ",
format(mydata$Margin_of_Error, big.mark = ",")
)

############################################################
# STEP 10: Build Color Palette
############################################################

qs <- unique(
  quantile(
    mydata$Estimate,
    probs = seq(0, 1, length.out = 6),
    na.rm = TRUE
  )
)

pal <- colorBin(
  palette = "Blues",
  domain = mydata$Estimate,
  bins = qs,
  pretty = FALSE
)

############################################################
# STEP 11: Create Choropleth Map
############################################################

AreaMap <- leaflet(mydata) |>
  addProviderTiles(
    providers$Esri.WorldTopoMap
  ) |>
  addPolygons(
    fillColor = ~pal(Estimate),
    fillOpacity = 0.5,
    color = "black",
    weight = 1,
    popup = ~popup
  ) |>
  addLegend(
    pal = pal,
    values = ~Estimate,
    title = "Estimate",
    labFormat = labelFormat(
      big.mark = ","
    )
  )

AreaMap

############################################################
# STEP 12: Create Plotly Dot Plot
############################################################

plotdf <- mydata |>
  st_drop_geometry() |>
  mutate(
    point_color = pal(Estimate),
    
    y_ordered = reorder(
      Area,
      Estimate
    ),
    
    hover_text = paste0(
      Area,
      "
Estimate: ",
      format(Estimate, big.mark = ","),
      "
Plus/Minus: ",
      format(Margin_of_Error, big.mark = ",")
    )
  )

DotPlot <- plot_ly(
  data = plotdf,
  x = ~Estimate,
  y = ~as.character(y_ordered),
  type = "scatter",
  mode = "markers",
  showlegend = FALSE,
  marker = list(
    color = ~point_color,
    size = 8,
    line = list(
      color = "rgba(120,120,120,0.9)",
      width = 0.5
    )
  ),
  error_x = list(
    type = "data",
    array = ~Margin_of_Error,
    arrayminus = ~Margin_of_Error,
    color = "rgba(0,0,0,0.65)",
    thickness = 1
  ),
  text = ~hover_text,
  hovertemplate = "%{text}"
) |>
  layout(
    title = list(
      text = paste0(
        "ACS Estimates by State",
        "
Error bars show ACS margins of error."
      )
    ),
    xaxis = list(
      title = "Estimate",
      tickformat = ",.0f",
      automargin = TRUE
    ),
    yaxis = list(
      title = "",
      automargin = TRUE,
      categoryorder = "array",
      categoryarray = levels(plotdf$y_ordered)
    ),
    margin = list(
      l = 200,
      r = 20,
      b = 60,
      t = 60,
      pad = 2
    )
  )

DotPlot

############################################################
# STEP 13: Export Dot Plot (Optional)
############################################################

saveWidget(
  widget = as_widget(DotPlot),
  file = "ACSGraph.html",
  selfcontained = TRUE
)

############################################################
# STEP 14: Export Map (Optional)
############################################################

saveWidget(
  widget = AreaMap,
  file = "ACSMap.html",
  selfcontained = TRUE
)

############################################################
# STEP 1: Load Packages
############################################################

if (!require("tidycensus")) install.packages("tidycensus")
if (!require("tidyverse")) install.packages("tidyverse")
if (!require("kableExtra")) install.packages("kableExtra")
if (!require("sf")) install.packages("sf")
if (!require("leaflet")) install.packages("leaflet")
if (!require("htmlwidgets")) install.packages("htmlwidgets")
if (!require("plotly")) install.packages("plotly")

library(tidycensus)
library(tidyverse)
library(kableExtra)
library(sf)
library(leaflet)
library(htmlwidgets)
library(plotly)

############################################################
# STEP 2: Verify Census API Key
############################################################

if (!nzchar(Sys.getenv("CENSUS_API_KEY"))) {
  
  stop(
    paste(
      "\nNo Census API key was found.",
      "\n\nTo obtain a free Census API key, visit:",
      "\nhttps://api.census.gov/data/key_signup.html",
      "\n\nThen run the one-time setup code at the",
      "\ntop of this script to store your key.",
      "\n"
    )
  )
  
}

############################################################
# STEP 3: Explore 2024 Codebooks and Choose an ACS Variable
############################################################

ProfileTables <- load_variables(2024, "acs5/profile")

# DetailedTables <- load_variables(2024, "acs5")
# SubjectTables  <- load_variables(2024, "acs5/subject")

# Replace the variable code below with the ACS variable
# you want to download.

my_variable <- "DP03_0062"

# Example variables:
#
# DP05_0001  = Total population
# DP03_0062  = Median Household Income
# DP03_0128P = Poverty Rate
# DP02_0068P = Bachelor's Degree or Higher Rate
# DP05_0018  = Median Age
# DP04_0046P = Homeownership Rate
# DP04_0134  = Median Gross Rent
# DP04_0142P = Pct. renters paying 35%+ HH income in rent (GRAPI)
# DP05_0090P = Hispanic or Latino Population Percentage
# DP02_0094P = Foreign-Born Percentage
# DP03_0009P = Unemployment Rate
# DP03_0025  = Average Commute Time

############################################################
# STEP 4: Download ACS Data
############################################################

# Common ACS geography options:
#
# "state"
# "county"
# "county subdivision"
# "tract"
# "block group"
# "zcta"
#
# Most of this script works for any ACS geography.
# Some geographies may require additional filtering
# or display adjustments.

mydata <- get_acs(
  geography = "county subdivision",
  state = "TN",
  year = 2024,
  survey = "acs5",
  variables = my_variable,
  geometry = TRUE
)

############################################################
# STEP 5: Keep Needed Variables and Geometry
############################################################

mydata <- mydata |>
  select(
    Area = NAME,
    Variable = variable,
    Estimate = estimate,
    Margin_of_Error = moe,
    geometry
  )

############################################################
# STEP 5B: Transform Coordinates for Leaflet
############################################################

mydata <- st_transform(
  mydata,
  crs = 4326
)

############################################################
# STEP 6: Optional Filtering
############################################################

# OPTION A: Filter for specific counties.
#
# Uncomment to keep only county subdivisions whose
# names contain the specified county names. To search
# for just a single county, use just the county name
# without the | character.

mydata <- mydata |>
  filter(
    str_detect(
      Area,
      "Hamilton"
    )
  )

############################################################
# STEP 7: Sort Results from Highest to Lowest
############################################################

mydata <- mydata |>
  arrange(desc(Estimate))

############################################################
# STEP 8: Display Results as a Table
############################################################

ResultsTable <- mydata |>
  st_drop_geometry() |>
  kbl(
    col.names = c(
      "Area",
      "Variable",
      "Estimate",
      "Margin of Error"
    ),
    format.args = list(big.mark = ",")
  ) |>
  kable_styling(
    bootstrap_options = c("striped", "hover"),
    full_width = FALSE
  )

ResultsTable

############################################################
# STEP 9: Create Pop-Up Content
############################################################

mydata$popup <- paste0(
  "", mydata$Area, "
",
"Estimate: ",
format(mydata$Estimate, big.mark = ","),
"
",
"Plus/Minus: ",
format(mydata$Margin_of_Error, big.mark = ",")
)

############################################################
# STEP 10: Build Color Palette
############################################################

qs <- unique(
  quantile(
    mydata$Estimate,
    probs = seq(0, 1, length.out = 6),
    na.rm = TRUE
  )
)

pal <- colorBin(
  palette = "Blues",
  domain = mydata$Estimate,
  bins = qs,
  pretty = FALSE
)

############################################################
# STEP 11: Create Choropleth Map
############################################################

HamiltonIncomeMap <- leaflet(mydata) |>
  addProviderTiles(
    providers$Esri.WorldTopoMap
  ) |>
  addPolygons(
    fillColor = ~pal(Estimate),
    fillOpacity = 0.5,
    color = "black",
    weight = 1,
    popup = ~popup
  ) |>
  addLegend(
    pal = pal,
    values = ~Estimate,
    title = "Estimate",
    labFormat = labelFormat(
      big.mark = ","
    )
  )

HamiltonIncomeMap

############################################################
# STEP 12: Create Plotly Dot Plot
############################################################

plotdf <- mydata |>
  st_drop_geometry() |>
  mutate(
    point_color = pal(Estimate),
    
    y_ordered = reorder(
      Area,
      Estimate
    ),
    
    hover_text = paste0(
      Area,
      "
Estimate: ",
      format(Estimate, big.mark = ","),
      "
Plus/Minus: ",
      format(Margin_of_Error, big.mark = ",")
    )
  )

DotPlot <- plot_ly(
  data = plotdf,
  x = ~Estimate,
  y = ~as.character(y_ordered),
  type = "scatter",
  mode = "markers",
  showlegend = FALSE,
  marker = list(
    color = ~point_color,
    size = 8,
    line = list(
      color = "rgba(120,120,120,0.9)",
      width = 0.5
    )
  ),
  error_x = list(
    type = "data",
    array = ~Margin_of_Error,
    arrayminus = ~Margin_of_Error,
    color = "rgba(0,0,0,0.65)",
    thickness = 1
  ),
  text = ~hover_text,
  hovertemplate = "%{text}"
) |>
  layout(
    title = list(
      text = paste0(
        "ACS Estimates by State",
        "
Error bars show ACS margins of error."
      )
    ),
    xaxis = list(
      title = "Estimate",
      tickformat = ",.0f",
      automargin = TRUE
    ),
    yaxis = list(
      title = "",
      automargin = TRUE,
      categoryorder = "array",
      categoryarray = levels(plotdf$y_ordered)
    ),
    margin = list(
      l = 200,
      r = 20,
      b = 60,
      t = 60,
      pad = 2
    )
  )

DotPlot

############################################################
# STEP 13: Export Dot Plot (Optional)
############################################################

saveWidget(
  widget = as_widget(DotPlot),
  file = "ACSGraph.html",
  selfcontained = TRUE
)

############################################################
# STEP 14: Export Map (Optional)
############################################################

saveWidget(
  widget = AreaMap,
  file = "ACSMap.html",
  selfcontained = TRUE
)

############################################################
# STEP 1: Load Packages
############################################################

if (!require("tidycensus")) install.packages("tidycensus")
if (!require("tidyverse")) install.packages("tidyverse")
if (!require("kableExtra")) install.packages("kableExtra")
if (!require("sf")) install.packages("sf")
if (!require("leaflet")) install.packages("leaflet")
if (!require("htmlwidgets")) install.packages("htmlwidgets")
if (!require("plotly")) install.packages("plotly")

library(tidycensus)
library(tidyverse)
library(kableExtra)
library(sf)
library(leaflet)
library(htmlwidgets)
library(plotly)

############################################################
# STEP 2: Verify Census API Key
############################################################

if (!nzchar(Sys.getenv("CENSUS_API_KEY"))) {
  
  stop(
    paste(
      "\nNo Census API key was found.",
      "\n\nTo obtain a free Census API key, visit:",
      "\nhttps://api.census.gov/data/key_signup.html",
      "\n\nThen run the one-time setup code at the",
      "\ntop of this script to store your key.",
      "\n"
    )
  )
  
}

############################################################
# STEP 3: Explore 2024 Codebooks and Choose an ACS Variable
############################################################

ProfileTables <- load_variables(2024, "acs5/profile")

# DetailedTables <- load_variables(2024, "acs5")
# SubjectTables  <- load_variables(2024, "acs5/subject")

# Replace the variable code below with the ACS variable
# you want to download.

my_variable <- "DP03_0062"

# Example variables:
#
# DP05_0001  = Total population
# DP03_0062  = Median Household Income
# DP03_0128P = Poverty Rate
# DP02_0068P = Bachelor's Degree or Higher Rate
# DP05_0018  = Median Age
# DP04_0046P = Homeownership Rate
# DP04_0134  = Median Gross Rent
# DP04_0142P = Pct. renters paying 35%+ HH income in rent (GRAPI)
# DP05_0090P = Hispanic or Latino Population Percentage
# DP02_0094P = Foreign-Born Percentage
# DP03_0009P = Unemployment Rate
# DP03_0025  = Average Commute Time

############################################################
# STEP 4: Download ACS Data
############################################################

# Common ACS geography options:
#
# "state"
# "county"
# "county subdivision"
# "tract"
# "block group"
# "zcta"
#
# Most of this script works for any ACS geography.
# Some geographies may require additional filtering
# or display adjustments.

mydata <- get_acs(
  geography = "county subdivision",
  state = "TN",
  year = 2024,
  survey = "acs5",
  variables = my_variable,
  geometry = TRUE
)

############################################################
# STEP 5: Keep Needed Variables and Geometry
############################################################

mydata <- mydata |>
  select(
    Area = NAME,
    Variable = variable,
    Estimate = estimate,
    Margin_of_Error = moe,
    geometry
  )

############################################################
# STEP 5B: Transform Coordinates for Leaflet
############################################################

mydata <- st_transform(
  mydata,
  crs = 4326
)

############################################################
# STEP 6: Optional Filtering
############################################################

# OPTION A: Filter for specific counties.
#
# Uncomment to keep only county subdivisions whose
# names contain the specified county names. To search
# for just a single county, use just the county name
# without the | character.

mydata <- mydata |>
  filter(
    str_detect(
      Area,
      "Hamilton"
    )
  )

############################################################
# STEP 7: Sort Results from Highest to Lowest
############################################################

mydata <- mydata |>
  arrange(desc(Estimate))

############################################################
# STEP 8: Display Results as a Table
############################################################

ResultsTable <- mydata |>
  st_drop_geometry() |>
  kbl(
    col.names = c(
      "Area",
      "Variable",
      "Estimate",
      "Margin of Error"
    ),
    format.args = list(big.mark = ",")
  ) |>
  kable_styling(
    bootstrap_options = c("striped", "hover"),
    full_width = FALSE
  )

ResultsTable

############################################################
# STEP 9: Create Pop-Up Content
############################################################

mydata$popup <- paste0(
  "", mydata$Area, "
",
"Estimate: ",
format(mydata$Estimate, big.mark = ","),
"
",
"Plus/Minus: ",
format(mydata$Margin_of_Error, big.mark = ",")
)

############################################################
# STEP 10: Build Color Palette
############################################################

qs <- unique(
  quantile(
    mydata$Estimate,
    probs = seq(0, 1, length.out = 6),
    na.rm = TRUE
  )
)

pal <- colorBin(
  palette = "Blues",
  domain = mydata$Estimate,
  bins = qs,
  pretty = FALSE
)

############################################################
# STEP 11: Create Choropleth Map
############################################################

HamiltonIncomeMap <- leaflet(mydata) |>
  addProviderTiles(
    providers$Esri.WorldTopoMap
  ) |>
  addPolygons(
    fillColor = ~pal(Estimate),
    fillOpacity = 0.5,
    color = "black",
    weight = 1,
    popup = ~popup
  ) |>
  addLegend(
    pal = pal,
    values = ~Estimate,
    title = "Estimate",
    labFormat = labelFormat(
      big.mark = ","
    )
  )

HamiltonIncomeMap

############################################################
# STEP 12: Create Plotly Dot Plot
############################################################

plotdf <- mydata |>
  st_drop_geometry() |>
  mutate(
    point_color = pal(Estimate),
    
    y_ordered = reorder(
      Area,
      Estimate
    ),
    
    hover_text = paste0(
      Area,
      "
Estimate: ",
      format(Estimate, big.mark = ","),
      "
Plus/Minus: ",
      format(Margin_of_Error, big.mark = ",")
    )
  )

DotPlot <- plot_ly(
  data = plotdf,
  x = ~Estimate,
  y = ~as.character(y_ordered),
  type = "scatter",
  mode = "markers",
  showlegend = FALSE,
  marker = list(
    color = ~point_color,
    size = 8,
    line = list(
      color = "rgba(120,120,120,0.9)",
      width = 0.5
    )
  ),
  error_x = list(
    type = "data",
    array = ~Margin_of_Error,
    arrayminus = ~Margin_of_Error,
    color = "rgba(0,0,0,0.65)",
    thickness = 1
  ),
  text = ~hover_text,
  hovertemplate = "%{text}"
) |>
  layout(
    title = list(
      text = paste0(
        "ACS Estimates by State",
        "
Error bars show ACS margins of error."
      )
    ),
    xaxis = list(
      title = "Estimate",
      tickformat = ",.0f",
      automargin = TRUE
    ),
    yaxis = list(
      title = "",
      automargin = TRUE,
      categoryorder = "array",
      categoryarray = levels(plotdf$y_ordered)
    ),
    margin = list(
      l = 200,
      r = 20,
      b = 60,
      t = 60,
      pad = 2
    )
  )

DotPlot

############################################################
# STEP 13: Export Dot Plot (Optional)
############################################################

saveWidget(
  widget = as_widget(DotPlot),
  file = "ACSGraph.html",
  selfcontained = TRUE
)

############################################################
# STEP 14: Export Map (Optional)
############################################################

saveWidget(
  widget = AreaMap,
  file = "ACSMap.html",
  selfcontained = TRUE
)