Labour Regression Models

Dataset

library(tidyverse)

## ── Attaching core tidyverse packages ──────────────────────── tidyverse 2.0.0 ──
## ✔ dplyr     1.1.1     ✔ readr     2.1.4
## ✔ forcats   1.0.0     ✔ stringr   1.5.0
## ✔ ggplot2   3.4.2     ✔ tibble    3.2.1
## ✔ lubridate 1.9.2     ✔ tidyr     1.3.0
## ✔ purrr     1.0.1     
## ── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
## ✖ dplyr::filter() masks stats::filter()
## ✖ dplyr::lag()    masks stats::lag()
## ℹ Use the conflicted package (<http://conflicted.r-lib.org/>) to force all conflicts to become errors

library(ggplot2)
library(readr)
library(dplyr)
library(stringr)
library(lmtest)

## Loading required package: zoo
## 
## Attaching package: 'zoo'
## 
## The following objects are masked from 'package:base':
## 
##     as.Date, as.Date.numeric

library(lfe)

## Loading required package: Matrix
## 
## Attaching package: 'Matrix'
## 
## The following objects are masked from 'package:tidyr':
## 
##     expand, pack, unpack
## 
## 
## Attaching package: 'lfe'
## 
## The following object is masked from 'package:lmtest':
## 
##     waldtest

df2 <- read.csv("Civil Labor Ballance.csv")
df2 <- df2 %>%
  mutate(
    DATE = as.Date(DATE),
    treat = as.numeric(STATE == "LA"),
    t = as.numeric(DATE > as.Date("2005-08-01")),
    year = as.numeric(format(DATE, "%Y"))
    )

## SUbset Data, because we want focus on before Feb 2010
df3 <- df2 %>%
    filter(DATE < as.Date("2010-02-01"))

## Makes Date Dummy
dloop <- data.frame(DATE = unique(df3$DATE))

dloop <- dloop %>%
  mutate(year = format(DATE,"%Y"),
         month=str_pad(string = format(DATE,"%m"),width = 2,side = "left",pad = 0))

## YEAR DUMMY
for (i in 1:nrow(dloop)) {
  df3[paste0("d.", dloop$year[i],".",dloop$month[i])] <- as.numeric(df3$DATE == unique(df3$DATE)[i])
}

Civil Labor

did.reg <- felm(log(Civil.Labor) ~ 
           treat:(d.2000.01+d.2000.02+d.2000.03+d.2000.04+d.2000.05+d.2000.06+d.2000.07+d.2000.08+d.2000.09+d.2000.10+d.2000.11+d.2000.12+
           d.2001.01+d.2001.02+d.2001.03+d.2001.04+d.2001.05+d.2001.06+d.2001.07+d.2001.08+d.2001.09+d.2001.10+d.2001.11+d.2001.12+
           d.2002.01+d.2002.02+d.2002.03+d.2002.04+d.2002.05+d.2002.06+d.2002.07+d.2002.08+d.2002.09+d.2002.10+d.2002.11+d.2002.12+
           d.2003.01+d.2003.02+d.2003.03+d.2003.04+d.2003.05+d.2003.06+d.2003.07+d.2003.08+d.2003.09+d.2003.10+d.2003.11+d.2003.12+
           d.2004.01+d.2004.02+d.2004.03+d.2004.04+d.2004.05+d.2004.06+d.2004.07+d.2004.08+d.2004.09+d.2004.10+d.2004.11+d.2004.12+
           d.2005.01+d.2005.02+d.2005.03+d.2005.04+d.2005.05+d.2005.06+d.2005.07+d.2005.09+d.2005.10+d.2005.11+d.2005.12+
           d.2006.01+d.2006.02+d.2006.03+d.2006.04+d.2006.05+d.2006.06+d.2006.07+d.2006.08+d.2006.09+d.2006.10+d.2006.11+d.2006.12+
           d.2007.01+d.2007.02+d.2007.03+d.2007.04+d.2007.05+d.2007.06+d.2007.07+d.2007.08+d.2007.09+d.2007.10+d.2007.11+d.2007.12+
           d.2008.01+d.2008.02+d.2008.03+d.2008.04+d.2008.05+d.2008.06+d.2008.07+d.2008.08+d.2008.09+d.2008.10+d.2008.11+d.2008.12+
           d.2009.01+d.2009.02+d.2009.03+d.2009.04+d.2009.05+d.2009.06+d.2009.07+d.2009.08+d.2009.09+d.2009.10+d.2009.11+d.2009.12+
           d.2010.01)| COUNTY + DATE | 0 | COUNTY,
              data = df3)

## Export Coeffiecent and test for Heteroscedasticity SE
coef_df <- coeftest(did.reg)

coef_df <- coef_df[,] %>%
    as_tibble() %>%
    mutate(variables = rownames(coef_df))

coef_df$date <- as.Date(paste0(coef_df$variables, ".01"),
    format = "treat:d.%Y.%m.%d")
 
event_df <- data.frame(date=seq.Date(as.Date("2000-01-01"),as.Date("2010-01-01"),by = "month"))
event_df <- left_join(event_df, coef_df[, c(1, 2, 6)])

## Joining with `by = join_by(date)`

event_df[is.na(event_df)] <- 0

colnames(event_df) <- c("time","coef","se")

## Makes Confident Interval
event_df$ci_upper <- event_df$coef + 1.96 * event_df$se
event_df$ci_lower <- event_df$coef - 1.96 * event_df$se
head(event_df,3)

##         time       coef          se   ci_upper    ci_lower
## 1 2000-01-01 0.03328925 0.008043848 0.04905519 0.017523304
## 2 2000-02-01 0.02356477 0.007633328 0.03852609 0.008603447
## 3 2000-03-01 0.02857348 0.007674360 0.04361523 0.013531738

dates <- seq.Date(as.Date("2000-01-01"),as.Date("2010-01-01"),by='month')
# Define start and end dates for plot
start_date <- as.Date("2004-08-01")
end_date <- as.Date("2010-01-01")

# Create sequence of dates for x-axis
dates_seq <- seq.Date(from = start_date, to = end_date, by = "month")

# Subset event_df to only include dates after or equal to start_date
event_df_sub <- event_df[dates >= start_date, ]

pl <- ggplot(event_df, aes(x = dates, y = coef)) +
  geom_point() +
  geom_line()+
  #geom_errorbar(aes(ymin = ci_lower, ymax = ci_upper), width = 0.2) +
  geom_ribbon(aes(ymin = ci_lower, ymax = ci_upper), alpha = 0.2) +
  geom_hline(yintercept = 0, linetype = "dashed", color = "red") +
  geom_vline(xintercept = as.Date("2005-08-01"), linetype = "dashed", color = "blue") +
  scale_x_date(date_labels = "%b %Y", limits = c(start_date, end_date), breaks = "1 year") +
  labs(title = "Civil Labor") +
  theme(axis.text.x = element_text(angle = 45, hjust = 1),
    plot.title = element_text(size = 14, face = "bold"),
    axis.title.x = element_blank(),
    axis.title.y = element_blank(),
    axis.text = element_text(size = 12),
    legend.text = element_text(size = 12),
    plot.caption = element_text(size = 12, face = c("italic", "bold")),
    axis.line = element_line(color = "black", size = 0.5),
    panel.background = element_blank(),
    panel.grid.major = element_blank(),
    panel.grid.minor = element_line(color = "black", size = 0.1))

## Warning: The `size` argument of `element_line()` is deprecated as of ggplot2 3.4.0.
## ℹ Please use the `linewidth` argument instead.
## This warning is displayed once every 8 hours.
## Call `lifecycle::last_lifecycle_warnings()` to see where this warning was
## generated.

## Warning: Vectorized input to `element_text()` is not officially supported.
## ℹ Results may be unexpected or may change in future versions of ggplot2.

pl

## Warning: Removed 55 rows containing missing values (`geom_point()`).

## Warning: Removed 55 rows containing missing values (`geom_line()`).

Employee Person

df4 <- read.csv("Merge Employee Person Full-Location Ballance.csv")
df4 <- df4 %>%
  mutate(
    DATE = as.Date(DATE)
  )

df4 <- df4 %>%
  mutate(
    treat = as.numeric(STATE == "LA"),
    t = as.numeric(DATE > as.Date("2005-08-01"))
  )

df4 <- df4 %>% 
  mutate(
    year = as.numeric(format(DATE, "%Y"))
  )

## SUbset Data, because we want focus on before Feb 2010
df5 <- df4 %>% filter(DATE <"2010-02-01")

## Makes Date Dummy
dloop<-data.frame(DATE=unique(df5$DATE))
dloop<-dloop%>%
  mutate(year=format(DATE,"%Y"),
         month=str_pad(string = format(DATE,"%m"),width = 2,side = "left",pad = 0))

## YEAR DUMMY
for (i in 1:nrow(dloop)) {
  df5[paste0("d.", dloop$year[i],".",dloop$month[i])] <- as.numeric(df5$DATE == unique(df5$DATE)[i])
}

## RUN MODEL
did.reg2 <- felm(log(Employee) ~ 
           treat:(d.2000.01+d.2000.02+d.2000.03+d.2000.04+d.2000.05+d.2000.06+d.2000.07+d.2000.08+d.2000.09+d.2000.10+d.2000.11+d.2000.12+
           d.2001.01+d.2001.02+d.2001.03+d.2001.04+d.2001.05+d.2001.06+d.2001.07+d.2001.08+d.2001.09+d.2001.10+d.2001.11+d.2001.12+
           d.2002.01+d.2002.02+d.2002.03+d.2002.04+d.2002.05+d.2002.06+d.2002.07+d.2002.08+d.2002.09+d.2002.10+d.2002.11+d.2002.12+
           d.2003.01+d.2003.02+d.2003.03+d.2003.04+d.2003.05+d.2003.06+d.2003.07+d.2003.08+d.2003.09+d.2003.10+d.2003.11+d.2003.12+
           d.2004.01+d.2004.02+d.2004.03+d.2004.04+d.2004.05+d.2004.06+d.2004.07+d.2004.08+d.2004.09+d.2004.10+d.2004.11+d.2004.12+
           d.2005.01+d.2005.02+d.2005.03+d.2005.04+d.2005.05+d.2005.06+d.2005.07+d.2005.09+d.2005.10+d.2005.11+d.2005.12+
           d.2006.01+d.2006.02+d.2006.03+d.2006.04+d.2006.05+d.2006.06+d.2006.07+d.2006.08+d.2006.09+d.2006.10+d.2006.11+d.2006.12+
           d.2007.01+d.2007.02+d.2007.03+d.2007.04+d.2007.05+d.2007.06+d.2007.07+d.2007.08+d.2007.09+d.2007.10+d.2007.11+d.2007.12+
           d.2008.01+d.2008.02+d.2008.03+d.2008.04+d.2008.05+d.2008.06+d.2008.07+d.2008.08+d.2008.09+d.2008.10+d.2008.11+d.2008.12+
           d.2009.01+d.2009.02+d.2009.03+d.2009.04+d.2009.05+d.2009.06+d.2009.07+d.2009.08+d.2009.09+d.2009.10+d.2009.11+d.2009.12+
           d.2010.01)| COUNTY + DATE | 0 | COUNTY,
              data = df5)

## Export Coeffiecent and test for Heteroscedasticity SE
coef_df <- coeftest(did.reg2)

coef_df <- coef_df[,] %>%
    as_tibble() %>%
    mutate(variables = rownames(coef_df))

coef_df$date <- as.Date(paste0(coef_df$variables, ".01"),
    format = "treat:d.%Y.%m.%d")
 
event_df <- data.frame(date=seq.Date(as.Date("2000-01-01"),as.Date("2010-01-01"),by = "month"))
event_df <- left_join(event_df, coef_df[, c(1, 2, 6)])

## Joining with `by = join_by(date)`

event_df[is.na(event_df)] <- 0

colnames(event_df)<-c("time","coef","se")

## Makes Confident Interval
event_df$ci_upper <- event_df$coef + 1.96 * event_df$se
event_df$ci_lower <- event_df$coef - 1.96 * event_df$se
head(event_df,3)

##         time       coef          se   ci_upper   ci_lower
## 1 2000-01-01 0.03699934 0.008874423 0.05439321 0.01960547
## 2 2000-02-01 0.03390130 0.008447374 0.05045815 0.01734445
## 3 2000-03-01 0.03369194 0.008599336 0.05054664 0.01683724

dates<-seq.Date(as.Date("2000-01-01"),as.Date("2010-01-01"),by='month')
# Define start and end dates for plot
start_date <- as.Date("2004-08-01")
end_date <- as.Date("2010-01-01")

# Create sequence of dates for x-axis
dates_seq <- seq.Date(from = start_date, to = end_date, by = "month")

# Subset event_df to only include dates after or equal to start_date
event_df_sub <- event_df[dates >= start_date, ]

pl2 <- ggplot(event_df, aes(x = dates, y = coef)) +
  geom_point() +
  geom_line()+
  #geom_errorbar(aes(ymin = ci_lower, ymax = ci_upper), width = 0.2) +
  geom_ribbon(aes(ymin = ci_lower, ymax = ci_upper), alpha = 0.2) +
  geom_hline(yintercept = 0, linetype = "dashed", color = "red") +
  geom_vline(xintercept = as.Date("2005-08-01"), linetype = "dashed", color = "blue") +
  scale_x_date(date_labels = "%b %Y", limits = c(start_date, end_date), breaks = "1 year") +
  labs(title = "Employee Person") +
  theme(axis.text.x = element_text(angle = 45, hjust = 1),
    plot.title = element_text(size = 14, face = "bold"),
    axis.title.x = element_blank(),
    axis.title.y = element_blank(),
    axis.text = element_text(size = 12),
    legend.text = element_text(size = 12),
    plot.caption = element_text(size = 12, face = c("italic", "bold")),
    axis.line = element_line(color = "black", size = 0.5),
    panel.background = element_blank(),
    panel.grid.major = element_blank(),
    panel.grid.minor = element_line(color = "black", size = 0.1))

## Warning: Vectorized input to `element_text()` is not officially supported.
## ℹ Results may be unexpected or may change in future versions of ggplot2.

pl2

## Warning: Removed 55 rows containing missing values (`geom_point()`).

## Warning: Removed 55 rows containing missing values (`geom_line()`).

Unemployment Rate

df6 <- read.csv("Merge UnemRate ballance V2 FULL-Location.csv")
df6 <- df6 %>% 
  mutate(
    DATE = as.Date(DATE)
  )

df6 <- df6 %>% 
  mutate(
    treat = as.numeric(STATE == "LA"),
    t = as.numeric(DATE > as.Date("2005-08-01"))
  )

df6 <- df6 %>% 
  mutate(
    year = as.numeric(format(DATE, "%Y"))
  )

## SUbset Data, because we want focus on before Feb 2010
df7 <- df6 %>% filter(DATE <"2010-02-01")

## Makes Date Dummy
dloop<-data.frame(DATE=unique(df7$DATE))
dloop<-dloop%>%
  mutate(year=format(DATE,"%Y"),
         month=str_pad(string = format(DATE,"%m"),width = 2,side = "left",pad = 0))

## YEAR DUMMY
for (i in 1:nrow(dloop)) {
  df7[paste0("d.", dloop$year[i],".",dloop$month[i])] <- as.numeric(df7$DATE == unique(df7$DATE)[i])
}

did.reg3 <- felm(log(Unemployee.Rate) ~ 
           treat:(d.2000.01+d.2000.02+d.2000.03+d.2000.04+d.2000.05+d.2000.06+d.2000.07+d.2000.08+d.2000.09+d.2000.10+d.2000.11+d.2000.12+
           d.2001.01+d.2001.02+d.2001.03+d.2001.04+d.2001.05+d.2001.06+d.2001.07+d.2001.08+d.2001.09+d.2001.10+d.2001.11+d.2001.12+
           d.2002.01+d.2002.02+d.2002.03+d.2002.04+d.2002.05+d.2002.06+d.2002.07+d.2002.08+d.2002.09+d.2002.10+d.2002.11+d.2002.12+
           d.2003.01+d.2003.02+d.2003.03+d.2003.04+d.2003.05+d.2003.06+d.2003.07+d.2003.08+d.2003.09+d.2003.10+d.2003.11+d.2003.12+
           d.2004.01+d.2004.02+d.2004.03+d.2004.04+d.2004.05+d.2004.06+d.2004.07+d.2004.08+d.2004.09+d.2004.10+d.2004.11+d.2004.12+
           d.2005.01+d.2005.02+d.2005.03+d.2005.04+d.2005.05+d.2005.06+d.2005.07+d.2005.09+d.2005.10+d.2005.11+d.2005.12+
           d.2006.01+d.2006.02+d.2006.03+d.2006.04+d.2006.05+d.2006.06+d.2006.07+d.2006.08+d.2006.09+d.2006.10+d.2006.11+d.2006.12+
           d.2007.01+d.2007.02+d.2007.03+d.2007.04+d.2007.05+d.2007.06+d.2007.07+d.2007.08+d.2007.09+d.2007.10+d.2007.11+d.2007.12+
           d.2008.01+d.2008.02+d.2008.03+d.2008.04+d.2008.05+d.2008.06+d.2008.07+d.2008.08+d.2008.09+d.2008.10+d.2008.11+d.2008.12+
           d.2009.01+d.2009.02+d.2009.03+d.2009.04+d.2009.05+d.2009.06+d.2009.07+d.2009.08+d.2009.09+d.2009.10+d.2009.11+d.2009.12+
           d.2010.01)| COUNTY + DATE | 0 | COUNTY,
              data = df7)

## Export Coeffiecent and test for Heteroscedasticity SE
coef_df <- coeftest(did.reg3)

coef_df <- coef_df[,] %>%
    as_tibble() %>%
    mutate(variables = rownames(coef_df))

coef_df$date <- as.Date(paste0(coef_df$variables, ".01"),
    format = "treat:d.%Y.%m.%d")
 
event_df <- data.frame(date=seq.Date(as.Date("2000-01-01"),as.Date("2010-01-01"),by = "month"))
event_df <- left_join(event_df, coef_df[, c(1, 2, 6)])

## Joining with `by = join_by(date)`

event_df[is.na(event_df)] <- 0

colnames(event_df)<-c("time","coef","se")

## Makes Confident Interval
event_df$ci_upper <- event_df$coef + 1.96 * event_df$se
event_df$ci_lower <- event_df$coef - 1.96 * event_df$se
head(event_df,3)

##         time        coef         se     ci_upper    ci_lower
## 1 2000-01-01 -0.03060959 0.02124344  0.011027556 -0.07224673
## 2 2000-02-01 -0.14777513 0.02316540 -0.102370943 -0.19317932
## 3 2000-03-01 -0.04987232 0.02244510 -0.005879936 -0.09386471

dates<-seq.Date(as.Date("2000-01-01"),as.Date("2010-01-01"),by='month')
# Define start and end dates for plot
start_date <- as.Date("2004-08-01")
end_date <- as.Date("2010-01-01")

# Create sequence of dates for x-axis
dates_seq <- seq.Date(from = start_date, to = end_date, by = "month")

# Subset event_df to only include dates after or equal to start_date
event_df_sub <- event_df[dates >= start_date, ]

pl3 <- ggplot(event_df, aes(x = dates, y = coef)) +
  geom_point() +
  geom_line()+
  #geom_errorbar(aes(ymin = ci_lower, ymax = ci_upper), width = 0.2) +
  geom_ribbon(aes(ymin = ci_lower, ymax = ci_upper), alpha = 0.2) +
  geom_hline(yintercept = 0, linetype = "dashed", color = "red") +
  geom_vline(xintercept = as.Date("2005-08-01"), linetype = "dashed", color = "blue") +
  scale_x_date(date_labels = "%b %Y", limits = c(start_date, end_date), breaks = "1 year") +
  labs(x = "Date", y = "Coefficient", title = "") +
  labs(title = "Unemployment Rate") +
  theme(axis.text.x = element_text(angle = 45, hjust = 1),
    plot.title = element_text(size = 14, face = "bold"),
    axis.title.x = element_blank(),
    axis.title.y = element_blank(),
    axis.text = element_text(size = 12),
    legend.text = element_text(size = 12),
    plot.caption = element_text(size = 12, face = c("italic", "bold")),
    axis.line = element_line(color = "black", size = 0.5),
    panel.background = element_blank(),
    panel.grid.major = element_blank(),
    panel.grid.minor = element_line(color = "black", size = 0.1))

## Warning: Vectorized input to `element_text()` is not officially supported.
## ℹ Results may be unexpected or may change in future versions of ggplot2.

pl3

## Warning: Removed 55 rows containing missing values (`geom_point()`).

## Warning: Removed 55 rows containing missing values (`geom_line()`).

– Plot of All Graphs –

ggplot() +
  annotation_custom(ggplotGrob(pl), xmin = 0, xmax = 0.5, ymin = 0.5, ymax = 1) +
  annotation_custom(ggplotGrob(pl2), xmin = 0.5, xmax = 1, ymin = 0.5, ymax = 1) +
  annotation_custom(ggplotGrob(pl3), xmin = 0.25, xmax = 0.75, ymin = 0, ymax = 0.5)+
  coord_cartesian(xlim = c(0, 1), ylim = c(0, 1), clip = "off") +
  theme_bw()

## Warning: Removed 55 rows containing missing values (`geom_point()`).

## Warning: Removed 55 rows containing missing values (`geom_line()`).

## Warning: Removed 55 rows containing missing values (`geom_point()`).

## Warning: Removed 55 rows containing missing values (`geom_line()`).

## Warning: Removed 55 rows containing missing values (`geom_point()`).

## Warning: Removed 55 rows containing missing values (`geom_line()`).