| Title: | Treatment Effect Estimation for Time-to-Event Data with Intercurrent Events |
|---|---|
| Description: | Analysis of treatment effects in clinical trials with time-to-event outcomes is complicated by intercurrent events. This package implements methods for estimating and inferring the cumulative incidence functions for time-to-event (TTE) outcomes with intercurrent events (ICE) under the five strategies outlined in the ICH E9 (R1) addendum, see Deng (2025) <doi:10.1002/sim.70091>. This package can be used for analyzing data from both randomized controlled trials and observational studies. In general, the data involve a primary outcome event and, potentially, an intercurrent event. Two data structures are allowed: competing risks, where only the time to the first event is recorded, and semicompeting risks, where the times to both the primary outcome event and intercurrent event (or censoring) are recorded. For estimation methods, users can choose nonparametric estimation (which does not use covariates) and semiparametrically efficient estimation. |
| Authors: | Yuhao Deng [aut], Yi Zhou [cre] |
| Maintainer: | Yi Zhou <[email protected]> |
| License: | GPL-3 |
| Version: | 1.1.4 |
| Built: | 2026-05-15 05:41:00 UTC |
| Source: | https://github.com/mephas/tteice |
This package aims to analyze treatment effects in clinical trials with time-to-event outcomes is complicated by intercurrent events. This package implements methods for estimating and inferring the cumulative incidence functions for time-to-event (TTE) outcomes with intercurrent events (ICE) under the five strategies outlined in the ICH E9 (R1) addendum, see Deng (2025) https://doi.org/10.1002/sim.70091. This package can be used for analyzing data from both randomized controlled trials and observational studies. In general, the data involve a primary outcome event and, potentially, an intercurrent event. Two data structures are allowed: competing risks, where only the time to the first event is recorded, and semicompeting risks, where the times to both the primary outcome event and intercurrent event (or censoring) are recorded. For estimation methods, nonparametric estimation (which does not use covariates) and semiparametrically efficient estimation are presented.
Main functions:
tteICE Using formula to fit cumulative incidence functions (CIFs) for competing/semicompeting risk time-to-event data with intercurrent events.
scr.tteICE Fit CIFs for semicompeting risk time-to-event data with intercurrent events.
surv.tteICE Fit CIFs for competing risk time-to-event with intercurrent events.
tteICEShiny Interactive Shiny app for the 'tteICE' package
Results output functions:
plot.tteICE Plot results from 'tteICE' objects.
print.tteICE Print a short summary of results from 'tteICE' objects
summary.tteICE Summarize results from 'tteICE' objects
predict.tteICE Predict risks for 'tteICE' objects at specific time
coef.tteICE Show the coefficient for 'tteICE' objects
bshaz.tteICE Extract the baseline hazards for 'tteICE' objects
zph.tteICE Perform a test for the proportional hazards assumption for the Cox models in 'tteICE' objects
Example data:
bmt Data from Section 1.3 of Klein and Moeschberger (1997)
Maintainer: Yi Zhou [email protected]
Authors:
Yuhao Deng [email protected]
Useful links:
Report bugs at https://github.com/mephas/tteICE/issues
The bmt data frame has 137 rows and 22 columns.
bmtbmt
This data frame contains the following columns:
Disease Group 1-ALL, 2-AML Low Risk, 3-AML High Risk
Time To Death Or On Study Time
Disease Free Survival Time (Time To Relapse, Death Or End Of Study)
Death Indicator 1-Dead 0-Alive
Relapse Indicator 1-Relapsed, 0-Disease Free
Disease Free Survival Indicator 1-Dead Or Relapsed, 0-Alive Disease Free)
Time To Acute Graft-Versus-Host Disease
Acute GVHD Indicator 1-Developed Acute GVHD 0-Never Developed Acute GVHD)
Time To Chronic Graft-Versus-Host Disease
Chronic GVHD Indicator 1-Developed Chronic GVHD 0-Never Developed Chronic GVHD
Time To Platelet Recovery
Platelet Recovery Indicator 1-Platelets Returned To Normal, 0-Platelets Never Returned to Normal
Patient Age In Years
Donor Age In Years
Patient Sex: 1-Male, 0-Female
Donor Sex: 1-Male, 0-Female
Patient CMV Status: 1-CMV Positive, 0-CMV Negative
Donor CMVStatus: 1-CMV Positive, 0-CMV Negative
Waiting Time to Transplant In Days
FAB: 1-FAB Grade 4 Or 5 and AML, 0-Otherwise
Hospital: 1-The Ohio State University, 2-Alferd, 3-St. Vincent, 4-Hahnemann
MTX used as a Graft-Versus-Host-Prophylactic: 1-Yes 0-No
Klein and Moeschberger (1997) Survival Analysis Techniques for Censored and Truncated Data, Springer.
data(bmt)data(bmt)
This function extracts the baseline cumulative hazards in the survival models
## S3 method for class 'tteICE' bshaz(x)## S3 method for class 'tteICE' bshaz(x)
x |
A fitted object returned by the function |
A data frame of baseline cumulative hazards in the working Kaplan-Meier or Cox models, stratified by treatment groups. The first column is time, the following columns are baseline cumulative hazards.
## load data data(bmt) bmt = transform(bmt, d4=d2+d3) A = as.numeric(bmt$group>1) X = as.matrix(bmt[,c('z1','z3','z5')]) bmt$A = A library(survival) fit = tteICE(Surv(t2, factor(d4))~A|z1+z3+z5, data=bmt, strategy="whileon", method='eff') bshaz(fit)## load data data(bmt) bmt = transform(bmt, d4=d2+d3) A = as.numeric(bmt$group>1) X = as.matrix(bmt[,c('z1','z3','z5')]) bmt$A = A library(survival) fit = tteICE(Surv(t2, factor(d4))~A|z1+z3+z5, data=bmt, strategy="whileon", method='eff') bshaz(fit)
This function extracts the coefficients in the Cox models
## S3 method for class 'tteICE' coef(object, ...)## S3 method for class 'tteICE' coef(object, ...)
object |
A fitted object returned by the function |
... |
Other arguments in function |
A list of coefficients of covariates in the working Cox models, stratified by treatment groups.
For the treatment policy strategy and composite variable strategy, only one Cox model is fit (for the primary
outcome event or the composite event). In these two strategies, coef1 is the coefficients in the treated
group, coef0 is the coefficients in the control group. For other strategies, Cox models are fitted for
each event (primary outcome event and intercurrent event). In these strategies, coef11 is the coefficients
for the primary outcome event in the treatment group, coef10 is the coefficients for the primary outcome
event in the control group, coef21 is the coefficients for the intercurrent event in the treated group,
coef20 is the coefficients for the intercurrent in the control group. If the nonparametric method is used,
the return is NULL.
## load data data(bmt) bmt = transform(bmt, d4=d2+d3) A = as.numeric(bmt$group>1) X = as.matrix(bmt[,c('z1','z3','z5')]) bmt$A = A library(survival) fit = tteICE(Surv(t2, factor(d4))~A|z1+z3+z5, data=bmt, strategy="whileon", method='eff') coef(fit)## load data data(bmt) bmt = transform(bmt, d4=d2+d3) A = as.numeric(bmt$group>1) X = as.matrix(bmt[,c('z1','z3','z5')]) bmt$A = A library(survival) fit = tteICE(Surv(t2, factor(d4))~A|z1+z3+z5, data=bmt, strategy="whileon", method='eff') coef(fit)
This function plots the estimated potential cumulative incidence functions or treatment effect curve with pointwise confidence intervals.
## S3 method for class 'tteICE' plot( x, type = c("ate", "inc")[1], decrease = FALSE, conf.int = 0.95, xlab = "Time", xlim = NULL, ylim = NULL, plot.configs = list(), ... )## S3 method for class 'tteICE' plot( x, type = c("ate", "inc")[1], decrease = FALSE, conf.int = 0.95, xlab = "Time", xlim = NULL, ylim = NULL, plot.configs = list(), ... )
x |
A fitted object returned by the function |
type |
Which plot to create: |
decrease |
|
conf.int |
#' Confidence level for the pointwise confidence intervals
If |
xlab |
Label for the x-axis. |
xlim |
A numeric vector of length 2 specifying the limits of the x-axis.
If |
ylim |
A numeric vector of length 2 giving the limits of the y-axis.
If |
plot.configs |
A named |
... |
Other arguments in function |
Plot the results from a tteICE object
plot_ate, plot_inc,
surv.tteICE, scr.tteICE,
tteICE
## load data data(bmt) bmt = transform(bmt, d4=d2+d3) A = as.numeric(bmt$group>1) bmt$A = A library(survival) ## simple model fitting and plotting fit1 = tteICE(Surv(t2, factor(d4))~A, data=bmt) plot(fit1, type="ate") plot(fit1, type="inc") ## plot cumulative incidence functions with p-values fit2 = surv.tteICE(A, bmt$t2, bmt$d4, "composite") plot(fit2, type="inc", decrease=TRUE, ylim=c(0,1), plot.configs=list(show.p.value=TRUE)) ## plot treatment effects for semicompeting risk data fit3 = scr.tteICE(A, bmt$t1, bmt$d1, bmt$t2, bmt$d2, "composite") plot(fit3, type="ate", ylim=c(-1,1), xlab="time", plot.configs=list(col="red"))## load data data(bmt) bmt = transform(bmt, d4=d2+d3) A = as.numeric(bmt$group>1) bmt$A = A library(survival) ## simple model fitting and plotting fit1 = tteICE(Surv(t2, factor(d4))~A, data=bmt) plot(fit1, type="ate") plot(fit1, type="inc") ## plot cumulative incidence functions with p-values fit2 = surv.tteICE(A, bmt$t2, bmt$d4, "composite") plot(fit2, type="inc", decrease=TRUE, ylim=c(0,1), plot.configs=list(show.p.value=TRUE)) ## plot treatment effects for semicompeting risk data fit3 = scr.tteICE(A, bmt$t1, bmt$d1, bmt$t2, bmt$d2, "composite") plot(fit3, type="ate", ylim=c(-1,1), xlab="time", plot.configs=list(col="red"))
This function predicts the potential cumulative incidence function and treatment effect at specific time points.
## S3 method for class 'tteICE' predict(object, timeset = NULL, ...)## S3 method for class 'tteICE' predict(object, timeset = NULL, ...)
object |
A fitted object returned by the function |
timeset |
Time at which to predict the risk.
If |
... |
Other arguments in function |
A matrix with each row being time points, potential cumulative incidences (under treated and under control), treatment effects, standard errors, and P-values.
predict a tteICE object. The meanings of each row are: time points, potential cumulative incidences (under treated and under control), treatment effects, standard errors, and P-values.
scr.tteICE, surv.tteICE, tteICE
surv.boot
surv.tteICE, scr.tteICE,
tteICE
## load data data(bmt) bmt = transform(bmt, d4=d2+d3) A = as.numeric(bmt$group>1) bmt$A = A X = as.matrix(bmt[,c('z1','z3','z5')]) ## predict results at specified time points ## model fitting using semicompeting risk data fit1 = scr.tteICE(A, bmt$t1, bmt$d1, bmt$t2, bmt$d2, "composite") predict(fit1, timeset=c(670,2000)) ## predict results without specifying any time points ## model fitting using competing risk data fit2 = surv.tteICE(A, bmt$t2, bmt$d4, "composite") predict(fit2) ## a simpler way library(survival) fit3 = tteICE(Surv(t2, factor(d4))~A|z1+z3+z5, data=bmt, strategy="composite", method='eff') predict(fit3, timeset=c(670,2000)) predict(fit3)## load data data(bmt) bmt = transform(bmt, d4=d2+d3) A = as.numeric(bmt$group>1) bmt$A = A X = as.matrix(bmt[,c('z1','z3','z5')]) ## predict results at specified time points ## model fitting using semicompeting risk data fit1 = scr.tteICE(A, bmt$t1, bmt$d1, bmt$t2, bmt$d2, "composite") predict(fit1, timeset=c(670,2000)) ## predict results without specifying any time points ## model fitting using competing risk data fit2 = surv.tteICE(A, bmt$t2, bmt$d4, "composite") predict(fit2) ## a simpler way library(survival) fit3 = tteICE(Surv(t2, factor(d4))~A|z1+z3+z5, data=bmt, strategy="composite", method='eff') predict(fit3, timeset=c(670,2000)) predict(fit3)
Print the summary of 'tteICE'
## S3 method for class 'summary.tteICE' print(x, digits = 3, ...)## S3 method for class 'summary.tteICE' print(x, digits = 3, ...)
x |
A fitted object returned by the function |
digits |
The digits of the results |
... |
Other arguments in function |
Print the summary of a tteICE object
This function summarizes the results
## S3 method for class 'tteICE' print(x, digits = 3, ...)## S3 method for class 'tteICE' print(x, digits = 3, ...)
x |
A fitted object returned by the function |
digits |
The digits of the results |
... |
Other arguments in function |
Print the summary of a tteICE object
surv.tteICE,
scr.tteICE,
tteICE
## load data data(bmt) bmt = transform(bmt, d4=d2+d3) A = as.numeric(bmt$group>1) bmt$A = A ## print the results fit1 = surv.tteICE(A, bmt$t2, bmt$d4, "composite") print(fit1) fit2 = scr.tteICE(A, bmt$t1, bmt$d1, bmt$t2, bmt$d2, "composite") print(fit2) library(survival) fit3 = tteICE(Surv(t2, factor(d4))~A|z1+z3+z5, data=bmt, strategy="composite", method='eff') print(fit3, digits=4)## load data data(bmt) bmt = transform(bmt, d4=d2+d3) A = as.numeric(bmt$group>1) bmt$A = A ## print the results fit1 = surv.tteICE(A, bmt$t2, bmt$d4, "composite") print(fit1) fit2 = scr.tteICE(A, bmt$t1, bmt$d1, bmt$t2, bmt$d2, "composite") print(fit2) library(survival) fit3 = tteICE(Surv(t2, factor(d4))~A|z1+z3+z5, data=bmt, strategy="composite", method='eff') print(fit3, digits=4)
This function estimates the potential cumulative incidence function for time-to event data under ICH E9 (R1) to address intercurrent events. The input data should be of a semicompeting risks structure.
scr.tteICE( A, Time, status, Time_int, status_int, strategy = "composite", cov1 = NULL, method = "np", weights = NULL, subset = NULL, na.rm = FALSE, nboot = 0, seed = 0 )scr.tteICE( A, Time, status, Time_int, status_int, strategy = "composite", cov1 = NULL, method = "np", weights = NULL, subset = NULL, na.rm = FALSE, nboot = 0, seed = 0 )
A |
Treatment indicator, 1 for treatment and 0 for control. |
Time |
Time to the primary (terminal) event. |
status |
Indicator of the primary (terminal) event, 1 for event and 0 for censoring. |
Time_int |
Time to the intercurrent event. |
status_int |
Indicator of the intercurrent event, 1 for event and 0 for censoring. |
strategy |
Strategy to address intercurrent events, |
cov1 |
Baseline covariates. |
method |
Estimation method, |
weights |
Weight for each subject. |
subset |
Subset, either numerical or logical. |
na.rm |
Whether to remove missing values. |
nboot |
Number of resamplings in the boostrapping method. If |
seed |
Seed for bootstrapping. |
Intercurrent events refer to the events occurring after treatment initiation of clinical trials that affect either the interpretation of or the existence of the measurements associated with the clinical question of interest. The International Conference on Harmonization (ICH) E9 (R1) addendum proposed five strategies to address intercurrent events, namely, treatment policy strategy, composite variable strategy, while on treatment strategy, hypothetical strategy, and principal stratum strategy. To answer a specific scientific question, a strategy with a particular estimand is chosen before the study design.
We adopt the potential outcomes framework that defines a causal estimand as the contrast between
functionals of potential outcomes. Consider a randomized controlled trial with individuals
randomly assigned to one of two treatment conditions, denoted by , where represents
the active treatment (a test drug) and represents the control (placebo). Assume that all
patients adhere to their treatment assignments and do not discontinue treatment. Associated with individual
are two potential time-to-event primary outcomes and ,
if any, which represent the time durations from treatment initiation to the primary outcome event under
two treatment assignments respectively. Let and denote the occurrence time of
potential intercurrent events, if any, under the two treatment assignments, respectively. Intercurrent
events are considered as absent if no post-treatment intercurrent events occur until the end of study.
We adopt the potential cumulative incidences under both treatment assignments as the target estimands. Potential cumulative incidences describe the probability of time-to-event outcomes occurring at each time point. We define the treatment effect as the contrast of two potential cumulative incidences. Cumulative incidences are model-free and collapsible, enjoying causal interpretations.
A list including the fitted object and input variables.
Deng, Y., Han, S., & Zhou, X. H. (2025). Inference for Cumulative Incidences and Treatment Effects in Randomized Controlled Trials With Time-to-Event Outcomes Under ICH E9 (R1). Statistics in Medicine. doi:10.1002/sim.70091
## load data data(bmt) bmt = transform(bmt, d4=d2+d3) A = as.numeric(bmt$group>1) X = as.matrix(bmt[,c('z1','z3','z5')]) ## Composite variable strategy, ## nonparametric estimation without covariates fit1 = scr.tteICE(A, bmt$t1, bmt$d1, bmt$t2, bmt$d2, "composite") ## Hypothetical strategy (natural effects), ## nonparametric estimation with inverse probability weighting fit2 = scr.tteICE(A, bmt$t1, bmt$d1, bmt$t2, bmt$d2, "natural", X, method='ipw') ## nonparametric estimation with weights as non-standardized inverse probability score ps = predict(glm(A ~ X, family='binomial'), type='response') w = A/ps + (1-A)/(1-ps) fit2 = scr.tteICE(A, bmt$t1, bmt$d1, bmt$t2, bmt$d2, "natural", weights=w) ## Hypothetical strategy (removing intercurrent events), ## semiparametrically efficient estimation with covariates fit3 = scr.tteICE(A, bmt$t1, bmt$d1, bmt$t2, bmt$d2, "removed", X, method='eff')## load data data(bmt) bmt = transform(bmt, d4=d2+d3) A = as.numeric(bmt$group>1) X = as.matrix(bmt[,c('z1','z3','z5')]) ## Composite variable strategy, ## nonparametric estimation without covariates fit1 = scr.tteICE(A, bmt$t1, bmt$d1, bmt$t2, bmt$d2, "composite") ## Hypothetical strategy (natural effects), ## nonparametric estimation with inverse probability weighting fit2 = scr.tteICE(A, bmt$t1, bmt$d1, bmt$t2, bmt$d2, "natural", X, method='ipw') ## nonparametric estimation with weights as non-standardized inverse probability score ps = predict(glm(A ~ X, family='binomial'), type='response') w = A/ps + (1-A)/(1-ps) fit2 = scr.tteICE(A, bmt$t1, bmt$d1, bmt$t2, bmt$d2, "natural", weights=w) ## Hypothetical strategy (removing intercurrent events), ## semiparametrically efficient estimation with covariates fit3 = scr.tteICE(A, bmt$t1, bmt$d1, bmt$t2, bmt$d2, "removed", X, method='eff')
This function summarizes the results
## S3 method for class 'tteICE' summary(object, digits = 3, ...)## S3 method for class 'tteICE' summary(object, digits = 3, ...)
object |
A fitted object returned by the function |
digits |
The digits of the results |
... |
Other arguments in function |
A list that consists of summaries of a tteICE object: data type, strategy, estimation method, maximum follow-up time, sample size, treated sample size, controlled sample size, p-value, and predicted risks at quartiles
surv.tteICE,
scr.tteICE,
tteICE,
print.tteICE
## load data data(bmt) bmt = transform(bmt, d4=d2+d3) A = as.numeric(bmt$group>1) bmt$A = A X = as.matrix(bmt[,c('z1','z3','z5')]) ## Composite variable strategy, ## nonparametric estimation without covariates fit1 = scr.tteICE(A, bmt$t1, bmt$d1, bmt$t2, bmt$d2, "composite") summary(fit1) fit2 = surv.tteICE(A, bmt$t2, bmt$d4, "composite") predict(fit2) library(survival) fit3 = tteICE(Surv(t2, factor(d4))~A|z1+z3+z5, data=bmt, strategy="composite", method='eff') summary(fit3)## load data data(bmt) bmt = transform(bmt, d4=d2+d3) A = as.numeric(bmt$group>1) bmt$A = A X = as.matrix(bmt[,c('z1','z3','z5')]) ## Composite variable strategy, ## nonparametric estimation without covariates fit1 = scr.tteICE(A, bmt$t1, bmt$d1, bmt$t2, bmt$d2, "composite") summary(fit1) fit2 = surv.tteICE(A, bmt$t2, bmt$d4, "composite") predict(fit2) library(survival) fit3 = tteICE(Surv(t2, factor(d4))~A|z1+z3+z5, data=bmt, strategy="composite", method='eff') summary(fit3)
This function estimates the hazard ratio for time-to event data under ICH E9 (R1) to address intercurrent events. Multiple strategies except the principal stratum strategy are allowed.
surv.HR( A, Time, cstatus, strategy = "composite", cov1 = NULL, conf.int = 0.95, weights = NULL, subset = NULL )surv.HR( A, Time, cstatus, strategy = "composite", cov1 = NULL, conf.int = 0.95, weights = NULL, subset = NULL )
A |
Treatment indicator, 1 for treatment and 0 for control. |
Time |
Time to event. |
cstatus |
Indicator of event, 1 for the primary event, 2 for the intercurrent event, 0 for censoring. |
strategy |
Strategy to address intercurrent events, |
cov1 |
Baseline covariates. |
conf.int |
Level of the confidence interval. |
weights |
Weight for each subject (not applied to the while on treatment strategy). |
subset |
Subset, either numerical or logical. |
For the treatment policy and hypothetical strategies, the hazard ratio (HR) is given by the Cox
regression regarding intercurrent events as censoring. For the composite variable strategy, the
hazard ratio is given by the Cox regression regarding the first occurrence of either intercurrent
event or primary event as the event of interest. For the while on treatment strategy, the hazard
ratio is given by the Fine-Gray subdistribution model. There is no existing method to estimate the
hazard ratio using principal stratum strategy.
The weakness of using hazard ratio to infer treatment effects is critical. First, the hazard ratio
relies on model specification. Second, the hazard ratio is not collapsible. Therefore, the hazard
ratio should only be treated as a descriptive or exploratory measure of the treatment effect.
A list including
Estimated log hazard ratio (logHR) of the treatment effect on the primary event.
Standard error of the estimated log hazard ratio (logHR).
Confidence interval of the hazard ratio (HR).
P value of the hazard ratio.
data(bmt) bmt = transform(bmt, d4=d2+d3) A = as.numeric(bmt$group>1) ## composite variable strategy fit = surv.HR(A, bmt$t2, bmt$d4, "composite") ## while on treatment strategy X = bmt[,c('z1','z3','z5')] fit = surv.HR(A, bmt$t2, bmt$d4, "whileon", cov1=X)data(bmt) bmt = transform(bmt, d4=d2+d3) A = as.numeric(bmt$group>1) ## composite variable strategy fit = surv.HR(A, bmt$t2, bmt$d4, "composite") ## while on treatment strategy X = bmt[,c('z1','z3','z5')] fit = surv.HR(A, bmt$t2, bmt$d4, "whileon", cov1=X)
This function estimates the potential cumulative incidence function for time-to event data under ICH E9 (R1) to address intercurrent events. The input data should be of a competing risks structure.
surv.tteICE( A, Time, cstatus, strategy = "composite", cov1 = NULL, method = "np", weights = NULL, subset = NULL, na.rm = FALSE, nboot = 0, seed = 0 )surv.tteICE( A, Time, cstatus, strategy = "composite", cov1 = NULL, method = "np", weights = NULL, subset = NULL, na.rm = FALSE, nboot = 0, seed = 0 )
A |
Treatment indicator, 1 for treatment and 0 for control. |
Time |
Time to event. |
cstatus |
Indicator of event, 1 for the primary event, 2 for the intercurrent event, 0 for censoring. |
strategy |
Strategy to address intercurrent events, |
cov1 |
Baseline covariates. |
method |
Estimation method, |
weights |
Weight for each subject. |
subset |
Subset, either numerical or logical. |
na.rm |
Whether to remove missing values. |
nboot |
Number of resamplings in the boostrapping method. If |
seed |
Seed for bootstrapping. |
Intercurrent events refer to the events occurring after treatment initiation of clinical trials that affect either the interpretation of or the existence of the measurements associated with the clinical question of interest. The International Conference on Harmonization (ICH) E9 (R1) addendum proposed five strategies to address intercurrent events, namely, treatment policy strategy, composite variable strategy, while on treatment strategy, hypothetical strategy, and principal stratum strategy. To answer a specific scientific question, a strategy with a particular estimand is chosen before the study design.
We adopt the potential outcomes framework that defines a causal estimand as the contrast between
functionals of potential outcomes. Consider a randomized controlled trial with individuals
randomly assigned to one of two treatment conditions, denoted by , where represents
the active treatment (a test drug) and represents the control (placebo). Assume that all
patients adhere to their treatment assignments and do not discontinue treatment. Associated with individual
are two potential time-to-event primary outcomes and ,
if any, which represent the time durations from treatment initiation to the primary outcome event under
two treatment assignments respectively. Let and denote the occurrence time of
potential intercurrent events, if any, under the two treatment assignments, respectively. Intercurrent
events are considered as absent if no post-treatment intercurrent events occur until the end of study.
We adopt the potential cumulative incidences under both treatment assignments as the target estimands. Potential cumulative incidences describe the probability of time-to-event outcomes occurring at each time point. We define the treatment effect as the contrast of two potential cumulative incidences. Cumulative incidences are model-free and collapsible, enjoying causal interpretations.
A list including the fitted object and input variables.
Deng, Y., Han, S., & Zhou, X. H. (2025). Inference for Cumulative Incidences and Treatment Effects in Randomized Controlled Trials With Time-to-Event Outcomes Under ICH E9 (R1). Statistics in Medicine. doi:10.1002/sim.70091
## load data data(bmt) bmt = transform(bmt, d4=d2+d3) A = as.numeric(bmt$group>1) X = as.matrix(bmt[,c('z1','z3','z5')]) ## Composite variable strategy, ## nonparametric estimation without covariates fit1 = surv.tteICE(A, bmt$t2, bmt$d4, "composite") ## Hypothetical strategy (natural effects), ## nonparametric estimation with inverse probability weighting fit2 = surv.tteICE(A, bmt$t2, bmt$d4, "natural", X, method='ipw') ## nonparametric estimation with weights as inverse propensity score ps = predict(glm(A ~ X, family='binomial'), type='response') w = A/ps + (1-A)/(1-ps) fit2 = surv.tteICE(A, bmt$t2, bmt$d4, "natural", weights=w) ## Hypothetical strategy (removing intercurrent events), ## semiparametrically efficient estimation with covariates fit3 = surv.tteICE(A, bmt$t2, bmt$d4, "removed", X, method='eff')## load data data(bmt) bmt = transform(bmt, d4=d2+d3) A = as.numeric(bmt$group>1) X = as.matrix(bmt[,c('z1','z3','z5')]) ## Composite variable strategy, ## nonparametric estimation without covariates fit1 = surv.tteICE(A, bmt$t2, bmt$d4, "composite") ## Hypothetical strategy (natural effects), ## nonparametric estimation with inverse probability weighting fit2 = surv.tteICE(A, bmt$t2, bmt$d4, "natural", X, method='ipw') ## nonparametric estimation with weights as inverse propensity score ps = predict(glm(A ~ X, family='binomial'), type='response') w = A/ps + (1-A)/(1-ps) fit2 = surv.tteICE(A, bmt$t2, bmt$d4, "natural", weights=w) ## Hypothetical strategy (removing intercurrent events), ## semiparametrically efficient estimation with covariates fit3 = surv.tteICE(A, bmt$t2, bmt$d4, "removed", X, method='eff')
This function estimates the potential cumulative incidence function for time-to event data under ICH E9 (R1) to address intercurrent events. The input data can be competing or semicompeting risks data structure.
tteICE( formula, add.scr = NULL, data, strategy = "composite", method = "np", weights = NULL, subset = NULL, na.rm = FALSE, nboot = 0, seed = 0 )tteICE( formula, add.scr = NULL, data, strategy = "composite", method = "np", weights = NULL, subset = NULL, na.rm = FALSE, nboot = 0, seed = 0 )
formula |
An object of class "formula" (or one that can be coerced to that class).
A symbolic description of the model to be fitted.
For example, |
add.scr |
Required for semicompeting data.
An object of class "Surv" (or one that can be coerced to that class).
For example, |
data |
Data or object coercible by as.data.frame to a data frame, containing the variables in the model. |
strategy |
Strategy to address intercurrent events, |
method |
Estimation method, |
weights |
Weight for each subject. |
subset |
Subset, either numerical or logical. |
na.rm |
Whether to remove missing values. |
nboot |
Number of resamplings in the boostrapping method. If |
seed |
Seed for bootstrapping. |
Intercurrent events refer to the events occurring after treatment initiation of clinical trials that affect either the interpretation of or the existence of the measurements associated with the clinical question of interest. The International Conference on Harmonization (ICH) E9 (R1) addendum proposed five strategies to address intercurrent events, namely, treatment policy strategy, composite variable strategy, while on treatment strategy, hypothetical strategy, and principal stratum strategy. To answer a specific scientific question, a strategy with a particular estimand is chosen before the study design.
We adopt the potential outcomes framework that defines a causal estimand as the contrast between
functionals of potential outcomes. Consider a randomized controlled trial with individuals
randomly assigned to one of two treatment conditions, denoted by , where represents
the active treatment (a test drug) and represents the control (placebo). Assume that all
patients adhere to their treatment assignments and do not discontinue treatment. Associated with individual
are two potential time-to-event primary outcomes and ,
if any, which represent the time durations from treatment initiation to the primary outcome event under
two treatment assignments respectively. Let and denote the occurrence time of
potential intercurrent events, if any, under the two treatment assignments, respectively. Intercurrent
events are considered as absent if no post-treatment intercurrent events occur until the end of study.
We adopt the potential cumulative incidences under both treatment assignments as the target estimands. Potential cumulative incidences describe the probability of time-to-event outcomes occurring at each time point. We define the treatment effect as the contrast of two potential cumulative incidences. Cumulative incidences are model-free and collapsible, enjoying causal interpretations.
The formula should be set in the following two ways.
When data take format of competing risk data, set the first argument formula = Surv(time, status) ~ treatment | covariate1+covariate2
or formula = Surv(time, status)~ A without any baseline covariates, where status is a factor variable with levels 0,1,2
(1 for the primary event, 2 for the intercurrent event, and 0 for censoring).
When data take the format of semicompeting risk data, set the first argument formula = Surv(time, status) ~ treatment | covariate1+covariate2
or formula = Surv(time, status) ~ A without any baseline covariates, where status=0,1 (1 for the primary event and 0 for censoring).
In addition, the second argument add.scr = ~ Surv(time.intercurrent, status.intercurrent) is required.
A list including the fitted object and input variables.
Deng, Y., Han, S., & Zhou, X. H. (2025). Inference for Cumulative Incidences and Treatment Effects in Randomized Controlled Trials With Time-to-Event Outcomes Under ICH E9 (R1). Statistics in Medicine. doi:10.1002/sim.70091
## load data data(bmt) bmt = transform(bmt, d4=d2+d3) A = as.numeric(bmt$group>1) X = as.matrix(bmt[,c('z1','z3','z5')]) bmt$A = A library(survival) ## Composite variable strategy, ## nonparametric estimation without covariates ## Composite variable strategy, ## nonparametric estimation without covariates ## model fitting for competing risk data without covariates fit1 = tteICE(Surv(t2, factor(d4)) ~ A, data=bmt, strategy="composite", method='np') print(fit1) ## model fitting for competing risk data without covariates ## with bootstrap confidence intervals fit.bt1 = tteICE(Surv(t2, factor(d4)) ~ A, data=bmt, strategy="composite", method='eff', nboot=20, seed=2) print(fit.bt1) ## model fitting for competing risk data with covariates fit2 = tteICE(Surv(t2, factor(d4)) ~ A | z1 + z3 + z5, data=bmt, strategy="composite", method='eff') print(fit2) ## model fitting for semicompeting risk data without covariates fitscr1 = tteICE(Surv(t1, d1) ~ A, ~Surv(t2, d2), data=bmt, strategy="composite", method='np') print(fitscr1) ## model fitting for semicompeting risk data without covariates fitscr2 = tteICE(Surv(t1, d1) ~ A | z1 + z3 + z5, ~Surv(t2, d2), data=bmt, strategy="composite", method='eff') print(fitscr2)## load data data(bmt) bmt = transform(bmt, d4=d2+d3) A = as.numeric(bmt$group>1) X = as.matrix(bmt[,c('z1','z3','z5')]) bmt$A = A library(survival) ## Composite variable strategy, ## nonparametric estimation without covariates ## Composite variable strategy, ## nonparametric estimation without covariates ## model fitting for competing risk data without covariates fit1 = tteICE(Surv(t2, factor(d4)) ~ A, data=bmt, strategy="composite", method='np') print(fit1) ## model fitting for competing risk data without covariates ## with bootstrap confidence intervals fit.bt1 = tteICE(Surv(t2, factor(d4)) ~ A, data=bmt, strategy="composite", method='eff', nboot=20, seed=2) print(fit.bt1) ## model fitting for competing risk data with covariates fit2 = tteICE(Surv(t2, factor(d4)) ~ A | z1 + z3 + z5, data=bmt, strategy="composite", method='eff') print(fit2) ## model fitting for semicompeting risk data without covariates fitscr1 = tteICE(Surv(t1, d1) ~ A, ~Surv(t2, d2), data=bmt, strategy="composite", method='np') print(fitscr1) ## model fitting for semicompeting risk data without covariates fitscr2 = tteICE(Surv(t1, d1) ~ A | z1 + z3 + z5, ~Surv(t2, d2), data=bmt, strategy="composite", method='eff') print(fitscr2)
This function opens the RShiny app for tteICE. RShiny application can be used for generating plots and basic analysis results. It provides a point-and-click interface, so users can obtain results without writing R code directly.
tteICEShiny()tteICEShiny()
Rshiny interface
if(interactive() && requireNamespace("shiny", quietly = TRUE)){ tteICEShiny() }if(interactive() && requireNamespace("shiny", quietly = TRUE)){ tteICEShiny() }
This function checks the proportional hazards assumption in the Cox model using Schoenfeld residuals. This function only return results for strategies based on efficient influence functions.
## S3 method for class 'tteICE' zph(x)## S3 method for class 'tteICE' zph(x)
x |
A fitted object returned by the function |
A list of P-values of testing the proportional hazards (PH) assumption in the working Cox models, for each
covariate and a global test, stratified by treatment groups.
For the treatment policy strategy and composite variable strategy, only one Cox model is fit (for the primary
outcome event or the composite event). In these two strategies, ph1 is the P-values in the treated
group, ph0 is the P-values in the control group. For other strategies, Cox models are fitted for
each event (primary outcome event and intercurrent event). In these strategies, ph11 is the P-values
for the primary outcome event in the treatment group, ph10 is the P-values for the primary outcome
event in the control group, ph21 is the P-values for the intercurrent event in the treated group,
ph20 is the P-values for the intercurrent in the control group. If the nonparametric method is used,
the return is NULL.
## load data data(bmt) bmt = transform(bmt, d4=d2+d3) A = as.numeric(bmt$group>1) X = as.matrix(bmt[,c('z1','z3','z5')]) bmt$A = A library(survival) fit = tteICE(Surv(t2, factor(d4))~A|z1+z3+z5, data=bmt, strategy="whileon", method='eff') print(fit$ph) zph(fit) plot(fit$ph$ph11) plot(fit$ph$ph10)## load data data(bmt) bmt = transform(bmt, d4=d2+d3) A = as.numeric(bmt$group>1) X = as.matrix(bmt[,c('z1','z3','z5')]) bmt$A = A library(survival) fit = tteICE(Surv(t2, factor(d4))~A|z1+z3+z5, data=bmt, strategy="whileon", method='eff') print(fit$ph) zph(fit) plot(fit$ph$ph11) plot(fit$ph$ph10)