rcure: Robust cure model

Description Usage Arguments Author(s) References Examples

Description

Fits robust cure model by incorporating a weakly informative prior distribution for uncure probability part in cure models

Usage

1
2
3
4
5
6
7
8
rcure(formula, cureform, offset.s = NULL, data, na.action = na.omit,
  model = c("aft", "ph"), link = "logit", Var = TRUE, emmax = 50,
  eps = 1e-07, nboot = 100, family = binomial(link = "logit"),
  method = "glm.fit", prior.mean = 0, prior.scale = NULL, prior.df = 1,
  prior.mean.for.intercept = 0, prior.scale.for.intercept = NULL,
  prior.df.for.intercept = 1, min.prior.scale = 1e-12, scaled = TRUE,
  n.iter = 100, Warning = TRUE, eva_model = NULL, cutpoint = c(0.1,
  0.25, 0.5, 0.75, 0.9))

Arguments

formula

a formula object for the survival part in cure model. left must be a survival object as returned by the Surv function

cureform

specifies the variables in the uncure probability part in cure model

offset.s

variable(s) with coefficient 1 in PH model or AFT model

data

a a data.frame

na.action

a missing-data filter function. By default na.action = na.omit

model

specifies survival part in cure model, "ph" or "aft"

link

specifies the link in incidence part. The "logit", "probit" or complementary loglog ("cloglog") links are available. By default link = "logit".

Var

If it is TRUE, the program returns Std.Error by bootstrap method. If set to False, the program only returns estimators of coefficients. By default, Var = TRUE

emmax

specifies the maximum iteration number. If the convergence criterion is not met, the EM iteration will be stopped after emmax iterations and the estimates will be based on the last maximum likelihood iteration. The default emmax = 100.

eps

sets the convergence criterion. The default is eps = 1e-7. The iterations are considered to be converged when the maximum relative change in the parameters and likelihood estimates between iterations is less than the value specified.

nboot

specifies the number of bootstrap sampling. The default nboot = 100

family

a description of the error distribution and link function to be used in the model. Default is binomial(link="logit")

method

the method to be used in fitting the glmbayes model. The default method "glm.fit" uses iteratively reweighted least squares (IWLS). The only current alternative is "model.frame" which returns the model frame and does no fitting

prior.mean

prior mean for the coefficients: default is 0. Can be a vector of length equal to the number of predictors (not counting the intercept, if any). If it is a scalar, it is expanded to the length of this vector.

prior.scale

prior scale for the coefficients: default is NULL; if is NULL, for a logit model, prior.scale is 2.5; for a probit model, prior scale is 2.5*1.6. Can be a vector of length equal to the number of predictors (not counting the intercept, if any). If it is a scalar, it is expanded to the length of this vector

prior.df

prior degrees of freedom for the coefficients. For t distribution: default is 1 (Cauchy). Set to Inf to get normal prior distributions. Can be a vector of length equal to the number of predictors (not counting the intercept, if any). If it is a scalar, it is expanded to the length of this vector

prior.mean.for.intercept

prior mean for the intercept: default is 0.

prior.scale.for.intercept

prior scale for the intercept: default is NULL; for a logit model, prior scale for intercept is 10; for probit model, prior scale for intercept is rescaled as 10*1.6.

prior.df.for.intercept

prior degrees of freedom for the intercept: default is 1.

min.prior.scale

Minimum prior scale for the coefficients: default is 1e-12.

scaled

scaled=TRUE, the scales for the prior distributions of the coefficients are determined as follows: For a predictor with only one value, we just use prior.scale. For a predictor with two values, we use prior.scale/range(x). For a predictor with more than two values, we use prior.scale/(2*sd(x)). If the response is Gaussian, prior.scale is also multiplied by 2 * sd(y). Default is TRUE

n.iter

integer giving the maximal number of bayesglm IWLS iterations, default is 100.

Warning

default is TRUE, which will show the error messages of not convergence and separation

eva_model

for Cox PH model, the default is "PH". For AFT model, it can be "PO".

cutpoint

the cut points for ROC to calculate AUC

...

further arguments passed to or from other methods

Author(s)

Xiaoxia Han

References

Cai, C., Zou, Y., Peng, Y., & Zhang, J. (2012). smcure: An R-Package for estimating semiparametric mixture cure models. Computer methods and programs in biomedicine, 108(3), 1255-1260.

Gelman, A., Jakulin, A., Pittau, M. G., & Su, Y. S. (2008). A weakly informative default prior distribution for logistic and other regression models. The Annals of Applied Statistics, 1360-1383.

Examples

 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
library(survival)
library(smcure)
library(arm)
data(e1684)

# fit PH robust cure model
pd <- rcure(Surv(FAILTIME,FAILCENS)~TRT+SEX+AGE,cureform=~TRT+SEX+AGE,
data=e1684,model="ph",Var =FALSE,
method = "glm.fit", prior.mean = 0, prior.scale = NULL, prior.df = 1,
prior.mean.for.intercept = 0, prior.scale.for.intercept = NULL,
prior.df.for.intercept = 1, min.prior.scale = 1e-12,
scaled = FALSE, n.iter = 100, Warning=F)
printrcure(pd,Var = FALSE, ROC=FALSE)
# plot predicted survival curves for male with median centered age by treatment groups
predm=predictrcure(pd,newX=cbind(c(1,0),c(0,0),c(0.579,0.579)),
newZ=cbind(c(1,0),c(0,0),c(0.579,0.579)),model="ph")
plotpredictrcure(predm,model="ph")

# just a test: this should be identical to classical cure model
pd2 <- rcure(Surv(FAILTIME,FAILCENS)~TRT+SEX+AGE,cureform=~TRT+SEX+AGE,
data=e1684,model="ph",Var = FALSE,
method = "glm.fit", prior.mean = 0, prior.scale = Inf, prior.df = Inf,
prior.mean.for.intercept = 0, prior.scale.for.intercept = Inf,
prior.df.for.intercept = Inf, Warning=F)
printrcure(pd2,Var = FALSE, ROC=FALSE)
pd3 <- smcure(Surv(FAILTIME,FAILCENS)~TRT+SEX+AGE,cureform=~TRT+SEX+AGE,
data=e1684,model="ph",Var = FALSE)

data(bmt)
# fit AFT robust cure model
bmtfit <- rcure(formula = Surv(Time, Status) ~ TRT, cureform = ~TRT,
data = bmt, model = "aft", Var = FALSE,
method = "glm.fit", prior.mean = 0, prior.scale = NULL, prior.df = 1,
prior.mean.for.intercept = 0, prior.scale.for.intercept = NULL,
prior.df.for.intercept = 1, min.prior.scale = 1e-12,
scaled = TRUE, n.iter = 100, Warning=F, eva_mode="PO")
printrcure(bmtfit,Var = FALSE, ROC=FALSE)
## plot predicted Survival curves by treatment groups
predbmt=predictrcure(bmtfit,newX=c(0,1),newZ=c(0,1),model="aft")
plotpredictrcure(predbmt,model="aft")

Example output

Loading required package: MASS
Loading required package: Matrix
Loading required package: lme4

arm (Version 1.10-1, built: 2018-4-12)

Working directory is /work/tmp

Program is running..be patient... done.
Call:
rcure(formula = Surv(FAILTIME, FAILCENS) ~ TRT + SEX + AGE, cureform = ~TRT + 
    SEX + AGE, data = e1684, model = "ph", Var = FALSE, method = "glm.fit", 
    prior.mean = 0, prior.scale = NULL, prior.df = 1, prior.mean.for.intercept = 0, 
    prior.scale.for.intercept = NULL, prior.df.for.intercept = 1, 
    min.prior.scale = 1e-12, scaled = FALSE, n.iter = 100, Warning = F)

Cure probability model:
               Estimate
(Intercept)  1.35020371
TRT         -0.56924133
SEX         -0.08326791
AGE          0.02021256


Failure time distribution model:
        Estimate
TRT -0.156862991
SEX  0.099745136
AGE -0.007646723


 Evaluation Metrics
      AUC         K         C 
0.6050845 0.5401964 0.5450531 
Program is running..be patient... done.
Call:
rcure(formula = Surv(FAILTIME, FAILCENS) ~ TRT + SEX + AGE, cureform = ~TRT + 
    SEX + AGE, data = e1684, model = "ph", Var = FALSE, method = "glm.fit", 
    prior.mean = 0, prior.scale = Inf, prior.df = Inf, prior.mean.for.intercept = 0, 
    prior.scale.for.intercept = Inf, prior.df.for.intercept = Inf, 
    Warning = F)

Cure probability model:
               Estimate
(Intercept)  1.36493298
TRT         -0.58847727
SEX         -0.08696490
AGE          0.02033857


Failure time distribution model:
        Estimate
TRT -0.153595097
SEX  0.099458470
AGE -0.007664013


 Evaluation Metrics
      AUC         K         C 
0.6050158 0.5399548 0.5448281 
Program is running..be patient... done.
Call:
smcure(formula = Surv(FAILTIME, FAILCENS) ~ TRT + SEX + AGE, 
    cureform = ~TRT + SEX + AGE, data = e1684, model = "ph", 
    Var = FALSE)

Cure probability model:
               Estimate
(Intercept)  1.36493298
TRT         -0.58847727
SEX         -0.08696490
AGE          0.02033857


Failure time distribution model:
        Estimate
TRT -0.153595097
SEX  0.099458470
AGE -0.007664013
Program is running..be patient... done.
Call:
rcure(formula = Surv(Time, Status) ~ TRT, cureform = ~TRT, data = bmt, 
    model = "aft", Var = FALSE, method = "glm.fit", prior.mean = 0, 
    prior.scale = NULL, prior.df = 1, prior.mean.for.intercept = 0, 
    prior.scale.for.intercept = NULL, prior.df.for.intercept = 1, 
    min.prior.scale = 1e-12, scaled = TRUE, n.iter = 100, Warning = F, 
    eva_model = "PO")

Cure probability model:
             Estimate
(Intercept) 1.0202106
TRT         0.3951315


Failure time distribution model:
              Estimate
(Intercept)  0.2101563
TRT         -0.3531250


 Evaluation Metrics
NULL

rcure documentation built on May 2, 2019, 7:01 a.m.

Related to rcure in rcure...