This is a R/Rcpp package BayesSurvive for Bayesian survival models with graph-structured selection priors for sparse identification of high-dimensional features predictive of survival (Madjar et al., 2021) and its extensions with the use of a fixed graph via a Markov Random Field (MRF) prior for capturing known structure of high-dimensional features, e.g. disease-specific pathways from the Kyoto Encyclopedia of Genes and Genomes (KEGG) database.
Install the latest released version from CRAN
install.packages("BayesSurvive")Install the latest development version from GitHub
#install.packages("remotes")
remotes::install_github("ocbe-uio/BayesSurvive")library("BayesSurvive")
# Load the example dataset
data("simData", package = "BayesSurvive")
dataset = list("X" = simData[[1]]$X,
"t" = simData[[1]]$time,
"di" = simData[[1]]$status)## Initial value: null model without covariates
initial = list("gamma.ini" = rep(0, ncol(dataset$X)))
# Prior parameters
hyperparPooled = list(
"c0" = 2, # prior of baseline hazard
"tau" = 0.0375, # sd (spike) for coefficient prior
"cb" = 20, # sd (slab) for coefficient prior
"pi.ga" = 0.02, # prior variable selection probability for standard Cox models
"a" = -4, # hyperparameter in MRF prior
"b" = 0.1, # hyperparameter in MRF prior
"G" = simData$G # hyperparameter in MRF prior
)
## run Bayesian Cox with graph-structured priors
fit <- BayesSurvive(survObj = dataset, model.type = "Pooled", MRF.G = TRUE,
hyperpar = hyperparPooled, initial = initial, nIter = 100)
## show posterior mean of coefficients and 95% credible intervals
library("GGally")
plot(fit) +
coord_flip() +
theme(axis.text.x = element_text(angle = 90, size = 7))
#plot(fit$output$beta.p[,1], type="l")
#fit$output$beta.margin
#fit$output$gamma.margin
#simData[[1]]$trueB
The function BayesSurvive::plotBrier() can show the
time-dependent Brier scores based on posterior mean of coefficients or
Bayesian model averaging.
plotBrier(fit, , survObj.new = dataset)
The integrated Brier score (IBS) can be obtained by the function
BayesSurvive::predict().
predict(fit, survObj.new = dataset)## IBS
## Null model 0.09147208
## Bayesian Cox model 0.03433363The function BayesSurvive::predict() can estimate the
survival probabilities and cumulative hazards.
predict(fit, survObj.new = dataset, type = c("cumhazard", "survival"))## observation times cumhazard survival
## <int> <num> <num> <num>
## 1: 1 3.3 2.11e-04 1.00e+00
## 2: 2 3.3 3.29e-01 7.20e-01
## 3: 3 3.3 2.06e-06 1.00e+00
## 4: 4 3.3 1.19e-02 9.88e-01
## 5: 5 3.3 5.36e-04 9.99e-01
## ---
## 9996: 96 9.5 2.67e+01 2.57e-12
## 9997: 97 9.5 1.08e+03 0.00e+00
## 9998: 98 9.5 2.23e+00 1.08e-01
## 9999: 99 9.5 3.72e+00 2.42e-02
## 10000: 100 9.5 3.37e+01 2.38e-15hyperparPooled <- append(hyperparPooled, list("lambda" = 3, "nu0" = 0.05, "nu1" = 5))
fit2 <- BayesSurvive(survObj = list(dataset), model.type = "Pooled", MRF.G = FALSE,
hyperpar = hyperparPooled, initial = initial, nIter = 10)# specify a fixed joint graph between two subgroups
hyperparPooled$G <- Matrix::bdiag(simData$G, simData$G)
dataset2 <- simData[1:2]
dataset2 <- lapply(dataset2, setNames, c("X", "t", "di", "X.unsc", "trueB"))
fit3 <- BayesSurvive(survObj = dataset2,
hyperpar = hyperparPooled, initial = initial,
model.type="CoxBVSSL", MRF.G = TRUE,
nIter = 10, burnin = 5)fit4 <- BayesSurvive(survObj = dataset2,
hyperpar = hyperparPooled, initial = initial,
model.type="CoxBVSSL", MRF.G = FALSE,
nIter = 3, burnin = 0)Katrin Madjar, Manuela Zucknick, Katja Ickstadt, Jörg Rahnenführer (2021). Combining heterogeneous subgroups with graph‐structured variable selection priors for Cox regression. BMC Bioinformatics, 22(1):586. DOI: 10.1186/s12859-021-04483-z.