tests/testthat/test-cic.R

In sccic: Synthetic Control Changes-in-Changes Estimator

test_that("cic returns correct structure", {
  set.seed(1)
  y_00 <- rnorm(100, 0, 1)
  y_01 <- rnorm(100, 0.5, 1)
  y_10 <- rnorm(50, 1, 1)
  y_11 <- rnorm(50, 2, 1)
  result <- cic(y_00, y_01, y_10, y_11)
  expect_s3_class(result, "cic")
  expect_true(!is.na(result$tau))
  expect_true(!is.na(result$se))
  expect_true(!is.na(result$pval))
  expect_equal(result$N, 300)
})

test_that("cic matches qte::CiC on workers comp data", {
  skip_if_not_installed("wooldridge")
  data("injury", package = "wooldridge")
  result <- cic(
    y_00 = injury$ldurat[injury$highearn == 0 & injury$afchnge == 0],
    y_01 = injury$ldurat[injury$highearn == 0 & injury$afchnge == 1],
    y_10 = injury$ldurat[injury$highearn == 1 & injury$afchnge == 0],
    y_11 = injury$ldurat[injury$highearn == 1 & injury$afchnge == 1]
  )
  # Known result from validation: tau = 0.0687
  expect_equal(round(result$tau, 2), 0.07, tolerance = 0.01)
})

test_that("cic with no treatment returns near-zero effect", {
  set.seed(42)
  y_00 <- rnorm(200, 0, 1)
  y_01 <- rnorm(200, 1, 1)  # time trend
  y_10 <- rnorm(200, 0, 1)  # same as control pre
  y_11 <- rnorm(200, 1, 1)  # same time trend, no treatment
  result <- cic(y_00, y_01, y_10, y_11)
  # With no treatment, tau should be near zero
  expect_lt(abs(result$tau), 0.5)
})

test_that("discrete CIC reproduces Athey-Imbens (2006) workers comp estimate", {
  skip_if_not_installed("wooldridge")
  data("injury", package = "wooldridge")
  result <- cic(
    y_00 = injury$ldurat[injury$highearn == 0 & injury$afchnge == 0],
    y_01 = injury$ldurat[injury$highearn == 0 & injury$afchnge == 1],
    y_10 = injury$ldurat[injury$highearn == 1 & injury$afchnge == 0],
    y_11 = injury$ldurat[injury$highearn == 1 & injury$afchnge == 1],
    discrete = TRUE, se = FALSE, boot = FALSE
  )
  expect_true(isTRUE(result$discrete))
  # A&I (2006) report 0.18 using discrete CIC on a subsample (N=5,624);
  # our full-sample (N=7,150) estimate is 0.184.
  expect_equal(round(result$tau, 2), 0.18, tolerance = 0.01)
  # Analytic SE should be NA for discrete = TRUE
  expect_true(is.na(result$se))
})

test_that("discrete CIC ecdf helpers are correct", {
  ec <- sccic:::make_ecdf(c(1, 2, 3, 4))
  # integral of Q(u) from 0 to 1: 0.25*1 + 0.25*2 + 0.25*3 + 0.25*4 = 2.5
  expect_equal(sccic:::integral_quantile(ec, 0, 1), 2.5)
  # ecdf_eval_left: Pr(Y < 2) = 1/4 for {1,2,3,4}
  expect_equal(sccic:::ecdf_eval_left(ec, 2), 0.25)
  # ecdf_eval_left at minimum: Pr(Y < 1) = 0
  expect_equal(sccic:::ecdf_eval_left(ec, 1), 0)
})

test_that("bootstrap SE is computed when requested", {
  set.seed(1)
  y_00 <- rnorm(50); y_01 <- rnorm(50)
  y_10 <- rnorm(30); y_11 <- rnorm(30)
  result <- cic(y_00, y_01, y_10, y_11, boot = TRUE, boot_iters = 100, seed = 42)
  expect_true(!is.na(result$boot_se))
  expect_gt(result$boot_se, 0)
})

# Helper: load and prepare Basque Country data
.basque_data <- function() {
  data("basque", package = "Synth")
  gdp <- reshape(basque[, c("regionno", "year", "gdpcap")],
                 idvar = "year", timevar = "regionno", direction = "wide")
  y_treated <- gdp[, "gdpcap.17"]
  donors <- as.matrix(gdp[, grep("gdpcap\\.", names(gdp))])
  donors <- donors[, !colnames(donors) %in% c("gdpcap.17", "gdpcap.1")]
  valid <- complete.cases(y_treated, donors)
  list(y_treated = y_treated[valid], donors = donors[valid, ])
}

test_that("sc_cic returns correct structure and matches Basque manuscript values", {
  skip_if_not_installed("Synth")
  d <- .basque_data()
  result <- sc_cic(d$y_treated, d$donors,
                   treatment_period = 16, boot = FALSE, seed = 42)
  expect_s3_class(result, "sc_cic")
  expect_s3_class(result, "cic")
  # Manuscript values: tau = -0.7943, pre_fit_rmse = 0.0493, 4 donors selected
  expect_equal(round(result$tau, 2),          -0.79, tolerance = 0.01)
  expect_equal(round(result$pre_fit_rmse, 3),  0.049, tolerance = 0.001)
  expect_equal(result$n_donors_eff, 4L)
  expect_gt(result$ks_pval, 0.10)  # distributional alignment OK
})

test_that("sensitivity_alpha returns expected structure and Basque alpha=1 value", {
  skip_if_not_installed("Synth")
  d <- .basque_data()
  sens <- sensitivity_alpha(d$y_treated, d$donors, treatment_period = 16)
  expect_s3_class(sens, "data.frame")
  expect_true(all(c("alpha", "tau_cic", "n_donors", "pre_rmse") %in% names(sens)))
  expect_equal(nrow(sens), 6L)
  # alpha=1 row should match the main sc_cic result
  expect_equal(round(sens$tau_cic[sens$alpha == 1], 2), -0.79, tolerance = 0.01)
})

test_that("loo_donors returns expected structure", {
  skip_if_not_installed("Synth")
  d <- .basque_data()
  # Use first 5 donors only to keep test fast
  loo <- loo_donors(d$y_treated, d$donors[, 1:5], treatment_period = 16)
  expect_s3_class(loo, "data.frame")
  expect_true(all(c("excluded", "tau_cic", "tau_did", "pre_rmse") %in% names(loo)))
  expect_equal(nrow(loo), 5L)
  expect_true(all(!is.na(loo$tau_cic)))
})

test_that("check_support returns TRUE for Basque (support condition satisfied)", {
  skip_if_not_installed("Synth")
  d <- .basque_data()
  result <- sc_cic(d$y_treated, d$donors,
                   treatment_period = 16, boot = FALSE, seed = 42)
  expect_true(check_support(result))
})