```{r setup, include=FALSE} knitr::opts_chunk$set(echo = TRUE) ``` ## Working in R These code downloads have been constructed as supplements to the full [Doing Economics projects](https://core-econ.org/doing-economics/). You’ll need to download the data before running the code that follows. ## Getting started in R For this project you will need the following packages: - `tidyverse`, to help with data manipulation - `mFilter`, to apply the Hodrick-Prescott (HP) filter You will also use the `ggplot2` package to produce accurate charts, but that comes as part of the `tidyverse` package. If you need to install any of these packages, run the following code: ``` r install.packages(c("tidyverse","mFilter"), repos="https://www.stats.bris.ac.uk/R/") ``` You can import these libraries now, or when they are used in the R walk-throughs below. ``` r library(tidyverse) ``` ### **R walk-through 1** Importing a `.csv` datafile into R and doing basic data cleaning We want to import the datafiles we have just downloaded into R, starting with the NBER US recession dates ('USRECQ.csv'). We start by setting our working directory using the `setwd` command. This command tells R where your codes and datafiles are stored. In the code below, replace 'YOURFILEPATH' with the full filepath that indicates the folder in which you have saved the code chunks file (not the 'data' subfolder, but its parent folder). Note that you have to use forward slashes ('/') rather than backslashes ('\\'). If you don’t know how to find the path to your working folder, see the [‘Technical Reference’ section](50-00-technical-reference.html#finding-your-path-on-mac). ``` r setwd('YOURFILEPATH') ``` You should now see the contents of your working directory in the bottom-right panel of RStudio (click the ‘Files’ tab). It should contain the 'data' subfolder and this '.rmd' file. Since our data is in csv format, we use the `read.csv` function to import the data into R, specifying that we want R to look for the file in the 'data' subfolder. We will call our file ‘US_recessions’. ``` r US_recessions <- read.csv("data/USRECQ.csv") ``` To check that the data has been imported correctly, you can use the `head` function to view the first six rows of the dataset, and confirm that they correspond to the columns in the csv file. ``` r head(US_recessions) ``` Before working with the data, we use the `str` function to check that the data is formatted correctly. ``` r str(US_recessions) ``` You can see that the variable `USRECQ` containing the quarterly recession dummy is an integer (which is correct as it is a binary dummy variable), while the `DATE` variable containing the quarterly dates is formatted as character (`chr`). We want to create an annual recession dummy because the other data is at the annual level (and we want to have consistent date variables when combining the various datasets). The first step is to create a numeric `Year` variable based on `DATE`, which can be done via the `substr` function (extracting the first four characters) inside `as.numeric` to convert the extracted character string to numeric (`num`) type. We then create an annual dummy variable called `Recession` that takes the value 1 for all observations within a year if one or more quarters of that year featured a recession (`USRECQ` = 1). ``` r US_recessions <- US_recessions %>% mutate(Year = as.numeric(substr(DATE,1,4))) %>% group_by(Year) %>% mutate(Recession = case_when(any(USRECQ == 1) ~ 1, TRUE ~ 0)) ``` Piping is a very useful technique in R for data analysis. If we have to perform a sequence of commands on the same object, then we can simplify the code by writing the sequence as a single command. We use the punctuation `%>%` to link the commands together. In this code above, we 1) created the numeric `Year` variable based on `DATE`, 2) grouped the data by `Year`, 3) created the annual `Recession` dummy using the `any` function within `case_when`, which sets the dummy to 1 if `USRECQ = 1` for any of the quarters within a year. Otherwise (indicated by `TRUE`) we set the dummy to 0. (For a more detailed introduction to piping, visit the University of Manchester’s [Econometric Computing Learning Resource](https://tinyco.re/5531433).) Next, we want to remove all variables from `US_recessions` except `Year` and `Recession`. Using the `distinct` function, this can be done in one step together with removing the duplicate observations (as we still have every year four times due to the initial data being at the quarterly level). We also want to create a new variable called `Country` indicating that this data is for the US (which will be helpful when merging it with the other datasets later). ``` r US_recessions <- US_recessions %>% distinct(Year, Recession) %>% mutate(Country ="United States") ``` ### **R walk-through 2** Recession dating based on quarterly GDP data We start by importing the quarterly GDP dataset and calling our file ‘GDP_quarterly’. ``` r GDP_quarterly <- read.csv("data/GDP_quarterly.csv") ``` The quarterly GDP values are stored in the `Value` column. We want to obtain a recession dummy for the non-US country, similar to the NBER one we downloaded for the US. In this example, we have chosen Germany in addition to the US. To compute the growth rate of GDP, we compute the percentage change from one quarter to the next. We only want to do this for the non-US country (Germany), so we `filter` for the German values only. We can use the `lag` function to take the previous row’s (i.e. quarter’s) value. ``` r GDP_quarterly <- GDP_quarterly %>% filter(Country == "Germany") %>% mutate(pct_change = (Value/lag(Value) - 1) * 100) ``` The NBER defines a recession as at least two consecutive quarters of negative GDP growth. We want to define a recession dummy variable that takes the value of 1 in all periods of all such instances, including the first quarter of negative growth. In all other periods (positive GDP growth or only one quarter of negative GDP growth), the dummy should be set to 0. ``` r GDP_quarterly <- GDP_quarterly %>% mutate(Recession_Germany_quarterly = case_when(pct_change < 0 & (lag(pct_change) < 0 | lead(pct_change) < 0) ~ 1, TRUE ~ 0)) ``` We used the `lag` and `lead` functions to check whether either the preceding or the following value (quarter) features a negative GDP growth. If a period has negative GDP growth (`pct_change < 0`) and either the lagging or leading period does as well, we set the recession dummy to 1 using the `case_when` function; otherwise (indicated by `TRUE`) we set the dummy to 0. The `TIME` variable contains dates as character strings (`chr`) in the format YYYY-QQ (year-quarter). As for the US above, we define a numeric `Year` variable, as we need the German recession period data on the annual level. We then convert the quarterly recession dummy `Recession_Germany_quarterly` to annual frequency in the same way as for the US (an annual dummy that takes the value 1 for all observations within a year if one or more quarters of that year featured a recession). The annual recession dummy is a new variable called `Recession_Germany`. We then keep only the `Country`, `Year`, and `Recession_Germany` variables and drop duplicate values (as we have four observations per year) in one step by using the `distinct` function. We can also rename the dataset `GDP_quarterly` to `Germany_recessions` in this step by assigning the dataset in the first line (using the `<-` operator). ``` r Germany_recessions <- GDP_quarterly %>% mutate(Year = as.numeric(substr(TIME,1,4))) %>% group_by(Year) %>% mutate(Recession_Germany = case_when(any(Recession_Germany_quarterly == 1) ~ 1, TRUE ~ 0)) %>% distinct(Country,Year,Recession_Germany) ``` To keep the workspace clean, we remove the `GDP_quarterly` dataset which is no longer needed. ``` r rm(GDP_quarterly) ``` ### **R walk-through 3** Reshaping data and further data cleaning Two of our datasets (`hours_worked` and `LFP_rates`) are currently in long format (the single variable `Value` contains values for multiple economic variables) and need to be reshaped (reformatted) to wide format (where each economic variable has its own column). For example in `hours_worked`, hours worked by all persons, females, and males should be their own variables (columns), rather than sharing the `SEX`, `Sex`, and `Value` columns. We first clean the `GDP_annual` and `labour_productivities` datasets, as these are already in wide format and do not need to be reshaped. We import them into R. ``` r GDP_annual <- read.csv("data/GDP_annual.csv") labour_productivities <- read.csv("data/labour_productivities.csv") ``` In `GDP_annual`, convert the `Value` variable to numeric. Remove all observations with missing values for `Value`. Rename the `Value` variable to `GDP`, the `Country.Name` variable to `Country`, and the `Time` variable to `Year`. Remove all columns except `Country`, `Year`, and `GDP`. ``` r GDP_annual <- GDP_annual %>% mutate(Value = as.numeric(Value)) %>% filter(!is.na(Value)) %>% rename(GDP = Value, Country = Country.Name, Year = Time) %>% select(Country,Year,GDP) ``` In `labour_productivities`, rename the `Time` variable to `Year` and the `Value` variable to `labour_productivity` and remove all columns except `Country`, `Year`, and `labour_productivity`. ``` r labour_productivities <- labour_productivities %>% select(Country,Time,Value) %>% rename(Year = Time, labour_productivity = Value) ``` We now import the two datasets (`hours_worked` and `LFP_rates`) in long format. These need to be shaped into wide format. ``` r hours_worked <- read.csv("data/hours_worked.csv") LFP_rates <- read.csv("data/LFP_rates.csv") ``` We can reshape the datasets with the `spread` function. The `key` input takes the variable that will be used to name the newly created variables (in our case the `Sex` column) while the `value` input contains the values that will be reshaped (`Value` in our data). Before we reshape, we need to remove the `SEX` column which also differs between `All persons`, `Women`, and `Men` (as `spread` can only take one `key` variable, which is `Subject`). This can be done using a minus sign inside the `select` function. ``` r hours_worked <- hours_worked %>% select(-SEX) %>% spread(key = Sex, value = Value) LFP_rates <- LFP_rates %>% select(-SEX) %>% spread(key = Sex, value = Value) ``` The datasets are now in wide format. However, as `spread` takes the (character string) values in the `key` variable as variable names for the newly generated variables, the new `All persons` variables have spaces in their variable name. While R allows variable names to have spaces, it is best practice to avoid this. Moreover, the new `women` and `men` variables are also not very descriptive in terms of what economic series those variables correspond to. Therefore, we rename the reshaped variables to the name of the data series (e.g. to `hours_worked_female`). Note that we have to put the ‘old’ variable names containing spaces (`All persons`) inside quotation marks when we rename them. We also rename the `Time` variable (containing the years as numeric values) to `Year`, consistent with the other datasets. Finally, we remove all columns except `Country`, `Year`, and the newly created (and renamed) variables. ``` r hours_worked <- hours_worked %>% rename(hours_worked_all = `All persons`, hours_worked_female = Women, hours_worked_male = Men, Year = Time) %>% select(Country,Year,hours_worked_all,hours_worked_female,hours_worked_male) LFP_rates <- LFP_rates %>% rename(LFP_rate_all = `All persons`, LFP_rate_female = Women, LFP_rate_male = Men, Year = Time) %>% select(Country,Year,LFP_rate_all,LFP_rate_female,LFP_rate_male) ``` ### **R walk-through 4** Merging multiple datasets into a single dataset To merge (combine) all six datasets into one, we first combine them in a list and then use the `reduce(full_join)` function, which automatically merges the datasets by the `Country` and `Year` variables. We save the new (merged) dataset as `dataset_main`. ``` r dataset_main <- list(GDP_annual,labour_productivities,hours_worked,LFP_rates,US_recessions,Germany_recessions) %>% reduce(full_join) ``` We replace the (missing) values of the `Recession` variable (containing the US recession dummy) with the values of the `Recession_Germany` variable, using the `coalesce` function which takes the first non-missing element of its two inputs. Then we delete the `Recession_Germany` variable. We also sort the data by `Country` and `Time` (in that order). We delete all pre-1960 data as the OECD data series are not available before then. ``` r dataset_main <- dataset_main %>% mutate(Recession = coalesce(Recession,Recession_Germany)) %>% select(-Recession_Germany) %>% arrange(Country,Year) %>% filter(Year >= 1960) ``` Finally, we delete all datasets except `dataset_main`. ``` r rm(GDP_annual,labour_productivities,hours_worked,LFP_rates,US_recessions,Germany_recessions) ``` ### **R walk-through 5** Creating line charts with `ggplot` and adding shaded vertical bars We now want to plot the female and male labour force participation (LFP) rates over time together in a line chart (one chart per country), including a legend. We use the `ggplot` function for this, which is part of the `tidyverse` package we loaded earlier. `ggplot` is one the most popular ways to plot data in R and comes with a lot of different chart types and customisation options. We start with the US: ``` r dataset_main %>% filter(Country == "United States") %>% ggplot(aes(Year)) + geom_line(aes(y = LFP_rate_female, colour = "LFP_rate_female")) + geom_line(aes(y = LFP_rate_male, colour = "LFP_rate_male")) + theme_bw() + # Add the vertical axis label labs(y = "Percentage") + # Scale horizontal axis tick marks in 10-year intervals scale_x_continuous(breaks = seq(1960, 2020, by = 10)) + # Add the title and legend labels ggtitle("Labour force participation rate by sex (United States)") + scale_colour_discrete(name = "", labels = c("Female", "Male")) + # Move the legend below the chart theme(legend.position="bottom") + # Remove the legend title theme(legend.title=element_blank()) ``` Next, we plot the same chart for Germany: ``` r dataset_main %>% filter(Country == "Germany") %>% ggplot(aes(Year)) + geom_line(aes(y = LFP_rate_female, colour = "LFP_rate_female")) + geom_line(aes(y = LFP_rate_male, colour = "LFP_rate_male")) + theme_bw() + # Add the vertical axis label labs(y = "Percentage") + # Scale horizontal axis tick marks in 10-year intervals scale_x_continuous(breaks = seq(1960, 2020, by = 10)) + # Add the title and legend labels ggtitle("Labour force participation rate by sex (Germany)") + scale_colour_discrete(name = "", labels = c("Female", "Male")) + # Move the legend below the chart theme(legend.position="bottom") + # Remove the legend title theme(legend.title=element_blank()) ``` We now add the shaded recession areas, using the `Recession` dummy, to produce the chart shown in the instructions. We use the `geom_rect` function to add the shaded areas for the recession periods, and the `subset` function to only plot the recession areas (where the recession dummy `Recession ==`1) for the US (`Country == "United States"`). ``` r dataset_main %>% filter(Country == "United States") %>% mutate(ymax = 90, ymin = 0) %>% ggplot(aes(Year)) + geom_line(aes(y = LFP_rate_female, colour = "LFP_rate_female")) + geom_line(aes(y = LFP_rate_male, colour = "LFP_rate_male")) + theme_bw() + # Add the vertical axis label labs(y = "Percentage") + # Scale horizontal axis tick marks in 10-year intervals scale_x_continuous(breaks = seq(1960, 2020, by = 10)) + # Add the title and legend labels ggtitle("Labour force participation rate by sex (United States)") + scale_colour_discrete(name = "", labels = c("Female", "Male")) + # Move the legend below the chart theme(legend.position="bottom") + # Remove the legend title theme(legend.title=element_blank()) + # Add recession periods as shaded areas geom_rect(data = subset(dataset_main, Recession == 1 & Country == "United States"), aes(ymin = -Inf, ymax = Inf, xmin = Year-0.5, xmax = Year+0.5), alpha = 0.2) ``` We also plot the chart with recession areas for Germany. ``` r dataset_main %>% filter(Country == "Germany") %>% ggplot(aes(Year)) + geom_rect(data = subset(dataset_main, Recession == 1 & Country == "Germany"), aes(ymin = -Inf, ymax = Inf, xmin = Year-0.5, xmax = Year+0.5), alpha = 0.2) + geom_line(aes(y = LFP_rate_female, colour = "LFP_rate_female")) + geom_line(aes(y = LFP_rate_male, colour = "LFP_rate_male")) + theme_bw() + # Add the vertical axis label labs(y = "Percentage") + # Scale horizontal axis tick marks in 10-year intervals scale_x_continuous(breaks = seq(1960, 2020, by = 10)) + # Add the title and legend labels ggtitle("Labour force participation rate by sex (Germany)") + scale_colour_discrete(name = "", labels = c("Female", "Male")) + # Move the legend below the chart theme(legend.position="bottom") + # Remove the legend title theme(legend.title=element_blank()) # Add recession periods as shaded areas geom_rect(data = subset(dataset_main, Recession == 1 & Country == "Germany"), aes(ymin = -Inf, ymax = Inf, xmin = Year-0.5, xmax = Year+0.5), alpha = 0.2) ``` ### **R walk-through 6** Using the Hodrick-Prescott (HP) filter and plotting cyclical components We use the `log` function to create variables containing the natural log of annual GDP and of hours worked (female, male, and all persons). ``` r dataset_main <- dataset_main %>% mutate(GDP_log = log(GDP), hours_worked_all_log = log(hours_worked_all), hours_worked_female_log = log(hours_worked_female), hours_worked_male_log = log(hours_worked_male)) ``` To apply the HP filter to the four log variables, we use the `hpfilter` function which is part of the `mFilter` package. We want to set the smoothing parameter lambda to 6.25 (a standard value for annual data like ours) and save both the trend and cyclical components. The `hpfilter` function outputs two vectors: `trend`, containing the trend component, and `cycle`, containing the cyclical component. To save the trend component as a new variable, we add `$trend` after calling the `hpfilter` function. Before applying the HP filter, we need to group by country (otherwise R would treat both countries’ log GDP as one data series) and filter out the missing values via `filter(!is.na(GDP_log))`. Note that this removes all observations which have missing values in `GDP_log`. We therefore merge the dataset with the HP-filtered values back to the main dataset. Finally, we remove all intermediate datasets (with their filenames all starting with `dataset_main_filter`) except for `dataset_main`, using `rm(list=ls(pattern="dataset_main_filter"))`. ``` r dataset_main_filter1 <- dataset_main %>% group_by(Country) %>% filter(!is.na(GDP_log)) %>% mutate(GDP_log_trend = hpfilter(GDP_log, freq=6.25)$trend, GDP_log_cycle = hpfilter(GDP_log, freq=6.25)$cycle) dataset_main_filter2 <- dataset_main %>% group_by(Country) %>% filter(!is.na(hours_worked_all_log)) %>% mutate(hours_worked_all_log_trend = hpfilter(hours_worked_all_log, freq=6.25)$trend, hours_worked_all_log_cycle = hpfilter(hours_worked_all_log, freq=6.25)$cycle) dataset_main_filter3 <- dataset_main %>% group_by(Country) %>% filter(!is.na(hours_worked_female_log)) %>% mutate(hours_worked_female_log_trend = hpfilter(hours_worked_female_log, freq=6.25)$trend, hours_worked_female_log_cycle = hpfilter(hours_worked_female_log, freq=6.25)$cycle) dataset_main_filter4 <- dataset_main %>% group_by(Country) %>% filter(!is.na(hours_worked_male_log)) %>% mutate(hours_worked_male_log_trend = hpfilter(hours_worked_male_log, freq=6.25)$trend, hours_worked_male_log_cycle = hpfilter(hours_worked_male_log, freq=6.25)$cycle) dataset_main <- list(dataset_main,dataset_main_filter1,dataset_main_filter2,dataset_main_filter3,dataset_main_filter4) %>% reduce(full_join) ``` ``` r rm(list=ls(pattern="dataset_main_filter")) ``` We plot the cyclical component of log GDP and of hours worked (by all persons) for each country, again adding the recession areas. ``` r dataset_main %>% filter(Country == "United States") %>% ggplot(aes(Year)) + geom_line(aes(y = GDP_log_cycle)) + theme_bw() + # Add the vertical axis label labs(y = "") + # Scale horizontal axis tick marks in 10-year intervals scale_x_continuous(breaks = seq(1960, 2020, by = 10)) + # Add the title and legend labels ggtitle("Cylical component of log GDP (United States)") + # Move the legend below the chart theme(legend.position="bottom") + # Remove the legend title theme(legend.title=element_blank()) + # Add recession periods as shaded areas geom_rect(data = subset(dataset_main, Recession == 1 & Country == "United States"), aes(ymin = -Inf, ymax = Inf, xmin = Year-0.5, xmax = Year+0.5), alpha = 0.2) ``` ``` r dataset_main %>% filter(Country == "United States") %>% ggplot(aes(Year)) + geom_line(aes(y = hours_worked_all_log_cycle)) + theme_bw() + # Add the vertical axis label labs(y = "") + # Scale horizontal axis tick marks in 10-year intervals scale_x_continuous(breaks = seq(1960, 2020, by = 10)) + # Add the title and legend labels ggtitle("Cylical component of average hours worked (United States)") + # Move the legend below the chart theme(legend.position="bottom") + # Remove the legend title theme(legend.title=element_blank()) + # Add recession periods as shaded areas geom_rect(data = subset(dataset_main, Recession == 1 & Country == "United States"), aes(ymin = -Inf, ymax = Inf, xmin = Year-0.5, xmax = Year+0.5), alpha = 0.2) ``` ``` r dataset_main %>% filter(Country == "Germany") %>% ggplot(aes(Year)) + geom_line(aes(y = GDP_log_cycle)) + theme_bw() + # Add the vertical axis label labs(y = "") + # Scale horizontal axis tick marks in 10-year intervals scale_x_continuous(breaks = seq(1960, 2020, by = 10)) + # Add the title and legend labels ggtitle("Cylical component of log GDP (Germany)") + # Move the legend below the chart theme(legend.position="bottom") + # Remove the legend title theme(legend.title=element_blank()) + # Add recession periods as shaded areas geom_rect(data = subset(dataset_main, Recession == 1 & Country == "Germany"), aes(ymin = -Inf, ymax = Inf, xmin = Year-0.5, xmax = Year+0.5), alpha = 0.2) ``` ``` r dataset_main %>% filter(Country == "Germany") %>% ggplot(aes(Year)) + geom_line(aes(y = hours_worked_all_log_cycle)) + theme_bw() + # Add the vertical axis label labs(y = "") + # Scale horizontal axis tick marks in 10-year intervals scale_x_continuous(breaks = seq(1960, 2020, by = 10)) + # Add the title and legend labels ggtitle("Cylical component of average hours worked (Germany)") + # Move the legend below the chart theme(legend.position="bottom") + # Remove the legend title theme(legend.title=element_blank()) + # Add recession periods as shaded areas geom_rect(data = subset(dataset_main, Recession == 1 & Country == "Germany"), aes(ymin = -Inf, ymax = Inf, xmin = Year-0.5, xmax = Year+0.5), alpha = 0.2) ``` ### **R walk-through 7** Computing cyclical properties of economic variables and plotting a grouped bar chart First, we study the cyclicality of average hours worked (of all persons). To do so, compute the correlation coefficient of the cyclical components of log hours worked (of all persons) with the cyclical component of log GDP (separately for each country). We use the `cor` function to compute the correlation, but need to add `use="complete.obs"` within the `cor` function to omit missing values in the computation. ``` r cyclical_correlation <- dataset_main %>% group_by(Country) %>% summarise(corr_hours_worked_GDP = cor(hours_worked_all_log_cycle, GDP_log_cycle, use="complete.obs")) print(cyclical_correlation) ``` Next, we compute the standard deviation of the cyclical components of log hours worked (of all persons, female, and male) separately for each country. We also compute the standard deviation of the cyclical component of log GDP. Then to obtain the relative standard deviation (also known as ‘cyclical volatility’) of the three `hours_worked` variables, we divide their standard deviation by the standard deviation of GDP. By only considering observations which are non-missing in all four variables (using `filter_at`), we ensure that the standard deviations are computed over the same time periods within each country. ``` r cyclical_volatility <- dataset_main %>% group_by(Country) %>% filter_at(vars(GDP,hours_worked_all_log_cycle,hours_worked_female_log_cycle,hours_worked_male_log_cycle),all_vars(!is.na(.))) %>% summarise(sd_total_hours_worked = sd(hours_worked_all_log_cycle), sd_female_hours_worked = sd(hours_worked_female_log_cycle), sd_male_hours_worked = sd(hours_worked_male_log_cycle), sd_GDP = sd(GDP_log_cycle)) %>% mutate(relative_sd_total_hours_worked = sd_total_hours_worked/sd_GDP, relative_sd_female_hours_worked = sd_female_hours_worked/sd_GDP, relative_sd_male_hours_worked = sd_male_hours_worked/sd_GDP) %>% select(Country,relative_sd_total_hours_worked,relative_sd_female_hours_worked,relative_sd_male_hours_worked) print(cyclical_volatility) ``` We now split the sample into a pre-global financial crisis (1983-2007) and a post-global financial crisis (2008-latest available) period, repeating the calculation of relative standard deviations, but only for female and male hours worked in the US. ``` r cyclical_volatility_pre_financial_crisis_US <- dataset_main %>% filter(Country == "United States" & Year >= 1983 & Year < 2008) %>% filter_at(vars(GDP,hours_worked_female_log_cycle,hours_worked_male_log_cycle),all_vars(!is.na(.))) %>% summarise(sd_female_hours_worked = sd(hours_worked_female_log_cycle), sd_male_hours_worked = sd(hours_worked_male_log_cycle), sd_GDP = sd(GDP_log_cycle)) %>% mutate(relative_sd_female_hours_worked = sd_female_hours_worked/sd_GDP, relative_sd_male_hours_worked = sd_male_hours_worked/sd_GDP) %>% select(relative_sd_female_hours_worked,relative_sd_male_hours_worked) %>% mutate(period = "1983-2007") cyclical_volatility_post_financial_crisis_US <- dataset_main %>% filter(Country == "United States" & Year >= 2008) %>% filter_at(vars(GDP,hours_worked_female_log_cycle,hours_worked_male_log_cycle),all_vars(!is.na(.))) %>% summarise(sd_female_hours_worked = sd(hours_worked_female_log_cycle), sd_male_hours_worked = sd(hours_worked_male_log_cycle), sd_GDP = sd(GDP_log_cycle)) %>% mutate(relative_sd_female_hours_worked = sd_female_hours_worked/sd_GDP, relative_sd_male_hours_worked = sd_male_hours_worked/sd_GDP) %>% select(relative_sd_female_hours_worked,relative_sd_male_hours_worked) %>% mutate(period = "2008--today") cyclical_volatility_comparison_US <- bind_rows(cyclical_volatility_pre_financial_crisis_US,cyclical_volatility_post_financial_crisis_US) %>% relocate(period) print(cyclical_volatility_comparison_US) ``` ``` r rm(cyclical_volatility_pre_financial_crisis_US,cyclical_volatility_post_financial_crisis_US) ``` Now we plot a grouped bar chart based on the `cyclical_volatility_comparison` table. First, we need to reshape the data into long format, with an extra column indicating sex. ``` r cyclical_volatility_comparison_US <- cyclical_volatility_comparison_US %>% gather(key = relative_sd_, value = relative_sd, relative_sd_female_hours_worked:relative_sd_male_hours_worked) %>% mutate(sex = case_when(relative_sd_ == "relative_sd_female_hours_worked" ~ "female", TRUE ~ "male")) %>% select(-relative_sd_) cyclical_volatility_comparison_US %>% ggplot() + aes(fill=sex, y=relative_sd, x=period) + geom_bar(position = "dodge", stat = "identity") + theme_bw() + # Add the vertical axis label labs(y = "Standard deviation of hours worked relative to GDP") + # Add the title and legend labels ggtitle("Cyclical volatility of hours worked by sex (United States)") + scale_fill_discrete(name = "", labels = c("Female", "Male")) + # Move the legend below the chart theme(legend.position="bottom") + # Remove the legend title and the horizontal axis label theme(legend.title=element_blank(), axis.title.x=element_blank()) ``` Then we create the same bar chart for Germany, using the same sample split as for the US: ``` r cyclical_volatility_pre_financial_crisis_Germany <- dataset_main %>% filter(Country == "Germany" & Year >= 1983 & Year < 2008) %>% filter_at(vars(GDP,hours_worked_female_log_cycle,hours_worked_male_log_cycle),all_vars(!is.na(.))) %>% summarise(sd_female_hours_worked = sd(hours_worked_female_log_cycle), sd_male_hours_worked = sd(hours_worked_male_log_cycle), sd_GDP = sd(GDP_log_cycle)) %>% mutate(relative_sd_female_hours_worked = sd_female_hours_worked/sd_GDP, relative_sd_male_hours_worked = sd_male_hours_worked/sd_GDP) %>% select(relative_sd_female_hours_worked,relative_sd_male_hours_worked) %>% mutate(period = "1983-2007") cyclical_volatility_post_financial_crisis_Germany <- dataset_main %>% filter(Country == "Germany" & Year >= 2008) %>% filter_at(vars(GDP,hours_worked_female_log_cycle,hours_worked_male_log_cycle),all_vars(!is.na(.))) %>% summarise(sd_female_hours_worked = sd(hours_worked_female_log_cycle), sd_male_hours_worked = sd(hours_worked_male_log_cycle), sd_GDP = sd(GDP_log_cycle)) %>% mutate(relative_sd_female_hours_worked = sd_female_hours_worked/sd_GDP, relative_sd_male_hours_worked = sd_male_hours_worked/sd_GDP) %>% select(relative_sd_female_hours_worked,relative_sd_male_hours_worked) %>% mutate(period = "2008--today") cyclical_volatility_comparison_Germany <- bind_rows(cyclical_volatility_pre_financial_crisis_Germany,cyclical_volatility_post_financial_crisis_Germany) %>% relocate(period) rm(cyclical_volatility_pre_financial_crisis_Germany,cyclical_volatility_post_financial_crisis_Germany) cyclical_volatility_comparison_Germany <- cyclical_volatility_comparison_Germany %>% gather(key = relative_sd_, value = relative_sd, relative_sd_female_hours_worked:relative_sd_male_hours_worked) %>% mutate(sex = case_when(relative_sd_ == "relative_sd_female_hours_worked" ~ "female", TRUE ~ "male")) %>% select(-relative_sd_) cyclical_volatility_comparison_Germany %>% ggplot() + aes(fill=sex, y=relative_sd, x=period) + geom_bar(position = "dodge", stat = "identity") + theme_bw() + # Add the vertical axis label labs(y = "Standard deviation of hours worked relative to GDP") + # Add the title and legend labels ggtitle("Cyclical volatility of hours worked by sex (Germany)") + scale_fill_discrete(name = "", labels = c("Female", "Male")) + # Move the legend below the chart theme(legend.position="bottom") + # Remove the legend title and the horizontal axis label theme(legend.title=element_blank(), axis.title.x=element_blank()) ``` To keep the workspace clean, delete all datasets except for `dataset_main`. ``` r rm(list=ls(pattern="cyclical")) ``` ### **R walk-through 8** Interpreting the trend component of the Hodrick-Prescott (HP) filter We apply the HP filter (`hpfilter`) to the log of average labour productivity (ALP, the variable `labour_productivity`) using the smoothing parameter lambda = 6.25, saving only the trend component. Before doing so, we group by country and filter out missing values in `labour_productivity`. Then we compute the annual growth rate (in %) of that trend component. ``` r dataset_main_filter <- dataset_main %>% group_by(Country) %>% filter(!is.na(labour_productivity)) %>% mutate(labour_productivity_trend = hpfilter(log(labour_productivity), freq=6.25)$trend, labour_productivity_trend_growth = (labour_productivity_trend/lag(labour_productivity_trend) - 1) * 100) dataset_main <- list(dataset_main,dataset_main_filter) %>% reduce(full_join) ``` ``` r rm(dataset_main_filter) ``` Plot the growth rate of the trend component of US productivity over time, separately for each country, adding in the recession periods as shaded areas. ``` r dataset_main %>% filter(Country == "United States") %>% ggplot(aes(Year)) + geom_line(aes(y = labour_productivity_trend_growth)) + theme_bw() + # Scale horizontal axis tick marks in 10-year intervals scale_x_continuous(breaks = seq(1960, 2020, by = 10)) + # Add the vertical axis label labs(y = "Growth rate in %") + # Add the title ggtitle("Trend productivity (United States)") + # Remove the horizontal axis label theme(axis.title.x=element_blank()) + # Add recession periods as shaded areas geom_rect(data = subset(dataset_main, Recession == 1 & Country == "United States"), aes(ymin = -Inf, ymax = Inf, xmin = Year-0.5, xmax = Year+0.5), alpha = 0.2) ``` Create the same chart for Germany. ``` r dataset_main %>% filter(Country == "Germany") %>% ggplot(aes(Year)) + geom_line(aes(y = labour_productivity_trend_growth)) + theme_bw() + # Scale horizontal axis tick marks in 10-year intervals scale_x_continuous(breaks = seq(1960, 2020, by = 10)) + # Add the vertical axis label labs(y = "Growth rate in %") + # Add the title ggtitle("Trend productivity (Germany)") + # Remove the horizontal axis label theme(axis.title.x=element_blank()) + # Add recession periods as shaded areas geom_rect(data = subset(dataset_main, Recession == 1 & Country == "Germany"), aes(ymin = -Inf, ymax = Inf, xmin = Year-0.5, xmax = Year+0.5), alpha = 0.2) ``` Apply the HP filter (`hpfilter`) to the log of female and male LFP rates with lambda = 6.25. To exclude the jumps in LFP due to changes in methodology, apply the filter only to years 1963 and onwards. Save the trend components, then compute their annual growth rates in percent. ``` r dataset_main_filter <- dataset_main %>% group_by(Country) %>% filter(!is.na(LFP_rate_female) & !is.na(LFP_rate_male) & Year >= 1963) %>% mutate(LFP_rate_female_trend = hpfilter(log(LFP_rate_female), freq=6.25)$trend, LFP_rate_male_trend = hpfilter(log(LFP_rate_male), freq=6.25)$trend, LFP_rate_female_trend_growth = (LFP_rate_female_trend/lag(LFP_rate_female_trend) - 1) * 100, LFP_rate_male_trend_growth = (LFP_rate_male_trend/lag(LFP_rate_male_trend) - 1) * 100) dataset_main <- list(dataset_main,dataset_main_filter) %>% reduce(full_join) ``` ``` r rm(dataset_main_filter) ``` Add the growth rates of the log LFP rates to the productivity charts. ``` r dataset_main %>% filter(Country == "United States") %>% ggplot(aes(Year)) + geom_line(aes(y = LFP_rate_female_trend_growth, colour = "labour_productivity_trend_growth")) + geom_line(aes(y = LFP_rate_male_trend_growth, colour = "LFP_rate_male_trend_growth")) + geom_line(aes(y = labour_productivity_trend_growth, colour = "LFP_rate_female_trend_growth")) + theme_bw() + # Scale horizontal axis tick marks in 10-year intervals scale_x_continuous(breaks = seq(1960, 2020, by = 10)) + # Add the vertical axis label labs(y = "Growth rate in %") + # Add the title and legend labels ggtitle("Trend productivity and labour force participation (United States)") + scale_colour_discrete(name = "", labels = c("Female LFP trend", "ALP trend", "Male LFP trend")) + # Move the legend below the chart theme(legend.position="bottom") + # Remove the horizontal axis label and the legend title theme(axis.title.x=element_blank(), legend.title=element_blank()) + # Add recession periods as shaded areas geom_rect(data = subset(dataset_main, Recession == 1 & Country == "United States"), aes(ymin = -Inf, ymax = Inf, xmin = Year-0.5, xmax = Year+0.5), alpha = 0.2) ``` ``` r dataset_main %>% filter(Country == "Germany") %>% ggplot(aes(Year)) + geom_line(aes(y = LFP_rate_female_trend_growth, colour = "LFP_rate_female_trend_growth")) + geom_line(aes(y = LFP_rate_male_trend_growth, colour = "LFP_rate_male_trend_growth")) + geom_line(aes(y = labour_productivity_trend_growth, colour = "labour_productivity_trend_growth")) + theme_bw() + # Scale horizontal axis tick marks in 10-year intervals scale_x_continuous(breaks = seq(1960, 2020, by = 10)) + # Add the vertical axis label labs(y = "Growth rate in %") + # Add the title and legend labels ggtitle("Trend productivity and labour force participation (Germany)") + scale_colour_discrete(name = "", labels = c("ALP trend","Female LFP trend", "Male LFP trend")) + # Move the legend below the chart theme(legend.position="bottom") + # Remove the horizontal axis label and the legend title theme(axis.title.x=element_blank(), legend.title=element_blank()) + # Add recession periods as shaded areas geom_rect(data = subset(dataset_main, Recession == 1 & Country == "Germany"), aes(ymin = -Inf, ymax = Inf, xmin = Year-0.5, xmax = Year+0.5), alpha = 0.2) ``` ### **R walk-through 9** Jobless recoveries Plot the log of the LFP rate of all persons over time (starting in 1963), adding in the recession periods as shaded areas. ``` r dataset_main %>% filter(Country == "United States" & Year >= 1963) %>% ggplot(aes(Year)) + geom_line(aes(y = log(LFP_rate_all))) + theme_bw() + # Scale horizontal axis tick marks in 10-year intervals scale_x_continuous(breaks = seq(1960, 2020, by = 10)) + # Add the vertical axis label labs(y = "Log value") + # Add the title and legend labels ggtitle("Labour force participation over time (United States)") + # Remove the horizontal axis label theme(axis.title.x=element_blank()) + # Add recession periods as shaded areas geom_rect(data = subset(dataset_main, Recession == 1 & Country == "United States"), aes(ymin = -Inf, ymax = Inf, xmin = Year-0.5, xmax = Year+0.5), alpha = 0.2) ``` Now redo this chart, split by sex. ``` r dataset_main %>% filter(Country == "United States" & Year >= 1963) %>% ggplot(aes(Year)) + geom_line(aes(y = log(LFP_rate_female))) + theme_bw() + # Scale horizontal axis tick marks in 10-year intervals scale_x_continuous(breaks = seq(1960, 2020, by = 10)) + # Add the vertical axis label labs(y = "Log value") + # Add the title and legend labels ggtitle("Female labour force participation over time (United States)") + # Remove the horizontal axis label theme(axis.title.x=element_blank()) + # Add recession periods as shaded areas geom_rect(data = subset(dataset_main, Recession == 1 & Country == "United States"), aes(ymin = -Inf, ymax = Inf, xmin = Year-0.5, xmax = Year+0.5), alpha = 0.2) ``` ``` r dataset_main %>% filter(Country == "United States" & Year >= 1963) %>% ggplot(aes(Year)) + geom_line(aes(y = log(LFP_rate_male))) + theme_bw() + # Scale horizontal axis tick marks in 10-year intervals scale_x_continuous(breaks = seq(1960, 2020, by = 10)) + # Add the vertical axis label labs(y = "Log value") + # Add the title and legend labels ggtitle("Male labour force participation over time (United States)") + # Remove the horizontal axis label theme(axis.title.x=element_blank()) + # Add recession periods as shaded areas geom_rect(data = subset(dataset_main, Recession == 1 & Country == "United States"), aes(ymin = -Inf, ymax = Inf, xmin = Year-0.5, xmax = Year+0.5), alpha = 0.2) ``` Redo all three charts for Germany. ``` r dataset_main %>% filter(Country == "Germany") %>% ggplot(aes(Year)) + geom_line(aes(y = log(LFP_rate_all))) + theme_bw() + # Scale horizontal axis tick marks in 10-year intervals scale_x_continuous(breaks = seq(1960, 2020, by = 10)) + # Add the vertical axis label labs(y = "Log value") + # Add the title and legend labels ggtitle("Labour force participation over time (Germany)") + # Remove the horizontal axis label theme(axis.title.x=element_blank()) + # Add recession periods as shaded areas geom_rect(data = subset(dataset_main, Recession == 1 & Country == "Germany"), aes(ymin = -Inf, ymax = Inf, xmin = Year-0.5, xmax = Year+0.5), alpha = 0.2) ``` ``` r dataset_main %>% filter(Country == "Germany") %>% ggplot(aes(Year)) + geom_line(aes(y = log(LFP_rate_female))) + theme_bw() + # Scale horizontal axis tick marks in 10-year intervals scale_x_continuous(breaks = seq(1960, 2020, by = 10)) + # Add the vertical axis label labs(y = "Log value") + # Add the title and legend labels ggtitle("Female labour force participation over time (Germany)") + # Remove the horizontal axis label theme(axis.title.x=element_blank()) + # Add recession periods as shaded areas geom_rect(data = subset(dataset_main, Recession == 1 & Country == "Germany"), aes(ymin = -Inf, ymax = Inf, xmin = Year-0.5, xmax = Year+0.5), alpha = 0.2) ``` ``` r dataset_main %>% filter(Country == "Germany") %>% ggplot(aes(Year)) + geom_line(aes(y = log(LFP_rate_male))) + theme_bw() + # Scale horizontal axis tick marks in 10-year intervals scale_x_continuous(breaks = seq(1960, 2020, by = 10)) + # Add the vertical axis label labs(y = "Log value") + # Add the title and legend labels ggtitle("Male labour force participation over time (Germany)") + # Remove the horizontal axis label theme(axis.title.x=element_blank()) + # Add recession periods as shaded areas geom_rect(data = subset(dataset_main, Recession == 1 & Country == "Germany"), aes(ymin = -Inf, ymax = Inf, xmin = Year-0.5, xmax = Year+0.5), alpha = 0.2) ```