Run a simple simulation based on the bounds. For each simulation, sample the set of counterfactual probabilities from a uniform distribution, translate into a multinomial distribution, and then compute the objective and the bounds in terms of the observable variables.
simulate_bounds(obj, bounds, nsim = 1000)
Object as returned by analyze_graph
Object as returned by optimize_effect
Number of simulation replicates
A data frame with columns: objective, bound.lower, bound.upper
b <- graph_from_literal(X -+ Y, Ur -+ X, Ur -+ Y)
V(b)$leftside <- c(0,0,0)
V(b)$latent <- c(0,0,1)
V(b)$nvals <- c(2,2,2)
E(b)$rlconnect <- E(b)$edge.monotone <- c(0, 0, 0)
obj <- analyze_graph(b, constraints = NULL, effectt = "p{Y(X = 1) = 1} - p{Y(X = 0) = 1}")
bounds <- optimize_effect(obj)
simulate_bounds(obj, bounds, nsim = 5)
#> objective bound.lower bound.upper
#> 1 -0.04861071 -0.5580883 0.4419117
#> 2 -0.00217234 -0.5494420 0.4505580
#> 3 -0.23536287 -0.6237436 0.3762564
#> 4 0.14257190 -0.3649807 0.6350193
#> 5 0.15133314 -0.4878743 0.5121257