crisk.bart {BART}  R Documentation 
Here we have implemented a simple and direct approach to utilize BART for competing risks that is very flexible, and is akin to discretetime survival analysis. Following the capabilities of BART, we allow for maximum flexibility in modeling the dependence of competing failure times on covariates. In particular, we do not impose proportional hazards.
To elaborate, consider data in the form: (s_i, \delta_i,
{x}_i)
where s_i
is the event time;
\delta_i
is an indicator distinguishing events,
\delta_i=h
due to cause h in {1, 2}
, from
rightcensoring, \delta_i=0
; {x}_i
is a vector of
covariates; and i=1, ..., N
indexes subjects.
We denote the K
distinct event/censoring times by
0<t_{(1)}<...<t_{(K)}<\infty
thus
taking t_{(j)}
to be the j^{th}
order statistic
among distinct observation times and, for convenience,
t_{(0)}=0
. Now consider event indicators for cause
h
: y_{hij}
for each subject i
at each distinct
time t_{(j)}
up to and including the subject's last
observation time s_i=t_{(n_i)}
with n_i=\arg \max_j
[t_{(j)}\leq s_i]
for cause 1, but only up to
n_iy_{1ij}
for cause 2.
We then denote by p_{hij}
the probability of an event at
time t_{(j)}
conditional on no previous event. We now write
the model for y_{hij}
as a nonparametric probit (or
logistic) regression of y_{hij}
on the time
t_{(j)}
and the covariates {x}_{hi}
, and then
utilize BART for binary responses. Specifically, y_{hij}\ =\
I[\delta_i=h] I[s_i=t_{(j)}],\ j=1, ..., n_iI[h=2]y_{1ij}
. Therefore, we have
p_{hij} = F(mu_{hij}),\ mu_{hij} = mu_h+f_h(t_{(j)},
{x}_{hi})
where
F
denotes the Normal (or Logistic) cdf. As in the binary response
case, f_h
is the sum of many tree models. Finally, based on
these probabilities, p_{hij}
, we can construct targets of
inference such as the cumulative incidence functions.
crisk.bart(x.train=matrix(0,0,0), y.train=NULL,
x.train2=x.train, y.train2=NULL,
times=NULL, delta=NULL, K=NULL,
x.test=matrix(0,0,0), x.test2=x.test, cond=NULL,
sparse=FALSE, theta=0, omega=1,
a=0.5, b=1, augment=FALSE,
rho=NULL, rho2=NULL,
xinfo=matrix(0,0,0), xinfo2=matrix(0,0,0),
usequants=FALSE,
rm.const=TRUE, type='pbart',
ntype=as.integer(
factor(type, levels=c('wbart', 'pbart', 'lbart'))),
k=2, power=2, base=0.95,
offset=NULL, offset2=NULL,
tau.num=c(NA, 3, 6)[ntype],
ntree=50, numcut=100, ndpost=1000, nskip=250,
keepevery = 10L,
printevery=100L,
id=NULL, ## crisk.bart only
seed=99, ## mc.crisk.bart only
mc.cores=2, ## mc.crisk.bart only
nice=19L ## mc.crisk.bart only
)
mc.crisk.bart(x.train=matrix(0,0,0), y.train=NULL,
x.train2=x.train, y.train2=NULL,
times=NULL, delta=NULL, K=NULL,
x.test=matrix(0,0,0), x.test2=x.test, cond=NULL,
sparse=FALSE, theta=0, omega=1,
a=0.5, b=1, augment=FALSE,
rho=NULL, rho2=NULL,
xinfo=matrix(0,0,0), xinfo2=matrix(0,0,0),
usequants=FALSE,
rm.const=TRUE, type='pbart',
ntype=as.integer(
factor(type, levels=c('wbart', 'pbart', 'lbart'))),
k=2, power=2, base=0.95,
offset=NULL, offset2=NULL,
tau.num=c(NA, 3, 6)[ntype],
ntree=50, numcut=100, ndpost=1000, nskip=250,
keepevery = 10L,
printevery=100L,
id=NULL, ## crisk.bart only
seed=99, ## mc.crisk.bart only
mc.cores=2, ## mc.crisk.bart only
nice=19L ## mc.crisk.bart only
)
x.train 
Covariates for training (in sample) data of cause 1. 
y.train 
Cause 1 binary response for training (in sample)
data. 
x.train2 
Covariates for training (in sample)
data of cause 2. Similar to 
y.train2 
Cause 2 binary response for training (in sample) data, i.e., failure
from any cause besides the cause of interest which is cause 1.
Similar to 
times 
The time of event or rightcensoring, 
delta 
The event indicator: 1 for cause 1, 2 for cause 2 and 0 is censored. 
K 
If provided, then coarsen 
x.test 
Covariates for test (out of sample) data of cause 1. 
x.test2 
Covariates for test (out of sample) data of cause 2.
Similar to 
cond 
A vector of indices for 
sparse 
Whether to perform variable selection based on a sparse Dirichlet prior; see Linero 2016. 
theta 
Set 
omega 
Set 
a 
Sparse parameter for 
b 
Sparse parameter for 
rho 
Sparse parameter: typically 
rho2 
Sparse parameter: typically 
augment 
Whether data augmentation is to be performed in sparse variable selection. 
xinfo 
You can provide the cutpoints to BART or let BART
choose them for you. To provide them, use the 
xinfo2 
Cause 2 cutpoints. 
usequants 
If 
rm.const 
Whether or not to remove constant variables. 
type 
Whether to employ probit BART via AlbertChib,

ntype 
The integer equivalent of 
k 
k is the number of prior standard deviations

power 
Power parameter for tree prior. 
base 
Base parameter for tree prior. 
offset 
Cause 1 binary offset. 
offset2 
Cause 2 binary offset. 
tau.num 
The numerator in the 
ntree 
The number of trees in the sum. 
numcut 
The number of possible values of cutpoints (see

ndpost 
The number of posterior draws returned. 
nskip 
Number of MCMC iterations to be treated as burn in. 
keepevery 
Every 
printevery 
As the MCMC runs, a message is printed every 
id 

seed 

mc.cores 

nice 

crisk.bart
returns an object of type criskbart
which is
essentially a list. Besides the items listed
below, the list has offset
, offset2
,
times
which are the unique times, K
which is the number of unique times, tx.train
and
tx.test
, if any.
yhat.train 
A matrix with 
yhat.test 
Same as 
surv.test 
test data fits for the survival function, 
surv.test.mean 
mean of 
prob.test 
The probability of suffering cause 1. 
prob.test2 
The probability of suffering cause 2. 
cif.test 
The cumulative incidence function of cause 1,

cif.test2 
The cumulative incidence function of cause 2,

cif.test.mean 
mean of 
cif.test2.mean 
mean of 
varcount 
a matrix with 
varcount2 
For each variable the total count of the number of times this variable is used for cause 2 in a tree decision rule is given. 
crisk.pre.bart
, predict.criskbart
,
mc.crisk.pwbart
, crisk2.bart
data(transplant)
pfit < survfit(Surv(futime, event) ~ abo, transplant)
# competing risks for type O
plot(pfit[4,], xscale=7, xmax=735, col=1:3, lwd=2, ylim=c(0, 1),
xlab='t (weeks)', ylab='AalenJohansen (AJ) CI(t)')
legend(450, .4, c("Death", "Transplant", "Withdrawal"), col=1:3, lwd=2)
## plot(pfit[4,], xscale=30.5, xmax=735, col=1:3, lwd=2, ylim=c(0, 1),
## xlab='t (months)', ylab='AalenJohansen (AJ) CI(t)')
## legend(450, .4, c("Death", "Transplant", "Withdrawal"), col=1:3, lwd=2)
delta < (as.numeric(transplant$event)1)
## recode so that delta=1 is cause of interest; delta=2 otherwise
delta[delta==1] < 4
delta[delta==2] < 1
delta[delta>1] < 2
table(delta, transplant$event)
times < pmax(1, ceiling(transplant$futime/7)) ## weeks
##times < pmax(1, ceiling(transplant$futime/30.5)) ## months
table(times)
typeO < 1*(transplant$abo=='O')
typeA < 1*(transplant$abo=='A')
typeB < 1*(transplant$abo=='B')
typeAB < 1*(transplant$abo=='AB')
table(typeA, typeO)
x.train < cbind(typeO, typeA, typeB, typeAB)
x.test < cbind(1, 0, 0, 0)
dimnames(x.test)[[2]] < dimnames(x.train)[[2]]
##test BART with token run to ensure installation works
set.seed(99)
post < crisk.bart(x.train=x.train, times=times, delta=delta,
x.test=x.test, nskip=1, ndpost=1, keepevery=1)
## Not run:
## run one long MCMC chain in one process
## set.seed(99)
## post < crisk.bart(x.train=x.train, times=times, delta=delta, x.test=x.test)
## in the interest of time, consider speeding it up by parallel processing
## run "mc.cores" number of shorter MCMC chains in parallel processes
post < mc.crisk.bart(x.train=x.train, times=times, delta=delta,
x.test=x.test, seed=99, mc.cores=8)
K < post$K
typeO.cif.mean < apply(post$cif.test, 2, mean)
typeO.cif.025 < apply(post$cif.test, 2, quantile, probs=0.025)
typeO.cif.975 < apply(post$cif.test, 2, quantile, probs=0.975)
plot(pfit[4,], xscale=7, xmax=735, col=1:3, lwd=2, ylim=c(0, 0.8),
xlab='t (weeks)', ylab='CI(t)')
points(c(0, post$times)*7, c(0, typeO.cif.mean), col=4, type='s', lwd=2)
points(c(0, post$times)*7, c(0, typeO.cif.025), col=4, type='s', lwd=2, lty=2)
points(c(0, post$times)*7, c(0, typeO.cif.975), col=4, type='s', lwd=2, lty=2)
legend(450, .4, c("Transplant(BART)", "Transplant(AJ)",
"Death(AJ)", "Withdrawal(AJ)"),
col=c(4, 2, 1, 3), lwd=2)
##dev.copy2pdf(file='../vignettes/figures/liverBART.pdf')
## plot(pfit[4,], xscale=30.5, xmax=735, col=1:3, lwd=2, ylim=c(0, 0.8),
## xlab='t (months)', ylab='CI(t)')
## points(c(0, post$times)*30.5, c(0, typeO.cif.mean), col=4, type='s', lwd=2)
## points(c(0, post$times)*30.5, c(0, typeO.cif.025), col=4, type='s', lwd=2, lty=2)
## points(c(0, post$times)*30.5, c(0, typeO.cif.975), col=4, type='s', lwd=2, lty=2)
## legend(450, .4, c("Transplant(BART)", "Transplant(AJ)",
## "Death(AJ)", "Withdrawal(AJ)"),
## col=c(4, 2, 1, 3), lwd=2)
## End(Not run)