# Seance 0: Identification strategy and independence
# Small artificial training data. No external data files required.
# ggplot2 is optional and used only for visualisation.

dat <- data.frame(
  worker = 1:12,
  D = c(1, 0, 1, 0, 1, 0, 1, 0, 1, 0, 1, 0),
  X = c(42, 41, 57, 56, 50, 49, 63, 62, 46, 47, 54, 55),
  Y = c(1, 0, 1, 0, 1, 1, 1, 0, 0, 0, 1, 1)
)

head(dat)

mean_X_treated <- mean(dat$X[dat$D == 1])
mean_X_control <- mean(dat$X[dat$D == 0])

mean_Y_treated <- mean(dat$Y[dat$D == 1])
mean_Y_control <- mean(dat$Y[dat$D == 0])

naive_difference <- mean_Y_treated - mean_Y_control

model <- lm(Y ~ D, data = dat)

cat("\nSeance 0: Identification strategy and independence\n")
cat("----------------------------------------------------\n")
cat("Mean X, D = 1: ", round(mean_X_treated, 2), "\n", sep = "")
cat("Mean X, D = 0: ", round(mean_X_control, 2), "\n", sep = "")
cat("Mean Y, D = 1: ", round(mean_Y_treated, 2), "\n", sep = "")
cat("Mean Y, D = 0: ", round(mean_Y_control, 2), "\n", sep = "")
cat("Naive treated-control difference: ",
    round(naive_difference, 2), "\n", sep = "")

cat("\nRegression output:\n")
print(summary(model))

cat("\nDataviz process:\n")
cat("Plot 1 compares baseline X by treatment status: balance is reassuring, but it is not proof of independence.\n")
cat("Plot 2 compares Y by treatment status: this is the observed contrast whose causal meaning depends on the identification strategy.\n")

# ggplot2 is used only for visualisation. The numerical exercise above works without it.
if (requireNamespace("ggplot2", quietly = TRUE)) {
  library(ggplot2)

  # Give the binary treatment a readable label for the x-axis and legend.
  dat$D_label <- ifelse(dat$D == 1, "D = 1", "D = 0")

  # Baseline balance plot: do treated and control workers look similar in X?
  plot_balance_x <- ggplot(dat, aes(x = D_label, y = X, colour = D_label)) +
    geom_jitter(width = 0.08, height = 0, size = 2, alpha = 0.75) +
    stat_summary(fun = mean, geom = "point", size = 4, colour = "black") +
    labs(
      title = "Baseline comparison by treatment status",
      x = "Treatment status",
      y = "Baseline covariate X"
    ) +
    theme_minimal() +
    theme(legend.position = "none")

  # Outcome plot: this visualises the same treated-control comparison as lm(Y ~ D).
  plot_outcome_y <- ggplot(dat, aes(x = D_label, y = Y, colour = D_label)) +
    geom_jitter(width = 0.08, height = 0.04, size = 2, alpha = 0.75) +
    stat_summary(fun = mean, geom = "point", size = 4, colour = "black") +
    labs(
      title = "Observed outcome comparison",
      x = "Treatment status",
      y = "Outcome Y"
    ) +
    theme_minimal() +
    theme(legend.position = "none")

  if (interactive()) {
    print(plot_balance_x)
    print(plot_outcome_y)
  } else {
    cat("Plots created as plot_balance_x and plot_outcome_y. Run the script in RStudio to display them.\n")
  }
} else {
  cat("ggplot2 is not installed. To draw the figures, run install.packages(\"ggplot2\") once, then rerun this script.\n")
}

cat("\nInterpretation prompt:\n")
cat("The coefficient on D is a treated-control comparison. It is causal only if the identification strategy makes D independent of Y(1) and Y(0).\n")