Comparing CIFs
Let’s demonstrate how to analyze competing risks using the cmprsk package in R. This package includes functions for estimating the cumulative incidence function (CIF) in the presence of competing risks.
Bone marrow transplant (BMT) data from the cmprsk package includes information for 35 leukemia cancer patients who underwent bone marrow transplantation. The data includes the following variables:
- dis: disease; 0 = ALL; 1 = AML (ALL, Acute lymphoblastic leukemia; AML, Acute myeloid leukemia)
- ftime: follow-up time
- status: 0 = censored (survival); 1 = relapse; 2 = Transplant-related mortality;
## 'data.frame': 35 obs. of 3 variables:
## $ dis : int 0 0 0 0 0 1 0 0 1 1 ...
## $ ftime : int 13 1 72 7 8 67 9 5 70 4 ...
## $ status: int 2 1 0 2 2 0 2 2 0 0 ...
## 'data.frame': 35 obs. of 3 variables:
## $ dis : Factor w/ 2 levels "ALL","AML": 1 1 1 1 1 2 1 1 2 2 ...
## $ ftime : int 13 1 72 7 8 67 9 5 70 4 ...
## $ status: Factor w/ 3 levels "Censored","Mortality",..: 3 2 1 3 3 1 3 3 1 1 ...
## Estimates and Variances:
## $est
## 20 40 60
## 1 Censored 0.1142857 0.2000000 0.2000000
## 1 Mortality 0.2571429 0.2571429 0.2571429
## 1 Relapse 0.3714286 0.4285714 0.4285714
##
## $var
## 20 40 60
## 1 Censored 0.002998057 0.004910750 0.004910750
## 1 Mortality 0.005636800 0.005636800 0.005636800
## 1 Relapse 0.006878377 0.007269181 0.007269181
## stat pv df
## Censored 5.9785107 0.014481226 1
## Mortality 0.9133497 0.339227192 1
## Relapse 9.4874094 0.002068867 1
Gray’s test for equality of CIFs across groups (ALL vs. AML) is split across different types of events, here (Censored, Mortality, and Relapse).
- For Censored, the p-value is small and suggest that there isa statistically significant difference in the CIF for censored events between the groups.
- For Mortality, the p-value is greater than 0.05, suggesting that there is no statistically significant difference in the CIF for mortality between the groups.
- For Relapse, the p-value is less than 0.05, indicating a statistically significant difference in the CIF for relapse between the groups.
- The above suggest that the groups are behaving differently in terms of the time to relapse and censored events, but not for mortality.
Competing Risks Regression
adopted from https://www.nature.com/articles/bmt2009359
Suppose that the BMT study was extended to include more participants and additional covariates and now includes:
- 177 observations
- Sex: gender of the individual
- D: disease; 0 = ALL; 1 = AML
- Phase: phase at transplant (Relapse, CR1, CR2, CR3)
- Age: age at the beginning of follow-up
- Status: 0 = censored; 1 = Transplant-related mortality; 2 = relapse
- Source: source of stem cells (BM+PB, PB)
- ftime: failure time
We are interested in modeling time to relapse in the presence of transplant-related death (competing event). We want to sutdy the effect on relapse of sex, disease type, phase at transplant, source of stem cells and age.
## 'data.frame': 177 obs. of 7 variables:
## $ Sex : Factor w/ 2 levels "F","M": 2 1 2 1 1 2 2 1 2 1 ...
## $ D : Factor w/ 2 levels "ALL","AML": 1 2 1 1 1 1 1 1 1 1 ...
## $ Phase : Factor w/ 4 levels "CR1","CR2","CR3",..: 4 2 3 2 2 4 1 1 1 4 ...
## $ Age : int 48 23 7 26 36 17 7 17 26 8 ...
## $ Status: int 2 1 0 2 2 2 0 2 0 1 ...
## $ Source: Factor w/ 2 levels "BM+PB","PB": 1 1 1 1 1 1 1 1 1 1 ...
## $ ftime : num 0.67 9.5 131.77 24.03 1.47 ...
## Competing Risks Regression
##
## Call:
## crr(ftime = data_bmtcrr$ftime, fstatus = data_bmtcrr$Status,
## cov1 = x, failcode = 1)
##
## coef exp(coef) se(coef) z p-value
## age -0.0185 0.982 0.0119 -1.554 0.1200
## sex_F -0.0352 0.965 0.2900 -0.122 0.9000
## D_AML -0.4723 0.624 0.3054 -1.547 0.1200
## P:CR1 -1.1018 0.332 0.3764 -2.927 0.0034
## P:CR2 -1.0200 0.361 0.3558 -2.867 0.0041
## P:CR3 -0.7314 0.481 0.5766 -1.268 0.2000
## source:PB 0.9211 2.512 0.5530 1.666 0.0960
##
## exp(coef) exp(-coef) 2.5% 97.5%
## age 0.982 1.019 0.959 1.005
## sex_F 0.965 1.036 0.547 1.704
## D_AML 0.624 1.604 0.343 1.134
## P:CR1 0.332 3.009 0.159 0.695
## P:CR2 0.361 2.773 0.180 0.724
## P:CR3 0.481 2.078 0.155 1.490
## source:PB 2.512 0.398 0.850 7.426
##
## Num. cases = 177
## Pseudo Log-likelihood = -267
## Pseudo likelihood ratio test = 24.4 on 7 df,
The first part of the output shows for each term in the design matrix the estimated coefficient \(\hat{\beta}_j\), the relative risk \(\mathrm{exp}(\hat{\beta}_j)\), the standard error, the z-value and the corresponding P-value for assessing significance.
Here, Sex is not significant, followed by Age and D (disease type), whereas Source is only marginally significant. Phase is a factor with relapse as baseline, so each P-value provides a test for the difference of each level with respect to the baseline.
An overall P-value for Phase can be obtained through the Wald test via aod R package:
## Wald test:
## ----------
##
## Chi-squared test:
## X2 = 14.0, df = 3, P(> X2) = 0.0029
The first argument to the function wald.test() is the estimated covariance matrix for the coefficients, followed by the vector of coefficients estimates, and the position of coefficients for which we want to assess significance. In our case, the P-value indicates that Phase is statistically significant.
The second part of the output for competing risks regression shows the relative risk for each term, and a 95% confidence interval. The relative risk or subdistribution hazard ratio for a categorical covariate is the ratio of subdistribution hazards for the actual group with respect to the baseline, with all other covariates being equal. If the covariate is continuous then the relative risk refers to the effect of a one unit increase in the covariate, with all other covariates being equal. In our data, exp(−0.0352)=0.965 is the relative risk of a woman with respect to a man, and exp(−0.0185)=0.982 is the relative risk for a 1 year increase in age.
The last part of the output shows the pseudo log-likelihood at maximum and the pseudo likelihood ratio test, that is, the difference in the objective function at the global null and at the final estimates. As this objective function is not a true likelihood, this test statistic is not asymptotically distributed as a χ2. As a consequence, model comparison based on likelihood ratio approach cannot be performed directly, but significance must be evaluated through simulations. However, a model selection criterion can be easily adopted as described in the following section.
For an example on model selection, model diagnostics or adding time-varying covariates in presence of competing events, please refer to the original tutorial article Scrucca, Santucci, and Aversa (2010) https://www.nature.com/articles/bmt2009359.
Scrucca, L, A Santucci, and F Aversa. 2010.
“Regression Modeling of Competing Risk Using r: An in Depth Guide for Clinicians.” Bone Marrow Transplantation 45: 1388–95.
https://doi.org/10.1038/bmt.2009.359.