riAFTBART {riAFTBART}

A flexible approach for causal inference with multiple treatments and clustered survival outcomes

Description

This function implements the random effect accelerated failure time BART (riAFT-BART) for causal inference with multiple treatments and clustered survival outcomes.

Usage

riAFTBART(
  M.burnin,
  M.keep,
  M.thin = 1,
  status,
  y,
  x,
  trt,
  cluster.id,
  verbose = FALSE,
  estimand = "ATE",
  reference_trt = NULL
)

Arguments

M.burnin

A numeric value indicating the number of MCMC iterations to be treated as burn in.

M.keep

A numeric value indicating the number of MCMC posterior draws after burn in.

M.thin

A numeric value indicating the thinning parameter.

status

A vector of event indicators: status = 1 indicates that the event was observed while status = 0 indicates the observation was right-censored.

y

A vector of follow-up times.

x

A dataframe or matrix, including all the covariates but not treatments with rows corresponding to observations and columns to variables.

trt

A numeric vector representing the treatment groups.

cluster.id

A vector of integers representing the clustering id. The cluster id should be an integer and start from 1.

verbose

A logical indicating whether to show the progress bar for riAFT-BART. The default is FALSE

estimand

A character string representing the type of causal estimand. Only "ATT" or "ATE" is allowed. When the estimand = "ATT", users also need to specify the reference treatment group by setting the reference_trt argument.

reference_trt

A numeric value indicating reference treatment group for ATT effect.

Value

A list of causal estimands in terms of log T between different treatment groups.

Examples

library(riAFTBART)
set.seed(20181223)
n = 5       # number of clusters
k = 50      # cluster size
N = n*k     # total sample size
cluster.id = rep(1:n, each=k)
tau.error = 0.8
b = stats::rnorm(n, 0, tau.error)
alpha = 2
beta1 = 1
beta2 = -1
sig.error = 0.5
censoring.rate = 0.02
x1 = stats::rnorm(N,0.5,1)
x2 = stats::rnorm(N,1.5,0.5)
trt.train = sample(c(1,2,3), N, prob = c(0.4,0.3,0.2), replace = TRUE)
trt.test = sample(c(1,2,3), N, prob = c(0.3,0.4,0.2), replace = TRUE)
error = stats::rnorm(N,0,sig.error)
logtime = alpha + beta1*x1 + beta2*x2 + b[cluster.id] + error
y = exp(logtime)
C = rexp(N, rate=censoring.rate) # censoring times
Y = pmin(y,C)
status = as.numeric(y<=C)
res_ate <- riAFTBART(M.burnin = 10, M.keep = 10, M.thin = 1, status = status,
                      y = Y, trt = trt.train,
                      x = cbind(x1,x2),
                      cluster.id = cluster.id, estimand = "ATE")