Nothing
context("Testing Darfur Example with fixest model")
# runs regression model
data("darfur")
# fixest using village as fixed effects
model <- fixest::feols(peacefactor ~ directlyharmed + age + farmer_dar + herder_dar +
pastvoted + hhsize_darfur + female|village, data = darfur)
# vanilla lm
model_lm <- lm(peacefactor ~ directlyharmed + age + farmer_dar + herder_dar +
pastvoted + hhsize_darfur + female + village, data = darfur)
# Model fixest using heteroskedastic se's
# -note adjusted estimates, partial R2, and the RV for point estimates should *always* use the "iid" SE.
# -this is because the SE is just a shortcut for computing the residual standard deviations.
# -Robust SEs should be used for inference (confidence intervals, RV_alpha),
# -but to do that properly we need to use the delta method.
# -this can be implemented later, for now sensemakr will simply ignore the "hetero" flag,
# -and provide a message to the user.
model_het <- fixest::feols(peacefactor ~ directlyharmed + age+ farmer_dar + herder_dar +
pastvoted + hhsize_darfur + female | village, data = darfur,
vcov = "hetero")
test_that(desc = "Identical Results lm and feols with village fixed effects",
{
darfur_out <- sensemakr(model, treatment = "directlyharmed", benchmark_covariates = "female", kd = 1:3)
darfur_out_lm <- sensemakr(model_lm, treatment = "directlyharmed", benchmark_covariates = "female", kd = 1:3)
plot1 <- plot(darfur_out)
plot1_lm <- plot(darfur_out_lm)
expect_equal(plot1, plot1_lm)
plot2 <- ovb_contour_plot(darfur_out)
plot2_lm <- ovb_contour_plot(darfur_out_lm)
expect_equal(plot2, plot2_lm)
plot3 <- ovb_extreme_plot(darfur_out)
plot3_lm <- ovb_extreme_plot(darfur_out_lm)
expect_equal(plot3, plot3_lm)
# info
expect_equal(darfur_out$info$formula,
peacefactor ~ directlyharmed + age + farmer_dar + herder_dar +
pastvoted + hhsize_darfur + female | village,
darfur_out_lm$info$formula,
peacefactor ~ directlyharmed + age + farmer_dar + herder_dar +
pastvoted + hhsize_darfur + female + village)
expect_equal(darfur_out$info$treatment, darfur_out_lm$info$treatment)
expect_equal(darfur_out$info$q, darfur_out_lm$info$q)
expect_equal(darfur_out$info$alpha, darfur_out_lm$info$alpha)
expect_equal(darfur_out$info$reduce, darfur_out_lm$info$reduce)
# sensitivity stats
expect_equal(darfur_out$sensitivity_stats$dof, darfur_out_lm$sensitivity_stats$dof)
expect_equal(darfur_out$sensitivity_stats$r2yd.x, darfur_out_lm$sensitivity_stats$r2yd.x, tolerance = 1e-5)
expect_equivalent(c(darfur_out$sensitivity_stats$rv_q), c(darfur_out_lm$sensitivity_stats$rv_q), tolerance = 1e-5)
expect_equivalent(c(darfur_out$sensitivity_stats$rv_qa), c(darfur_out_lm$sensitivity_stats$rv_qa), tolerance = 1e-5)
# bounds
expect_equivalent(darfur_out$bounds, darfur_out_lm$bounds)
# out1 <- capture.output(darfur_out)
# out1_lm <- capture.output(darfur_out_lm)
# expect_equal(out1, out1_lm)
# out2 <- capture.output(summary(darfur_out))
# out2_lm <- capture.output(summary(darfur_out_lm))
# expect_equal(out2, out2_lm)
out3 <- capture.output(ovb_minimal_reporting(darfur_out))
out3_lm <- capture.output(ovb_minimal_reporting(darfur_out_lm))
expect_equal(out3, out3_lm)
darfur_out2 <- sensemakr(model, treatment = "directlyharmed")
darfur_out2_lm <- sensemakr(model_lm, treatment = "directlyharmed")
plot1 <- plot(darfur_out2)
plot1_lm <- plot(darfur_out2)
expect_equal(plot1, plot1_lm)
plot2 <- ovb_contour_plot(darfur_out2)
plot2_lm <- ovb_contour_plot(darfur_out2_lm)
expect_equal(plot2, plot2_lm)
plot3 <- ovb_extreme_plot(darfur_out2)
plot3_lm <- ovb_extreme_plot(darfur_out2_lm)
expect_equal(plot3, plot3_lm)
out3 <- capture.output(ovb_minimal_reporting(darfur_out2))
out3_lm <- capture.output(ovb_minimal_reporting(darfur_out2_lm))
expect_equal(out3, out3_lm)
})
test_that(desc = "Identifical results with using different vcov (we should ignore vcov for adjusted estimates)" ,
{
darfur_out <- sensemakr(model, treatment = "directlyharmed", benchmark_covariates = "female", kd = 1:3)
darfur_out_het <- sensemakr(model_het, treatment = "directlyharmed", benchmark_covariates = "female", kd = 1:3)
plot1 <- plot(darfur_out)
plot1_lm <- plot(darfur_out_het)
expect_equal(plot1, plot1_lm)
plot2 <- ovb_contour_plot(darfur_out)
plot2_lm <- ovb_contour_plot(darfur_out_het)
expect_equal(plot2, plot2_lm)
plot3 <- ovb_extreme_plot(darfur_out)
plot3_lm <- ovb_extreme_plot(darfur_out_het)
expect_equal(plot3, plot3_lm)
# info
expect_equal(darfur_out$info$formula,
peacefactor ~ directlyharmed + age + farmer_dar + herder_dar +
pastvoted + hhsize_darfur + female | village,
darfur_out_het$info$formula,
peacefactor ~ directlyharmed + age + farmer_dar + herder_dar +
pastvoted + hhsize_darfur + female + village)
expect_equal(darfur_out$info$treatment, darfur_out_het$info$treatment)
expect_equal(darfur_out$info$q, darfur_out_het$info$q)
expect_equal(darfur_out$info$alpha, darfur_out_het$info$alpha)
expect_equal(darfur_out$info$reduce, darfur_out_het$info$reduce)
# sensitivity stats
expect_equal(darfur_out$sensitivity_stats$dof, darfur_out_het$sensitivity_stats$dof)
expect_equal(darfur_out$sensitivity_stats$r2yd.x, darfur_out_het$sensitivity_stats$r2yd.x, tolerance = 1e-5)
expect_equivalent(c(darfur_out$sensitivity_stats$rv_q), c(darfur_out_het$sensitivity_stats$rv_q), tolerance = 1e-5)
expect_equivalent(c(darfur_out$sensitivity_stats$rv_qa), c(darfur_out_het$sensitivity_stats$rv_qa), tolerance = 1e-5)
# bounds
expect_equivalent(darfur_out$bounds, darfur_out_het$bounds)
out1 <- capture.output(darfur_out)
out1_lm <- capture.output(darfur_out_het)
expect_equal(out1, out1_lm)
out2 <- capture.output(summary(darfur_out))
out2_lm <- capture.output(summary(darfur_out_het))
expect_equal(out2, out2_lm)
out3 <- capture.output(ovb_minimal_reporting(darfur_out))
out3_lm <- capture.output(ovb_minimal_reporting(darfur_out_het))
expect_equal(out3, out3_lm)
darfur_out2 <- sensemakr(model, treatment = "directlyharmed")
darfur_out_het2 <- sensemakr(model_het, treatment = "directlyharmed")
plot1 <- plot(darfur_out2)
plot1_lm <- plot(darfur_out2)
expect_equal(plot1, plot1_lm)
plot2 <- ovb_contour_plot(darfur_out2)
plot2_lm <- ovb_contour_plot(darfur_out_het2)
expect_equal(plot2, plot2_lm)
plot3 <- ovb_extreme_plot(darfur_out2)
plot3_lm <- ovb_extreme_plot(darfur_out_het2)
expect_equal(plot3, plot3_lm)
out3 <- capture.output(ovb_minimal_reporting(darfur_out2))
out3_lm <- capture.output(ovb_minimal_reporting(darfur_out_het2))
expect_equal(out3, out3_lm)
})
test_that("testing darfur manual bounds",
{
sense.out <- sensemakr(model, treatment = "directlyharmed", benchmark_covariates = "female", r2dz.x = .1)
bounds.check <- sense.out$bounds
to_check <- bounds.check$adjusted_se[1]
true_check <- adjusted_se(model, treatment = "directlyharmed", r2dz.x = .1, r2yz.dx = .1)
expect_equal(to_check, unname(true_check))
})
test_that(desc = "testing darfur sensemakr with formula",
{
darfur_out <- sensemakr(formula = peacefactor ~ directlyharmed + age + farmer_dar + herder_dar +
pastvoted + hhsize_darfur + female + village, data = darfur, method = "feols", treatment = "directlyharmed", benchmark_covariates = "female", kd = 1:3)
# info
expect_equal(darfur_out$info$formula, peacefactor ~ directlyharmed + age + farmer_dar + herder_dar +
pastvoted + hhsize_darfur + female + village)
expect_equal(darfur_out$info$treatment, "directlyharmed")
expect_equal(darfur_out$info$q, 1)
expect_equal(darfur_out$info$alpha, 0.05)
expect_equal(darfur_out$info$reduce, TRUE)
# sensitivity stats
expect_equal(darfur_out$sensitivity_stats$dof, 783)
expect_equal(darfur_out$sensitivity_stats$r2yd.x, 0.02187, tolerance = 1e-5)
expect_equivalent(c(darfur_out$sensitivity_stats$rv_q), 0.13878, tolerance = 1e-5)
expect_equivalent(c(darfur_out$sensitivity_stats$rv_qa), 0.07626, tolerance = 1e-5)
# bounds
check_bounds <- structure(list(bound_label = c("1x female", "2x female", "3x female"),
r2dz.x = c(0.00916428667504862, 0.0183285733500972, 0.0274928600251459),
r2yz.dx = c(0.12464092303637, 0.249324064199975, 0.374050471038094),
treatment = rep("directlyharmed", 3),
adjusted_estimate = c(0.0752202712144491, 0.0529151723844518, 0.0303960234641548),
adjusted_se = c(0.0218733277437572, 0.0203500620779637, 0.0186700648170924),
adjusted_t = c(3.43890386024675, 2.60024623913809, 1.62806202131271),
adjusted_lower_CI = c(0.032282966, 0.012968035,-0.006253282),
adjusted_upper_CI = c(0.11815758, 0.09286231, 0.06704533)),
.Names = c("bound_label", "r2dz.x", "r2yz.dx", "treatment",
"adjusted_estimate", "adjusted_se", "adjusted_t",
"adjusted_lower_CI", "adjusted_upper_CI"),
row.names = c(NA, -3L), class = c("ovb_bounds", "data.frame"))
expect_equivalent(darfur_out$bounds, check_bounds)
})
test_that(desc = "testing darfur sensemakr manually",
{
model.treat <- fixest::feols(directlyharmed ~ age + farmer_dar + herder_dar +
pastvoted + hhsize_darfur + female | village, data = darfur)
darfur_out <- sensemakr(estimate = 0.09731582,
se = 0.02325654,
dof = 783,
treatment = "directlyharmed",
benchmark_covariates = "female",
r2dxj.x = partial_r2(model.treat, covariates = "female"),
r2yxj.dx = partial_r2(model, covariates = "female"),
kd = 1:3)
plot(darfur_out)
plot(darfur_out, type = "extreme")
plot(darfur_out, sensitivity.of = "t-value")
add_bound_to_contour(0.3,0.3, "test")
# info
expect_equal(darfur_out$info$formula, "Data provided manually")
expect_equal(darfur_out$info$treatment, "directlyharmed")
expect_equal(darfur_out$info$q, 1)
expect_equal(darfur_out$info$alpha, 0.05)
expect_equal(darfur_out$info$reduce, TRUE)
expect_equal(darfur_out$sensitivity_stats$dof, 783)
expect_equal(darfur_out$sensitivity_stats$r2yd.x, 0.02187, tolerance = 1e-5)
expect_equivalent(c(darfur_out$sensitivity_stats$rv_q), 0.13878, tolerance = 1e-5)
expect_equivalent(c(darfur_out$sensitivity_stats$rv_qa), 0.07626, tolerance = 1e-5)
# bounds
check_bounds <-
structure(
list(
bound_label = c("1x female", "2x female", "3x female"),
r2dz.x = c(0.00916428667504862, 0.0183285733500972, 0.0274928600251459),
r2yz.dx = c(0.12464092303637, 0.249324064199975, 0.374050471038094),
adjusted_estimate = c(0.0752202698486415, 0.0529151689180575,
0.0303960178770157),
adjusted_se = c(0.0218733298036818, 0.0203500639944344,
0.0186700665753491),
adjusted_t = c(3.43890347394571, 2.60024582392121,
1.62806156873318),
adjusted_lower_CI = c(0.0322829603180086,
0.0129680276030601, -0.00625329133645187),
adjusted_upper_CI = c(0.118157579379274,
0.092862310233055, 0.0670453270904833)
),
row.names = c(NA, -3L),
class = "data.frame"
)
expect_equivalent(as.data.frame(darfur_out$bounds), as.data.frame(check_bounds), tolerance = 1e-5)
})
test_that(desc = "testing darfur sensitivity stats",
{
# checks RV
## RV q = 1
rv <- robustness_value(model, covariates = "directlyharmed", alpha = 1)
expect_equivalent(c(rv), 0.138776, tolerance = 1e-5)
expect_equivalent(attributes(rv)$q, 1)
expect_equivalent(attributes(rv)$names, "directlyharmed")
## RV q = 1, alpha = 0.05
rv <- robustness_value(model, covariates = "directlyharmed", q = 1, alpha = 0.05)
expect_equivalent(c(rv), 0.07625797, tolerance = 1e-5)
expect_equivalent(attributes(rv)$q, 1)
expect_equivalent(attributes(rv)$alpha, 0.05)
expect_equivalent(attributes(rv)$names, "directlyharmed")
# checks partial R2
r2 <- partial_r2(model, covariates = "directlyharmed")
expect_equivalent(r2, 0.02187309, tolerance = 1e-5)
# checks partial f2
f2 <- partial_f2(model, covariates = "directlyharmed")
expect_equivalent(f2, 0.02236222, tolerance = 1e-5)
# sensitivity stats
sens_stats <- sensitivity_stats(model, treatment = "directlyharmed")
expect_equivalent(sens_stats$treatment, "directlyharmed")
expect_equivalent(sens_stats$estimate, 0.09731582, tolerance = 1e5)
expect_equivalent(sens_stats$se, 0.02325654, tolerance = 1e5)
expect_equivalent(sens_stats$t_statistic, 4.18445, tolerance = 1e5)
expect_equivalent(sens_stats$r2yd.x, 0.02187309, tolerance = 1e5)
expect_equivalent(sens_stats$rv_q , 0.1387764, tolerance = 1e5)
expect_equivalent(sens_stats$rv_qa , 0.07625797, tolerance = 1e5)
expect_equivalent(sens_stats$f2yd.x , 0.02236222, tolerance = 1e5)
expect_equivalent(sens_stats$dof , 783, tolerance = 1e5)
expect_equivalent(group_partial_r2(model, covariates = "directlyharmed"), partial_r2(model, covariates = "directlyharmed"), tolerance = 1e5)
expect_error(group_partial_r2(model))
expect_equal(group_partial_r2(model, covariates = c("directlyharmed", "female")), 0.1350435, tolerance = 1e-5)
})
test_that(desc = "testing darfur adjusted estimates",
{
should_be_zero <- adjusted_estimate(model, treatment = "directlyharmed", r2yz.dx = 1, r2dz.x = partial_r2(model, covariates = "directlyharmed"))
expect_equivalent(should_be_zero, 0)
rv <- robustness_value(model, covariates = "directlyharmed", alpha = 1)
should_be_zero <- adjusted_estimate(model, treatment = "directlyharmed", r2yz.dx = rv, r2dz.x = rv)
expect_equivalent(should_be_zero, 0)
rv <- robustness_value(model, covariates = "directlyharmed", alpha = 0.05)
thr <- qt(0.975, df = 783 - 1)
should_be_1.96 <- adjusted_t(model, treatment = "directlyharmed", r2yz.dx = rv, r2dz.x = rv)
expect_equivalent(c(should_be_1.96), thr)
should_be_estimate <- bias(model, treatment = "directlyharmed", r2yz.dx = 1, r2dz.x = partial_r2(model, covariates = "directlyharmed"))
expect_equivalent(should_be_estimate, coef(model)["directlyharmed"])
rv <- robustness_value(model, covariates = "directlyharmed", alpha = 1)
should_be_estimate <- bias(model, treatment = "directlyharmed", r2yz.dx = rv, r2dz.x = rv)
expect_equivalent(should_be_estimate, coef(model)["directlyharmed"])
rv <- robustness_value(model, covariates = "directlyharmed", q = 0.5, alpha = 1)
should_be_half_estimate <- bias(model, treatment = "directlyharmed", r2yz.dx = rv, r2dz.x = rv)
expect_equivalent(should_be_half_estimate, coef(model)["directlyharmed"]/2)
})
test_that(desc = "testing darfur plots",
{
# testing countour output with internal functions
## point estimate
contour_out <- ovb_contour_plot(model, treatment = "directlyharmed",
benchmark_covariates = "female",kd = 1:3)
add_bound_to_contour(model, treatment = "directlyharmed",
benchmark_covariates = "age", kd = 10)
add_bound_to_contour(model, treatment = "directlyharmed",
benchmark_covariates = "age", kd = 200, ky = 20)
adj_est <- adjusted_estimate(model, treatment = "directlyharmed",
r2dz.x = contour_out$r2dz.x[10],
r2yz.dx = contour_out$r2yz.dx[1])
expect_equivalent(adj_est, contour_out$value[10])
bounds <- ovb_bounds(model,
treatment = "directlyharmed",
benchmark_covariates = "female",
kd = 1:3)
expect_equivalent(contour_out$bounds, bounds[c(2,3,1)])
## t-value
contour_out <- ovb_contour_plot(model, treatment = "directlyharmed",
benchmark_covariates = "female",kd = 1:3,
sensitivity.of = "t-value")
add_bound_to_contour(model, treatment = "directlyharmed", benchmark_covariates = "age",
kd = 200, ky = 10, sensitivity.of = "t-value")
adj_t <- adjusted_t(model, treatment = "directlyharmed",
r2dz.x = contour_out$r2dz.x[10],
r2yz.dx = contour_out$r2yz.dx[1])
expect_equivalent(adj_t, contour_out$value[10])
bounds <- ovb_bounds(model,
treatment = "directlyharmed",
benchmark_covariates = "female",
kd = 1:3)
expect_equivalent(contour_out$bounds, bounds[c(2,3,1)])
# tests bounds numerically
# tests bounds numerically
expect_equivalent(bounds$adjusted_estimate, c(0.0752202712144491, 0.0529151723844518, 0.0303960234641548))
expect_equivalent(bounds$adjusted_t, c(3.43890386024675, 2.60024623913809, 1.62806202131271))
# test extreme scenario plot
extreme_out <- ovb_extreme_plot(model, treatment = "directlyharmed", kd = 1:3)
adj_est <- adjusted_estimate(model, treatment = "directlyharmed",
r2yz.dx = 1,
r2dz.x = extreme_out$scenario_r2yz.dx_1$r2dz.x[5])
expect_equivalent(adj_est, extreme_out$scenario_r2yz.dx_1$adjusted_estimate[5])
})
test_that("testing darfur different q",
{
darfur_out <- sensemakr(model, treatment = "directlyharmed", benchmark_covariates = "female", q = 2, kd = 1:3)
rvq <- darfur_out$sensitivity_stats$rv_q
rvqa <- darfur_out$sensitivity_stats$rv_qa
expect_equivalent(rvq, robustness_value(model, covariates = "directlyharmed", q = 2, alpha = 1))
expect_equivalent(rvqa, robustness_value(model, covariates = "directlyharmed", q = 2,alpha = 0.05))
})
model5 <- fixest::feols(peacefactor ~ directlyharmed + age + farmer_dar + herder_dar +
pastvoted + hhsize_darfur + female + village, data = darfur)
test_that("Darfur group benchmarks",
{
village <- grep(pattern = "village", names(coef(model5)), value = T)
village_lm <- grep(pattern = "village", names(coef(model_lm)), value = T)
sensitivity.lm <- sensemakr(model_lm, treatment = "directlyharmed",
benchmark_covariates = list(village = village_lm),
kd = 0.3)
sensitivity <- sensemakr(model5, treatment = "directlyharmed",
benchmark_covariates = list(village = village),
kd = 0.3)
expect_equal(sensitivity.lm$bounds, sensitivity$bounds)
r2y <- group_partial_r2(model5, covariates = village)
treat.model <- update(model5, directlyharmed ~ .-directlyharmed)
r2d <- group_partial_r2(treat.model, covariates = village)
bounds.check <- ovb_partial_r2_bound(r2dxj.x = r2d, r2yxj.dx = r2y, kd = 0.3)
bounds <- sensitivity$bounds
expect_equal(bounds$r2dz.x, bounds.check$r2dz.x)
expect_equal(bounds$r2yz.dx, bounds.check$r2yz.dx)
})
Any scripts or data that you put into this service are public.
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.