curelps {blapsr} | R Documentation |
Promotion time cure model with Laplace P-splines.
Description
Fits a promotion time cure model with the Laplace-P-spline methodology. The routine can be applied to survival data for which a plateau is observed in the Kaplan-Meier curve. In this case, the follow-up period is considered to be sufficiently long to intrinsically account for long-term survivors and hence a cured fraction. The user can separately specify the model covariates influencing the cure probability (long-term survival) and the population hazard dynamics (short-term survival).
Usage
curelps(formula, data, K = 30, penorder = 2, tmax = NULL, constr = 6)
Arguments
formula |
A formula object where the ~ operator separates the response
from the covariates. In a promotion time cure model, it takes the form
response ~ covariates, where response is a survival object
returned by the |
data |
Optional. A data frame to match the variable names provided in
|
K |
A positive integer specifying the number of cubic B-spline
functions in the basis. Default is |
penorder |
The penalty order used on finite differences of the coefficients of contiguous B-splines. Can be either 2 for a second-order penalty (the default) or 3 for a third-order penalty. |
tmax |
A user-specified value for the upper bound of the B-spline basis. The default is NULL, so that the B-spline basis is specified in the interval [0, tup], where tup is the upper bound of the follow-up times. It is required that tmax > tup. |
constr |
Constraint imposed on last B-spline coefficient (default is 6). |
Details
The log-baseline hazard is modeled as a linear combination of
K
cubic B-splines as obtained from cubicbs
. A
robust Gamma prior is imposed on the roughness penalty parameter.
A grid-based approach is used to explore the posterior penalty space and
the resulting quadrature points serve to compute the approximate (joint)
posterior of the latent field vector. Point and set estimates of latent
field elements are obtained from a finite mixture of Gaussian densities.
The routine centers the columns of the covariate matrix around their mean
value for numerical stability. See print.curelps
for a
detailed explanation on the output printed by the curelps
function.
Value
An object of class curelps
containing various components
from the promotion time cure model fit. Details can be found in
curelps.object
. Estimates on the baseline survival,
population survival (for a chosen covariate profile) and cure probability
can be obtained with the plot.curelps
and
curelps.extract
routines.
Author(s)
Oswaldo Gressani oswaldo_gressani@hotmail.fr.
References
Cox, D.R. (1972). Regression models and life-tables. Journal of the Royal Statistical Society: Series B (Methodological) 34(2): 187-202.
Bremhorst, V. and Lambert, P. (2016). Flexible estimation in cure survival models using Bayesian P-splines. Computational Statistical & Data Analysis 93: 270-284.
Gressani, O. and Lambert, P. (2018). Fast Bayesian inference using Laplace approximations in a flexible promotion time cure model based on P-splines. Computational Statistical & Data Analysis 124: 151-167.
Lambert, P. and Bremhorst, V. (2019). Estimation and identification issues in the promotion time cure model when the same covariates influence long- and short-term survival. Biometrical Journal 61(2): 275-289.
See Also
curelps.object
, curelps.extract
,
plot.curelps
, print.curelps
,
Surv
.
Examples
## Fit a promotion time cure model on malignant melanoma data
data(melanoma)
medthick <- median(melanoma$thickness)
# Kaplan-Meier estimate to check the existence of a plateau
KapMeier <- survfit(Surv(time,status) ~ 1, data = melanoma)
plot(KapMeier, mark.time = TRUE, mark = 4, xlab = "Time (in years)")
# Fit with curelps
fit <- curelps(Surv(time , status) ~ lt(thickness + ulcer) +
st(thickness + ulcer), data = melanoma, K = 40)
fit
# Cure prediction for median thickness and absence of ulceration
curelps.extract(fit, time = c(2, 4 ,6, 8), curvetype = "probacure",
cred.int = 0.90, covar.profile = c(medthick, 0, medthick, 0))
# Plot of baseline and population survival functions
opar <- par(no.readonly = TRUE)
par(mfrow = c(1, 2))
# Baseline survival
plot(fit, curvetype = "baseline", plot.cred = FALSE, ylim = c(0,1))
# Population survival
plot(fit, curvetype = "population", covar.profile = c(medthick, 0, medthick, 0),
plot.cred = FALSE, ylim = c(0,1))
par(opar)