Estimating causal `policyFX` with `clusteredinterference`

  collapse = TRUE,
  comment = "#>"

First, load the clusteredinterference package


Now load a quick data example that's included in the package



Estimation is implemented with the policyFX() function:

suppressWarnings(RNGversion("3.5.0")) ## For backwards compatibility 
causal_fx <- policyFX(
  data = toy_data,
  formula = Outcome | Treatment ~ Age + Distance + (1 | Cluster_ID) | Cluster_ID,
  alphas = c(.3, .5), 
  k_samps = 1

The policyFX() function outputs a "policyFX" object, which works well with a few methods, including:


Necessary arguments


A data.frame. At present, tibbles are coerced back to standard data.frames. I also recommend against using factors in the columns.


A numeric vector of probabilities corresponding the the policies of interest. Each must be between 0 and 1.


The user must specify the number of sum-sampled vectors for estimating the counterfactual probabilities (nuisance parameters). It is recommended to choose k_samps <=5. To avoid the sub-sampling approximation and use all possible vectors, set k_samps=0.


The formula may be the trickiest, and it has plenty of information. It provides:

outcome | treatment ~ predictors and random intercept | clustering specification

Note that the middle section is passed to glmer() to fit the mixed effects model, so this is how to specify the modeling formula.

Treatment ~ Age + Distance + (1 | Cluster_ID)

See below for the model output.

Formal arguments with defaults

verbose = FALSE

root_options = NULL

This is for rootSolve::multiroot() used in the point estimation procedure. E.g., this will be passed in:

root_options = list(atol = 1e-7) 


This is for lme4::glmer(). The default in glmer() is nAGQ=1, which indicates a Laplace approximation to the log-likelihood. Instead, in this package the default is nAGQ=2, which indicates that n=2 Adaptive Gaussian Quadrature points will be used. This is slightly slower but is a more accurate calculation. In limited testing, it seems that nAGQ=2 was almost as accurate as higher values, so 2 was chosen as the default. See their documentation for more details.

return_matrices = FALSE

If TRUE, this will return the "bread" and "meat" matrices in the variance calculations. The default is FALSE as these matrices can be quite large.

Dots tricks

In the event you're only interested in a subset of contrasts, you can pass a customized grid into the function.

my_grid <- makeTargetGrid(alphas = (3:8)/20, small_grid = TRUE) 

This can be particularly useful for plotting, as you can "turn off" the variance estimates

my_grid$estVar <- FALSE

This is available through the dots argument. Note that when supplying a custom target_grid, it's not necessary to specify the alphas argument, as that is taken directly from target_grid.

causal_fx2 <- policyFX(
  data = toy_data,
  formula = Outcome | Treatment ~ Age + Distance + (1 | Cluster_ID) | Cluster_ID,
  # alphas = c(.3, .5), 
  target_grid = my_grid,
  k_samps = 5,
  verbose = FALSE,
  root_options = list(atol=1e-4)

print(causal_fx, nrows = 9)


The tidy dataframe estimates can be easily used for plotting:

plotdat <- causal_fx2$estimates[causal_fx2$estimates$estimand_type=="mu",]
plot(x = plotdat$alpha1, y = plotdat$estimate, main = "Estimated Population Means")

Treatment model

As mentioned above, the treatment model is specified via the formula argument. For example, compare:

# Returns the specified formula, coerced to a Formula object
# causal_fx$model is a glmerMod S4 object
lme4::getME(causal_fx$model, c("beta", "theta"))

Try the clusteredinterference package in your browser

Any scripts or data that you put into this service are public.

clusteredinterference documentation built on May 1, 2019, 9:26 p.m.