geecure {geecure}

Marginal proportional hazards mixture cure model with generalzied estimating equations

Description

Fit the marginal proportional hazards mixture cure (PHMC) model with the generalized estimating equations (GEE). GEE approach is generalized to the marginal PHMC model through the expectation-solution (ES) algorithm to account for the correlation among the cure statuses and the dependence among the failure times of uncured patients to improve the estimation efficiency.

Usage

geecure(formula, cureform, data, id, model = c("para", "semi"),  
        corstr = c("independence", "exchangeable"), Var = TRUE, stdz = FALSE, 
        boots = FALSE, nboot = 100, esmax = 100, eps = 1e-06)

Arguments

formula

a formula expression, of the form response ~ predictors. The response is a Surv object with right censoring. It is used to specify the covariate effects on the failure time of uncured subjects. See the documentation for survreg, Surv for details. The expression to the right of the "~" specifies the effect of covariates on the failure time of uncured patients.

cureform

a formula expression, of the form cureform ~ predictors. It is used to specify the effects of covariates on the cure rate. A covariate may be used in both formula and cureform.

data

a data frame in which to interpret the variables named in the formula and the cureform.

id

a vector which identifies the clusters. The length of id should be the same as the number of observations.

model

specifies your model, it can be para which represents the parametric PHMC model with two-parameter Weibull baseline survival function, or semi which represents the semiparametric PHMC model.

corstr

a character string specifying the correlation structure. The following are permitted: independence and exchangeable.

Var

If it is TRUE, the program returns Std.Error by the sandwich method. By default, Var = TRUE.

stdz

If it is TRUE, all the covariates in the formula and cureform are standardized. By default, stdz = FALSE.

boots

If it is TRUE, the program returns Std.Error by the bootstrap method. By default, boots = FALSE.

nboot

the number of bootstrap samples. The default is nboot = 100.

esmax

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

eps

tolerance for convergence. The default is eps = 1e-6. Iteration stops once the relative change in deviance is less than eps.

Details

The marginal PHMC model is considered in this function. For cure rate, a logistic regression model is employed and the probability of being cured is given by (1+\exp(\gamma Z))^{(-1)}. For uncured subject, the failure time is modeled by either the parametric proportional hazards model with Weibull baseline distributions or the semiparametric proportional hazards model. A covariate can be used either in formula or in cureform or in both. The model parameters are estimated by the expectation-solution (ES) algorithm and the standard error estimates are obtained from sandwich variance formula based on Niu and Peng (2014) and Niu et al. (2018).

Value

An object of class geecure is returned. It can be examined by print.geecure().

References

Niu, Y. and Peng, Y. (2014) Marginal regression analysis of clustered failure time data with a cure fraction. Journal of Multivariate Analysis, 123: 129-142.

Niu, Y., Song, L., Liu, Y, and Peng, Y. (2018) Modeling clustered long-term survivors using marginal mixture cure model. Biometrical Journal, doi: 10.1002/bjmj.201700114.

Examples

# Be patient, the following examples may take several minites on a faster computer.
# Example 1. Fit the marginal parametric PHMC model for the smoking cessation data.
data(smoking)
smoking$Time <- ifelse(smoking$Relapse == 0, smoking$Timept1, 
                       (smoking$Timept1 + smoking$Timept2)/2)

# plot the KM survival curve of smoking cessation data
plot(survfit(Surv(Time, Relapse) ~ SexF + (SI.UC), data = smoking), 
      ylab = "Survival function", xlab = "Years", ylim = c(0.5, 1), 
      xlim = c(0, 6), lty = 1:4, col = 1:4)
legend(0.5, 0.63, c("SI/Male", "SI/Female", "UC/Female", "UC/Male"), cex = 1, 
        col = c(2, 4, 3, 1), lty = c(2, 4, 3, 1))

geesmokingind <- geecure(Surv(Time, Relapse) ~ SexF + Duration + SI.UC + F10Cigs + 
      SexF * SI.UC, cureform = ~ SexF + Duration + SI.UC + F10Cigs + SexF * SI.UC, 
      data = smoking,  model = "para", id = smoking$Zip, corstr = "independence")

geesmokingexch <- geecure(Surv(Time, Relapse) ~ SexF + Duration + SI.UC + F10Cigs + 
      SexF * SI.UC, cureform = ~ SexF + Duration + SI.UC + F10Cigs + SexF * SI.UC, 
      data = smoking,  model = "para", id = smoking$Zip, corstr = "exchangeable")

# Example 2. Fit the marginal semiparametric PHMC model for the bmt data.
data(bmt)
bmt$g <- factor(bmt$g, label = c("ALL", "AML low risk", "AML high risk"))
bmt$Z8 <- factor(bmt$Z8, label = c("Otherwise", "FAB"))

# plot the KM survival curve of bmt data 
plot(survfit(Surv(T2, d3) ~ 1, data = bmt), xlab = "Days", ylab = "Survival Probability", 
     cex.lab = 1.7, cex.axis = 2, cex.main = 1.7, mark.time = TRUE)

geebmtind <- geecure(Surv(T2, d3) ~ factor(g) + Z8, cureform = ~ factor(g) + Z8, 
             data = bmt, model = "semi", id = bmt$Z9, corstr= "independence")           

geebmtexch <- geecure(Surv(T2, d3) ~ factor(g) + Z8, cureform = ~ factor(g) + Z8, 
              data = bmt, model = "semi", id = bmt$Z9, corstr= "exchangeable", 
              stdz = TRUE, boots = TRUE)

# Example 3. Fit the marginal semiparametric PHMC model for the tonsil data.
data(tonsil)
tonsil<-tonsil[-c(141,136,159),]
tonsil$Sex <- ifelse(tonsil$Sex == 1, 0, 1)
tonsil$Cond <- ifelse(tonsil$Cond == 1, 0, 1)
tonsil$T <- ifelse(tonsil$T < 4, 0, 1)

# plot the KM survival curve of tonsil data 
plot(survfit(Surv(Time, Status) ~ 1, data = tonsil), xlab = "Days", ylab = "Survival 
     Probability", cex.lab = 1.7, cex.axis = 2, cex.main = 1.7, mark.time = TRUE)

geetonsilind <- geecure2(Surv(Time, Status) ~ Sex + factor(Grade) + Age + Cond + T, 
                cureform = ~ Sex + factor(Grade) + Age + Cond + T, data = tonsil, 
                id = tonsil$Inst, corstr = "independence")

geetonsilexch <- geecure2(Surv(Time, Status) ~ Sex + factor(Grade) +Age + Cond + T, 
                 cureform = ~ Sex + factor(Grade) + Age + Cond + T, data = tonsil, 
                id = tonsil$Inst, corstr = "exchangeable",  stdz = TRUE, Var = FALSE)