· 6 years ago · Mar 06, 2020, 03:40 PM
1####################################################################
2### Graph US Unemployment data from FRED
3### - <Mathew.Binkley@Vanderbilt.edu>
4####################################################################
5
6### We need the following packages for this example.
7packages <- c("fredr", "lubridate", "fredr", "ggplot2", "ggfortify", "seasonal",
8 "tidyverse", "ggthemes", "tsibble", "dplyr", "magrittr", "broom")
9
10### Install packages if needed, then load them quietly
11new_packages <- packages[!(packages %in% installed.packages()[, "Package"])]
12if (length(new_packages)) install.packages(new_packages, quiet = TRUE)
13invisible(lapply(packages, "library", quietly = TRUE,
14 character.only = TRUE, warn.conflicts = FALSE))
15
16### Set my FRED API key to access the FRED database.
17### You may request an API key at:
18### https://research.stlouisfed.org/useraccount/apikeys
19api_key_fred <- "GO_GET_YOUR_OWN"
20fredr_set_key(api_key_fred)
21
22### We only want to plot the last ~3 years of data
23date_start <- as.Date("2017-01-01")
24date_end <- as.Date(now())
25
26
27####################################################################
28### Fetch unemployment data from FRED
29####################################################################
30data <- fredr(series_id = "UNRATENSA", frequency = "m") %>%
31 as_tsibble(index = "date")
32
33date <- data %>% pull("date")
34values <- data %>% pull("value")
35
36
37####################################################################
38### Deseasonalize using X13-ARIMA-SEATS
39####################################################################
40values_ts <- values %>% ts(frequency = 12,
41 start = c(year(min(date)), month(min(date))))
42
43x13arima <- seas(values_ts, estimate.maxiter = 10000, x11 = "")
44x13arima_fit <- x13arima$data %>% as_tibble()
45x13arima_trend <- x13arima_fit$final
46x13arima_seasonal <- x13arima_fit$seasonal
47
48data_df <- data.frame(date = date,
49 values = values,
50 trend = x13arima_trend,
51 seasonal = x13arima_seasonal)
52
53data_subset <- data_df %>% filter(date >= date_start & date <= date_end)
54
55####################################################################
56### Graph: Deseasonalized unemployment compared to raw data
57####################################################################
58c1 <- "U.S. Bureau of Labor Statistics, Unemployment Rate in US [UNRATENSA]\n"
59c2 <- "Retrieved from FRED, Federal Reserve Bank of St. Louis;\n"
60c3 <- "https://fred.stlouisfed.org/series/UNRATENSA\n"
61c4 <- paste("Data retrieved on", format(Sys.time(), "%B %d, %Y at %I:%M %p %Z"))
62caption <- paste(c1, c2, c3, c4)
63
64title <- paste("US Unemployment Rate\nbased on unemployment data from",
65 month.name[month(last(date))], year(last(date)))
66
67s1 <- "Raw unemployment rate in grey, "
68s2 <- "deseasonalized unemployment rate in gold\n"
69s3 <- "Seasonal signal in dotted blue\n"
70
71subtitle <- paste(s1, s2, s3, sep = "")
72xlab <- "Year"
73ylab <- "Unemployment Rate"
74caption <- paste(c1, c2, c3, c4)
75
76p <- ggplot(data = data_subset, aes(x = date, y = values)) +
77
78 theme_economist() +
79
80 geom_point(color = "grey60", shape = 1) +
81
82 geom_line(data = data_subset, aes(x = date, y = values),
83 size = 1.3, color = "black", alpha = 0.1) +
84
85 geom_line(data = data_subset, aes(x = date, y = trend),
86 size = 1.3, color = "goldenrod2", alpha = 0.8) +
87
88 geom_line(data = data_subset,
89 aes(x = date,
90 y = (max(trend) + min(trend)) / 1.5 + seasonal),
91 size = 1.3, color = "steelblue2",
92 alpha = 0.5, linetype = "dotted") +
93
94 labs(title = title, subtitle = subtitle, caption = caption,
95 x = xlab, y = ylab)
96print(p)