Dataset Clean

civil <- read_csv("Civil Labor Ballance.csv")
civil2 <- civil %>%
  mutate(FIPS = sprintf("%05d", FIPS),
         DATE = as.Date(DATE),
         treat = as.numeric(STATE == "LA"),
         #treat = as.numeric(STATE == "LA" & DATE >= as.Date("2005-08-01")), For SDID
         t = as.numeric(DATE > as.Date("2005-08-01")),
         year = as.numeric(format(DATE, "%Y"))) %>%
  group_by(FIPS) %>%
  filter(!any(is.na(Civil.Labor))) %>%
  ungroup()

## SUbset Data, because we want focus on before 2011
civil3 <- civil2 %>%
  filter(DATE <= as.Date("2010-12-01") & !is.na(Civil.Labor))

Ploting

I skip, to the final ploting and subseting Fips for Felm. Here I simulate subset the FIPS using sampling function as much 50 FIPS, I also manually add FIPS in LA, because what makes error is because none of LA county are selected in cnty_sdid$FIPS, it will makes treat in felm value all 0,

I also change a bit for Fix effect from County to FIPS, since County may have duplicated.

## Makes Date Dummy
#cnty_sdid <- map_sdid %>%
#  distinct(FIPS, COUNTY)

#df4 <- civil3 %>%
#  filter(FIPS %in% cnty_sdid$FIPS) %>%
#  filter(DATE < as.Date("2010-02-01"))

## For Sampling Code ~ MY REVISE,
##### NOTE : The reason that cause, error in your code, is because none of LA county are selected in cnty_sdid$FIPS, 
##### since it comes from weight of SDID and SC, why because county in LA State are the treated one.
cnty_sdid<-sample(unique(civil3$FIPS),size = 50,replace = F)
cnty_sdid[51]<-"22001" ## here I manually ADD LA Treatment 

df4 <- civil3 %>%
  filter(FIPS %in% cnty_sdid) %>%
  filter(DATE < as.Date("2010-02-01"))

View(df4)

dloop <- data.frame(DATE = unique(df4$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)) {
  df4[paste0("d.", dloop$year[i],".",dloop$month[i])] <- as.numeric(df4$DATE == unique(df4$DATE)[i])
}

## RUN MODEL
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)| FIPS + DATE | 0 | FIPS,
                data = df4)

## 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)])
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.04515872 0.01447151 0.07352288 0.016794553
## 2 2000-02-01 0.03988815 0.01653763 0.07230191 0.007474389
## 3 2000-03-01 0.04593513 0.01627898 0.07784193 0.014028331
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))
pl