Exploring the Seasonality of Unemployment in Northern Michigan
Gerrit Anderson
2024-03-12
Introduction
In this report, I explore BLS LAUS data - which reports unemployment statistics on a county level for the United States over a 14 month period. The data was downloaded from the latest LAUS release, which was March 1st, 2024.
For this lesson, I am specifically interested in exploring the seasonality of unemployment in a selection of Northern Michigan Counties. The counties I am interested in are Antrim, Benzie, Charlevoix, Emmet, Grand Traverse, Kalkaska, Leelanau, Manistee, Missaukee, and Wexford. These ten counties are part of a larger economic region call “Prosperity Region 2”. You can view the Michigan Prosperity Regions on the map below.
Load Packages
library(tidyverse)
library(here)
library(ggbeeswarm)
library(RColorBrewer)
library(lubridate)Selecting Prosperity Region 2 from the data.
#read in the LAUS data that was cleaned and manipulated last week, and select counties in Prosperity Region 2
plotLAUS <- read.csv(here("data", "testLAUS.csv")) %>%
filter(CountyName_StateAbbreviation %in% c("Emmet County, MI", "Charlevoix County, MI", "Antrim County, MI",
"Kalkaska County, MI", "Missaukee County, MI", "Wexford County, MI",
"Grand Traverse County, MI", "Manistee County, MI", "Benzie County, MI",
"Leelanau County, MI" ))
#view the data
glimpse(plotLAUS)## Rows: 126
## Columns: 14
## $ X <int> 1241, 1246, 1255, 1259, 1271, 1276, 1282,…
## $ LAUSCode <chr> "CN2601900000000", "CN2602900000000", "CN…
## $ StateFIPSCode <int> 26, 26, 26, 26, 26, 26, 26, 26, 26, 26, 2…
## $ CountyFIPSCode <int> 19, 29, 47, 55, 79, 89, 101, 113, 165, 19…
## $ CountyName_StateAbbreviation <chr> "Benzie County, MI", "Charlevoix County, …
## $ County <chr> "Benzie County", "Charlevoix County", "Em…
## $ State <chr> " MI", " MI", " MI", " MI", " MI", " MI",…
## $ Period <chr> "Nov-22", "Nov-22", "Nov-22", "Nov-22", "…
## $ LaborForce <chr> "8,586 ", "12,253 ", "16,064 …
## $ Employed <chr> "8,201 ", "11,731 ", "15,263 …
## $ Unemployed <chr> "385 ", "522 ", "801 ", "1,57…
## $ UnemploymentRate <dbl> 4.5, 4.3, 5.0, 3.2, 5.4, 3.6, 5.3, 4.1, 4…
## $ UR_all <lgl> TRUE, TRUE, TRUE, FALSE, TRUE, FALSE, TRU…
## $ UR_State <lgl> FALSE, FALSE, TRUE, FALSE, TRUE, FALSE, T…Data Preperation
Some additional data cleaning had to be done to prepare the columns for plotting.
#convert the period column to a date format
plotLAUS$Period <- my(plotLAUS$Period)
#remove extra spaces from the LaborForce Column
plotLAUS$LaborForce <- trimws(plotLAUS$LaborForce, "r")
#change the LaborForce Column to numeric, removing the column
plotLAUS$LaborForce <- as.numeric(gsub(",","",plotLAUS$LaborForce))
#remove extra spaces from the Unemployed Column
plotLAUS$Unemployed <- trimws(plotLAUS$Unemployed, "r")
#change the Unemployed Column to numeric, removing the column
plotLAUS$Unemployed <- as.numeric(gsub(",","",plotLAUS$Unemployed))Analysis
First, lets just look at how each of these counties unemployment rates look for the most recent period (22-Nov to Dec-23p (p = preliminary data)).
#plot county unemployment rate by period by county
plotLAUS %>%
na.omit() %>%
ggplot(aes(x = Period, y = UnemploymentRate, colour = CountyName_StateAbbreviation)) +
geom_point() +
geom_smooth() +
facet_wrap(~CountyName_StateAbbreviation) +
coord_flip() +
labs(title = "Unemployment Rates",
subtitle = "for selected periods",
caption = "data from https://www.bls.gov/web/metro/laucntycur14.txt",
x = "Period",
y = "Unemployment Rate")
This data is not seasonally adjusted, so it is expected that there would be some peaks and troughs. However, as you can see from the plot, there are some extremely high unemployment rates from December to March for most of the selected counties. An unemployment rate that is doubling from summer to winter is unusual, so there must be something else occurring in the data. Let’s see if the labor force for each county is also drastically changing between the summer and winter months.
#plot county labor force by period by county
plotLAUS %>%
na.omit() %>%
group_by(Period, CountyName_StateAbbreviation) %>%
ggplot(aes(x = Period, y = LaborForce, fill = CountyName_StateAbbreviation)) +
geom_col() +
geom_smooth() +
facet_wrap(~CountyName_StateAbbreviation) +
coord_flip() +
labs(title = "Labor Force",
subtitle = "for selected periods",
caption = "data from https://www.bls.gov/web/metro/laucntycur14.txt",
x = "Period",
y = "Labor Force")
It seems that not all of the counties show a strong seasonality in their labor force (with the main exception being Grand Traverse County). Let’s just take a look at how labor force and unemployment (not unemployment rate), influence each other.
#plot relationship between unemployed and labor force
plotLAUS %>%
na.omit() %>%
ggplot(aes(x = Unemployed, y = LaborForce, color = UnemploymentRate)) +
geom_point() +
geom_smooth() +
theme_classic() +
scale_color_distiller(palette = "RdYlBu") +
labs(title = "Unemployed vs. Labor Force",
subtitle = "for selected periods",
caption = "data from https://www.bls.gov/web/metro/laucntycur14.txt",
x = "Unemployed",
y = "Labor Force")
This plot should make sense to anyone who has worked with unemployment data before. A high number of unemployed people, with a low number of people in the overall labor force, means that there will be a high unemployment rate. While this seems obvious, it is still helpful to see this relationship in the selected data. Moreover, this can also help us see the raw numbers behind those extreme unemployment rates.
Addtional Analysis
From just the LAUS data we have identified that there is some seasonality to the unemployment rates in the selected counties, and we have taken a closer look at the relationship between unemployment and the labor force. However, we have still not fully identified what could be driving this seasonality in unemployment. For this, lets take a look at some industry data from the Bureau of Labor Statistics: Quarterly Census of Employment and Wages. Specifically, we are interested in the establishment location quotient (LQ) of the Amusement, Gambling, and Recreation Industries (NAICS code 713). A LQ is an analytical statistic that measures a region’s industrial specialization relative to a larger geographic unit. In this case the region is each county and the larger geographic unit is the United States. An LQ greater than 1.2 indicates an industry specialization that is stronger than the larger geographic unit.
#Use BLS API example to pull LQs for selection area
#713 Amusement, Gambling, and Recreation Industries
qcewGetIndustryData <- function (year, qtr, industry) {
url <- "http://data.bls.gov/cew/data/api/YEAR/QTR/industry/INDUSTRY.csv"
url <- sub("YEAR", year, url, ignore.case=FALSE)
url <- sub("QTR", tolower(qtr), url, ignore.case=FALSE)
url <- sub("INDUSTRY", industry, url, ignore.case=FALSE)
read.csv(url, header = TRUE, sep = ",", quote="\"", dec=".", na.strings=" ", skip=0)
}#pulling relevant LQs for the past 3 available years for the selected counties
Recreation2020 <- qcewGetIndustryData("2020", "a", "713") %>%
filter(area_fips %in% c("26009", "26019", "26029", "26047", "26055", "26079", "26089", "26101", "26113", "26165"))
Recreation2021 <- qcewGetIndustryData("2021", "a", "713") %>%
filter(area_fips %in% c("26009", "26019", "26029", "26047", "26055", "26079", "26089", "26101", "26113", "26165"))
Recreation2022 <- qcewGetIndustryData("2022", "a", "713") %>%
filter(area_fips %in% c("26009", "26019", "26029", "26047", "26055", "26079", "26089", "26101", "26113", "26165"))
#append each year into one dataframe
RecreationTotal <- rbind(Recreation2020, Recreation2021, Recreation2022)
#check the data types for each column
glimpse(RecreationTotal)## Rows: 30
## Columns: 38
## $ area_fips <chr> "26009", "26019", "26029", "26047", "…
## $ own_code <int> 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5…
## $ industry_code <int> 713, 713, 713, 713, 713, 713, 713, 71…
## $ agglvl_code <int> 75, 75, 75, 75, 75, 75, 75, 75, 75, 7…
## $ size_code <int> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0…
## $ year <int> 2020, 2020, 2020, 2020, 2020, 2020, 2…
## $ qtr <chr> "A", "A", "A", "A", "A", "A", "A", "A…
## $ disclosure_code <chr> "", "", "", "", "", "", "", "", "N", …
## $ annual_avg_estabs <int> 9, 8, 27, 31, 45, 5, 20, 9, 3, 13, 11…
## $ annual_avg_emplvl <int> 139, 74, 213, 405, 375, 27, 100, 172,…
## $ total_annual_wages <dbl> 4057194, 2242111, 7650777, 14786044, …
## $ taxable_annual_wages <dbl> 1776497, 969095, 2423014, 5073476, 42…
## $ annual_contributions <int> 177006, 84959, 139633, 291319, 242066…
## $ annual_avg_wkly_wage <int> 560, 586, 691, 702, 503, 433, 574, 64…
## $ avg_annual_pay <int> 29119, 30470, 35947, 36494, 26155, 22…
## $ lq_disclosure_code <chr> "", "", "", "", "", "", "", "", "N", …
## $ lq_annual_avg_estabs <dbl> 1.83, 1.95, 3.90, 2.46, 1.50, 1.75, 3…
## $ lq_annual_avg_emplvl <dbl> 3.17, 2.06, 2.56, 2.63, 0.90, 0.76, 1…
## $ lq_total_annual_wages <dbl> 6.06, 4.01, 4.64, 5.11, 1.12, 0.79, 3…
## $ lq_taxable_annual_wages <dbl> 4.97, 3.04, 3.35, 3.85, 1.30, 0.85, 2…
## $ lq_annual_contributions <dbl> 10.97, 6.68, 5.37, 5.71, 2.38, 1.81, …
## $ lq_annual_avg_wkly_wage <dbl> 1.91, 1.96, 1.81, 1.94, 1.25, 1.03, 1…
## $ lq_avg_annual_pay <dbl> 1.91, 1.96, 1.82, 1.94, 1.25, 1.03, 1…
## $ oty_disclosure_code <chr> "", "", "", "", "", "", "", "", "N", …
## $ oty_annual_avg_estabs_chg <int> -1, 1, 0, 0, 0, 1, 0, 1, 0, 0, 2, 1, …
## $ oty_annual_avg_estabs_pct_chg <dbl> -10.0, 14.3, 0.0, 0.0, 0.0, 25.0, 0.0…
## $ oty_annual_avg_emplvl_chg <int> -16, -20, -45, -107, -183, -2, -29, -…
## $ oty_annual_avg_emplvl_pct_chg <dbl> -10.3, -21.3, -17.4, -20.9, -32.8, -6…
## $ oty_total_annual_wages_chg <dbl> -179876, -201005, -650811, -1481365, …
## $ oty_total_annual_wages_pct_chg <dbl> -4.2, -8.2, -7.8, -9.1, -20.6, -8.4, …
## $ oty_taxable_annual_wages_chg <dbl> -129112, -176020, -593231, -1151168, …
## $ oty_taxable_annual_wages_pct_chg <dbl> -6.8, -15.4, -19.7, -18.5, -16.8, -30…
## $ oty_annual_contributions_chg <int> -4632, -16715, -41601, -93444, -51626…
## $ oty_annual_contributions_pct_chg <dbl> -2.6, -16.4, -23.0, -24.3, -17.6, -32…
## $ oty_annual_avg_wkly_wage_chg <int> 36, 85, 71, 91, 77, -3, 113, 96, 0, 6…
## $ oty_annual_avg_wkly_wage_pct_chg <dbl> 6.9, 17.0, 11.5, 14.9, 18.1, -0.7, 24…
## $ oty_avg_annual_pay_chg <int> 1856, 4410, 3729, 4748, 4005, -163, 5…
## $ oty_avg_annual_pay_pct_chg <dbl> 6.8, 16.9, 11.6, 15.0, 18.1, -0.7, 24…#replace county fips codes with county names
RecreationTotal2 <- RecreationTotal %>%
mutate(CountyName = case_when(
endsWith(area_fips, "26009") ~ "Antrim County, MI",
endsWith(area_fips, "26019") ~ "Benzie County, MI",
endsWith(area_fips, "26029") ~ "Charlevoix County, MI",
endsWith(area_fips, "26047") ~ "Emmet County, MI",
endsWith(area_fips, "26055") ~ "Grand Traverse County, MI",
endsWith(area_fips, "26079") ~ "Kalkaska County, MI",
endsWith(area_fips, "26089") ~ "Leelanau County, MI",
endsWith(area_fips, "26101") ~ "Manistee County, MI",
endsWith(area_fips, "26113") ~ "Missaukee County, MI",
endsWith(area_fips, "26165") ~ "Wexford County, MI"
))Now that we have the LQ data prepared, lets plot each county’s LQ (if it is greater than 1.2) over the past 3 years available (2020 to 2022).
#plot LQs higher than national average (>1.2) in recreation for the selected counties and periods
RecreationTotal2 %>%
na.omit() %>%
filter(lq_annual_avg_estabs > 1.2) %>%
group_by(CountyName) %>%
ggplot(aes(x = year, y = lq_annual_avg_estabs, fill = CountyName)) +
scale_x_reverse() +
geom_col() +
geom_smooth() +
facet_wrap(~CountyName) +
coord_flip() +
labs(title = "Recreation Establishments LQs",
subtitle = "for selected periods and counties",
caption = "data from http://data.bls.gov/cew/data/api/YEAR/QTR/industry/INDUSTRY.csv",
x = "Period",
y = "LQ")
As expected, each of the selected counties have very high LQs for recreational establishments. This likely means that a large portion of the labor force in these counties are employed in seasonal jobs that end during the winter months, leaving them unemployed for that period (which means a higher unemployment rate).
ggsave(here("output", "RecreationEstablishmentLQ.png"))