Business Cycle Statistics

Loading packages, defining colors and using data

We first clear the workspace using rm(list = ls())and then include all packages we need. If a package is missing in your R distribution (which is quite likely initially), just use install.packages("package_name") with the respective package name to install it on your system. If you execute the code in the file install_packages.R, then all necessary packages will be installed into your R distribution. If the variable export_graphs is set to TRUE, then the graphs will be exported as pdf-files. In addition, we define a set of colors here to make graphs look more beautiful. Furthermore, we define some axis breaks for graphs to look nicer.

rm(list = ls())
library(zoo)
library(xts)
library(dynlm)
library(reshape2)
library(base)
library(ggplot2)
library(grid)
library(scales)
library(stringr)
library(tidyverse)
library(fredr)
library(formattable)

# should graphs be exported to pdf
export_pdf <- FALSE

# define some colors
mygreen <- "#00BA38"
myblue  <- "#619CFF"
myred   <- "#F8766D"

# define plotting breaks
xbreaks <- c(seq(from = as.Date("1950-01-01"), to = as.Date("2020-01-01"),by = "10 years"))

Quarterly data from the FRED database

Now we want download several data series from the FRED database at once. To this end, we’ve written a function get_fred_data that allows us to do this. This function is stored in a separate file named chap06_func_get_fred_data.R. The function get_fred_data receives as inputs two arrays that need to have the same length as well as three scalar values. The array series contains the ID of each series that should be downloaded. You can find the ID on the respective series website, in the case of GDP per capita this is this one.. The array names contains the respective names under which the series should be stored in our final data structure. Finally, we need to pass the starting date (start_date), the ending date (end_date), and the frequency (freq) to the function. Note that the starting and end date as well as the frequency need to be identical for all series that should be downloaded. Once we call the function get_fred_data with these arguments, it is going to return a data structure that contains all the data we’ve asked for.

# Load function to get FRED data
source("chap06_func_get_fred_data.R")

# US recession dates
selection <- data.frame(series = "USRECQP", names = "RECESSION")

# nominal time series
selection <- rbind(selection, c("GDP", "GDP"))
selection <- rbind(selection, c("PCEC", "CONS"))
selection <- rbind(selection, c("PCDG", "CONS_D"))
selection <- rbind(selection, c("PCND", "CONS_N"))
selection <- rbind(selection, c("PCESV", "CONS_S"))
selection <- rbind(selection, c("GPDI", "INV"))
selection <- rbind(selection, c("PRFI", "INV_R"))
selection <- rbind(selection, c("PNFI", "INV_N"))
selection <- rbind(selection, c("CBI", "INV_I"))
selection <- rbind(selection, c("GCE", "GOV"))
selection <- rbind(selection, c("EXPGS", "EXP"))
selection <- rbind(selection, c("IMPGS", "IMP"))

# price deflators
selection <- rbind(selection, c("GDPDEF", "P_GDP"))
selection <- rbind(selection, c("DPCERD3Q086SBEA", "P_CONS"))
selection <- rbind(selection, c("DDURRD3Q086SBEA", "P_CONS_D"))
selection <- rbind(selection, c("DNDGRD3Q086SBEA", "P_CONS_N"))
selection <- rbind(selection, c("DSERRD3Q086SBEA", "P_CONS_S"))
selection <- rbind(selection, c("A006RD3Q086SBEA", "P_INV"))
selection <- rbind(selection, c("A011RD3Q086SBEA", "P_INV_R"))
selection <- rbind(selection, c("A008RD3Q086SBEA", "P_INV_N"))
selection <- rbind(selection, c("A371RD3Q086SBEA", "P_INV_I"))
selection <- rbind(selection, c("A822RD3Q086SBEA", "P_GOV"))
selection <- rbind(selection, c("A020RD3Q086SBEA", "P_EXP"))
selection <- rbind(selection, c("A021RD3Q086SBEA", "P_IMP"))

# additional series
selection <- rbind(selection, c("CE16OV", "EMPL"))
selection <- rbind(selection, c("UNRATE", "UNEMP"))
selection <- rbind(selection, c("AWHMAN", "HOURS"))
selection <- rbind(selection, c("PRS85006092", "LPROD"))
selection <- rbind(selection, c("PRS85006151", "REALCOMP"))
selection <- rbind(selection, c("TB3MS", "NOMRATE"))

# pull all data series from FRED
fred <- get_fred_data(selection$series, selection$names, "1948-01-01", "2023-12-31", "q")

## [1] "Loading series GDP"
## [1] "Loading series PCEC"
## [1] "Loading series PCDG"
## [1] "Loading series PCND"
## [1] "Loading series PCESV"
## [1] "Loading series GPDI"
## [1] "Loading series PRFI"
## [1] "Loading series PNFI"
## [1] "Loading series CBI"
## [1] "Loading series GCE"
## [1] "Loading series EXPGS"
## [1] "Loading series IMPGS"
## [1] "Loading series GDPDEF"
## [1] "Loading series DPCERD3Q086SBEA"
## [1] "Loading series DDURRD3Q086SBEA"
## [1] "Loading series DNDGRD3Q086SBEA"
## [1] "Loading series DSERRD3Q086SBEA"
## [1] "Loading series A006RD3Q086SBEA"
## [1] "Loading series A011RD3Q086SBEA"
## [1] "Loading series A008RD3Q086SBEA"
## [1] "Loading series A371RD3Q086SBEA"
## [1] "Loading series A822RD3Q086SBEA"
## [1] "Loading series A020RD3Q086SBEA"
## [1] "Loading series A021RD3Q086SBEA"
## [1] "Loading series CE16OV"
## [1] "Loading series UNRATE"
## [1] "Loading series AWHMAN"
## [1] "Loading series PRS85006092"
## [1] "Loading series PRS85006151"
## [1] "Loading series TB3MS"

The series for GDP and its components we’ve downloaded are nominal series. Therefore we have to deflate them into real series before we can use them. Deflators, however, are usually calculated separately for each component of GDP. Therefore, the sum of the deflated GDP components usually doesn’t add up exactly to deflated GDP anymore. We therefore adjust each series proportionally in order to make them consistent with aggregate real GDP. This involves calculating the ratio between real GDP and the sum of its components in the variable factor. We then multiply each individual series with this factor to make our data consistent. The argument that inflated series don’t usually add up well also holds for real total consumption and real total investment. We therefore construct these aggregates by adding up the different real components of consumption and investment, respectively. Finally, we calculate calculate net exports in both nominal and real terms.

# calculate real variables by deflating nominal ones
for(var in c("GDP", "CONS_D", "CONS_N", "CONS_S", "INV_R", "INV_N", "INV_I", "GOV", "EXP", "IMP")) {
  fred[paste("R_", var, sep="")] = fred[var]/fred[paste("P_", var, sep="")]*100
}

# calculate adjusted real GDP (because of different deflators that don't finally add up correctly)
factor <- fred$R_GDP/(fred$R_CONS_D + fred$R_CONS_N + fred$R_CONS_S +fred$R_INV_R + fred$R_INV_N + fred$R_INV_I + fred$R_GOV + fred$R_EXP - fred$R_IMP)

# scale series
for(var in c("CONS_D", "CONS_N", "CONS_S", "INV_R", "INV_N", "INV_I", "GOV", "EXP", "IMP")) {
  fred[paste("R_", var, sep="")] = fred[paste("R_", var, sep="")]*factor
}

# calculate aggregate real consumption and investment
fred["R_CONS"] = fred["R_CONS_D"] + fred["R_CONS_N"] + fred["R_CONS_S"]
fred["R_INV"] = fred["R_INV_R"] + fred["R_INV_N"] + fred["R_INV_I"]

# determine net exports
fred$XM   <- fred$EXP   - fred$IMP
fred$R_XM <- fred$R_EXP - fred$R_IMP

Calculate Statistics Table

The next step is to determine how recessions are different compared to normal times. We start with calculating normal GDP growth in y_normal, which is just the average growth rate of GDP over the entire time horizon of our data. We can do this using the diff function, which just calculates the forward difference for a series, i.e. \(Y_{t+1} - Y_t\). We then just have to take the mean of all the calculated GDP growth rates. Next, we only select data from recessions into the variable y_rec to calculate GDP growth rates again, but only in recession times. After this exercise, y_normal contains the normal GDP growth rate and y_rec the GDP growth rates in recessions only.

# calculate regression beginnings and ends
fred$REC_DATES <- 0
fred$REC_DATES[2:length(fred$REC_DATES)] <- diff(fred$RECESSION)

# calculate normal GDP growth (in quarterly terms)
y_normal <- mean(diff(fred$R_GDP)/fred$R_GDP[1:length(fred$R_GDP)-1])

# GDP growth in recessions relative to normal growth
y_rec <- vector();
for(t in 1:length(fred$R_GDP)) {
  if(fred$RECESSION[t] == 1 & fred$REC_DATES[t] != -1) {
    y_rec <- append(y_rec, (fred$R_GDP[t+1] - fred$R_GDP[t])/fred$R_GDP[t])
  }
}

# calculate mean across recessions
y_rec <- mean(y_rec)

Next, we calculate the share of each component in GDP as well as its contribution to GDP growth in recessions relative to normal times. The share of each component in total GDP is stored in dat[i, 1]. To calculate the contribution of each variable to GDP growth in recessions relative to normal times is a bit more difficult. We first calculate the first difference of each component of GDP (for aggregate consumption this for example is \(C_{t+1} - C_t\)) and relate it to the respective value of GDP \(Y_t\). We do this both for the entire time series of data (v_normal) as well as for the recession dates only (v_rec). We finally relate the difference in v_rec and v_normal to the difference in y_rec and y_normal to get to each component’s contribution in GDP growth.

# here statistics are stored
dat <- array(dim=c(10,2))

# iterate over relevant variables
i <- 1
for(var in c("CONS", "CONS_D", "CONS_N", "CONS_S", "INV", "INV_R", "INV_N", "INV_I", "GOV", "XM")) {
  
  # get mean share in GDP
  dat[i, 1] <- mean(unlist(fred[var]/fred["GDP"]))
  
  # get real variable name
  rvar <- paste("R_", var, sep="")
  
  # calculate normal share in GDP growth of variable
  v_normal <- mean(diff(unlist(fred[rvar]))/unlist(fred[1:length(fred$R_GDP)-1,]["R_GDP"]))
  
  # variable share of GDP growth in recessions
  v_rec <- vector();
  for(t in 1:length(fred$R_GDP)) {
    if(fred$RECESSION[t] == 1 & fred$REC_DATES[t] != -1) {
      v_rec <- append(v_rec, unlist(fred[t+1,][rvar] - fred[t,][rvar])/fred$R_GDP[t])
    }
  }
  
  # calculate mean across recessions
  v_rec <- mean(v_rec)
  
  # calculate share in GDP fall relative to normal growth
  dat[i, 2] <- (v_rec - v_normal)/(y_rec - y_normal)
  
  # increment i
  i = i + 1
}

Generate Table

From the data we calculated above, we then create a simple table using the formattable statement.

table_data <- data.frame("Component" = 
                          c("Consumption", " - Durables", " - Nondurables", " - Services",
                            "Investment", " - Residential", " - Fixed Nonresidential", " - Inventories",
                            "Government Purchases", "Net Exports"))
table_data <- cbind(table_data, "Average share in GDP" = percent(dat[, 1]))
table_data <- cbind(table_data, "Average share in fall in GDP in recessions relative to normal growth" = percent(dat[, 2]))

# output table
formattable(table_data,
            align= c("l", "c", "c"))

Component	Average share in GDP	Average share in fall in GDP in recessions relative to normal growth
Consumption	63.57%	45.54%
Durables	8.39%	6.31%
Nondurables	18.89%	10.98%
Services	36.30%	28.25%
Investment	17.39%	52.59%
Residential	4.61%	10.76%
Fixed Nonresidential	12.25%	19.73%
Inventories	0.53%	22.09%
Government Purchases	20.46%	7.87%
Net Exports	-1.42%	-6.00%

Calculate Recession Start and End Dates

Next, we determine the recession start and end dates. This is simply done to be able to shade recessions in the following graphs.

# date at beginning and end of recessions
rstart <- fred$date[fred$REC_DATES == 1]
rend   <- fred$date[fred$REC_DATES == -1]
rec_df <- data.frame(rstart, rend)

Employment

Now we want to plot a couple of important series and study their properties during recessions. With start out with percentage changes in aggregate employment. The growth rates of employment can again be calculated using the diff function.

# calculate employment changes
fred$D_EMPL <- NA
fred$D_EMPL[1:(length(fred$EMPL)-1)] <- diff(fred$EMPL)/fred$EMPL[1:(length(fred$EMPL)-1)]*100

# set minimum and maximum
ymin <- -14
ymax <- 8
dy   <- (abs(ymin)+abs(ymax))/2

# create plot
myplot <- ggplot(data = fred) + 
  geom_rect(data=rec_df, aes(xmin=rstart, xmax=rend, ymin=ymin-0.05*dy, ymax=ymax+0.05*dy), fill='gray', alpha=0.5) + 
  geom_hline(yintercept=0, color="gray", linewidth=0.5) +
  geom_line(aes(x=date, y=D_EMPL), color="darkblue", linewidth=0.5) +
  scale_x_date(breaks = xbreaks, date_labels = "%Y", expand=c(0, 0)) +
  scale_y_continuous(breaks=seq(ymin, ymax, 2), expand=c(0, 0)) +  
  labs(x = "Year t",
       y = "Employment Change (in %)") +
  theme_bw()

# print the plot
print(myplot)

Unemployment

Next, we take a look at unemployment rates.

# set minimum and maximum
ymin <- 0
ymax <- 15
dy   <- (abs(ymin)+abs(ymax))/2

# create plot
myplot <- ggplot(data = fred) + 
  geom_rect(data=rec_df, aes(xmin=rstart, xmax=rend, ymin=ymin-0.05*dy, ymax=ymax+0.05*dy), fill='gray', alpha=0.5) + 
  geom_line(aes(x=date, y=UNEMP), color="darkblue", linewidth=0.5) +
  scale_x_date(breaks = xbreaks, date_labels = "%Y", expand=c(0, 0)) +
  scale_y_continuous(breaks=seq(ymin, ymax, 2.5), expand=c(0, 0)) +  
  labs(x = "Year t",
       y = "Unemployment Rate (quarterly, in %)") +
  theme_bw()

# print the plot
print(myplot)

Average weekly hours of production workers in manufacturing

This series plots average weekly hours of production workers in manufacturing.

# set minimum and maximum
ymin <- 36
ymax <- 44
dy   <- (abs(ymin)+abs(ymax))/2

# create plot
myplot <- ggplot(data = fred) + 
  geom_rect(data=rec_df, aes(xmin=rstart, xmax=rend, ymin=ymin-0.05*dy, ymax=ymax+0.05*dy), fill='gray', alpha=0.5) + 
  geom_line(aes(x=date, y=HOURS), color="darkblue", linewidth=0.5) +
  scale_x_date(breaks = xbreaks, date_labels = "%Y", expand=c(0, 0)) +
  scale_y_continuous(breaks=seq(ymin, ymax, 2), expand=c(0, 0)) +  
  labs(x = "Year t",
       y = "Average Weekly Hours, Prod. Workers") +
  theme_bw()

# print the plot
print(myplot)

Change in output per hour, non-farm business

This series contains changes in output per hour in non-farm businesses.

# set minimum and maximum
ymin <- -8
ymax <- 18
dy   <- (abs(ymin)+abs(ymax))/2

# create plot
myplot <- ggplot(data = fred) + 
  geom_rect(data=rec_df, aes(xmin=rstart, xmax=rend, ymin=ymin-0.05*dy, ymax=ymax+0.05*dy), fill='gray', alpha=0.5) + 
  geom_hline(yintercept=0, color="gray", linewidth=0.5) +
  geom_line(aes(x=date, y=LPROD), color="darkblue", linewidth=0.5) +
  scale_x_date(breaks = xbreaks, date_labels = "%Y", expand=c(0, 0)) +
  scale_y_continuous(breaks=seq(ymin, ymax, 4), expand=c(0, 0)) +  
  labs(x = "Year t",
       y = "Change in Output per Hour (in %)") +
  theme_bw()

# print the plot
print(myplot)

Change in real compensation per hour, non-farm business

Again for non-farm businesses, but this time changes in real compensation per hour.

# set minimum and maximum
ymin <- -4
ymax <- 10
dy   <- (abs(ymin)+abs(ymax))/2

# create plot
myplot <- ggplot(data = fred) + 
  geom_rect(data=rec_df, aes(xmin=rstart, xmax=rend, ymin=ymin-0.05*dy, ymax=ymax+0.05*dy), fill='gray', alpha=0.5) + 
  geom_hline(yintercept=0, color="gray", linewidth=0.5) +
  geom_line(aes(x=date, y=REALCOMP), color="darkblue", linewidth=0.5) +
  scale_x_date(breaks = xbreaks, date_labels = "%Y", expand=c(0, 0)) +
  scale_y_continuous(breaks=seq(ymin, ymax, 4), expand=c(0, 0)) +  
  labs(x = "Year t",
       y = "Change in Real Compensation per Hour (in %)") +
  theme_bw()

# print the plot
print(myplot)

Inflation

Finally, we take a look at prices and interest rates. This is the inflation rate.

# calculate annualized inflation from GDP deflator
fred$INFL <- NA
fred$INFL[5:(length(fred$P_GDP))] <- diff(fred$P_GDP, lag = 4)/fred$P_GDP[1:(length(fred$P_GDP)-4)]*100

# set minimum and maximum
ymin <- -2
ymax <- 12
dy   <- (abs(ymin)+abs(ymax))/2

# create plot
myplot <- ggplot(data = fred) + 
  geom_rect(data=rec_df, aes(xmin=rstart, xmax=rend, ymin=ymin-0.05*dy, ymax=ymax+0.05*dy), fill='gray', alpha=0.5) + 
  geom_hline(yintercept=0, color="gray", linewidth=0.5) +
  geom_line(aes(x=date, y=INFL), color="darkblue", linewidth=0.5) +
  scale_x_date(breaks = xbreaks, date_labels = "%Y", expand=c(0, 0)) +
  scale_y_continuous(breaks=seq(ymin, ymax, 2), expand=c(0, 0)) +  
  labs(x = "Year t",
       y = "Inflation Relative to Previous Year (in %)") +
  theme_bw()

# print the plot
print(myplot)

Nominal interest rate on 3-month treasuries

An this is the nominal interest rate on 3-month treasuries.

# set minimum and maximum
ymin <- 0
ymax <- 16
dy   <- (abs(ymin)+abs(ymax))/2

# create plot
myplot <- ggplot(data = fred) + 
  geom_rect(data=rec_df, aes(xmin=rstart, xmax=rend, ymin=ymin-0.05*dy, ymax=ymax+0.05*dy), fill='gray', alpha=0.5) + 
  geom_line(aes(x=date, y=NOMRATE), color="darkblue", linewidth=0.5) +
  scale_x_date(breaks = xbreaks, date_labels = "%Y", expand=c(0, 0)) +
  scale_y_continuous(breaks=seq(ymin, ymax, 2), expand=c(0, 0)) +  
  labs(x = "Year t",
       y = "Nominal Interest on 3-month Treasuries (in %)") +
  theme_bw()

# print the plot
print(myplot)

Ex post real interest rate on 3-month treasuries

The difference between the nominal rate and inflation leads us to the ex-post real rate.

# set minimum and maximum
ymin <- -8
ymax <- 8
dy   <- (abs(ymin)+abs(ymax))/2

# create plot
myplot <- ggplot(data = fred) + 
  geom_rect(data=rec_df, aes(xmin=rstart, xmax=rend, ymin=ymin-0.05*dy, ymax=ymax+0.05*dy), fill='gray', alpha=0.5) + 
  geom_hline(yintercept=0, color="gray", linewidth=0.5) +
  geom_line(aes(x=date, y=NOMRATE-INFL), color="darkblue", linewidth=0.5) +
  scale_x_date(breaks = xbreaks, date_labels = "%Y", expand=c(0, 0)) +
  scale_y_continuous(breaks=seq(ymin, ymax, 2), expand=c(0, 0)) +  
  labs(x = "Year t",
       y = "Real Interest on 3-month Treasuries (in %)") +
  theme_bw()

# print the plot
print(myplot)