Figures
Abstract
Background
Evolution of metastatic melanoma (MM) under B-RAF inhibitors (BRAFi) is unpredictable, but anticipation is crucial for therapeutic decision. Kinetics changes in metastatic growth are driven by molecular and immune events, and thus we hypothesized that they convey relevant information for decision making.
Patients and methods
We used a retrospective cohort of 37 MM patients treated by BRAFi only with at least 2 close CT-scans available before BRAFi, as a model to study kinetics of metastatic growth before, under and after BRAFi. All metastases (mets) were individually measured at each CT-scan. From these measurements, different measures of growth kinetics of each met and total tumor volume were computed at different time points. A historical cohort permitted to build a reference model for the expected spontaneous disease kinetics without BRAFi. All variables were included in Cox and multistate regression models for survival, to select best candidates for predicting overall survival.
Results
Before starting BRAFi, fast kinetics and moreover a wide range of kinetics (fast and slow growing mets in a same patient) were pejorative markers. At the first assessment after BRAFi introduction, high heterogeneity of kinetics predicted short survival, and added independent information over RECIST progression in multivariate analysis. Metastatic growth rates after BRAFi discontinuation was usually not faster than before BRAFi introduction, but they were often more heterogeneous than before.
Citation: Hartung N, Huynh CT-K, Gaudy-Marqueste C, Flavian A, Malissen N, Richard-Lallemand M-A, et al. (2017) Study of metastatic kinetics in metastatic melanoma treated with B-RAF inhibitors: Introducing mathematical modelling of kinetics into the therapeutic decision. PLoS ONE 12(5): e0176080. https://doi.org/10.1371/journal.pone.0176080
Editor: Suzie Chen, Rutgers University, UNITED STATES
Received: December 1, 2016; Accepted: April 5, 2017; Published: May 4, 2017
Copyright: © 2017 Hartung et al. This is an open access article distributed under the terms of the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original author and source are credited.
Data Availability: All relevant data are within the paper and its Supporting Information files.
Funding: This work was funded by a grant from F. Hoffmann-La Roche Pharmaceutical. The funders had no role in study design, data collection and analysis, decision to publish, or preparation of the manuscript.
Competing interests: The authors received funding from F. Hoffmann-La Roche Pharmaceutical, a commercial company, for this study. There are no patents, products in development or marketed products to declare. This does not alter our adherence to all the PLOS ONE policies on sharing data and materials.
Introduction
The course of a metastatic melanoma (MM) is currently unpredictable since aggressiveness depends on a network of variables related to tumour and host reaction [1, 2]. B-RAF inhibitors (BRAFi), and MEKi, have led to a major improvement on survival in B-RAF mutated patients [3–7]. However, resistance mechanisms are mostly unpredictable [8–10] and heterogeneity of resistance mechanisms within the same individual [11–14] makes the monitoring difficult, even if liquid biopsies are under development [15]. In addition, resistance may not only be due to molecular events and immune changes may interfere [16].
Growth kinetics, measuring the change in tumor load over time, may be one of the best ways to characterize disease scenarios for therapeutic trials [17]. In a historical cohort of MM patients treated before the era of new treatments, we have shown that initial kinetics of metastases (mets) measured by 2 successive CT-scans is highly predictive for survival [18]. It has also been shown that the homogeneity of response under BRAFi had a prognostic impact on survival [19].
We hypothesized that kinetics of changes in tumor load before and under BRAFi was per se reflecting molecular, genetic, and immune mechanisms driving the disease, and could be easier to monitor than a huge number of biomarkers.
Our objective was to show that monitoring of metastatic disease kinetics under targeted therapy is a source of relevant predictive information, which mathematical modelling could use to anticipate events for decision-making.
Materials and methods
Study populations
BRAFi treated population.
To have the simplest approach for modelling, we retrospectively selected from the cohort of MM patients treated in our department (Dermatology and skin cancer department, La Timone Hospital Marseille, France), a series of patients with the following inclusion criteria: stage IIIC or IV AJCC [20], BRAF V600E/K mutation, treatment with BRAFi monotherapy only, and at least two whole-body CT-scans available before BRAFi treatment and at least one CT-scan after BRAFi treatment, all performed on the same machine (in one of the radiology department of our institution), with the same procedure, at most three months apart. Data collection was performed between June 2014 and April 2015.
Historical cohort.
To estimate the natural kinetics of the metastatic disease, we used a historical cohort of patients who never received any treatment with a demonstrated impact on survival [18] retrospectively selected with the following inclusion criteria: stage IV MM patients treated in our Institution between September 2007 and October 2011 who had two total body computed tomography (CT) scans on the same machine with the same procedure within a maximum of 4 months period after first distant metastases diagnosis, and who meanwhile received either no treatment or only monochemotherapy with dacarbazine or fotemustine or vaccines. Data collection was performed between June and September 2012.
Assessment of metastatic volumes and kinetics
Volumes of mets were computed using the two native axial measurements and the third measurement from a coronal reconstruction (General Electric Medical Systems, Advantage Workstation 4·4), and assuming an ellipsoidal shape. Each met of which any diameter exceeded 1cm at some time during the follow-up was measured (Fig 1) in all patients at all CT-scans along the follow-up. In this pilot study, measurements were manual and required many hours per patient, which led to limit the sample of patients.
Different variables characterizing metastatic growth kinetics were computed and tested as prognostic markers for survival. Growth kinetics were evaluated as exponential growth rate, linear slope, and relative volume change (RVC). As a measure of the overall kinetic tendency in an individual, kinetics of the total tumor volume and mean kinetics of all mets were calculated. In addition, range and standard deviation of individual metastatic growth kinetics were computed as a measure of homogeneity vs heterogeneity of kinetics between different mets of a given individual. Tumor size-based kinetic indices were computed before and after start of BRAFi therapy (Table 1).
All kinetic measures were computed both on the global volume and as a mean measure on the individual metastases. Spread of kinetics was evaluated as standard deviation and range of measures computed on the individual metastases, thus representing a measure of inter-metastasis heterogeneity.
Modelling the natural spontaneous course of the disease
We used the historical cohort to predict, from initial kinetics, what would have been the natural course of the metastatic disease. An exponential model based on pre-treatment tumor volumes was used to calculate a theoretical total metastatic volume at a time t after first kinetics measure: . We assumed that a patient would die when V(t) reached a critical metastatic burden Vlethal, the most realistic value of which was computed from the historical cohort. The best variable to predict spontaneous survival was found to be the inverse of the time to reach the critical metastatic burden
, thereafter called “lethal burden risk score”. More complex tumor growth models were not considered because the risk score had to be derived from only two volume measurements per patient.
Assessing the predictive value of kinetics indices
Kinetics variables were included in multivariate models with best recognized pejorative predictive markers in patients treated by BRAFi, i.e. RECIST (Response Evaluation Criteria In Solid Tumors version 1.1) progression [21] appearance of new mets under BRAFi, and mixed vs uniform initial response [19]. Since RECIST evaluation was not always available, the baseline lesions with largest volumes were considered as theoretical RECIST targets, up to a maximum of five in total and two per organ. A doubling of any non-target lesion volume was considered an unequivocal progression.
Statistical analyses
Overall survival from start of BRAFi therapy was taken as the clinical endpoint. Patient survival was censored on 16/02/2015 and left-truncated at the time where inclusion criteria were met, i.e. at the first evaluation after start of BRAFi. Impact of covariates available at entering the group at-risk on patient survival was evaluated in Cox proportional hazards (PH) models. To assess whether the tested kinetic indices are independent prognostic markers with respect to the above-mentioned other clinical criteria, they were also evaluated together in multiple explanatory covariate Cox PH models. Furthermore, to account for the appearance of brain metastases and change in RECIST status, multistate models with transition-specific covariates were used [22].
All regression analyses were carried out in R software version 3·1·3 [23] using packages “survival” for Cox and multistate regression [24, 25], and “mstate” for multistate modelling pre and postprocessing [26–28]. Continuous covariates ranging over several orders of magnitude were log-transformed prior to inclusion. The adequacy of functional form of a covariate was evaluated graphically via martingale residuals vs. covariate plots [29]. Likelihood ratios tests were used for model selection, and p-values below 0.05 were considered to be statistically significant.
A survival tree was generated using the “rpart” R package version 4·1–9 with default options [30]. In the algorithm, a two-step procedure is used: in the first step, a tree is built up by recursively partitioning patients into binary subgroups, maximizing the relative risk at each node. Then, to avoid overfitting, the tree is pruned back using a cross-validation criterion, a method that mimics prediction of a future patient not used in building the model. The local full likelihood model allows to deal with censored observation times [31].
Results
Study populations
BRAFi-treated population.
37 patients fulfilled the inclusion criteria. 27 were treated by vemurafenib only, 5 by dabrafenib and 5 successively by vemurafenib and dabrafenib. Population characteristics were as follows: 20 m/17 f, median age 54 (range 20–84 years), AJCC III C in 1 patient, IV M1a in 3 patients, IV M1b in 11 patients and IV M1C in 22 patients. At the end of follow-up, 12 patients were alive on BRAFi, 4 alive with another treatment after BRAFi, 7 died after BRAFi discontinuation, and 14 died during BRAFi therapy.
Historical cohort.
109 patients of the historical cohort fulfilled the inclusion criteria. Population characteristics are described elsewhere [18] but were similar as in the BRAFi-treated population. In the above-mentioned “lethal burden risk score”, the optimal critical value Vlethal = 1200 cm3 was determined from the historical cohort.
Cox regression
A high “lethal burden risk score”, which was the most significant predictors of overall survival (OS) in the historical cohort (p<0·0001) was shown to be also predictor of a poor prognosis in the BRAFi cohort, although less strongly associated (p = 0·04), confirming that initial kinetics before BRAFi is influencing survival under BRAFi. However, the hazard ratio (HR) of a unit increase in risk score tends to belower in the BRAFi cohort, than in the historical cohort with no active treatment: HR = 1·6 (95% CI 1·0–2·3) vs 2·5 (95% CI 1·9–3·4, respectively (p = 0·07 in a test for same HR in both cohorts).
Heterogeneity of initial kinetics expressed as range of RVC before BRAFi, as well as heterogeneity of kinetics in response to BRAFi were significantly predicting poor survival (p = 0·02, and p = 0·005, respectively). Standard deviations were not considered further to describe heterogeneity, due to the high correlation with range of RVC.
As expected, RECIST progressive disease (PD), appearance of new lesions, and a “mixed response” [19] at 1st disease assessment after BRAFi were significantly associated with a poorer outcome (RECIST PD p = 0·004, new lesions p = 0·01, mixed response p = 0·005). Since these three criteria were highly correlated, only RECIST PD, the most significant variable, was retained for further evaluation.
All predictors retained in single explanatory covariate Cox regression, namely “lethal burden risk score”, range of RVC at baseline, and range of RVC at initial response were tested together with RECIST at 1st treatment assessment in a multiple explanatory covariate Cox model. Only range of RVC, hereafter referred to as “response heterogeneity” significantly contributed to the predictive power of the model (p = 0·05). High response heterogeneity was a predictor of poor outcome among RECIST progressors, but not significantly linked to survival in RECIST responders.
Multistate modelling
To evaluate the impact of dynamic changes of the disease status likely to have an impact on survival, two multistate models were considered, one for brain mets and one for RECIST progression (see Fig 2). Brain-metastatic status was defined by presence or absence of brain mets and by distinguishing the cause of death (linked to or independent from brain mets). Both RECIST PD at 1st assessment and response heterogeneity were associated to the risk of dying for other reasons than brain mets, but not to the risk of developing brain mets, or the risk of dying from brain mets.
hx: Transition hazards. Asterisks (*) indicate transitions on which covariates are modelled.
Compared to the standard Cox model, the level of significance of “response heterogeneity” increased both independently of RECIST PD at 1st assessment (p = 0·001 in multistate vs p = 0·005 in Cox model) and also when considering RECIST PD and response heterogeneity jointly (p = 0·03 in multistate vs p = 0·05 in Cox model). A proportional baseline hazard assumption for the different causes of death was tested but not supported by the data [26].
We tested whether response heterogeneity could be also prognostic at a RECIST PD occurring later in the course of the diseases. Considering a multistate model for RECIST progression, response heterogeneity was shown to be a prognostic marker for survival from RECIST PD (p = 0·02), but no link was found between time to RECIST PD and survival from RECIST PD.
Comparing metastatic kinetics after BRAFi with pre-treatment kinetics
In the 11 patients with CT-scans performed after BRAFi discontinuation, empirical exponential growth rates computed from two volume measurements just before treatment start and two measurements just after treatment were compared to assess the impact of BRAFi on metastatic growth kinetics. Exponential rates of global burden and mean exponential rates were decreased in 9 and 7 patients, stable in 0 and 2, increased in 2 and 2, respectively. Heterogeneity of kinetics between the different mets of a given individual (range of individual growth rates for each met), was decreased in 5 patients, stable in 2, but increased in 5.
Discussion
Although manual measurement is extremely time-consuming until automatic measures are developed, this study is the proof of concept that a mathematical modelling using kinetics data could be helpful in predicting disease course under BRAFi. Several variables seem important: those characterizing initial kinetics, but more importantly those characterizing chaos in kinetics either before BRAFi or at 1st assessment under BRAFi.
Initial kinetics of MM disease are known to be a good reflect of tumor aggressiveness [18], which can be measured by 2 CT-scans at least 1 month interval before treatment. Ethical considerations about potential loss of time during this interval are not a real issue since usual delay before treatment (eg, validation of mutational status, surgery, pretreatment assessment…) is often longer than a month. A variable like “Lethal burden risk score” representing the predicted spontaneous survival according to initial kinetics, confirms that a high initial kinetics is associated with poor prognosis even in patients treated by BRAFi. This score however shows a trend for a lower impact on the outcome in BRAFi treated patients than in the patients of the historical cohort, suggesting that BRAFi can change the natural history of the disease. It should be noticed that due to the limited cohort size no distinction was made between the two different BRAFi molecules which seem to have similar response profile and response rate.
This study underlines the major prognostic role of heterogeneity of the kinetics as expressed by range of RVC or other similar variables. It is true for heterogeneity of initial kinetics, which may reflect the diversity of molecular pathways at work before treatment, thus increasing the probability of primary resistance in some mets. It is also true for heterogeneity of kinetics at 1st response under BRAFi, which may reflect the diversity of secondary resistance mechanisms in different mets. This evaluation of heterogeneity based on a comprehensive assessment of all mets is supporting the results of a previous study measuring changes in a sample of mets showing that a “mixed” initial response with BRAFi was associated to poorer outcome than a “uniform” response [19]. High heterogeneity in the kinetics of response under BRAFi was still able to add an independent pejorative information to RECIST progression (Cox and multistate models), although initial RECIST progression was one of the most powerful predictor of short survival for any treatment of MM. Furthermore, if we focus on the prediction of survival after RECIST progression, response heterogeneity between mets at progression was a better pejorative marker than a short time to RECIST progression. Unexpectedly, deaths due to brain a mets did not correlate with response heterogeneity perhaps because brain lesions were treated by Gamma-knife in most of our patients.
It is a common perception among clinicians that, when a patient is escaping the B-RAF blockade, the disease is even faster than before treatment. Although potentially biased since cases with early deaths did not permit post-treatment evaluation, our data are not supporting this perception. However, they show more heterogeneity in kinetics after treatment than before, which may be accounted by the fact that BRAFi may promote different resistance mechanisms in the different mets. Describing the whole time course of metastatic growth under BRAFi with a mathematical model could potentially alleviate this bias, but was beyond the scope of this manuscript.
Monitoring of kinetics can be used as a predictor at 2 levels. As a simple clinical indicator, it is more or less supporting the instinct of the clinician experienced with BRAFi: 1- fast metastatic growth before treatment, especially if there is a clear disorder in growth among the different mets is probably not an ideal situation for BRAFi, and 2-a response to BRAFi with contrasting kinetics in different mets should prompt clinician to switch treatment early if possible. Another level of application is to build up a real predictive model to facilitate therapeutic decision. An example of survival tree has selected two factors (Fig 3): expected spontaneous survival deduced from initial kinetics before BRAFi and 2) heterogeneity of kinetics at 1st assessment of response under BRAFi expressed as range of RVC. This model requires a comprehensive measure of all mets. Such a model could obviously be improved if it was generated from a larger sample of patients, and validated in another sample. This is out of reach for manual measurements but automatization or computer assistance seems quite feasible.
Using a partitioning algorithm, 2 kinetics-related factors were selected from the four variables considered in the multiple explanatory covariate Cox model described in the Results. The derived classification algorithm first uses the expected survival deduced from initial kinetics and then heterogeneity of kinetics at 1st assessment expressed as range of RVC.
This pilot study shows that kinetics monitoring is a useful source of information in MM patients treated with BRAFi, which may also apply to BRAFi-MEKi combinations. Although most of research is focusing on molecular and immunological biomarkers for the selection and adaptation of the treatment in a given patient, mathematical modelling based on usual imaging systems deserves more interest.
Supporting information
S1 File. Dataset specification.
Specification of the variables in Supporting Datasets S1 and S2 Datasets.
https://doi.org/10.1371/journal.pone.0176080.s001
(PDF)
Author Contributions
- Conceptualization: JJG CG NH FH CH.
- Data curation: CH AF NH.
- Formal analysis: NH FH.
- Funding acquisition: JJG CG.
- Investigation: CG CH NM AF.
- Methodology: JJG FH NH.
- Project administration: JJG CG FH.
- Resources: JJG CG NM FH NH MAR.
- Software: NH FH.
- Supervision: JJG FH.
- Validation: JJG CG FH NH.
- Visualization: CG NH JJG CH.
- Writing – original draft: NH CG FH JJG NM MAR.
- Writing – review & editing: CG JJG NH FH.
References
- 1. Ma X.H., Piao S., Wang D., McAfee Q.W., Nathanson K. L., Lum J. J., et al. Measurements of tumor cell autophagy predict invasiveness, resistance to chemotherapy, and survival in melanoma. Clinical Cancer Research 2011;17:3478–3489. pmid:21325076
- 2. Qendro V., Lundgren D. H., Rezaul K., Mahony F., Ferrell N., Bi A., et al. Large-scale proteomic characterization of melanoma expressed proteins reveals nestin and vimentin as biomarkers that can potentially distinguish melanoma subtypes. Journal of Proteome Research 2014;13: 5031–5040. pmid:25322343
- 3. Larkin J., Ascierto P. A., Dréno B., Atkinson V., Liszkay G., Maio M., et al. Combined Vemurafenib and Cobimetinib in BRAF -Mutated Melanoma. New England Journal of Medicine 2014;371:1867–1876. pmid:25265494
- 4. Larkin J., Del Vecchio M., Ascierto P. A., Krajsova I., Schachter J., Neyns B., et al. Vemurafenib in patients with BRAFV600 mutated metastatic melanoma: An open-label, multicentre, safety study. The Lancet Oncology 2014; 15:436–444. pmid:24582505
- 5. McArthur G. A., Chapman P. B., Robert C., Larkin J., Haanen J. B., Dummer R., et al. Safety and efficacy of vemurafenib in BRAFV600E and BRAFV600K mutation-positive melanoma (BRIM-3): extended follow-up of a phase 3, randomised, open-label study. The Lancet Oncology 2014 15:323–332. pmid:24508103
- 6. Robert C., Karaszewska B., Schachter J., Rutkowski P., Mackiewicz A., Stroiakovski D., et al. Improved Overall Survival in Melanoma with Combined Dabrafenib and Trametinib. New England Journal of Medicine 2015;372:30–39. pmid:25399551
- 7. Sosman J. A., Kim K. B., Schuchter L., Gonzalez R., Pavlick A. C., Weber J. S., et al. Survival in BRAF V600-mutant advanced melanoma treated with vemurafenib. The New England Journal of Medicine 2012;366:707–14. pmid:22356324
- 8. Bucheit A. D., & Davies M. A. Emerging insights into resistance to BRAF inhibitors in melanoma. Biochemical Pharmacology, 2014;87:381–389. pmid:24291778
- 9. Fallahi-Sichani M., Moerke N. J., Niepel M., Zhang T., Gray N. S., & Sorger P. K. Systematic analysis of BRAF V 600 E melanomas reveals a role for JNK / c-Jun pathway in adaptive resistance to drug-induced apoptosis. Molecular Systems Biology 2015;11:797. pmid:25814555
- 10. Long G. V, Fung C., Menzies A. M., Pupo G. M., Carlino M. S., Hyman J., et al. Increased MAPK reactivation in early resistance to dabrafenib/trametinib combination therapy of BRAF-mutant metastatic melanoma. Nature Communications 2014;5: 5694 pmid:25452114
- 11. Hoogstraat M., Gadellaa-van Hooijdonk C. G., Ubink I., Besselink N. J., Pieterse M., Veldhuis W., et al. Detailed imaging and genetic analysis reveal a secondary BRAF(L505H) resistance mutation and extensive intrapatient heterogeneity in metastatic BRAF mutant melanoma patients treated with vemurafenib. Pigment Cell & Melanoma Research 2015;28:318–323.
- 12. Kemper K., Krijgsman O., Cornelissen-Steijger P., Shahrabi A., Weeber F., Song J. Y., et al. Intra- and inter-tumor heterogeneity in a vemurafenib-resistant melanoma patient and derived xenografts. EMBO Molecular Medicine 2015;7, e201404914.
- 13. Monsma D. J., Cherba D. M., Eugster E. E., Dylewski D. L., Davidson P. T., Peterson C. A., et al. Melanoma patient derived xenografts acquire distinct Vemurafenib resistance mechanisms. American Journal of Cancer Research 2015;5:1507–1518. pmid:26101714
- 14. Rizos H., Menzies A. M., Pupo G. M., Carlino M. S., Fung C., Hyman J., et al. BRAF Inhibitor Resistance Mechanisms in Metastatic Melanoma: Spectrum and Clinical Impact. Clinical Cancer Research 2014;20:1965–1977. pmid:24463458
- 15. Girotti M. R., Gremel G., Lee R., Galvani E., Rothwell D., Viros A., et al. Application of sequencing, liquid biopsies and patient-derived xenografts for personalized medicine in melanoma. Cancer Discovery 2015; http://doi.org/10.1158/2159-8290.CD-15-1336
- 16. Smith M. P., Sanchez-Laorden B., O’Brien K., Brunton H., Ferguson J., Young H., et al. The Immune Microenvironment Confers Resistance to MAPK Pathway Inhibitors through Macrophage-Derived TNFalpha. Cancer Discovery 2014;4:1214–1229. pmid:25256614
- 17. Grob J.-J., Long G. V., Schadendorf D., Flaherty K. Disease kinetics for decision-making in advanced melanoma: a call for scenario-driven strategy trials. Lancet Oncol 2015;16:e522–6. pmid:26433825
- 18. Gaudy-Marqueste C., Archier E., Grob A., Durieux O., Loundou A., Richard M.-A., et al. Initial metastatic kinetics is the best prognostic indicator in stage IV metastatic melanoma. European Journal of Cancer 2014;50:1120–4. pmid:24440089
- 19. Menzies A. M., Haydu L. E., Carlino M. S., Azer M. W. F., Carr P. J. A., Kefford R. F., et al. Inter- and intra-patient heterogeneity of response and progression to targeted therapy in metastatic melanoma. PloS One 2014;9, e85004. pmid:24400126
- 20. Balch C. M., Gershenwald J. E., Soong S.-J., Thompson J. F., Atkins M. B., Byrd D. R., et al. Final Version of 2009 AJCC Melanoma Staging and Classification. Journal of Clinical Oncology 2009;27:6199–6206. pmid:19917835
- 21. Eisenhauer E. A., Therasse P., Bogaerts J., Schwartz L. H., Sargent D., Ford R., et al. New response evaluation criteria in solid tumours: Revised RECIST guideline (version 1.1). European Journal of Cancer 2009;45:228–247. pmid:19097774
- 22. Andersen P. K., & Keiding N. Multi-state models for event history analysis. Statistical Methods in Medical Research 2002, 11(2), 91–115. http://doi.org/10.1191/0962280202SM276ra pmid:12040698
- 23.
R Core Team. (2015). R: A language and environment for statistical computing. R Foundation for Statistical Computing, Vienna, Austria. http://doi.org/ISBN3-900051-07-0
- 24.
Therneau T. M., & Grambsch P. M. (2000). Modeling Survival Data: Extending the Cox Model (1st ed.). Springer-Verlag New York. http://doi.org/10.1007/978-1-4757-3294-8
- 25.
Therneau, T. M. (2015). A package for Survival Analysis in S. version 2.38. Retrieved from http://cran.r-project.org/package=survival
- 26. de Wreede L. C., Fiocco M., & Putter H. The mstate package for estimation and prediction in non- and semi-parametric multi-state and competing risks models. Computer Methods and Programs in Biomedicine 2010, 99(3), 261–274. http://doi.org/10.1016/j.cmpb.2010.01.001 pmid:20227129
- 27. de Wreede L. C., Fiocco M., & Putter H. mstate: An R Package for the Analysis of Competing Risks and Multi-State Models. Journal of Statistical Software 2011, 38(7), 1–30. http://www.jstatsoft.org/v38/i07
- 28. LeBlanc M., & Crowley J. Relative risk trees for censored survival data. Biometrics 1992, 48(2), 411–415. http://doi.org/10.2307/2532300 pmid:1637970
- 29. Putter H., Fiocco M., & Geskus R. B. Tutorial in biostatistics: competing risks and multi-state models. Statistics in Medicine 2007, 26(11), 2389–430. http://doi.org/10.1002/sim.2712 pmid:17031868
- 30. Therneau T. M., Grambsch P. M., & Fleming T. R. (1990). Martingale-based residuals for survival models. Biometrika. http://doi.org/10.1093/biomet/77.1.147
- 31. Therneau T., Atkinson B., & Ripley B. (2015). rpart: Recursive Partitioning and Regression Trees. R package version. 4.1–9. Retrieved from http://cran.r-project.org/package=rpart