Description Usage Arguments Details Value Author(s) Examples
This function performs MetropolisHastings exploration of the posterior distribution for Bayesian discrete time survival models. The models may have competing risks, semiparametric covariate effects (like arbitrary baseline hazards), and random effects shared across repeated events and correlated across competing risks.
1 
x 
a matrix of positive whole numbers or a dataframe with positive whole or factor columns specifying the values of the (discretized) predictors at each persontime point; the (i,j) entry of 
y 
a nonnegative whole number vector or matrix indicating event occurrence at each persontime point; with only one event type, the i^th entry of the vector 
priors 
a list with one element for each predictor variable (column of 
S 
the number of MCMC iterations to perform postburnin (postburnin iterations are stored in the object returned by the 
B 
the number of burnin MCMC iterations to perform (burnin iterations are not stored) 
n 
a vector of positive whole numbers with length equal to the number of persontime observations; the i^th entry of 
K 
a positive whole number vector whose m^th entry is the number of distinct values that the m^th discretized predictor may assume; if 
store_re 
if 
The data provided to the brea_mcmc
function is specified at the persontime level: there is one row in x
and y
for each discrete time point each person or thing was at risk for event occurrence. All predictors in x
must be encoded as factors or their corresponding integer codes. The underlying type of predictor is specified in the priors
argument, which is a list with one element for each predictor variable which specifies both the type of that predictor and the prior distribution to use. The allowed predictor types are:
"cat"
for categorical variables. The first element of the prior specification list is the string "cat"
, and the second element is a positive conditional prior standard deviation parameter.
"gmrf"
for underlying continuous predictors; continuous predictors should be discretized with cut()
before being included into x
; Gaussian Markov random field (GMRF) priors are then used to smooth the effects of adjacent categories of the discretized continuous predictor. The first element of the prior specification list is the string "gmrf"
, the second element is a prior degrees of freedom for the scaled inverse chisquared prior for the random walk increment variance, and the third element is a prior scale for the scaled inverse chisquared.
"re"
for variables (like subject id numbers) that represent random effects. The first element of the prior specification list is the string "re"
, the second element is a prior degrees of freedom for an inverse Wishart prior for the random effects covariance, and the third element is a prior scale matrix for the random effects covariance.
"zero"
for predictors that are not used in the current MCMC run. This is provided as a convenient way to exclude predictors from certain runs. The first and only element of the prior specification list is the string "zero"
.
A list providing the values from each of the S
postburnin MCMC iterations for the intercept, the other linear predictor parameters, standard deviations of the GMRF random walk increments, and covariances of the random effects:
b_0_s 
an 
b_m_s 
a list whose m^th element is the 
s_m_s 
a list whose m^th element is, in the case of a GMRF prior, the 
b_m_a 
a list whose m^th element is the length 
Adam King
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53  # leukemia remission times data from Table 1 of "Regression Models and
# LifeTables" (Cox, 1972)
# times of event occurrence or right censoring:
time < c(6,6,6,6,7,9,10,10,11,13,16,17,19,20,22,23,25,32,32,34,35,
1,1,2,2,3,4,4,5,5,8,8,8,8,11,11,12,12,15,17,22,23)
# event/censoring indicator (1 for event occurrence, 0 for right censoring):
event < c(0,1,1,1,1,0,0,1,0,1,1,0,0,0,1,1,0,0,0,0,0,
1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1)
# treatment group indicator (1 for treatment, 2 for control):
treat < c(rep(1,21),rep(2,21))
# total number of persontime observations among all subjects:
N < sum(time)
# create and fill in persontime predictor and event variables:
x < matrix(0,N,2) # first column is time, second column is treatment
y < numeric(N)
next_row < 1 # next available row in the persontime variables
for (i in seq_along(time)) {
rows < next_row:(next_row + time[i]  1) # observations for subject i
x[rows,1] < seq_len(time[i]) # time is integers 1,2,...,time[i]
x[rows,2] < rep(treat[i],time[i]) # treatment group is constant
y[rows] < c(rep(0,time[i]  1),event[i]) # outcome is 0's until event
next_row < next_row + time[i] # increment the next available row pointer
}
# group the timeatrisk variable into 3week intervals:
x[,1] < cut(x[,1],seq(0,36,3),labels=FALSE)
# use GMRF prior for time, and categorical prior for treatment group:
priors < list(list("gmrf",3,.01),list("cat",4))
# run the default of 100 burnin iterations followed by 1,000 stored iterations:
fit < brea_mcmc(x,y,priors)
# examine the structure of the returned posterior samples and acceptance counts:
str(fit)
# posterior samples of the treatment and control group parameters:
b_treatment < fit$b_m_s[[2]][1,1,]
b_control < fit$b_m_s[[2]][1,2,]
# posterior sample of treatment effect on linear predictor scale:
d < b_control  b_treatment
# posterior mean, median, and sd of treatment effect on linear predictor scale:
mean(d); median(d); sd(d)
# posterior mean and median hazard ratios:
mean(exp(d)); median(exp(d))

Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.