cfd.FUN {cfdecomp}

Flexible Function Decomposition: decompose any function that returns a vector

Description

Flexible Function Decomposition: decompose any function that returns a vector

Usage

cfd.FUN(
  formula.y,
  formula.m,
  mediator,
  group,
  data,
  family.y = "binomial",
  family.m = "binomial",
  bs.size = 250,
  mc.size = 50,
  FUN.y = mean,
  alpha = 0.05,
  cluster.sample = FALSE,
  cluster.name = NA,
  cluster.mrows = FALSE,
  sample.resid.y = FALSE,
  sample.resid.m = FALSE,
  print.iteration = FALSE,
  ...
)

Arguments

formula.y

the formula for the multivariable model (see glm) for the outcome Y.

formula.m

the formula for the multivariable model (see glm) for the mediator M.

mediator

the column name of the mediator M.

group

column name of a factor variable containing the group identifier.

data

a data frame containing the variables in the model.

family.y

a description of the error distribution to be used in the model, see family for details. For the outcome variable any member of the glm family can be used.

family.m

a description of the error distribution to be used in the model, see family for details. For the mediator, currently gaussian, binomial and poisson are supported.

bs.size

the number of bootstrap iterations to be performed.

mc.size

the number of Monte Carlo iterations to be performed (more = more MC error reduction).

FUN.y

a function to compute the statistics which can be applied to all data subsets, this function should return a vector and should be ran on pred_y (simulated y values in the natural course or counterfactual) and optional additional columns.

alpha

the alpha level used to construct confidence intervals (0.05 = 95 percent confidence interval).

cluster.sample

set to TRUE if data are clustered in the long format (i.e. multiple rows per individual or other cluster).

cluster.name

the name (as a character) of the column containing the cluster identifiers.

cluster.mrows

for the mediator model, only allows 1 observation per mediator so that the mediator model is not weighted by number of observations. e.g. set to TRUE if the mediator is time constant in longitudinal analysis of long format data.

sample.resid.y

sample.resid if the outcome is Gaussian, should the simulation sample from the residuals of the linear regression model of the outcome to approximate the empirical distribution of the outcome in the simulation (Monte Carlo integration) (if so, set to TRUE), or should it sample from a Gaussian distribution with the standard deviation of the outcome? If the true distribution of the continuous outcome is not very Gaussian, the former may be preferred.

sample.resid.m

sample.resid if the mediator is Gaussian, should the simulation sample from the residuals of the linear regression model of the mediator to approximate the empirical distribution of the mediator in the simulation (Monte Carlo integration) (if so, set to TRUE), or should it sample from a Gaussian distribution with the standard deviation of the mediator? If the true distribution of the continuous mediator is not very Gaussian, the former may be preferred.

print.iteration

print the bootstrap iteration

...

further arguments passed to or used by methods.

Value

out_nc_m returns the mean level of the mediator under the natural course, which is a value that should be close to the empirically observed value of the mediator for each group. out_nc_quantile provides the alpha/2 and 1-alpha/2 bootstrap quantiles for this mean (AKA bootstrap percentile confidence intervals). out_nc_y provides the output of the function fed into FUN.y for each bootstrap iteration, with out_nc_quantile_y providing the alpha/2 and 1-alpha/2 bootstrap quantiles of that output. Similarly, out_cf_m, out_cf_quantile_m,out_cf_y, and out_cf_quantile_y provide the corresponding values for the counterfactual scenario where the mediators of the groups are equalized. mediation and mediation_quantile are not provided for this function, so should be calculated by the user based on the output. mc_conv_info_m and mc_conv_info_y provide information that can help determine the number of Monte Carlo and Bootstrap iterations needed to achieve stability. See the Examples for more information.

Examples

set.seed(100)
# the decomposition functions in our package are computationally intensive
# to make the example run quick, I perform it on a subsample (n=250) of the data:
cfd.example.sample <- cfd.example.data[sample(250),]
# define some function (here one that calculates the mean from the data)
# such a function already exists, but this is to demonstrate how to do it for one that
# will be implemented in cfd.FUN:
mean.fun <- function(data,yname) {
x <- data
return(mean(x[,yname],na.rm=TRUE))
}
# test if the function works on normal data:
mean.fun(cfd.example.sample,yname="med.pois")
# then enter it into cfd.FUN and run:
mean.results <- cfd.FUN(formula.y='out.gauss ~ SES + med.gauss + med.binom + age',
                          formula.m='med.gauss ~ SES + age',
                          mediator='med.gauss',
                          group='SES',
                          data=cfd.example.sample,
                          family.y='gaussian',
                          family.m='gaussian',
                          FUN.y=mean.fun,
                          bs.size=15,
                          mc.size=5,
                          alpha=0.05,
                          print.iteration=TRUE,
                          yname="pred_y")
# more advanced code demonstrating how to do this with a function that calculates
# the age-adjusted rate ratio and life expectancy will hopefully soon be available
# in a publication.
#' @import stats utils