R/gbiqq_designer.R

In lilymedina/gbiqqtools: Tools for generalized BIQQ

#' Integrate gbiqq with DeclareDesign
#'
#' Define a model, estimand, and answer strategy.
#'
#' @param model A model object generated by \code{make_model()}.
#' @param restrictions A list of character vectors specifying nodal types to be kept or removed from the model. Use \code{get_nodal_types} to see syntax.
#' @param priors A numeric vector. Dirichlet hyperparameters.
#' @param parameters A numeric vector. True parameters.
#' @param inquiry A named list of strings. Names are estimand labels and values are strings containing estimand definitions. See `Examples` for correct format.
#' @param data_strat should be of the form "list(n_obs = n_obs, vars = list("X", "Y"), probs = list(...), n = NULL, subsets = list(...))
#' @param answer_strat priors used at the analysis stage, if different
#' @import DeclareDesign
#' @importFrom dplyr %>%
#' @import Rcpp
#' @export
#' @examples
#' require("DeclareDesign")
#' my_design <- gbiqq_designer(
#'   model = make_model("X -> Y"),
#'   parameters = c(.5, .5, .1, .7, .1, .1),  # Truth
#'   data_strat = list(n_obs = 5, vars = list(NULL), probs = list(NULL),
#'                     ns = NULL,
#'      subsets = list(NULL)),
#'   inquiry = list(ATE = "Y[X=1] - Y[X=0]"))
#'
#' \dontrun{draw_data(my_design)}
#' \dontrun{draw_estimands(my_design)}
#' \dontrun{draw_estimates(my_design)}
#' \dontrun{diagnose_design(my_design, sims = 2)}
#'
#' my_design <- gbiqq_designer(
#' model = make_model("X -> M -> Y"),
#' inquiry = list(ATE = "Y[X=1] - Y[X=0]"),
#' data_strat = list(n_obs = 5,
#' vars = list(c("X", "Y"), "M"),
#' probs = list(1, .5),
#' n = NULL,
#' subsets = list(NULL, "X==1 & Y==0")),
#' answer_strat  = NULL
#' )
#'
#' df <- draw_data(my_design)
#' draw_estimands(my_design)
#' \dontrun{
#' get_estimates(my_design, df)
#' }

gbiqq_designer <- function(
	model = make_model("X -> Y"),
	restrictions = NULL,
	priors = "uniform",
	parameters = NULL,        # True parameters
	inquiry = list(ATE = "Y[X=1] - Y[X=0]"),
	data_strat = NULL,
	answer_strat = NULL
) {

	answer_model <- model <- 	model %>%
		set_restrictions(restrictions) %>%
		set_priors(prior_distribution = priors) %>%
		set_prior_distribution()

	if(!is.null(answer_strat)) answer_model <- set_priors(prior_distribution = answer_strat) %>% # Need to be generalized
			set_prior_distribution()


	if(is.null(parameters)) {message("No true parameters provided; parameters drawn from prior on each run")}

	data_step <- declare_population(data =
																		data_strategy(
																			model,
																			parameters = parameters,
																			n_obs =      data_strat$n_obs,
																			vars =       data_strat$vars,
																			probs =      data_strat$probs,
																			n =          data_strat$n,
																			subsets =    data_strat$subsets))

	# Estimand given parameters
	estimand <- declare_estimand(handler = function(data) {
		value <- query_model(model,
													 using = "parameters",
													 parameters = parameters,
													 queries = inquiry)
		data.frame(estimand_label = names(inquiry),
							 estimand = value$mean,
							 stringsAsFactors = FALSE)})

	# Estimator runs gbiqq assuming answer-strategy model
	estimate <- declare_estimator(handler = function(data) {
		updated <- gbiqq(model = answer_model,  data = data)
		value   <- gbiqq::query_model(updated, using = "posteriors", queries = inquiry)
		data.frame(estimate_label = paste0("est_", names(inquiry)),
							 estimand = names(inquiry),
							 estimate = value$mean,
							 sd_estimate = value$sd, stringsAsFactors = FALSE)
	})

	# Declare design
	data_step + estimand + estimate

}