|Home | About | Journals | Submit | Contact Us | Français|
Early detection and effective treatments have dramatically improved breast cancer survivorship, yet the risk of relapse persists even 15 years after the initial diagnosis. It is important to identify prognostic factors for late breast cancer events. The authors investigated time-varying effects of tumor characteristics on breast-cancer-free survival using data on 3,088 breast cancer survivors from 4 US states who participated in a randomized dietary intervention trial in 1995–2006, with maximum follow-up through 15 years (median, 9 years). A piecewise constant penalized spline approach incorporating time-varying coefficients was adopted, allowing for deviations from the proportional hazards assumption. This method is more flexible than standard approaches, provides direct estimates of hazard ratios across time intervals, and is computationally tractable. Having a stage II or III tumor was associated with a 3-fold higher hazard of breast cancer than having a stage I tumor during the first 2.5 years after diagnosis; this hazard ratio decreased to 2.1 after 7.7 years, but higher tumor stage remained a significant risk factor. Similar diminishing effects were found for poorly differentiated tumors. Interestingly, having a positive estrogen receptor status was protective up to 4 years after diagnosis but detrimental after 7.7 years (hazard ratio = 1.5). These results emphasize the importance of careful statistical modeling allowing for possibly time-dependent effects in long-term survivorship studies.
Breast cancer is a heterogeneous disease comprising multiple subtypes that exhibit varying treatment response and survival rates. Known prognostic factors include tumor size, number of lymph nodes affected, histologic grade, and estrogen receptor (ER) and progesterone receptor (PR) status (1). Only the smallest tumors are considered “cured” after primary treatment, and recurrence risk remains elevated for at least 15 years postdiagnosis (2). Since early detection initiatives have increased survivorship (3), treatment focus has shifted towards characteristics of the original tumor that are predictive of survival, and the time frame has changed from short-term to 15 or more years (4). With long-term follow-up, the possibility of time-dependent risk should be considered, necessitating more complicated statistical modeling than standard approaches based on Cox proportional hazards models.
In several predominantly lymph-node-positive clinical cohorts, women with ER-positive (ER+) tumors had greater recurrence-free survival than those with ER-negative (ER−) tumors for the first 2 years after diagnosis, but this effect reversed (5) or diminished (6) at later times. Dignam et al. (7) reported similar results and noted that the hazard rate for ER− patients declined dramatically after the initial 2-year period, even in patients receiving surgery alone. Other studies have suggested that S-phase fraction (8), tumor size (8), tumor grade (9), and nodal status (9) have time-dependent effects on prognosis.
The issue of proportional hazards violation in Cox models has generated considerable interest in the statistics literature. Many authors recommend using stratified analyses when the proportional hazards assumption is not met (10). Although the stratified approach is simple to implement, it is unsatisfactory when the hypothesis of interest involves a variable that violates the proportional hazards assumption. Other suggestions include graphical tests (10), change-point models (11), adding a time-by-covariate interaction term in the Cox model (12, 13), or using fractional polynomials (14). Gray's penalized-spline approach (5) provides flexibility by modeling effects of prognostic variables as smooth functions over time and thus does not require a proportional hazards assumption. Most importantly, this method partitions the time scale into subintervals and provides direct estimates of hazard ratios (with 95% confidence intervals) for each subinterval. These hazard ratio estimates are useful for identifying time periods in which risk estimates change, which could inform decisions on treatment and disease management.
In this paper, we utilize Gray's penalized splines (5) to make time-dependent inferences regarding tumor characteristics using a large database from the Women's Healthy Eating and Living (WHEL) Study (15, 16). This study enrolled 3,088 breast cancer patients diagnosed between 1991 and 2000 (median year, 1997) and had follow-up data on 96% of participants through June 2006. Participants could enroll in the WHEL Study up to 4 years after diagnosis. Although the study is not representative of women whose breast cancer recurred in the first years following diagnosis, its detailed long-term follow-up data provide a unique opportunity for examining factors associated with “late” breast cancer events. A total of 518 WHEL Study participants had a confirmed second breast cancer event (16).
The WHEL Study was a multisite randomized trial that tested the efficacy of a diet high in fruit, vegetables, and fiber and low in fat in reducing new breast cancer events and improving survival among 3,088 women previously diagnosed with early-stage breast cancer. Details on the study design, protocol, and overall effect of the dietary intervention on the outcome have been previously reported (15, 16). Briefly, between 1995 and 2000, the WHEL Study recruited participants from 7 clinical sites in California, Oregon, Arizona, and Texas according to the following eligibility criteria: 1) having been diagnosed within the past 4 years with primary operable invasive stage I (≥1 cm), stage II, stage IIIA, or stage IIIC breast carcinoma (17); 2) having been aged 18–70 years at diagnosis; 3) having been treated with axillary dissection and total mastectomy or lumpectomy followed by primary breast radiation; 4) not being scheduled for or currently undergoing chemotherapy; 5) having had no evidence of recurrent or new breast cancer since completion of the initial treatment; and 6) having had no other invasive cancer in the past 10 years. Participants were randomized to either an intensive dietary intervention arm or a comparison arm (15).
The institutional review boards of all participating institutions approved the procedures for the WHEL Study, and written informed consent was obtained from all study participants before enrollment.
Data on the original breast cancer (diagnosis date, tumor stage, tumor grade, ER/PR status, and treatment) were obtained and confirmed via medical records. Reported breast cancer recurrences or new primary breast cancers diagnosed after study enrollment were independently adjudicated by 2 breast oncologists. Recurrent cancers were classified as local/regional or distant metastasis; ductal carcinoma in situ was not counted as a study outcome. The breast cancer event-free interval was defined as the time from diagnosis to the development of an additional breast cancer. Follow-up time was censored at the time of death (if not from breast cancer), at the last documented staff contact date, or at study completion (June 1, 2006).
To examine time-varying effects of covariates, we fitted Gray's spline model (5), which is specified as
where hi(t) is the hazard function for subject i, β(t)′Zi(t) = ∑βj(t)Zij(t), and Zij(t) equals the value of the jth covariate for subject i at time t. Here we focus on time-fixed covariates, so for each i and j, Zij(t) = Zij. In the Cox proportional hazards model, βj(t) = βj—that is, covariate effects are assumed to be independent of time. In Gray's approach, βj(t) = ∑θjkBjk(t), where Bjk(t) is a basis of spline functions, representing possibly time-varying effects for covariate Zj, and k indexes knots (i.e., join-points of spline sections). We used piecewise constant splines and varied the number of knots between 4 and 10, with knots placed to yield equal numbers of events in each time interval. Parameters were estimated from penalized partial likelihoods, and hypothesis tests used penalized versions of Wald statistics (5).
We conducted preliminary univariate analyses evaluating the effects of demographic factors (age at diagnosis, education, ethnicity), tumor characteristics (stage, grade, type, ER/PR status), treatment (radiation, chemotherapy, use of antiestrogen therapy), clinical site, and randomization arm (intervention vs. comparison) on breast-cancer-free survival in separate spline models, allowing for time-varying effects for each covariate. Variables that were significantly associated with the outcome at the 0.05 significance level or which violated the proportional hazards assumption at the 0.1 level were then included in multiple regression models; we adopted a liberal criterion for nonproportionality, because a primary focus of this analysis was to elucidate time-varying effects of prognostic factors.
We then developed multiple regression models using univariate results. All models were stratified by time from diagnosis to study entry using four 1-year time intervals. To compare standard approaches with the spline method, we fitted Cox models with time-fixed coefficients and Gray's spline model with time-varying coefficients. Cox models were fitted under 3 scenarios: 1) the “usual” approach encompassing the entire follow-up period; 2) follow-up truncated at ≤5 years since diagnosis (“early model”); and 3) follow-up conditional on surviving more than 5 years from diagnosis (“late model”). The early and late models arbitrarily used a cutpoint of 5 years and provided a crude assessment of changes in risks associated with prognostic factors.
Table 1 gives demographic, medical, and treatment characteristics of the study sample. WHEL women had an average age of 51 years (standard deviation, 8.8; range, 26–71) at diagnosis and were primarily Caucasian (85%) and college graduates (54%). The prevalence of stage I tumors was 38.6%; 15.7% of tumors displayed well-differentiated pathology; 63.1% of tumors were ER+ and PR+; and 85.9% of tumors were of the ductal (with or without lobular) type. Regarding treatment, 52% of WHEL women underwent a mastectomy (the remainder underwent lumpectomy with radiation therapy), 70% received adjuvant chemotherapy, and 68% received antiestrogen therapy (66% took tamoxifen, 1% took raloxifene, 0.1% took anastrazole/arimidex, and for 1% the type of antiestrogen medication prescribed was unknown). Median time from diagnosis was 8.99 years (range, 0.79–15.01).
In the univariate analysis, age at diagnosis, education, tumor stage, tumor type, tumor grade, ER/PR status, receipt of chemotherapy, receipt of antiestrogens, and clinical site were significantly associated with outcomes (P < 0.05), while age at diagnosis, tumor stage, tumor type, tumor grade, ER/PR status, and antiestrogen therapy exhibited non-proportional-hazards effects (P < 0.1). These results guided the development of multiple regression models (Tables 2 and and33).
In Cox models (Table 2), stage II or III tumors were associated with the highest risk of breast cancer events in the early model; this effect appeared to diminish with time, although it was still a major risk factor in the late model. Interestingly, although the finding was not statistically significant, having an ER+ tumor was protective against early disease but was associated with poorer prognosis more than 5 years postdiagnosis. Similar time-related inversions of risks were apparent for tumor type; that is, having a ductal tumor was associated with higher risk for earlier disease in our sample but not for risk more than 5 years postdiagnosis. As expected, the log hazard ratios in the overall model were a weighted average of the log hazard ratios from the early and late models. Older age at diagnosis, having completed college, and having received antiestrogen therapy exhibited similar protective effects across all 3 models.
The Cox models in Table 2 involved stratification on each year from diagnosis to study entry, with different baseline hazard functions being computed for these 4 strata; thus, varying recurrence risks were modeled by year from diagnosis. Additionally, we fitted delayed-entry Cox models (10) to further examine possible biases due to variability in time from diagnosis to study entry. The results (not presented) were almost identical to those of Table 2, with changes of less than 1% in the delayed-entry model hazard ratios, suggesting that this source of bias was not a factor in our analyses.
Time-varying coefficient models (5) further explicate the Cox model findings by formally testing the proportional hazards assumption. They permit more flexibility by allowing the estimation of hazard ratios across shorter time intervals. Table 3 and Figure 1 present results from the time-varying coefficient model with 4 knots and 2 degrees of freedom. Hazard ratios for education, age at diagnosis, PR status, and receipt of antiestrogen therapy were constant over time and similar to those for the overall Cox model (Table 2). Tumor stage, ER status, and tumor grade violated the proportional hazards assumption. Higher stage was associated with a 3-fold greater hazard of breast cancer events than a stage I tumor during the first 2.5 years after diagnosis; this effect diminished to a hazard ratio of 2.1 after 7.7 years, but higher stage remained a significant risk factor. Similarly, poorly differentiated tumors were associated with higher risk (hazard ratios ranging from 1.25 to 1.5) of a breast cancer event up to 5.5 years postdiagnosis, with no apparent effect after 7.7 years. Interestingly, although it was not statistically significant at the 5% level in each time interval, having ER+ status was protective up to 4 years after diagnosis; this trend was reversed in later years, with ER+ status exhibiting a hazard ratio of 1.5 after 7.7 years.
We also fitted a time-varying coefficient model with 10 knots. The results (not shown) were similar to those of Table 3 for ER status and stage, with ER+ status conferring striking protection early (for <1.7 years since diagnosis, hazard ratio = 0.63, 95% confidence interval: 0.4, 0.98) and detrimental effects later (for >9.5 years since diagnosis, hazard ratio = 2.01, 95% confidence interval: 1.22, 3,30). Interestingly, in the 10-knot model, tumor type violated the proportional hazards assumption, with hazard ratios increasing up to 5 years postdiagnosis and decreasing thereafter for ductal tumors versus nonductal tumors.
In additional sensitivity analyses, we included other covariates (e.g., ethnicity, receipt of radiation or chemotherapy, clinical site, and diet group (intervention vs. comparison)) in the model. These covariates were excluded from the final parsimonious models because either they were not significantly associated with outcomes or including them did not meaningfully alter the hazard ratios for other factors in the model.
It is important to contrast the Cox model findings (Table 2) with those from the spline approach (Table 3). At first glance, both give similar results—namely, that ER status and stage exhibit time-dependent effects on disease-free survival. However, the spline approach affords many advantages. First, it formally tests the proportional hazards assumption and identifies variables that violate it. Second, the spline approach does not arbitrarily choose a time point of 5 years to delineate early versus late recurrence. Instead, this method partitions the time axis on the basis of knot locations and estimates hazard ratios for each prognostic factor on each time subinterval. For instance, using splines, WHEL data suggest that the risk of second breast cancer events is 50% higher for poorly differentiated tumors in the first 2.5 years after diagnosis, whereas ER+ tumors appear to be associated with higher risk after 7.7 years. Thus, Gray's method allows for the possibility that time periods in which hazard ratios change are different for different prognostic factors, rather than imposing an a priori cutpoint of 5 years. Third, this approach also provides finer estimates of hazard ratios (with 95% confidence intervals) for each time subinterval, whereas the early and late models simply provide “averaged” hazard ratio estimates for time intervals of ≤5 years versus >5 years. In particular, the lack of a significant effect for poorly differentiated tumors in the early model (Table 2) is probably due to this factor's only having an impact on very early breast cancer events (<2.5 years), so that using a 5-year time scale nullifies the risk. Finally, most analyses of such data would simply use the Cox model for the entire follow-up period (i.e., our “overall” model), which would lead to potentially incorrect conclusions regarding the prognostic effects of ER status.
Modeling the time course of prognostic factors in cancer survival is important from a clinical standpoint, yet it poses statistical challenges. Much clinical and computational research has focused on this issue. Illustrating his analytic method, Gray (5) showed that ER+ tumors conferred significantly lower recurrence risk in the first 2 years postdiagnosis than ER− tumors, but beyond 4 years the risk profiles reversed, with ER+ patients having greater recurrence risk. Similar results on time-varying effects of ER+ tumors have been reported by others (6–8). Hilsenbeck et al. (8) also found that S-phase fraction and tumor size had diminishing detrimental effects on disease-free survival among 2,873 breast cancer patients whose tumors were diagnosed and resected between 1978 and 1992, with follow-up ranging from 1 month to 207 months (median, 44 months). Additionally, using a parametric Cox model with exponential decay terms to model time-varying effects among 2,299 breast cancer cases followed through 1998, Warwick et al. (9) concluded that higher tumor grade and positive lymph node status were associated with large risks of breast cancer death during the first 5 years postdiagnosis, whereas these risks were significantly reduced after 15 or 20 years. Finally, Brewster et al. (18) reported that stage II–III, grade 1, and hormone-receptor-positive tumors conferred higher risks of late recurrence (>5 years) among 2,838 breast cancer patients who received neoadjuvant systemic therapy between 1985 and 2001.
In the current analysis, we used the sample of WHEL Study women diagnosed and treated for stage I–III breast cancer between 1991 and 2000, with up to 15 years of follow-up from diagnosis. Our results regarding time-varying effects of cancer stage, ER status, and tumor grade on breast cancer outcomes are similar to previously published findings (6–9, 18). Since treatment options for breast cancer have evolved and improved dramatically over the last 20 years, it is important to evaluate the effects of prognostic factors on survival in multiple cohorts, diagnosed in different time periods, with long-term follow-up.
Importantly, our statistical approach differed from that of most previous studies (6–9, 18). We utilized splines (5), which are more flexible than the semiparametric methods of Warwick et al. (9) and Brewster et al. (18). The choice of a particular time interval from diagnosis (e.g., 5 or 10 years) is arbitrary from a statistical and clinical standpoint, since the time courses of different prognostic factors are likely to exhibit different patterns. The time-varying coefficient model we adopted addresses this issue by estimating hazard ratios across multiple time intervals defined by knot locations. For instance, our data suggest that higher stage is associated with more than 2-fold risk of additional breast cancer events across all time intervals, whereas poorly differentiated tumors appear to increase risk by 25%–50% only during the first 5 years postdiagnosis. It is important not to overinterpret the time intervals in our analysis but instead use them as a rough guide for gauging when prognostic effects might be changing.
Hess et al. (6), Dignam et al. (7), and Hilsenbeck et al. (8) adopted graphical approaches and tests based on smoothed, scaled Schoenfeld residuals. These methods, which include as special cases the addition of time-by-covariate interaction terms to the Cox model (10, 12, 19), provide excellent visual aids and in some cases statistical tests for assessing the proportional hazards assumption, yet these methods tend to be sensitive to transformations of the time scale (10). In addition, natural cubic splines (6, 8) were used for obtaining smoothed estimates of time-varying hazard ratios. A variety of spline-based methods have been proposed, ranging from entirely nonparametric models with knots at all data points (20) to regression splines without penalty functions (21). We prefer Gray's approach (5), which is low-dimensional, thus avoiding the computational complexity of fully nonparametric models, and incorporates penalty terms for the amount of smoothing and hence is not as sensitive to the number and location of knots as regression splines. The parametric penalized Gray method is best suited to our application, since it is more flexible than standard approaches yet is computationally tractable.
The consistent finding that ER+ tumors confer protection early but could be associated with greater risk in later years is of particular interest. The WHEL Study sample, which recruited women up to 4 years postdiagnosis, possibly underrepresents ER− cancers with early events. Nevertheless, our results on the association between ER status and late recurrences should be valid. This association may reflect the clinical practice of initiating tamoxifen therapy at the completion of primary treatment(s) and continuing therapy for 5 years. The discontinuation of antiestrogen therapy may have resulted in an increased risk of “late” disease among ER+ women. It could also be postulated that tumors are not 100% ER+, and administering antiestrogen therapy while controlling ER+ recurrent disease may result in a greater risk of ER− recurrent disease. Furthermore, future studies might examine time to recurrence in ER+ women treated with tamoxifen who were not able to metabolize the drug extensively because of genetic factors (e.g., the cytochrome P450 2D6 (CYP2D6) gene). The role of ER status in breast cancer mortality also needs examination. Having an ER+ tumor suggests biologic susceptibility to breast cancer and hence may increase recurrence risk throughout the woman's life, yet ER+ tumors may not necessarily confer additional risk of mortality. Regardless, this prognostic factor should be carefully modeled in long-term follow-up studies; not allowing for time-varying effects would lead to a possibly null ER+ effect (canceling the early vs. late risks), as observed in the overall Cox model (Table 2).
The statistical literature is replete with competing analytic methods for modeling time-varying effects in survival studies. Other modeling paradigms such as accelerated failure time models (22), additive hazard models, Bayesian methods (23, 24), reduced rank methods (13), or Cox models with fractional polynomial time transformations (14) might be viable alternatives to the approach we adopted. However, we did not pursue these alternate approaches, since our goal was to demonstrate general uses of the popular Cox model even when standard modeling assumptions fail. Comparing goodness of fit among competing modeling methods for time-varying effects was beyond the scope of this work.
A limitation of our study is the time frame of data collection. Because the initial breast cancer occurred between 1991 and 2000, we did not have information on molecular subtypes of the original tumor. Knowledge of human epidermal growth factor receptor 2 (Her2) status would permit further classification of ER+ tumors as Her2− (luminal A) or Her2+ (luminal B), which may refine predictive value (1, 25). It is also possible that risks associated with certain subtypes change over time (26). However, not including this finer classification of ER status in the models should not nullify the consistent findings regarding the time-varying effects of other tumor characteristics such as stage and grade. At study entry, WHEL participants were within 4 years of their primary breast cancer diagnosis, had completed initial treatment, and had no evidence of recurrent disease. Thus, the sample was left-truncated (10) and biased towards “healthy survivors.” This selection bias could have caused us to underestimate the effects of prognostic factors in the models during the first 4 years after diagnosis but would not have affected risk estimates beyond this initial time frame. The spline approach, although flexible, is limited by the diminishing number of breast cancer events as time since diagnosis increases. With 4 knots, we would have approximately 100 events in each time interval, an arguably small number with which to precisely estimate hazard ratios and avoid problems with overfitting the data. Thus, it is vitally important that these results be replicated in large breast cancer cohorts with long-term follow-up.
In conclusion, this analysis reiterates the importance of careful statistical modeling of prognostic factors in survival studies. Precise risk estimates of the effects of primary tumor characteristics on late breast cancer events are needed to ensure adequate statistical power when designing long-term cancer survivorship studies. Treatment options and strategies for reducing long-term risks associated with ER+ and higher-stage tumors need to be further investigated. In future studies, investigators using molecular subtyping derived from gene-based classification schemes might apply flexible statistical models such as those used herein to improve predictive value and identify subsets of patients who are at high risk for early recurrence as opposed to late recurrence.
Author affiliations: Rebecca and John Moores UCSD Cancer Center, School of Medicine, University of California, San Diego, La Jolla, California (Loki Natarajan, Minya Pu, Barbara A. Parker, Shirley W. Flatt, Lisa Madlensky, Wael K. Al-Delaimy, Nazmus Saquib, John P. Pierce); Arizona Cancer Center, University of Arizona, Tucson, Arizona (Cynthia A. Thomson); Division of Research, Kaiser Permanente Northern California, Oakland, California (Bette J. Caan); M. D. Anderson Cancer Center, University of Texas, Houston, Texas (Richard A. Hajek); and Department of Public Health Sciences, School of Medicine, University of California, Davis, Davis, California (Ellen B. Gold).
The Women's Healthy Eating and Living (WHEL) Study was initiated with the support of the Walton Family Foundation and was continued with funding from the National Cancer Institute (grant CA 69375). Some of the data were collected from General Clinical Research Centers under National Institutes of Health grants M01-RR00070, M01-RR00079, and M01-RR00827.
WHEL Study Coordinating Center: Cancer Prevention and Control Program, Rebecca and John Moores UCSD Cancer Center, School of Medicine, University of California, San Diego, La Jolla, California (Dr. John P. Pierce, Susan Faerber, Dr. Barbara A. Parker, Dr. Loki Natarajan, Dr. Cheryl L. Rock, Vicky A. Newman, Shirley W. Flatt, Sheila Kealey, Dr. Linda Wasserman, Dr. Wayne A. Bardwell, Dr. Lisa Madlensky, and Dr. Wael Al-Delaimy).
WHEL Study clinical sites: Center for Health Research–Portland, Portland, Oregon (Drs. Njeri Karanja and Mark U. Rarick); Kaiser Permanente Northern California, Oakland, California (Drs. Bette J. Caan and Lou Fehrenbacher); Stanford Prevention Research Center, Department of Medicine, School of Medicine, Stanford University, Stanford, California (Drs. Marcia L. Stefanick and Robert Carlson); Arizona Cancer Center, University of Arizona Health Sciences Center, Tucson, Arizona (Drs. Cynthia Thomson, James Warneke, and Cheryl Ritenbaugh); Department of Public Health Sciences, School of Medicine, University of California, Davis, Davis, California (Drs. Ellen B. Gold and Sidney Scudder); Rebecca and John Moores UCSD Cancer Center, School of Medicine, University of California, San Diego, La Jolla, California (Drs. Kathryn A. Hollenbach and Vicky Jones); and M. D. Anderson Cancer Center, University of Texas, Houston, Texas (Drs. Lovell A. Jones, Richard Hajek, and Richard Theriault).
Conflict of interest: none declared.