#!/usr/bin/env Rscript
# 02_event_study_did.R
# Estimate group-time ATT and event-study aggregations for the matched unique origin-device-month panel.
# Requires: install.packages(c('did','data.table','ggplot2'))

suppressPackageStartupMessages({
  library(data.table)
  library(did)
  library(ggplot2)
})

args <- commandArgs(trailingOnly = TRUE)
if (length(args) < 4) {
  stop('Usage: Rscript 02_event_study_did.R matched_panel.csv outcome_col group_col output_prefix')
}

panel_path <- args[[1]]
outcome <- args[[2]]
group_col <- args[[3]]
out_prefix <- args[[4]]

dt <- fread(panel_path)
if (!group_col %in% names(dt)) stop(paste('group_col not found:', group_col))
if (!outcome %in% names(dt)) stop(paste('outcome not found:', outcome))

# Unique unit and time indexes. Device is part of the unit when present.
if ('device' %in% names(dt)) {
  dt[, id_unit := paste(origin, device, sep = '||')]
} else {
  dt[, id_unit := origin]
}
dt[, id := as.integer(factor(id_unit))]
dt[, time_index := as.integer(substr(as.character(yyyymm),1,4))*12 + as.integer(substr(as.character(yyyymm),5,6))]

g_to_idx <- function(x) {
  fifelse(is.na(x) | x == 0, 0L,
          as.integer(substr(as.character(x),1,4))*12 + as.integer(substr(as.character(x),5,6)))
}
dt[, gname := g_to_idx(get(group_col))]
dt <- unique(dt, by = c('id', 'time_index'))
dt <- dt[!is.na(get(outcome))]

att <- att_gt(
  yname = outcome,
  tname = 'time_index',
  idname = 'id',
  gname = 'gname',
  data = dt,
  control_group = 'notyettreated',
  clustervars = 'id',
  allow_unbalanced_panel = TRUE,
  bstrap = TRUE,
  biters = 999
)

dyn <- aggte(att, type = 'dynamic', min_e = -6, max_e = 12)
res <- data.table(event_time = dyn$egt, att = dyn$att.egt, se = dyn$se.egt)
res[, ci_low := att - 1.96*se]
res[, ci_high := att + 1.96*se]
fwrite(res, paste0(out_prefix, '_dynamic_att.csv'))

png(paste0(out_prefix, '_dynamic_att.png'), width = 1200, height = 800, res = 150)
print(ggplot(res, aes(event_time, att)) +
  geom_hline(yintercept = 0, linetype = 'dashed') +
  geom_vline(xintercept = -0.5, linetype = 'dotted') +
  geom_ribbon(aes(ymin = ci_low, ymax = ci_high), alpha = 0.2) +
  geom_line() + geom_point() +
  labs(x = 'Event time in months', y = paste('ATT on', outcome), title = paste('Dynamic ATT:', outcome)) +
  theme_minimal())
dev.off()

print(summary(att))
print(summary(dyn))