Robust Experimental Design

Overview

robustExpDesign designs experiments that will be combined with observational or other external evidence. Because the external evidence may be misspecified, the design balances statistical precision against worst-case bias exposure.

A typical workflow is:

  1. specify a target sensitivity vector omega, or matrix Omega;
  2. provide covariance and experimental-variance inputs;
  3. specify costs and feasibility restrictions;
  4. solve a proportional-regret design problem;
  5. inspect experiment selection, estimator weights, sample allocation, and regret components.

The package implements methods from Epanomeritakis and Viviano (2026), Learning What to Learn: Experimental Design when Combining Experimental with Observational Evidence, manuscript dated August 4, 2026.

Installation and optional solvers

Matrix and rlang are required dependencies. Optional features use:

install.packages(c("quadprog", "CVXR", "ggplot2"))
library(robustExpDesign)

Baseline proportional-regret design

Consider three parameters. Sigma_obs is the covariance matrix of the external estimates, and v2[j] / n_j is the variance of experimental estimate j when sample size n_j is assigned to that experiment.

Sigma_obs <- matrix(
  c(
    0.08, 0.01, 0.00,
    0.01, 0.05, 0.01,
    0.00, 0.01, 0.09
  ),
  nrow = 3,
  byrow = TRUE
)

v2 <- c(direct = 1.00, income = 1.40, wage = 1.10)
omega <- c(direct = 0.40, income = 1.00, wage = -0.70)

fit <- solve_minimax_design(
  Sigma_obs = Sigma_obs,
  v2 = v2,
  n_total = 500,
  omega = omega,
  costs = c(1, 1, 1.5),
  bias_weights = c(1, 1, 1),
  h = 2,
  min_experiments = 1,
  solver = "quadprog"
)

The objective is the maximum of normalized variance and bias regret. The main outputs are:

fit$selected
fit$x_opt
fit$gamma_opt
fit$s_opt
fit$a_exp_opt
fit$a_obs_opt
fit$n_opt
fit$alpha_ratio
fit$beta_ratio
fit$regret

Evaluating a supplied design

evaluate_design() computes allocation, variance, and bias for a selection and weighting rule supplied by the user.

evaluated <- evaluate_design(
  Sigma_obs = Sigma_obs,
  v2 = v2,
  n_total = 500,
  omega = omega,
  x = c(1, 1, 0),
  gamma = c(0.8, 0.6, 0),
  costs = c(1, 1, 1.5)
)

evaluated$n_opt
evaluated$alpha
evaluated$beta

Weighted ambiguity sets

bias_weights implements coordinate-specific bounds

|b_j| <= bias_weights[j] * B.

Zero weights mark coordinates treated as not misspecified.

fit_weighted <- solve_minimax_design(
  Sigma_obs = Sigma_obs,
  v2 = v2,
  n_total = 500,
  omega = omega,
  costs = c(1, 1, 1.5),
  bias_weights = c(1.0, 0.5, 2.0),
  h = 2,
  min_experiments = 1,
  solver = "quadprog"
)

Feasibility restrictions

The feasible menu can be restricted with experiment-count bounds, eligibility indicators, explicit feasible sets, or linear constraints on the selection vector.

selection_constraints <- list(
  A = matrix(c(1, 1, 0), nrow = 1),
  sense = "<=",
  rhs = 1
)

fit_restricted <- solve_minimax_design(
  Sigma_obs = Sigma_obs,
  v2 = v2,
  n_total = 500,
  omega = omega,
  costs = c(1, 1, 1.5),
  h = 2,
  selection_constraints = selection_constraints,
  min_experiments = 1,
  solver = "quadprog"
)

Positive n_min imposes a minimum allocation on each selected experiment and currently requires Gurobi. The floor is applied to candidate designs and oracle benchmarks.

fit_floor <- solve_minimax_design(
  Sigma_obs = Sigma_obs,
  v2 = v2,
  n_total = 500,
  omega = omega,
  costs = c(1, 1, 1.5),
  h = 2,
  min_experiments = 1,
  n_min = 50,
  solver = "gurobi"
)

Budget sweeps

sweep_n_total() evaluates a grid of budgets and experiment-count limits.

sweep <- sweep_n_total(
  Sigma_obs = Sigma_obs,
  v2 = v2,
  omega = omega,
  n_grid = c(200, 300, 400, 500, 750, 1000),
  h_values = c(1, 2),
  costs = c(1, 1, 1.5),
  experiment_names = c("Direct", "Income", "Wage"),
  include_variance_design = TRUE,
  make_plots = TRUE,
  solver = "quadprog"
)

head(sweep$data$arms)
head(sweep$data$regret)
sweep$plots$allocation
sweep$plots$regret

Moment selection and GMM weights

For general moment-loading problems, provide a moment Jacobian Lambda, moment covariance Sigma, and target sensitivity omega or Omega.

With optimize_W = FALSE, the function evaluates a supplied weighting matrix.

Lambda <- matrix(
  c(
    1.0, 0.0,
    0.0, 1.0,
    1.0, 0.5,
    0.5, 1.0
  ),
  nrow = 4,
  byrow = TRUE
)

Sigma_mom <- diag(c(0.05, 0.06, 0.10, 0.12))
Omega <- matrix(c(0.4, 1.0), nrow = 1)

fixed_moments <- solve_moment_design(
  Lambda = Lambda,
  Sigma = Sigma_mom,
  Omega = Omega,
  W = diag(4),
  biased_moments = c(3, 4),
  norm = "l2",
  optimize_W = FALSE
)

With optimize_W = TRUE, CVXR optimizes the equivalent linear GMM estimator over the supplied moment sets or masks.

optimized_moments <- solve_moment_design(
  Lambda = Lambda,
  Sigma = Sigma_mom,
  Omega = Omega,
  biased_moments = c(3, 4),
  bias_weights = c(0, 0, 1, 1.5),
  norm = "linf",
  moment_sets = list(c(1, 2, 3), c(1, 2, 4), c(1, 2, 3, 4)),
  optimize_W = TRUE,
  cvxr_solver = "CLARABEL"
)

Audience-regret designs

solve_audience_regret_design() solves the Appendix C.2 finite-grid objective. The scalarization grid represents

lambda = B^2 / (1 + B^2)

and must include 0 and 1. The current implementation requires Gurobi.

fit_audience <- solve_audience_regret_design(
  Sigma_obs = Sigma_obs,
  v2 = v2,
  n_total = 500,
  omega = omega,
  costs = c(1, 1, 1.5),
  bias_weights = c(1, 1, 1),
  h = 2,
  min_experiments = 1
)

fit_audience$r_opt
fit_audience$risk_by_lambda