--- title: "Exercise template for 'Post-earnings announcement drift'" author: Your name format: html bibliography: book.bib --- ```{r} #| include: false #| message: false library(dplyr, warn.conflicts = FALSE) library(tidyr) library(stringr) library(lubridate) library(purrr) library(DBI) library(dbplyr) library(ggplot2) library(farr) library(modelsummary) library(knitr) library(furrr) ``` ```{r} #| label: database-connection #| include: false #| cache: false db <- dbConnect(RPostgres::Postgres(), bigint = "integer", check_interrupts = TRUE) funda <- tbl(db, Id(schema = "comp", table = "funda")) fundq <- tbl(db, Id(schema = "comp", table = "fundq")) company <- tbl(db, Id(schema = "comp", table = "company")) ccmxpf_lnkhist <- tbl(db, Id(schema = "crsp", table = "ccmxpf_lnkhist")) funda_mod <- funda |> filter(indfmt == "INDL", datafmt == "STD", consol == "C", popsrc == "D") fundq_mod <- fundq |> filter(indfmt == "INDL", datafmt == "STD", consol == "C", popsrc == "D") ``` ## Fiscal years {#sec-fyear} ```{r} #| label: fiscal-data #| include: false #| cache: true #| cache.vars: [fyear_data] fyear_data <- funda_mod |> distinct(fyear, datadate) |> mutate(month = month(datadate), year = year(datadate)) |> filter(!is.na(fyear)) |> collect() ``` ```{r} #| echo: true fms <- list(lm(fyear ~ factor(month) + year - 1, data = fyear_data), lm(fyear ~ month + year, data = fyear_data)) ``` ```{r} #| include: false fyear_plot_data <- fyear_data |> mutate(pred_1 = predict(fms[[1]]), pred_2 = predict(fms[[2]]), resid_1 = pred_1 - fyear, resid_2 = pred_2 - fyear) ``` ```{r} #| label: tbl-fyear #| tbl-cap: Regression of fiscal year on month indicators #| echo: false #| tbl-pos: H modelsummary(fms, estimate = "{estimate}", statistic = NULL, gof_map = "r.squared") ``` ```{r} #| echo: true #| fig-pos: H #| label: fig-plot-sample #| fig-height: 3 #| fig-cap: Plot of fyear against datadate for sample values #| fig-alt: "Plot of fyear against datadate. Plot shows scatterplot of actual data plus predicted values from two models. One model does a relatively poor job and the other model perfectly predicts fyear given datadate." fyear_plot_data |> filter(year %in% c(2001, 2002)) |> select(datadate, fyear, pred_1, pred_2) |> distinct() |> arrange(datadate) |> ggplot(aes(x = datadate)) + geom_point(aes(y = fyear, color = "fyear")) + geom_line(aes(y = pred_1, color = "pred_1")) + geom_line(aes(y = pred_2, color = "pred_2")) + scale_x_date(date_breaks = "1 month") + theme(axis.text.x = element_text(angle = 90), legend.position = "inside", legend.position.inside = c(.85, .35)) ``` ### Exercises 1. What is different between the first and second models in `fms`? What is `factor()` doing here? 2. What does the inclusion of `- 1` do in the first model in `fms`? Would the omission of `- 1` affect the fit of that model? Would it affect the interpretability of results? Would the inclusion of `- 1` affect the fit of the second model in `fms`? Would it affect the interpretability of results? 3. Does @fig-plot-sample help you understand what's going on? Why did we focus on a relatively short period in @fig-plot-sample? (*Hint*: What happens if you remove the line `filter(year %in% c(2001, 2002)) |>` from the code?) 4. Using `year()` and `month()`, add some code along the lines of `mutate(fyear_calc = ...)` to *calculate* `fyear`. Use this code to create a function `fiscal_year(datadate)`. Check that you match `fyear` in each case. ## Quarterly data {#sec-qtrs} ```{r} #| label: quarterly-data #| include: false #| cache: true #| cache.vars: [fundq_probs, link_table, ni_annual, ni_qtrly, ni_merged, plot_data] fundq_probs <- fundq_mod |> group_by(gvkey, datadate) |> filter(n() > 1) |> ungroup() |> select(gvkey, datadate, fyearq, fqtr, fyr, rdq) |> arrange(gvkey, datadate) |> collect() link_table <- fundq_mod |> rename(datadateq = datadate) |> select(gvkey:fyr) |> mutate(year = if_else(fyr <= 5L, fyearq + 1L, fyearq)) |> mutate(month = lpad(as.character(fyr), 2L, "0")) |> mutate(datadate = as.Date(str_c(year, month, '01', sep = "-"))) |> mutate(datadate = as.Date(datadate + months(1) - days(1))) |> select(-month, -year, -fqtr) |> collect() ni_annual <- funda_mod |> select(gvkey, datadate, fyr, sale, ni) |> collect() ni_qtrly <- fundq_mod |> select(gvkey, datadate, fyr, fqtr, saleq, niq, ibq) |> rename(datadateq = datadate) |> collect() ni_merged <- ni_annual |> inner_join(link_table, by = c("gvkey", "datadate", "fyr")) |> inner_join(ni_qtrly, by = c("gvkey", "fyr", "datadateq")) plot_data <- ni_merged |> mutate(decade = str_c(as.integer(floor(fyearq / 10) * 10), "s")) |> filter(!is.na(fqtr), fyearq < 2020) |> group_by(decade, fqtr) |> summarize(prop_ni = sum(niq, na.rm = TRUE)/ sum(ni, na.rm = TRUE), prop_sale = sum(saleq, na.rm = TRUE)/ sum(sale, na.rm = TRUE), .groups = "drop") |> mutate(fqtr = factor(fqtr)) |> pivot_longer(cols = c(prop_ni, prop_sale), names_to = "metric", values_to = "value") ``` ### Exercises 1. Pick a couple of `gvkey` values from `fundq_probs`. Is it possible to construct a "clean" sequence of consecutive quarterly earnings announcements for each of these firms? (Here "clean" means that, at the very least, each quarter shows up just once in the series.) What challenges does one face in this task? 2. The code below produces @fig-ni-annual. From @fig-ni-annual, it seems that Q2 has been the most profitable on average over the last three decades, while in all decades, Q4 has seen the most sales. Can you speculate as to what might explain these facts? ```{r} #| echo: true #| fig-pos: H #| label: fig-ni-annual #| fig-cap: Sales and net income by quarter over decades #| fig-alt: "Plot shows the average proportion of sales and net income by quarter for each decade from the 1960s to the 2010s. Net income is generally highest in the second quarter. Sales are generally highest in the fourth quarter and next highest in the second quarter." plot_data |> ggplot(aes(x = fqtr, y = value, fill = fqtr)) + geom_bar(stat = "identity") + facet_grid(metric ~ decade) + theme(legend.position = "none") ``` 3. Create another plot using data in `ni_merged` that you think might be interesting. (Feel free to add variables to `ni_annual` or `ni_qtrly` before merging.) ## Time-series properties of earnings ```{r} #| label: forecast-data #| include: false #| cache: true #| cache.vars: [n_qtrs, n_firms, focus_years, reg_data, reg_data_fos, ccm_link] n_qtrs <- 20 n_firms <- 70 focus_years <- c(1974L, 2019L) # See Table 1 of Foster (1977) for SICs sic2s <- as.character(c(29, 49, 28, 35, 32, 33, 37, 20, 26, 10, 36, 59)) companies <- company |> mutate(sic2 = str_sub(sic, 1L, 2L)) |> filter(sic2 %in% sic2s) |> select(gvkey, sic2) fundq_local <- fundq_mod |> semi_join(companies, by = "gvkey") |> filter(saleq > 0 & !is.na(saleq)) |> select(gvkey, datadate, fyr, fqtr, fyearq, rdq, niq, saleq, ibq) |> rename(datadateq = datadate) |> collect() firm_years <- funda_mod |> select(gvkey, datadate) |> collect() merged_data <- firm_years |> inner_join(link_table, by = c("gvkey", "datadate")) |> inner_join(fundq_local, by = c("gvkey", "datadateq", "fyearq", "fyr")) qtr_num <- merged_data |> group_by(gvkey, datadate) |> count(name = "num_quarters") |> ungroup() regular_fyears <- firm_years |> inner_join(qtr_num, by = c("gvkey", "datadate")) |> group_by(gvkey) |> arrange(gvkey, datadate) |> mutate(fyear_length = datadate - lag(datadate)) |> ungroup() |> mutate(regular_year = num_quarters == 4 & (is.na(fyear_length) | fyear_length %in% c(365, 366))) |> filter(regular_year) |> select(gvkey, datadate) reg_data <- merged_data |> semi_join(companies, copy = TRUE, by = "gvkey") |> semi_join(regular_fyears, by = c("gvkey", "datadate")) |> select(gvkey, datadateq, fyearq, rdq, niq, saleq) |> group_by(gvkey) |> arrange(datadateq) |> mutate(sale_lag_1 = lag(saleq, 1L), sale_lag_4 = lag(saleq, 4L), sale_lag_5 = lag(saleq, 5L), sale_diff = saleq - sale_lag_1, sale_seas_diff = saleq - sale_lag_4, lag_sale_seas_diff = lag(sale_seas_diff, 1L), ni_lag_1 = lag(niq, 1L), ni_lag_4 = lag(niq, 4L), ni_lag_5 = lag(niq, 5L), ni_diff = niq - ni_lag_1, ni_seas_diff = niq - ni_lag_4, lag_ni_seas_diff = lag(ni_seas_diff, 1L)) |> ungroup() reg_data_fos <- merged_data |> semi_join(regular_fyears, by = c("gvkey", "datadate")) |> select(gvkey, datadateq, fyearq, rdq, ibq) |> group_by(gvkey) |> arrange(datadateq) |> mutate(ib_lag_4 = lag(ibq, 4L), ib_seas_diff = ibq - ib_lag_4, lag_ib_seas_diff = lag(ib_seas_diff, 1L), qtr = as.numeric(quarter(datadateq, with_year = TRUE))) |> ungroup() ccm_link <- ccmxpf_lnkhist |> filter(linktype %in% c("LC", "LU", "LS"), linkprim %in% c("C", "P")) |> mutate(linkenddt = coalesce(linkenddt, max(linkenddt, na.rm = TRUE))) |> rename(permno = lpermno) |> collect() ``` ```{r} #| label: fit-model #| include: false fit_model <- function(gvkey, datadateq) { firm_data <- reg_data |> filter(gvkey == !!gvkey) train_data <- firm_data |> filter(datadateq < !!datadateq) |> top_n(n_qtrs, datadateq) if (nrow(train_data) < n_qtrs) return(NULL) test_data <- firm_data |> filter(datadateq == !!datadateq) # Estimate models 2 & 4 model_24 <- train_data |> group_by(gvkey) |> summarize(sale_diff = mean(sale_diff, na.rm = TRUE), ni_diff = mean(ni_diff, na.rm = TRUE), sale_seas_diff = mean(sale_seas_diff, na.rm = TRUE), ni_seas_diff = mean(ni_seas_diff, na.rm = TRUE)) # Fit models 1, 2, 3 & 4 df_model_1234 <- test_data |> # We drop these variables because we will replace them with # their means from model_24 select(-sale_diff, -ni_diff, -sale_seas_diff, -ni_seas_diff) |> inner_join(model_24, by = "gvkey") |> mutate(ni_m1 = ni_lag_4, sale_m1 = sale_lag_4, ni_m2 = ni_lag_4 + ni_seas_diff, sale_m2 = sale_lag_4 + sale_seas_diff, ni_m3 = ni_lag_1, sale_m3 = sale_lag_1, ni_m4 = ni_lag_1 + ni_diff, sale_m4 = sale_lag_1 + sale_diff) # Fit model 5 sale_fm5 <- tryCatch(lm(sale_seas_diff ~ lag_sale_seas_diff, data = train_data, model = FALSE), error = function(e) NULL) ni_fm5 <- tryCatch(lm(ni_seas_diff ~ lag_ni_seas_diff, data = train_data, model = FALSE), error = function(e) NULL) # Fit model 6 sale_fm6 <- tryCatch(lm(saleq ~ sale_lag_1 + sale_lag_4 + sale_lag_5, data = train_data, model = FALSE), error = function(e) NULL) ni_fm6 <- tryCatch(lm(niq ~ ni_lag_1 + ni_lag_4 + ni_lag_5, data = train_data, model = FALSE), error = function(e) NULL) if (!is.null(sale_fm5) & !is.null(ni_fm5)) { results <- df_model_1234 |> mutate(ni_m5 = ni_lag_4 + predict(ni_fm5, newdata = test_data)) |> mutate(sale_m5 = sale_lag_4 + predict(sale_fm5, newdata = test_data)) |> mutate(ni_m6 = predict(ni_fm6, newdata = test_data)) |> mutate(sale_m6 = predict(sale_fm6, newdata = test_data))|> select(gvkey, datadateq, fyearq, niq, saleq, matches("(ni|sale)_m[0-9]")) |> pivot_longer(cols = ni_m1:sale_m6, names_to = "item", values_to = "predicted") |> mutate(actual = if_else(str_detect(item, "^ni"), niq, saleq), abe = abs(predicted - actual) / predicted, se = abe^2 * sign(predicted)) |> separate(item, into = c("item", "model"), sep = "_m") |> select(-niq, -saleq) results } } ``` ```{r} #| include: false top_firms <- reg_data |> filter(fyearq %in% focus_years) |> group_by(gvkey, fyearq) |> summarize(total_sales = sum(saleq), .groups = "drop") |> group_by(fyearq) |> arrange(desc(total_sales)) |> mutate(rank = row_number()) |> filter(rank <= n_firms) test_years <- reg_data |> semi_join(top_firms, by = c("gvkey", "fyearq")) |> select(gvkey, datadateq) ``` ```{r} #| label: forecast-results #| include: false #| warning: false #| cache: true #| cache.vars: [forecast_results] #| dependson: forecast-data #| cache.extra: !expr digest::digest(list(formals(fit_model), body(fit_model), test_years)) forecast_results <- pmap(test_years, fit_model) |> list_rbind() ``` ```{r} #| include: false fix_outliers <- function(x) { if_else(x < 0 | x > 1, 1, x) } ``` ```{r} #| include: false results_summ <- forecast_results |> mutate(abe = fix_outliers(abe), se = fix_outliers(se)) |> group_by(fyearq, item, model) |> summarize(mabe = mean(abe, na.rm=TRUE), mse = mean(se, na.rm=TRUE), .groups = "drop") ``` ```{r} #| echo: false #| fig-pos: H #| label: fig-table-3 #| fig-cap: Plot of results like Table 3 of Foster (1977) #| fig-alt: "Plots of results analogous to those in Table 3 of Foster (1977). The relative performance of the models differs across periods. Error rates are generally much higher for net income models than for sales models." results_summ |> pivot_longer(cols = mabe:mse, names_to = "metric", values_to = "value") |> ggplot(aes(x = model, y = value, fill = model)) + geom_bar(stat = "identity") + facet_grid(fyearq ~ item + metric) + theme(legend.position = "none") ``` ### Exercises 1. How do the denominators in our calculations of `abe` and `se` differ from those in the analogous calculations in @Foster:1977wy? Is one approach to the denominator more correct than the other? 2. What does `fix_outliers()` do? Does @Foster:1977wy do anything to address outliers? If so, how does the approach in @Foster:1977wy compare to that of `fix_outliers()`? Do you agree with the approach taken in `fix_outliers()`? What would you do differently? 3. How do the results in @fig-table-3 compare with those in @Foster:1977wy? 4. What do you make of the significantly "worse" performance of models predicting `ni` than those predicting `sale`? Does this imply that `ni` is simply more difficult to forecast? Can you suggest an alternative approach to measuring performance that might place these models on a more "level playing field"? ## Post-earnings announcement drift The supplied baseline uses 1974--1986. The full forecasting and return extraction can take several minutes on a first run. ```{r} #| label: pead-functions #| include: false pead_firms <- split(reg_data_fos, reg_data_fos$gvkey) fit_model_fos <- function(gvkey, quarter) { n_qtrs <- 24 min_qtrs_fos <- 16 min_qtrs <- 10 firm_data <- pead_firms[[gvkey]] train_data <- firm_data |> filter(qtr < !!quarter) |> top_n(n_qtrs, datadateq) test_data <- firm_data |> filter(qtr == !!quarter) if (nrow(train_data) < min_qtrs) return(NULL) if (nrow(train_data) >= min_qtrs_fos) { # Fit model 5 ib_fm <- tryCatch(lm(ib_seas_diff ~ lag_ib_seas_diff, data = train_data, na.action = na.exclude, model = FALSE), error = function(e) NULL) } else { ib_fm <- NULL } if (!is.null(ib_fm)) { train_results <- train_data |> mutate(fib = ib_lag_4 + predict(ib_fm)) } else { train_results <- train_data |> mutate(fib = ib_lag_4) } denom_m2 <- train_results |> mutate(fe = ibq - fib) |> pull() |> sd() if (is.null(ib_fm)) { results <- test_data |> mutate(fib = ib_lag_4) } else { results <- test_data |> mutate(fib = ib_lag_4 + predict(ib_fm, newdata = test_data)) } results |> mutate(fe1 = (ibq - fib) / abs(ibq), fe2 = (ibq - fib) / denom_m2) } get_results <- function(quarter) { gvkeys <- reg_data_fos |> filter(qtr == quarter) |> distinct(gvkey) |> pull() map2(gvkeys, quarter, fit_model_fos) |> list_rbind() } get_deciles <- function(x) { breaks <- quantile(x, probs = seq(0, 1, .1), na.rm = TRUE) breaks[1] <- -Inf breaks[length(breaks)] <- Inf list(breaks) } ``` ```{r} #| label: pead-analysis #| include: false run_pead <- function(start_year, end_year) { quarters <- reg_data_fos |> filter(qtr >= start_year, qtr < end_year + 1) |> distinct(qtr) |> arrange(qtr) |> pull() # All models are deterministic. Limit workers to avoid large memory copies. plan(multisession, workers = 2) results <- future_map(quarters, get_results) |> list_rbind() plan(sequential) decile_cuts <- results |> group_by(qtr) |> summarize(fe1_deciles = get_deciles(fe1), fe2_deciles = get_deciles(fe2), .groups = "drop") |> arrange(qtr) |> mutate(fe1_deciles_lag = lag(fe1_deciles), fe2_deciles_lag = lag(fe2_deciles)) results_deciles <- results |> inner_join(decile_cuts, by = "qtr") |> # The first quarter has no prior-quarter breakpoints. filter(lengths(fe1_deciles_lag) > 1, lengths(fe2_deciles_lag) > 1) |> rowwise() |> mutate(fe1_decile = cut(fe1, fe1_deciles_lag, labels = FALSE, include.lowest = TRUE), fe2_decile = cut(fe2, fe2_deciles_lag, labels = FALSE, include.lowest = TRUE)) |> filter(!is.na(fe1_decile) | !is.na(fe2_decile)) |> ungroup() |> select(-matches("^fe[12]_deciles")) event_links <- results_deciles |> filter(!is.na(rdq)) |> distinct(gvkey, rdq) |> inner_join(ccm_link, join_by(gvkey, between(rdq, linkdt, linkenddt))) |> distinct(gvkey, rdq, permno) rets <- event_links |> distinct(permno, rdq) |> get_event_rets(db, event_date = "rdq", win_start = -60, win_end = 60) # Return only the summaries used below, rather than a large nested return table. plot_data <- results_deciles |> filter(!is.na(rdq)) |> group_by(gvkey, rdq) |> filter(datadateq == max(datadateq)) |> ungroup() |> mutate(decile = fe2_decile) |> inner_join(event_links, by = c("gvkey", "rdq")) |> inner_join(rets, by = c("rdq", "permno")) |> group_by(decile, relative_td) |> summarize(ar = mean(ret - decret, na.rm = TRUE), .groups = "drop") list(plot_data = plot_data, forecasts = nrow(results), classified = nrow(results_deciles), events = nrow(event_links), return_rows = nrow(rets), duplicate_returns = nrow(rets) - nrow(distinct(rets, permno, rdq, date)), first_quarter = min(quarters)) } plot_pead <- function(data, start_day = 0, end_day = 60) { data |> filter(between(relative_td, start_day, end_day), !is.na(decile)) |> mutate(decile = factor(decile)) |> group_by(decile) |> arrange(relative_td) |> mutate(ar = if_else(relative_td == min(relative_td), 0, ar), car = cumsum(ar), label = if_else(relative_td == max(relative_td), as.character(decile), NA)) |> ggplot(aes(relative_td, car, group = decile, color = decile)) + geom_line() + geom_label(aes(label = label), na.rm = TRUE) + theme(legend.position = "none") } ``` ```{r} #| label: pead-baseline #| include: false #| warning: false #| cache: true #| cache.vars: [pead_base] #| dependson: forecast-data #| cache.extra: !expr digest::digest(lapply(list(fit_model_fos, get_results, get_deciles, run_pead), function(f) list(formals(f), body(f)))) pead_base <- run_pead(1974, 1986) ``` ```{r} #| label: fig-pre #| fig-cap: Pre-announcement returns, 1974--1986 #| echo: false #| fig-pos: H plot_pead(pead_base$plot_data, -60, 0) ``` ```{r} #| label: fig-post #| fig-cap: Post-announcement returns, 1974--1986 #| echo: false #| fig-pos: H plot_pead(pead_base$plot_data) ``` ### Discussion questions 1. Both @Bernard:1989uu and @Ball:1968ub were addressing issues with "conventional wisdom" at their respective times. How had conventional wisdom changed in the years between 1968 and 1989? 2. Evaluate the introduction of @Bernard:1989uu. How clear is the research question to you from reading this? How would you compare this introduction with other papers you have read in this course or elsewhere? 3. How persuasive do you find @Bernard:1989uu to be? What evidence and arguments do you find persuasive or not persuasive? (Answering this question requires you to read the paper fairly closely.) 4. The analysis underlying @fig-post considers the 13-year period from 1974 to 1986. What changes would you need to make to the code to run the analysis for the 13-year period from 2007 to 2019? (If you choose to make this change and run the code, what do you notice about the profile of returns in the post-announcement period? Does it seem necessary to make an additional tweak to the code to address this?) 5. Considering a single stock, what trading strategy is implicit in calculating `ar` as `ret - decret`? 6. In calculating mean returns by `decile` and `relative_td` (i.e., first using `group_by(decile, relative_td)` and then calculating `ar` by aggregating `mean(ret - decret, na.rm = TRUE))`, are we making assumptions about the trading strategy? What issues are created by this trading strategy? Can you suggest an alternative trading strategy? What changes to the code would be needed to implement this alternative? 7. Is it appropriate to *add* returns to get cumulative abnormal returns as is done in `car = cumsum(ar)`? What would be one alternative approach? ```{r} #| include: false #| cache: false dbDisconnect(db) ```