Martingale Residuals Research Articles

Association between Split Selection Instability and Predictive Error in Survival Trees

To evaluate split selection instability in six survival tree algorithms and its relationship with predictive error by means of a bootstrap study. We study the following algorithms: logrank statistic with multivariate p-value adjustment without pruning (LR), Kaplan-Meier distance of survival curves (KM), martingale residuals (MR), Poisson regression for censored data (PR), within-node impurity (WI), and exponential log-likelihood loss (XL). With the exception of LR, initial trees are pruned by using split-complexity, and final trees are selected by means of cross-validation. We employ a real dataset from a clinical study of patients with gallbladder stones. The predictive error is evaluated using the integrated Brier score for censored data. The relationship between split selection instability and predictive error is evaluated by means of box-percentile plots, covariate and cutpoint selection entropy, and cutpoint selection coefficients of variation, respectively, in the root node. We found a positive association between covariate selection instability and predictive error in the root node. LR yields the lowest predictive error, while KM and MR yield the highest predictive error. The predictive error of survival trees is related to split selection instability. Based on the low predictive error of LR, we recommend the use of this algorithm for the construction of survival trees. Unpruned survival trees with multivariate p-value adjustment can perform equally well compared to pruned trees. The analysis of split selection instability can be used to communicate the results of tree-based analyses to clinicians and to support the application of survival trees.

Proposal for revision of the TNM classification system for renal cell carcinoma

The current study defined an optimal tumor size breakpoint to stratify localized renal cell carcinoma (RCC) into groups with significantly different cancer-related outcomes and proposed a revision of the TNM classification system. The authors analyzed the data from 1138 patients who had undergone partial or radical nephrectomy for localized RCC at 7 European urologic centers. The optimal pathologic size breakpoint was calculated using the martingale residuals from a Cox proportional hazards regression model. The mean follow-up time was 87 months. The scatterplot of tumor size versus expected risk of death per patient suggested that an interval of 5-6 cm was appropriate. A total of 720 (63.3%) and 418 (36.7%) patients had tumors measuring < or = 5.5-cm and tumors measuring > 5.5-cm, respectively. Significant cancer-specific survival differences between the two groups of patients were reported in the series by all the centers participating in the study. On univariate analysis, the other variables found to be associated with cancer-specific survival were the patient's age, symptomatic tumor presentation, and the Fuhrman nuclear grade. On multivariate analysis, the pathologic stage of the primary tumor defined according to the 5.5-cm breakpoint was found to be an independent predictor of cancer-specific survival, as well as age, mode of presentation, and nuclear grade. According to the multivariate analysis, the authors clustered patients into 3 groups with statistically significant outcome differences: 1) patients with < or = 5.5-cm incidentally detected RCC; 2) patients with < or = 5.5-cm symptomatic RCC; and 3) patients with > 5.5-cm RCC. This cancer-related outcome stratification was valid regardless of the patient's age. The 5.5-cm breakpoint was found to be the optimal tumor size breakpoint with which to stratify patients with organ-confined RCC. The study supported the upgrade of the TNM classification system according to this breakpoint.

PT2 CLASSIFICATION FOR RENAL CELL CARCINOMA. CAN ITS ACCURACY BE IMPROVED?

The 2002 tumor classification for renal cell carcinoma (RCC) classifies pT2 tumors as more than 7 cm in greatest dimension, limited to the kidney. In this study we determined whether a size cutoff point exists within pT2 tumors and whether such subclassification would further improve the accuracy of the current tumor classification. We studied 544 patients with unilateral, sporadic pT2 RCC treated with radical nephrectomy or nephron sparing surgery between 1970 and 2000. The association of tumor size with death from RCC was examined using martingale residuals from a Cox proportional hazards regression model to determine the optimal size cutoff point. There were 204 deaths from RCC a median of 3.8 years following nephrectomy. Univariately tumor size was significantly associated with death from RCC (risk ratio 1.08, 95% CI 1.04 to 1.13, p <0.001). A scatterplot of tumor size vs expected risk of death per patient suggested that a cutoff point between 9 and 10 cm was appropriate. When adjusted for regional lymph node involvement and distant metastases, the 10 cm cutoff point performed better than the 9 cm point (risk ratio 1.42, 95% CI 1.07 to 1.90, p = 0.017 vs 1.22, 95% 0.86 to 1.72, p = 0.268). Therefore, we propose using a 10 cm cutoff point to subclassify patients into pT2a and pT2b. Our data suggest that the prognostic accuracy of the 2002 pT2 tumor classification can be further improved by subclassifying patients with tumors greater than 7 and less than 10 cm into a pT2a category, and those with tumors 10 cm or greater into a pT2b category.

Testing association of a pathway with survival using gene expression data

A recent surge of interest in survival as the primary clinical endpoint of microarray studies has called for an extension of the Global Test methodology to survival. We present a score test for association of the expression profile of one or more groups of genes with a (possibly censored) survival time. Groups of genes may be pathways, areas of the genome, clusters from a cluster analysis or all genes on a chip. The test allows one to test hypotheses about the influence of these groups of genes on survival directly, without the intermediary of single gene testing. The test is based on the Cox proportional hazards model and is calculated using martingale residuals. It is possible to adjust the test for the presence of covariates. We also present a diagnostic graph to assist in the interpretation of the test result, visualizing the influence of genes. The test is applied to a tumor dataset, revealing pathways from the gene ontology database that are associated with survival of patients. The Global Test for survival has been incorporated into the R-package globaltest (version 3.0), available at http://www.bioconductor.org

Responder identification in clinical trials with censored data

A new technique is presented for identification of positive and negative responders to a new treatment which was compared to a classical treatment (or placebo) in a randomized clinical trial with respect to survival time. This method was primarily developed for trials in which the two treatment arms do not differ in overall survival. It checks in a systematic manner if certain subgroups, described by so-called predictive factors, do show difference in survival due to the new treatment. The method relies on a good prognostic model built on one arm treated with placebo or the classical therapy. It employs the martingale residuals of the prognostic model in a stabilized bump hunting procedure, which finds groups of responders in the new treatment group with large positive or large negative residuals. The results of a simulation study are presented, comparing the performance of the new method to that of the Cox-PH model with treatment–covariate interactions. In a simulation study on average the proposed method recognizes in 90% the correct positive responder group and in 99% the correct negative responder group. The method is to be used in explorative analysis for hypothesis generation. The results are to be validated in future studies.

Cytokines in pancreatic carcinoma

Cytokines have been implicated in diverse processes that are relevant to pancreatic carcinoma, including cachexia, asthenia, and tumor growth. The objective of this study was to examine the association between serum levels of proinflammatory and antiinflammatory or angiogenic cytokines and the outcomes of patients with pancreatic carcinoma. Serum cytokine levels were measured by enzyme-linked immunosorbent assay from 51 patients with pancreatic carcinoma and from 48-62 healthy volunteers. Cytokine levels were compared with disease manifestations and overall survival. Circulating levels of vascular endothelial growth factor, tumor necrosis factor alpha, interleukin-1alpha (IL-1alpha), and IL-1beta were not elevated significantly in patients with pancreatic carcinoma, but levels of IL-6, IL-8, IL-10, and IL-1 receptor antagonist (IL-1RA) were elevated significantly (P <0.05). Cytokine levels were dichotomized based on an analysis of null Martingale residuals. Patients who had IL-6 levels > 5.2 pg/mL or IL-10 levels >9.8 pg/mL had significantly worse survival compared with patients who had lower IL-6 or IL-10 levels (P <0.05). IL-8 levels were not associated with survival differences. Patients who had IL-1RA levels <159 pg/mL had significantly worse survival compared with patients who had higher IL-1RA levels (P <0.05). Higher IL-6, IL-10, and IL-8 levels were associated with poor performance status and/or weight loss. In multivariate analysis, only T4 tumors and high IL-6 levels were selected as independent prognostic factors for poor survival. Circulating levels of several cytokines were high in patients with pancreatic carcinoma, and their association with weight loss and poor performance status suggested that they may be involved in these disease manifestations. Furthermore, serum cytokine levels, in particular IL-6, may be a useful prognostic marker.

Alternative Tree-Structured Survival Analysis Based on Variance of Survival Time

Tree-structured survival analysis (TSSA) is a popular alternative to the Cox proportional hazards regression in medical research of survival data. Several methods for constructing a tree of different survival profiles have been developed, including TSSA based on log-rank statistics, martingale residuals, Lp Wasserstein metrics between Kaplan-Meier survival curves, and a method based on a weighted average of the within-node impurity of the death indicator and the within-node loss function of follow-up times. Lu and others used variance of restricted mean lifetimes as an index of degree of separation (DOS) to measure the efficiency in separations of survival profiles by a classification method. Like tree-based regression analysis that uses variance as a criterion for node partition and pruning, the variance of restricted mean lifetimes between different groups can be an alternative index to log-rank test statistics in construction of survival trees. In this article, the authors explore the use of DOS in TSSA. They propose an algorithm similar to the least square regression tree for survival analysis based on the variance of the restricted mean lifetimes. They apply the proposed method to prospective cohort data from the Study of Osteoporotic Fracture that motivated the research and then compare their classification rule to those rules based on the conventional TSSA mentioned above. A limited simulation study suggests that the proposed algorithm is a competitive alternative to the log-rank or martingale residual-based TSSA approaches.

A test for informative censoring in clustered survival data.

Frailty models are frequently used to analyse clustered survival data. The assumption of non-informative censoring is commonly used by these models, even though it may not be true in many situations. This article proposes a test for this assumption. It uses the estimated correlation between two types of martingale residuals, one from a model for failure and the other from a model for censoring. It distinguishes two types of censoring, namely withdrawal and the end of the study. Simulation studies show that the proposed test works well under various scenarios. For illustration, the test is applied to a data set for kidney disease patients from multiple dialysis centres.

Dose-response for biochemical control among high-risk prostate cancer patients after external beam radiotherapy

Introduction The literature on dose–response characteristics of high-risk prostate cancer has been scarce in this era, when these patients are treated with hormone therapy along with radiotherapy. In this study, we estimated the dose–response of prostate-specific antigen (PSA) control probability in high-risk prostate cancer patients treated with radiotherapy alone. Methods and materials The data set contains information on 363 high-risk prostate cancer patients who were treated with external beam radiotherapy without hormonal treatment between February 1987 and September 1998. These patients have one or more of the following adverse prognostic features: digital rectal examination stage ≥cT3, PSA >20 ng/mL, and biopsy Gleason score ≥8. These patients had biopsy-proven adenocarcinoma of prostate and were staged according to the 1992 AJCC staging system that was based on digital rectal examination. The logistic model was fitted to the data at various time points after treatment, and the dose–response parameters were estimated. Results The dose required to have 50% tumor control, TCD50 (95% confidence interval), for high-risk patients is 75.5 (range: 70.7–80.2) Gy. The γ50 (95% confidence interval) is 1.7 (range: 0.7–2.7) around 75.5 Gy. Recursive partitioning analysis based on the null Martingale residuals identifies two subgroups within the high-risk group. The TCD50 estimates of the two subgroups (PSA ≤ vs. >20 ng/mL) differ by 15 Gy at 5 years. There is a dose response in both subgroups. Conclusion We recognize that this study has the usual limitations of a retrospective study that includes treatment policy change that spanned a long time frame. However, our data strongly suggest a benefit of dose escalation for all the patients in the entire high-risk group. There is a steep dose response in PSA control probability around a modern dose of 78 Gy. A 5-Gy dose increase beyond 78 Gy may improve PSA control by about 10%.

Comparison of tree-based methods for prognostic stratification of survival data

Tree-based methods can be used to generate rules for prognostic classification of patients that are expressed as logical combinations of covariate values. Several splitting algorithms have been proposed for generating trees from survival data. However, the choice of an appropriate algorithm is difficult and may also depend on clinical considerations. By means of a prognostic study of patients with gallbladder stones and of a simulation study, we compare the following splitting algorithms: log-rank statistic adjusted for measurement scale with (AP) and without (AU) pruning, exponential log-likelihood loss (EP), Kaplan–Meier (KP) distance of survival curves, unadjusted log-rank statistic (LP), martingale residuals (MP), and node impurity (ZP). With the exception of the AU algorithm (based on a Bonferroni-adjusted p-value driven stopping rule), trees are pruned using the measure of split-complexity, and optimally-sized trees are selected using cross-validation. The integrated Brier score is used for the evaluation of predictive models. According to the results of our simulation study and of the clinical example, we conclude that the AU, AP, EP, and LP algorithm may yield superior predictive accuracy. The choice among these four algorithms may be based on the required parsimonity and on medical considerations.

Functional form of the effect of the numbers of axillary nodes on survival in early breast cancer.

The change in survival in function of the numbers of involved and uninvolved axillary nodes in early breast cancer - i.e. the functional form - was investigated to search for prognostic cutoffs and to assess if ratio-based characterization of node involvement is a significant prognostic factor or not. Women aged 40-69, diagnosed in 1988-1997 with T1-T2 invasive breast carcinoma, who underwent axillary dissection, are selected from the SEER public database. The method determines the functional form by applying smoothed plots to the martingale residuals obtained from a proportional hazards model. The results on 55,267 selected patients find that the ratio of involved nodes on examined nodes, in a multivariate model that takes into account known prognostic factors (age, race, tumor size, topography, histology, grade, hormone receptors), is associated with a relative mortality hazard of 1.012 (95% confidence interval 1.010-1.014; relative increase of mortality of 1.2% for each 1% increase in the percentage of involved nodes). The functional form for the number of uninvolved nodes shows that the relative mortality hazard initially steeply decreases and then tends to level off beyond 5-10 uninvolved nodes. For the number of involved nodes, the relative mortality hazard continues to increase with each involved node without any obvious cutpoint. Even when the number of involved nodes is already large, each additional involved node increases the relative mortality hazard by at least 1.3%.

Is there a minimum number of lymph nodes that should be histologically assessed for a reliable nodal staging of T3N0M0 colorectal carcinomas?

Because of the existing controversy, we searched for a cutoff value for the number of lymph nodes (LNs) to be examined in order to establish a reliable node-negative stage in colorectal carcinomas (CRCs). From the SEER database, 8,574 T3N0M0 first, single, histologically confirmed, surgically treated CRCs, with at least 1 LN examined histologically, were considered. As a first approach, the relationships between number of examined LNs and 5- and 10-year overall survival (OS) rates, computed by the Kaplan-Meier method, were assessed. Next, multivariate analysis was performed; a proportional hazards model was fitted to the data and used to obtain a smoothed plot of the martingale residuals vs. the number of negative LNs. Both OS rates displayed an improvement with an increase of number of LNs examined. The smoothed plot of the martingale residuals against the number of negative LNs was reasonably linear. Both approaches suggest that there is no cutoff value for the number of LNs to be examined for an adequate nodal staging; for a reliable pN0 staging, as many LNs should be assessed as possible. However, qualitative features of lymph nodes (e.g., those identified by sentinel lymphadenectomy) may alter this recommendation.

Practical application of residuals from survival models in quantitative trait linkage analysis.

A number of familial diseases have an age-of-onset component, which can be considered as a censored quantitative trait. However, few software resources are available for the use of time-to-event endpoints in linkage analysis. The purpose of this analysis was to examine the use of martingale residuals from Cox survival models as quantitative traits for familial diseases with variable age at onset. We used these residuals as quantitative traits in variance components linkage scans for the 50 replicates of the general population simulated data for chromosomes 6 and 7. The region on chromosome 6 containing markers D06G034 and D06G035 demonstrated evidence for linkage, consistent with the underlying genetic model. This analysis demonstrates the applicability of using martingale residuals as a quantitative trait in linkage analyses of diseases that depend on age of onset.

Variance components analysis for genetic linkage of time to onset for disease.

We compared two variance components methods for detecting genes that influence time to onset for a complex disease using simulated data. We first divided the extended families into nuclear families. The first method fitted variance components to the martingale residuals, which were obtained from first fitting a proportional hazards model to the time to onset data for the trait, allowing for the quantitative traits Q1-Q5, sex, age, and the environmental factor. The second method treated time to onset among the affected individuals as a quantitative trait adjusting for the same factors as in the first method. Power of these analyses were similar for either approach. However, we found an excess of false-positive results when fitting the martingale residual model or the affected-only model to identify genetic factors linked to chromosome 6. Applying a power transformation to the martingale residuals decreased the type I error rate and increased the power of tests for genetic linkage. We also found that robust variance correction lead to test with a slightly lower type I error rate, perhaps because the robust variance correction adjusts for the fact that we did not specifically model the effects of the mitochondrial factor in our analysis.

A nonparametric dynamic additive regression model for longitudinal data

In this work we study additive dynamic regression models for longitudinal data. These models provide a flexible and nonparametric method for investigating the time-dynamics of longitudinal data. The methodology is aimed at data where measurements are recorded at random time points. We model the conditional mean of responses given the full internal history and possibly time-varying covariates. We derive the asymptotic distribution for a new nonparametric least squares estimator of the cumulative time-varying regression functions. Based on the asymptotic results, confidence bands may be computed and inference about time-varying coefficients may be drawn. We propose two estimators of the cumulative regression function. One estimator that involves smoothing and one that does not. The latter, however, has twice the variance as the smoothing based estimator. Goodness of fit of the model is considered using martingale residuals. Finally, we also discuss how partly-conditional mean models in which the mean of the response is regressed onto selected time-varying covariates may be analysed in the same framework. We apply the methods to longitudinal data on height development for cystic fibrosis patients.

Standardized martingale residuals applied to grouped left truncated observations of dementia cases.

The use of martingale residuals have been proposed for model checking and also to get a non-parametric estimate of the effect of an explanatory variable. We apply this approach to an epidemiological problem which presents two characteristics: the data are left truncated due to delayed entry in the cohort; the data are grouped into geographical units (parishes). This grouping suggests a natural way of smoothing the graph of residuals which is to compute the sum of the residuals for each parish. It is also natural to present a graph with standardized residuals. We derive the variances of the estimated residuals for left truncated data which allows computing the standardized residuals. This method is applied to the study of dementia in a cohort of old people, and to the possible effect of the concentration of aluminum and silica in drinking water on the risk of developing dementia.

Model choice and influential cases for survival studies.

Deletion diagnostics are developed for identifying observations that influence the estimates of regression parameters and the mixture parameter in the families of relative risk functions for failure time data. The diagnostic for the regression parameters is a generalization of Cain and Lange's (1984, Biometrics 40, 493-499) measure of individual influence. The generalizations of martingale residuals, Schoenfeld's partial residuals (1982, Biometrika 69, 239-241), and score residuals by Therneau, Grambsch, and Fleming (1990, Biometrika 77, 147-160) are also obtained. The influence of some observations on regression parameters can be drastically modified as the mixture parameter changes, even for a very small change. In addition, adding or deleting some observations might result in choosing different models. The diagnostics are applied to a family proposed by Guerrero and Johnson (1982, Biometrika 69, 309-314). One illustrative example is presented.

A Goodness-of-Fit Test for Cox's Proportional Hazards Model Based on Martingale Residuals

A goodness-of-fit test for Cox's proportional hazards model is proposed on the basis of martingale residuals. The test is derived as the score test for the presence of an extra random effect in the exponential part of the hazard function. Assumptions about the distribution of the random effects are not made; only the correlation matrix must be specified. This can be done on the basis of an existing division of individuals into groups, such as families, or on the basis of the distance between observed covariates. The variance of the test statistic is derived using a counting process approach. Corrections for the estimation of parameters are available. The test is demonstrated in a simulated data set as well as in a real data set on ovarian cancer.

Two-sample goodness-of-fit tests for additive risk models with censored observations

SUMMARY The additive risk model assumes that the hazard function associated with a set of covariates is the sum of the baseline hazard function and the regression function of covariates. We propose two different test procedures for checking the adequacy of two-sample additive risk models for randomly censored observations. One is based on the martingale residuals and the other on the difference between weighted estimators of the excess risk. The test statistics are shown to be asymptotically normal under appropriate regularity conditions and consistent under any model misspecifications. Finally, two real examples are provided, along with results of a simulation study.

Tests of Independence for Bivariate Survival Data

We propose two test statistics based on the covariance process of the martingale residuals for testing independence of bivariate survival data. The first test statistic takes the supremum over time of the absolute value of the covariance process, and the second test statistic is a time-weighted summary of the process. We derive asymptotic properties of the two test statistics under the null hypothesis of independence. In addition, we derive the asymptotic distribution of the weighted test and construct optimal weights for contiguous alternatives to independence. Through simulations, we compare the performance of the proposed tests and the inner product of the Savage scores statistics of Clayton and Cuzick (1985, Journal of the Royal Statistical Society, Series A 148, 82-108). These demonstrate that the supremum test is generally more powerful with comparatively little power loss relative to their test when Clayton's family alternative holds, and the weighted test is more powerful when the weight is appropriately chosen.