This script reproduces a stylized version of the treatment effects from the paper “Early Impacts of M365 Copilot” using simulated data. It aligns closely with the original methodology—implementing fixed effects for worker, time, and firm-by-month, and applying Newey-West standard errors. The code calibrates intent-to-treat (ITT) effects to match…

replicate_m365_copilot_event_study.R

library(tidyverse)
library(fixest)
library(ggtext)

set.seed(123)
n_workers <- 6000
n_firms <- 56
rel_months <- -1:6 # –1 = pre-rollout month, 1…6 = months after rollout

# –– 2) Calibrate “true” effects to paper’s ITT estimates -------------------

# Month 1  ITT = –12.30 min
# Month 6  ITT = –31.16 min
effects <- tibble(
  rel_month = rel_months,
  effect = case_when(
    rel_month < 1 ~ 0,
    TRUE ~ -12.30 + (rel_month - 1) * ((-31.16) - (-12.30)) / 5
  )
)

# –– 3) Build the panel with firm, worker, and “calendar” month --------------

df_event <- expand_grid(
  worker_id = 1:n_workers,
  rel_month = rel_months
) %>%
  mutate(
    firm_id     = sample(1:n_firms, n_workers, replace = TRUE)[worker_id],
    treat_group = sample(0:1, n_workers, replace = TRUE)[worker_id]
  ) %>%
  left_join(effects, by = "rel_month") %>%
  mutate(
    # apply effect only for treated workers
    treatment_effect = if_else(treat_group == 1, effect, 0),
    base_time        = 170.16, # control‐group mean from paper
    epsilon          = rnorm(n(), 0, 10),
    y                = base_time + treatment_effect + epsilon,
    year_month       = rel_month, # stand‐in for calendar time FE
    rel_month_factor = factor(rel_month)
  )

# –– 4) Estimate event‐study DiD with firm×month & Newey–West SEs ------------

model_event <- feols(
  y ~ i(rel_month_factor, treat_group, ref = "-1") |
    worker_id + year_month + firm_id:year_month,
  data = df_event,
  panel.id = ~ worker_id + year_month,
  vcov = NW(lag = 2) # 2-month lag for Newey–West
)

# –– 5) Tidy up coefficients -------------------------------------------------

event_df <- broom::tidy(model_event, conf.int = TRUE) %>%
  filter(str_detect(term, "^rel_month_factor")) %>%
  mutate(
    month    = as.integer(str_extract(term, "-?\\d+")),
    estimate = estimate,
    lower    = conf.low,
    upper    = conf.high
  )

# –– 6) Plot ------------------------------------------------------------------

ggplot(event_df, aes(x = month, y = estimate)) +
  geom_line(color = "steelblue", size = 1) +
  geom_point(color = "steelblue", size = 2) +
  geom_ribbon(aes(ymin = lower, ymax = upper), fill = "steelblue", alpha = 0.2) +
  geom_vline(xintercept = 0, linetype = "dashed", color = "gray40") +
  geom_hline(yintercept = 0, linetype = "dotted") +
  geom_richtext(
    data = subset(event_df, month %in% c(1, 6)),
    aes(label = paste0(round(estimate, 1), " min")),
    fill = "white",
    label.color = NA,
    size = 5,
    hjust = c(0, 1),
    nudge_y = c(2, -2),
    family = "sans"
  ) +
  labs(
    title = "Copilot Cuts Corporate Email Time by Over Half an Hour in Six Months",
    subtitle = "A randomized trial of more than 6,000 workers across 56 firms shows that simply assigning a Microsoft<br>Copilot license reduced weekly email reading by 12.3 minutes in Month 1 and grew to 31.2 minutes by<br>Month 6.",
    x = "Months Since Copilot Rollout (0 = Pre)",
    y = "Estimated Treatment Effect (minutes)",
    caption = paste(
      "Assumptions not made in the original paper:",
      "1) Linearly interpolate the ITT effect from –12.3 min at Month 1 to –31.16 min at Month 6;",
      "2) Assume treatment compliance jumps to 100% starting Month 1 (no gradual uptake);",
      "3) Effects are homogeneous across all workers and firms;",
      "4) Calendar‐time FE is set equal to rel_month for simplicity;",
      "5) Use a 2-month Newey–West lag (paper used 9-week lag).",
      sep = "\n"
    )
  ) +
  theme_minimal(base_family = "sans", base_size = 14) +
  theme(
    plot.title.position = "plot",
    plot.title          = element_text(hjust = 0, size = 20, face = "bold"),
    plot.subtitle       = element_markdown(hjust = 0, size = 15, margin = margin(b = 10)),
    plot.caption.position = "plot",
    plot.caption        = element_text(hjust = 0, size = 10, margin = margin(t = 10)),
    legend.position     = "top",
    legend.title        = element_blank(),
    plot.margin         = margin(t = 25, r = 25, b = 15, l = 25),
    axis.title.x        = element_text(size = 14, face = "bold"),
    axis.title.y        = element_text(size = 14, face = "bold"),
    axis.text           = element_text(size = 12),
    panel.grid.minor    = element_blank(),
    panel.grid.major    = element_line(color = "gray90")
  ) +
  guides(
    color = guide_colorbar(
      title.position = "top",
      title.hjust    = 0.5,
      barwidth       = unit(15, "lines"),
      barheight      = unit(0.5, "lines")
    )
  ) +
  coord_cartesian(expand = FALSE, clip = "off")