wgtmed {twangMediation}

Weighted mediation analysis.

Description

Estimate causal mediation mechanism of a treatment using propensity score weighting.

Usage

wgtmed(
  formula.med,
  data,
  a_treatment,
  y_outcome = NULL,
  med_interact = NULL,
  total_effect_wts = NULL,
  total_effect_ps = NULL,
  total_effect_stop_rule = NULL,
  method = "ps",
  sampw = NULL,
  ps_n.trees = 10000,
  ps_interaction.depth = 3,
  ps_shrinkage = 0.01,
  ps_bag.fraction = 1,
  ps_n.minobsinnode = 10,
  ps_perm.test.iters = 0,
  ps_verbose = FALSE,
  ps_stop.method = c("ks.mean", "ks.max"),
  ps_version = "gbm",
  ps_ks.exact = NULL,
  ps_n.keep = 1,
  ps_n.grid = 25,
  ps_cv.folds = 10,
  ps_keep.data = FALSE
)

Arguments

formula.med

A object of class formula relating the mediatior(s) to the covariates (potential confounding variables).

data

A dataset of class data.frame that includes the treatment indicator, mediator(s), and covariates.

a_treatment

The (character) name of the treatment variable, which must be dichotomous (0, 1).

y_outcome

The (character) name of the outcome variable, y. If this is not provided, then no effects will be calculated and a warning will be raised. Default : NULL.

med_interact

The (character) vector of names of variables specified on the right-hand side of formula.med that consist of crossproducts or interactions between a covariate and the mediator. See the tutorial for details on these variables.

total_effect_wts

A vector of total effect weights, which if left NULL then total_effect_ps must be supplied. Default : NULL.

total_effect_ps

A ps object that contains the total effect weights,

total_effect_stop_rule

The stopping rule (ks.mean, ks.max, es.mean, es.max) for the total effect weights, which only needs to be specified if total_effect_ps is provided. Default : NULL.

method

The method for getting weights ("ps", "logistic", or "crossval"). Default : "ps".

sampw

Optional sampling weights Default : NULL.

ps_n.trees

Number of gbm iterations passed on to gbm. Default: 10000.

ps_interaction.depth

A positive integer denoting the tree depth used in gradient boosting. Default: 3.

ps_shrinkage

A numeric value between 0 and 1 denoting the learning rate. See gbm for more details. Default: 0.01.

ps_bag.fraction

A numeric value between 0 and 1 denoting the fraction of the observations randomly selected in each iteration of the gradient boosting algorithm to propose the next tree. See gbm for more details. Default: 1.0.

ps_n.minobsinnode

An integer specifying the minimum number of observations in the terminal nodes of the trees used in the gradient boosting. See gbm for more details. Default: 10.

ps_perm.test.iters

A non-negative integer giving the number of iterations of the permutation test for the KS statistic. If perm.test.iters=0 then the function returns an analytic approximation to the p-value. Setting perm.test.iters=200 will yield precision to within 3% if the true p-value is 0.05. Use perm.test.iters=500 to be within 2%. Default: 0.

ps_verbose

If TRUE, lots of information will be printed to monitor the the progress of the fitting. Default: FALSE.

ps_stop.method

A method or methods of measuring and summarizing balance across pretreatment variables. Current options are ks.mean, ks.max, es.mean, and es.max. ks refers to the Kolmogorov-Smirnov statistic and es refers to standardized effect size. These are summarized across the pretreatment variables by either the maximum (.max) or the mean (.mean). Default: c("ks.mean", "ks.max").

ps_version

"gbm", "xgboost", or "legacy", indicating which version of the twang package to use.

• "gbm" uses gradient boosting from the gbm package.

• "xgboost" uses gradient boosting from the xgboost package.

• "legacy" uses the prior implementation of the ps function.

ps_ks.exact

NULL or a logical indicating whether the Kolmogorov-Smirnov p-value should be based on an approximation of exact distribution from an unweighted two-sample Kolmogorov-Smirnov test. If NULL, the approximation based on the exact distribution is computed if the product of the effective sample sizes is less than 10,000. Otherwise, an approximation based on the asymptotic distribution is used. Warning: setting ks.exact = TRUE will add substantial computation time for larger sample sizes. Default: NULL.

ps_n.keep

A numeric variable indicating the algorithm should only consider every n.keep-th iteration of the propensity score model and optimize balance over this set instead of all iterations. Default : 1.

ps_n.grid

A numeric variable that sets the grid size for an initial search of the region most likely to minimize the stop.method. A value of n.grid=50 uses a 50 point grid from 1:n.trees. It finds the minimum, say at grid point 35. It then looks for the actual minimum between grid points 34 and 36.If specified with n.keep>1, n.grid corresponds to a grid of points on the kept iterations as defined by n.keep. Default: 25.

ps_cv.folds

A numeric variable that sets the number of cross-validation folds if using method='crossval'. Default: 10.

ps_keep.data

A logical variable that determines if the dataset should be saved in the resulting ps model objects. Default: FALSE.

Details

For users comfortable with ps, any options prefaced with ps_ are passed directly to the ps() function. Model A is used to estimate Pr(A=1 | X) where X is the vector of background covariates specified in formula.med. If method equals "ps" model A is fit using the twang ps function with estimand= "ATE". If method equals "logistic" then model A is fit using logistic regression. If method equals "crossval" then gbm using cross-validation is used to estimate model A. Because X might include variables not used to estimate the user-provided total effect weights, model A is fit rather than using the user-provided total effect weights to derive Pr(A | X). If the user uses the same set of variables to estimate their provided total effect weights as they enter in the wgtmed function to estimate the cross-world weights and the user uses the same estimation method and arguments as specified in the wgtmed function, then the estimated model A will match the model the user used to obtain the provided total effect weights.

Value

mediation object The mediation object includes the following:

• model_a The model A ps() results.

• model_m1 The model M1 ps() results.

• model_m0 The model M0 ps() results.

• data The data set used to compute models

• stopping_methods The stopping methods passed to stop.method.

• datestamp The date when the analysis was run.

• For each stop.method, a list with the following:

  • TE The total effect.

  • NDE_0 The natural direct effect, holding the mediator constant at 0.

  • NIE_1 The natural indirect effect, holding the exposure constant at 1.

  • NDE_1 The natural direct effect, holding the mediator constant at 1.

  • NIE_0 The natural indirect effect, holding the exposure constant at 0.

  • expected_treatment0_mediator0 E(Y(0, M(0)))

  • expected_treatment1_mediator1 E(Y(1, M(1)))

  • expected_treatment1_mediator0 E(Y(1, M(0)))

  • expected_treatment0_mediator1 E(Y(0, M(1)))

• dx.wts A list with information for checking covariate balance of for each estimated effect. Elements are TE, NIE1, NDE0, NIE0, NDE1, with results of twang dx.wts for the covariates when weighted by weights used in the estimating the effect.

See Also

ps

Examples

data("tMdat")

## tMdat is small simulated data set included in twangMediation for 
## demonstrating the functions. See ?tMdat for details

head(tMdat)

## The tMdat data contains the following variables:
##   w1, w2, w3 -- Simulatad covariates
##   A   -- Simulated dichotomous exposure indicator
##   M   -- Simulated discrete mediator (11 values)
##   Y   -- Simulated continuous outcome
##   te.wgt -- Estimated inverse probability weight, estimated using 
##             GBM via the twang ps function

fit.es.max <- wgtmed(M ~ w1 + w2 + w3,
                      data = tMdat,
                      a_treatment = "A",
                      y_outcome = "Y",
                      total_effect_wts = tMdat$te.wgt,
                      method = "ps",
                      ps_n.trees=1500,
                      ps_shrinkage=0.01,
                      ps_stop.method=c("es.max")
                      )

fit.es.max