Time-varying effect in older patients with early-stage breast cancer: a model considering the competing risks based on a time scale

Background: Patients with early-stage breast cancer may have a higher risk of dying from other diseases, making a competing risks model more appropriate. Considering subdistribution hazard ratio, which is used often, limited to model assumptions and clinical interpretation, we aimed to quantify the effects of prognostic factors by an absolute indicator, the difference in restricted mean time lost (RMTL), which is more intuitive. Additionally, prognostic factors of breast cancer may have dynamic effects (time-varying effects) in long-term follow-up. However, existing competing risks regression models only provide a static view of covariate effects, leading to a distorted assessment of the prognostic factor.

Methods: To address this issue, we proposed a dynamic effect RMTL regression that can explore the between-group cumulative difference in mean life lost over a period of time and obtain the real-time effect by the speed of accumulation, as well as personalized predictions on a time scale.

Results: A simulation validated the accuracy of the coefficient estimates in the proposed regression. Applying this model to an older early-stage breast cancer cohort, it was found that 1) the protective effects of positive estrogen receptor and chemotherapy decreased over time; 2) the protective effect of breast-conserving surgery increased over time; and 3) the deleterious effects of stage T2, stage N2, and histologic grade II cancer increased over time. Moreover, from the view of prediction, the mean C-index in external validation reached 0.78.

Conclusion: Dynamic effect RMTL regression can analyze both dynamic cumulative effects and real-time effects of covariates, providing a more comprehensive prognosis and better prediction when competing risks exist.

1 Background

Older patients with early-stage breast cancer (particularly when comorbidities are advanced) tend to die from other diseases. That is, the number of deaths from non-breast cancer is large (1). Patients may experience a variety of outcomes: death from breast cancer, death from heart disease and so on. Assuming the outcome we are interested in is death from breast cancer (event of interest), then we hope to observe the time from the start of follow-up to the occurrence of the event of interest, but this will not be observed for the patient dying from heart disease (competing event). In this case, traditional single-endpoint survival analysis, such as the Cox proportional hazards model, only considers the event of interest and treats patients who die of heart disease as censored simply. However, this doesn’t meet the non-informative censoring hypothesis. That is, the risk of dying from breast cancer among women who have already died of heart disease needs to be the same as that among women who remain in follow-up. However, patients who have already died of heart disease will not die of breast cancer again. So it is improper to simply treat patients who experience competing events as censored, which will overestimate the cumulative incidence of the event of interest and result in bias (2–4). To solve the situation where multiple outcomes compete with each other, we should consider competing risks models.

In traditional multivariate analysis of competing risks, cause-specific Cox regression and Fine-Gray regression are often used, and the corresponding effect sizes are the cause-specific hazard ratio (cHR) and subdistribution hazard ratio (sHR). However, both the cHR and sHR are relative indicators, defined as the ratio of the hazard function. It is difficult for clinicians to interpret them as intuitive clinical benefits and communicate with patients (5, 6). For example, when the sHR of estrogen receptor (ER) is 0.43, that is, the risk of death in the ER-positive group is 0.43 times that in the ER-negative group, and because the baseline hazard is unknown generally, the absolute risk of death of two groups cannot be known. In addition, both cause-specific Cox regression and Fine-Gray regression models need to satisfy the proportional hazards assumption. From a clinical perspective, clinicians or patients are more interested in direct (absolute) effect sizes on a time scale, e.g., how long will I live? How long will surgery extend my life expectancy? As Blagoev (6) points out, “While a hazard ratio has some value, for the clinician caring for a patient and, more importantly, the patient, it does not convey benefit in terms that are meaningful—how much longer will the patient live or live without experiencing disease progression.” Therefore, restricted mean time lost (RMTL) has been proposed as an alternative measure to the hazard ratio (7–12). The RMTL is the area under the cause-specific cumulative incidence function (CIF) over a period of time (from 0 to a restricted time point $\tau$), which can be interpreted as the life expectancy lost due to a specific cause in this period of time. Compared with the hazard ratio, the interpretation of RMTL is more intuitive, giving the life lost of each group due to death from breast cancer over a period of time and measuring the effect of a factor by the difference in RMTL between groups. For example, Figure 1 shows areas under the CIF for patients in the ER-positive group and ER-negative group. During the 10.5 years ($\tau =10.5$), the mean life lost due to death from breast cancer was 0.9 ($S_0$) years for patients in the ER-positive group and 2.4 ($S_0+S_1$) years for patients in the ER-negative group, and ER-negative patients lost an additional 1.5 (2.4–0.9) years of life on average. At the same time, RMTL does not need to meet the proportional hazards assumption.

Figure 1

Figure 1 Cumulative incidence curves for death from breast cancer in the ER-positive group and ER-negative group. $S_0,S_1$ correspond to the blue and red areas, respectively.

The existing multivariate analysis of RMTL includes two methods: regression based on the pseudo-value method (11) and regression based on inverse probability of censoring weighting(IPCW) (12). The regressions mentioned above only concern the cumulative effects of prognostic factors during a $\tau$-year follow-up, which are constant values (static or time-fixed effects). However, the real-time effects of many covariates (e.g.,: covariates with time-varying effects) vary, for example, the real-time effect of chemotherapy tends to decrease with increasing follow-up time. It has been documented that the effects of age, histological grade, and ER status on the survival of patients with breast cancer change over time (13, 14), so it may not be comprehensive to fit only the static effects of prognostic factors.

Given covariates with time-varying effects in the field of competing risks, we proposed a dynamic effect RMTL regression model. Monte Carlo simulation was used to assess the accuracy of the coefficient estimates of the model. At the same time, we applied this model to older patients with early-stage breast cancer in the Surveillance, Epidemiology, and End Results (SEER) database to explore the dynamic cumulative effects and real-time effects of prognostic factors, as well as to establish a prediction model to predict mean life lost due to death from breast cancer over a period of time among patients. We hope to guide doctors to better determine the prognosis of patients, select better therapeutic regimens, and improve the survival time of patients.

2 Methods

2.1 Model construction

Let $T$ be the time to event and C be the censoring time so that the observed time is $U=\min (T,C)$. An event indicator ε equals 1, 2, or 0 when the observed outcome is an event of interest, a competing event, or censoring, respectively. At the same time, let ${\mathit{\text{Z}}}^{\ast }=(1,Z)$ denote the $n\times (p+1)$ matrix of covariates allowing an intercept term. Thus, for patient $i(i=1,\dots ,n)$, the observed data include $\left\{U_i,{\epsilon }_i,Z_i^{\ast }\right\}$.

Let a continuous variable $l(0\le l\le \tau \text{, }\tau \le t_{\max })$ be the pre-specified end time of follow-up, where $\tau$ is the pre-specified maximum follow-up time, and t_max is the natural maximum follow-up time of data. J time points l_j are selected from 0 to $\tau$ in ascending order and recorded as $(l_1,l_2,\dots ,l_J)$.

For the need of the method, we advance the end time of follow-up from t_max to l_j. Correspondingly, each patient’s survival outcome will change at different pre-specified end times of follow-up. When a patient experienced the event of interest or the competing event before $l_j$, $\epsilon (l_j)$ is equal to 1 or 2; in other cases, $\epsilon (l_j)$ is equal to 0. $\epsilon (l_j)$ denotes the survival outcome after restraint. Let $T(l_j)=\min (T,l_j)$ and U(l_j) = min(U,l_j) be the event time and observed time after constraint, respectively. For patient $i$, the observed data consist of $\left\{U_i(l_j),{\epsilon }_i(l_j),Z_i^{\ast }\right\}$.

A regression model is developed to assess the dynamic effects of covariates in RMTL:

$$g\left\{{\mu }_1(l|Z_i^{\ast })\right\}=g\left\{E[(l-T_i(l))\times I({\epsilon }_i(l)=1)]\right\}=Z_i^{\ast }\beta (l)$$

with link function $g(·)$ using the identity function. ${\mu }_1(l|Z_i^{\ast })$ is the life expectancy lost due to the event of interest of patient i during the l-year follow-up. Regression coefficients $\beta (l)=({\beta }_1(l),{\beta }_2(l),\dots ,{\beta }_{p+1}(l))$, where ${\beta }_k(l)$ is defined as $({\beta }_{k0},{\beta }_{k1},{\beta }_{k2})\times (1,l,l^2)$.

Because $Z_{ik}^{\ast }\times {\beta }_k(l)=({\beta }_{k0},{\beta }_{k1},{\beta }_{k2})\times (Z_{ik}^{\ast },l\times Z_{ik}^{\ast },l^2\times Z_{ik}^{\ast })={\beta }_k\times Z_{ik}^{\ast }(l)$, the model can be rewritten as

$$g\left\{{\mu }_1(l|Z_i^{\ast })\right\}=Z_i^{\ast }(l)\beta$$

Regression coefficients are estimated by solving the estimating equation

$$\Phi (\beta )={\displaystyle \sum _{j=1}^J\sum _{i=1}^{\text{n}}I({\epsilon }_i(l_{\text{j}})\ne 0)Z_i^{\ast }(l_{\text{j}})\left\{(l_{\text{j}}-T_i(l_{\text{j}}))\times I({\epsilon }_i(l_{\text{j}})=1)-g^{-1}(Z_i^{\ast }(l_{\text{j}})\beta )\right\}}=0$$

The model assumes the life lost of individuals with the event of interest is $(l_{\text{j}}-T_i(l_{\text{j}}))\times I({\epsilon }_i(l_{\text{j}})=1)=l_{\text{j}}-T_i(l_{\text{j}})$; the life lost of individuals with the competing event is $(l_{\text{j}}-T_i(l_{\text{j}}))\times I({\epsilon }_i(l_{\text{j}})=1)=0$; and censored observations have $I({\epsilon }_i(l_{\text{j}})\ne 0)=0$, which means they do not contribute to the estimating equation. $E(\Phi (\beta ))\ne 0$ in the presence of censoring. However, when applying IPCW to the estimating equation, its expectation is 0 (15). Therefore, the estimating equation changes to

$$\Phi (\beta )={\displaystyle \sum _{j=1}^J\sum _{i=1}^{\text{n}}\frac{I({\epsilon }_i(l_{\text{j}})\ne 0)}{\widehat{G}(T_i(l_{\text{j}}),l_{\text{j}})}{Z}_i^{\ast }(l_{\text{j}})\left\{(l_{\text{j}}-T_i(l_{\text{j}}))\times I({\epsilon }_i(l_{\text{j}})=1)-g^{-1}(Z_i^{\ast }(l_{\text{j}})\beta )\right\}}=0.$$

The fitted data are actually obtained by stacking $J$ datasets. The j-th dataset is the risk set with l_j-year follow-up (as shown in Figure 2). $\widehat{G}(t,l_j)$ is the Kaplan-Meier estimator of the non-censoring distribution in the $j$-th dataset.

Figure 2

Figure 2 Composition of the stacked dataset.

We treat the inverse probability censoring weight as a fixed value rather than a random variable (16). Thus, the variance in the regression coefficients will not consider the variation brought by the weight. Therefore, we have $\sqrt{n}(\widehat{\beta }-\beta )\sim N(0,A^{-1}BA)$, $\widehat{V}ar(\widehat{\beta })={\widehat{A}}^{-1}\widehat{B}{\widehat{A}}^{-1}$, where

$$\widehat{A}=E\left[\sum _{j=1}^J{Z}_i^{\ast }{(l_j)}^{\otimes 2}h(Z_i^{\ast }(l_{\text{j}})\widehat{\beta })\right]$$

$$\widehat{\mathit{B}}=E\left[\sum _{j=1}^J{\epsilon }_{ij}{(\widehat{\mathit{\beta }})}^{\otimes 2}\right]$$

$${\epsilon }_{\text{ij}}(\widehat{\beta })=\frac{I({\epsilon }_i(l_{\text{j}})\ne 0)}{\widehat{G}(T_i(l_{\text{j}}),l_{\text{j}})}{Z}_i^{\ast }(l_{\text{j}})\left\{(l_{\text{j}}-T_i(l_{\text{j}}))\times I({\epsilon }_i(l_{\text{j}})=1)-g^{-1}(Z_i^{\ast }(l_{\text{j}})\widehat{\beta })\right\}$$

where $a^{\otimes 2}=a{a}^T$ and $h(x)=\partial g^{-1}(x)/\partial x$.

Then, we estimate regression coefficients by a generalized estimating equation, thereby correcting for data correlation.

2.2 Simulation designs

Next, we assessed the performance of the estimation of dynamic-effect RMTL regression by a simulation. We used the mean bias, mean relative bias, root mean squared error, relative standard error, and empirical coverage rate as evaluation indicators.

2.2.1 Data generation

First, we generated two independent variables Z = (Z₁,Z₂), which were generated by an independent Bernoulli distribution. We let the subdistribution hazard function for the event of interest follow a Gompertz distribution, ${\lambda }_1(t|Z)={\gamma }_z\exp ({\rho }_{z}t)$, where ${\gamma }_z$ and ${\rho }_z$ were set according to four strata of $(Z_1,Z_2)$.

We defined the CIF for the event of interest as $F_1(t|Z)=P(T\le t,\epsilon =1|Z)=1-\exp \left\{-{\int }_0^t{\lambda }_1(s|Z)ds\right\}$ and the CIF for the competing event as $F_2(t|Z)=P(T\le t,\epsilon =2|Z)=\exp ({\gamma }_z/{\rho }_z)\{1-\exp (-t)\}$. Survival outcome was generated by Bernoulli distribution, $P(\epsilon =1|Z)=1-\exp ({\gamma }_z/{\rho }_z)$. Thus, the conditional CIF for the event of interest and the competing event were $P(T\le t|\epsilon =1,Z)=[1-\exp \left\{-{\int }_0^t{\lambda }_1(s|Z)ds\right\}]/[1-\exp ({\gamma }_z/{\rho }_z)]$ and $P(T\le t|\epsilon =2,Z)=1-\exp (-t)$, respectively. Next, we used the inverse method to generate event time $T$. Finally, we generated right censoring and determined the final observed survival outcome.

We calculated the sHRs of independent variables based on the subdistribution hazard function: $sHR(Z_1,t)={\gamma }_{(1,0)}\exp ({\rho }_{(1,0)}t)/{\gamma }_{(0,0)}\exp ({\rho }_{(0,0)}t)$, $sHR(Z_2,t)={\gamma }_{(0,1)}\exp ({\rho }_{(0,1)}t)/{\gamma }_{(0,0)}\exp ({\rho }_{(0,0)}t)$. We let ${\rho }_{(1,0)}\ne {\rho }_{(0,0)},{\rho }_{(0,1)}\ne {\rho }_{(0,0)}$; then, the sHRs changed over time, and the proportional subdistribution hazards assumption of two variables was not satisfied.

2.2.2 Parameters, scenarios, and true values of simulation

Consider $({\gamma }_{(0,0)},{\gamma }_{(0,1)},{\gamma }_{(1,0)},{\gamma }_{(1,1)})=(2.88,1.95,2.29,1.55)$, $({\rho }_{(0,0)},{\rho }_{(0,1)},{\rho }_{(1,0)},{\rho }_{(1,1)})=(-1.7,-1.4,-2.9,-2.8)$ (12). The range of $l$ was between the 10th percentile and 95th percentile of the time of patients with the event of interest in the simulated data. Figure 3 shows that the sHRs of the independent variables changed over time. $Z_1$ was a protective factor in the early period and a risk factor in the late period, while the protective effect of Z₂ increased over time.

Figure 3

Figure 3 Subdistribution hazard ratios (sHRs) of two independent variables for the event of interest in the simulation.

Twelve scenarios were simulated considering varying sample sizes (250, 500, 1000), proportions of exposure (both $P(Z_1=1)$ and $P(Z_2=1)$ equal to 0.25 or 0.5), and censoring rates (0.1 or 0.25). Each scenario was simulated 2000 times.

By integrating the CIF, we obtained the true value of the RMTL of each group at different Z, ${\mu }_1(l|Z)={\int }_0^l{F}_1(t|Z)dt$. Moreover, the true values of regression coefficients were obtained by the difference in RMTL between groups: the baseline is ${\beta }_0(l)(={\mu }_1(l|z_1=0,z_2=0))$; the cumulative effect of Z₁ is ${\beta }_1(l)(={\mu }_1(l|z_1=1,z_2=0)-{\mu }_1(l|z_1=0,z_2=0))$; and the cumulative effect of Z₂ is ${\beta }_2(l)(={\mu }_1(l|z_1=0,z_2=1)-{\mu }_1(l|z_1=0,z_2=0))$. Taking $l=(0.75,1,1.5)$, then $({\beta }_0(0.75),{\beta }_0(1),{\beta }_0(1.5))=(0.368,0.550,0.937)$, $({\beta }_1(0.75),{\beta }_1(1),{\beta }_1(1.5))=(-0.093,-0.15,-0.284)$, $({\beta }_2(0.75),{\beta }_2(1),{\beta }_2(1.5))=(-0.074,-0.100,-0.146)$. We found that the regression coefficient ${\beta }_k(l)$ was a cumulative quantity, and its absolute value increased as l increased.

3 Results

3.1 Simulation results

Tables 1, 2 demonstrate the accuracy of $\widehat{\beta }(l)$ at different l. In all cases, the mean relative bias was small, in which the mean relative bias of ${\beta }_0(l)$ was less than 2%; the relative standard error was approximately 1; and the coverage rate was approximately 95%. Because the absolute value of true value of the regression coefficient increased with increasing l, it was reasonable that the mean bias increased with increasing l. The simulation showed that the estimation of dynamic effect RMTL regression was accurate.

Table 1

Table 1 Performance of dynamic-effect RMTL regression in the simulation when the proportion of exposure is 0.25.

Table 2

Table 2 Performance of dynamic-effect RMTL regression in the simulation when the proportion of exposure is 0.5.

With the increase in sample size, the mean bias, mean relative bias, and root mean squared error were more likely to decrease. Moreover, different censoring rates had little effect on the mean relative bias and root mean squared error.

3.2 Model application

In this study, we extracted data from the SEER database for older patients with early-stage breast cancer.

Covariates included race, age, marriage, T stage, N stage, histological grade, estrogen receptor (ER) status, progesterone receptor (PR) status, breast surgery, axillary surgery, chemotherapy, and radiotherapy (17). Breast surgery includes mastectomy or breast-conserving surgery (BCS); while axillary surgery includes axillary lymph node dissection (ALND) and lymph node biopsy (SLNB).

We used 3892 patients diagnosed from 2000 to 2012 as a training set and another 1561 patients diagnosed from 2013 to 2015 as an externally validated set. Details of data collection and variables can be found in the Supplementary Material.

There were 769 deaths from breast cancer and 998 deaths from non-breast cancer in the training set, giving an approximately 55% censoring rate. The follow-up time ranged from 0.17 to 18.92 years, with a median of 8 years.

Table 3 shows the result of the static effect RMTL regression ($\tau =\mathbf{10.5}$ years) (12). The regression coefficient β indicates a cumulative difference in mean life lost during the 10.5-year follow-up due to death from breast cancer between groups of the prognostic factor. For example, ${\beta }_{ER}=-0.638$ showed that patients in the ER-positive group died of breast cancer 0.638 years later than those in the ER-negative group during the 10.5 years, so ER positivity was a protective factor. In general, a prognostic factor was protective when β was negative and deleterious when β was positive. Table 3 shows that patients with ER positivity, PR positivity, breast-conserving surgery (relative to mastectomy), and chemotherapy had a better prognosis, while patients with older age, higher T stage, higher N stage, and higher histological grade had a worse prognosis. Race, marriage, axillary surgery, and radiation therapy had no statistical significance on survival.

Table 3

Table 3 Regression coefficients of static-effect RMTL regression (τ = 10.5 years).

The significant covariates did not meet the proportional subdistribution hazards assumption, indicating that dynamic effects might exist. However, static effect RMTL regression only gives the static cumulative effect, and the real-time effect in the cumulative process cannot be known. Therefore, we fitted the proposed dynamic effect RMTL regression and used a backward stepwise approach to screen covariates. As a result, race, marriage, axillary surgery, and radiation were screened out.

Table 4 shows the results of the dynamic effect RMTL regression. Because the cumulative effect of the k-th prognostic factor was assumed by ${\beta }_k(l)=({\beta }_{k0},{\beta }_{k1},{\beta }_{k2})\times (1,l,l^2)$, which would be screened by the stepwise method, the regression coefficients would include at least one of ${\beta }_{k0},{\beta }_{k1},{\beta }_{k2}$. For example, ${\beta }_{ER}(l)=0.291-0.141l+0.005l^2$, which was dynamic and varied with $l$. In the case of, ${\beta }_{ER}(4.5)\approx -0.25(=0.291-0.141\times 4.5+0.005\times {4.5}^2)$, which means ER-negative patients lost an additional 0.25 years of life on average during the 4.5 years follow-up. Figures 4A–J shows the regression coefficients of different prognostic factors in different $l$. For example, Figure 4G shows ${\beta }_{ER}(l)$ (solid black line), with breast cancer deaths occurring an average of 0.25 years later (y-axis ${\beta }_{ER}(4.5)\approx -0.25$) in ER-positive patients than in ER-negative patients during the 4.5-year follow-up (x-axis l = 4.5); breast cancer deaths occurring an average of 0.42 years later (y-axis ${\beta }_{ER}(6.5)\approx -0.42$) in ER-positive patients than in ER-negative patients during the 6.5 years (x-axis $l=6.5$).

Table 4

Table 4 Regression coefficients of dynamic-effect RMTL regression (2.5 years ≤ l ≤ 10.5 years).

Figure 4

Figure 4 The curves of the regression coefficient changing over time. The panels (A–J) represent different variables. The black solid line represents the regression coefficient β(l), the black dashed line represents the 95% confidence interval of β(l), and the red dotted line is an auxiliary line (straight line between two endpoints of β(l)), which is used to judge whether β(l) is a curve. In (G*, H*), the blue dashed line corresponds to the right coordinate and is the absolute value of the slope of β(l).

In the dynamic effect RMTL regression, in addition to obtaining dynamic cumulative effects β(l), real-time effects of prognostic factors can also be obtained by the speed of accumulation of β(l) at different moments. We obtained real-time effects by the absolute value of the slope of the curve of β(l), which represents the speed of accumulation. Figures 4G*, 4H* show the regression coefficients (black solid lines) and the speed of accumulation (blue dashed lines) of ER and PR, respectively. In Figure 4G*, the speed of decline of β(l) was decreasing, and the speed of decline was 0.097 when l = 4.5, which is in units of the difference in life lost between the positive group and negative group in the 1-year follow-up; the speed of decline was 0.077 when $l=6.5$. Therefore, the real-time effect of ER decreased with time. In Figure 4H*, the speed of decline of β(l) was constant, and the speeds of decline were both 0.027 when l = 4.5 and l = 6.5. Therefore, the real-time effect of PR remained unchanged with time.

We added an auxiliary line (red dotted line), which is the line between the two endpoints of the regression coefficient curve β(l) in Figures 4A–J, to determine whether the real-time effect of the prognostic factor changed. In general, 1) when β(l) coincided with the auxiliary line (Figures 4A, D, F, H), the real-time effect was unchanged; 2) when β(l) decreased with increasing l (Figures 4G, I, J), if the regression coefficient curve was below the auxiliary line, that is, the real-time effect decreased (Figures 4G, J), and conversely, the real-time effect increased (Figure 4I); and 3) when β(l) increased with increasing l (Figures 4B, C, E), the regression coefficient curve above the auxiliary line corresponded to a decrease in the real-time effect, and conversely, it corresponded to an increase in the real-time effect (Figures 4B, C, E). Therefore, it was concluded that the real-time effects of age, stage N3 (relative to stage N1), histological grade III&IV (relative to grade II), and PR positivity were unchanged; the real-time effects of ER positivity and chemotherapy decreased; and the real-time effects of T2 (relative to T1), N2 (relative to N1), histological grade II (relative to grade I), and breast-conserving surgery increased.

In addition to exploring the dynamic cumulative effects and real-time effects of prognostic factors, another role of dynamic effect RMTL regression is providing personalized prediction for patients. Three patients were selected (see Table 5 for details). Figure 5 shows the predicted RMTL during the l-year follow-up of each patient, and Table 5 also shows the predicted RMTL during the 5-year and 10-year follow-up. In the case of patient A, the predicted mean life lost due to death from breast cancer was 1.5 years in the 5-year follow-up; in the decade of follow-up, the predicted value was 4.2 years. Patients B and A differed only in the choice of treatment. Compared with patient A, patient B received breast-conserving surgery and chemotherapy, and his predicted RMTL was less than that of patient A; that is, breast-conserving surgery and chemotherapy could prolong the survival time of older patients with early-stage breast cancer. Patient C differed from patient B in N stage and histological grade, and because patient C had lower N stage and histological grade, his predicted RMTL was lower than that of patient B.

Table 5

Table 5 The definition of three example patients.

Figure 5

Figure 5 Predicted trajectories of RMTL for different patients. The predicted mean life lost of patient A due to death from breast cancer was 1.5 years in the 5-year follow-up; in the decade of follow-up, the corresponding predicted value was 4.2 years.

In addition, the accuracy of prediction was evaluated by an external validation set. Figure 6 shows the C-index and relative prediction error when the pre-specified end time of follow-up was different (18). The mean C-index was 0.78, indicating good discrimination of the model, and the relative prediction error was within 10%.

Figure 6

Figure 6 C-index (A) and relative prediction error (B) at different end times of follow-up. The C-index refers to the accuracy of the model in predicting the sequence of occurrence of death from breast cancer in the $l$-year follow-up. Relative prediction error is the proportion of prediction error to length of follow-up.

The prediction formula can be seen in Table 4, and the prediction model has been converted into a web-based prediction tool available on the web at https://m92imi-oscar-0.shinyapps.io/newapp/.

4 Discussion

When the effect of a prognostic factor on competing events is large, we should use a competing risks approach; otherwise, the estimate of the effect of this factor on the event of interest will be biased greatly (19). In our data, the sHRs of age and chemotherapy on death from non-breast cancer (the competing event) were 2.486 (95% CI: 2.181 to 2.834) and 0.627 (95% CI: 0.545 to 0.722), respectively. Moreover, the number of those who experienced the competing event accounted for 26% of the total sample size and 56% of the total number of events, so it is necessary to consider competing risks in these data.

In the static effect RMTL regression, it only gives the cumulative effect during the $\tau$-year follow-up, and it is impossible to know the real-time effect in the cumulative process. In particular, this result is incomplete for covariates with time-varying effects. Additionally, for patients who have been followed up for some time, the cumulative effect from 0 to $\tau$ years is no longer applicable. In contrast, the dynamic effect RMTL regression can not only obtain the dynamic cumulative effect in the l-year follow-up but also explore the real-time effect. The real-time effect can help doctors and patients to have a better understanding of the prognosis of breast cancer. For example, the real-time effect of ER positivity decreased, which means its protective effect is larger in the first period and smaller in the later period, suggesting that estrogen therapy should be used as early as possible; the real-time effect of breast-conserving surgery increased, which means its protective effect is larger in the later period, suggesting that the effect of breast-conserving surgery is delayed.

Regarding the prognostic analysis of death from breast cancer, Yao used Cox regression and cause-specific Cox regression to analyze the difference in the effects of prognostic factors on breast cancer in men and women (20), and Xu used Fine-Gray regression to develop a prediction model for patients with inflammatory breast cancer (21). However, none of these studies considered the potential time-varying effects of prognostic factors. Moreover, some studies analyzed the time-varying effects of prognostic factors (13, 14, 17, 22), but these were the results of single-endpoint survival analysis and did not consider the impact of competing events, which may result in competing bias.

In this paper, both competing risks and time-varying effects were considered for the first time, and the real-time effects of the following prognostic factors were found to be different from the previous single-endpoint analysis results. First, in a single-endpoint analysis of breast cancer, the risk effect of stage N2 relative to stage N1 decreased over time (17). In contrast, we found that stage N2 was also a risk factor, but the real-time effect increased over time (Figure 4C). Second, previous single-endpoint studies have shown that the deleterious effect of histological grade II relative to grade I decreased over time (13, 22). However, we found that the deleterious effect of histological grade II increased over time (Figure 4E). Third, in previous single-endpoint studies, ER positivity was a protective factor in the early period and a deleterious factor in the late period (13, 22, 23). This was different from our results, which showed that the protective effect of ER positivity decreased over time (Figure 4G). Fourth, in terms of treatment, we found that patients with breast-conserving surgery had a better prognosis than those with mastectomy (Figure 4I). This is consistent with Kim’s study and a meta-analysis, which showed that patients who underwent breast-conserving surgery had a higher overall survival rate than those who underwent mastectomy (24, 25). However, we further discovered that the protective effect of breast-conserving surgery increased over time (Figure 4I). Finally, chemotherapy was the protective factor, and its real-time effect decreased (Figure 4J). This is similar to Rakovich’s study, which found that chemotherapy after breast-conserving surgery in patients with ductal carcinoma in situ reduced the risk of early local recurrence but not the risk of late recurrence (26).

Finally, the final dynamic RMTL model was constructed with the full dataset (see Web Supplementary Table 2 in Supplementary Material), and the result was similar to that constructed with the training set (Table 4).

Time-varying covariate and covariate with time-varying effect are two different types of data, which requires different statistical methods to analyze (27). Time-varying covariate means the value of a covariate changes over time, which needs methods related to longitudinal data to analyze. While covariate with time-varying effect means the effect on the outcome is time-varying (28). Meanwhile, covariates do not meet the proportional subdistribution hazards assumption, tending to have time-varying effect in the competing risks. Because time-varying effect is difficult to identify, we often ignore it. And then biased estimates will be obtained, and the significant effect occurring only in part of the follow-up period will be missed (29). Among the two types of covariates, this paper focuses on the latter and proposes an extended RMTL regression model to depict time-varying effects, which also can be used in single-endpoint survival data. The extension for time-varying covariates will be the focus of our future research.

There are still some shortcomings in this study. First, the model uses IPCW. It should be noted that there are very few patients remaining at-risk at the end of follow-up, which may lead to large and unstable weights. 2) The life lost is the time lost due to death from breast cancer over a period of time (the l-year follow-up) rather than the reduction in total life in the traditional sense. 3) The HER2 status is also an important prognostic factor for breast cancer. Due to the SEER database only beginning to record HER2 status in 2010, we have chosen not to include this variable in our analysis.

5 Conclusion

To explore the potential time-varying effects of prognostic factors under competing risks survival data, we develop a dynamic effect RMTL regression to model the stacked dataset by generalized estimating equation and IPCW technique. The simulation of regression coefficients and external validation of prediction demonstrate that dynamic effect RMTL regression is accurate in both prognosis and prediction when competing risks exist. The new model can explore dynamic cumulative effects and real-time effects of prognostic factors on a time scale, which gives clinical researchers a more comprehensive understanding of the progression of breast cancer. Moreover, time-scale-based individual prediction also allows physicians and patients to more intuitively determine the disease and choose the best treatment.

Data availability statement

The data analyzed in this study is subject to the following licenses/restrictions: The data were available upon request to the SEER website (www.seer.cancer.gov). Requests to access these datasets should be directed to www.seer.cancer.gov.

Ethics statement

Ethical approval was not required for the studies involving humans because the data was from Surveillance, Epidemiology, and End Results (SEER) database. The studies were conducted in accordance with the local legislation and institutional requirements. Written informed consent for participation was not required from the participants or the participants’ legal guardians/next of kin in accordance with the national legislation and institutional requirements because the data was from Surveillance, Epidemiology, and End Results (SEER) database.

Author contributions

ZY: Methodology, Data curation, Formal analysis, Writing – original draft, Writing – review & editing. ZL: Writing – review & editing. CZ: Writing – review & editing. DZ: Visualization, Writing – review & editing. YH: Conceptualization, Writing – review & editing. ZC: Conceptualization, Methodology, Funding acquisition, Supervision, Writing – review & editing.

Funding

The author(s) declare financial support was received for the research, authorship, and/or publication of this article. This work was supported by the National Natural Science Foundation of China (grant numbers 82173622, 81903411, 81673268) and the Guangdong Basic and Applied Basic Research Foundation (grant number 2022A115011525, 2024A1515011402).

Acknowledgments

The authors thank National Cancer Institute (NCI) for providing the surveillance, epidemiology, and end results (SEER) database.

Conflict of interest

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Publisher’s note

All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.

Author disclaimer

Any perspectives, results, or conclusions found in this paper are those of the authors.

Supplementary material

The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fonc.2024.1352111/full#supplementary-material

Abbreviations

cHR, cause-specific hazard ratio; sHR, subdistribution hazard ratio; RMTL, restricted mean time lost; CIF, cause-specific cumulative incidence function; IPCW, inverse probability of censoring weighting; SEER, Surveillance, Epidemiology, and End Results; PR, progesterone receptor status; ER, breast-conserving surgery; ALND, axillary lymph node dissection; SLNB, lymph node biopsy; 95%CI, 95% confidence interval.

References

1. Freedman RA, Keating NL, Lin NU, Winer EP, Vaz-Luis I, Lii J, et al. Breast cancer-specific survival by age: Worse outcomes for the oldest patients. Cancer. (2018) 124:2184–91. doi: 10.1002/cncr.31308

PubMed Abstract | CrossRef Full Text | Google Scholar

2. Coemans M, Verbeke G, Döhler B, Süsal C, Naesens M. Bias by censoring for competing events in survival analysis. BMJ. (2022) 378:e071349. doi: 10.1136/bmj-2022-071349

PubMed Abstract | CrossRef Full Text | Google Scholar

3. Schuster NA, Hoogendijk EO, Kok AAL, Twisk JWR, Heymans MW. Ignoring competing events in the analysis of survival data may lead to biased results: a nonmathematical illustration of competing risk analysis. J Clin Epidemiol. (2020) 122:42–8. doi: 10.1016/j.jclinepi.2020.03.004

PubMed Abstract | CrossRef Full Text | Google Scholar

4. Ramspek CL, Teece L, Snell KIE, Evans M, Riley RD, van Smeden M, et al. Lessons learnt when accounting for competing events in the external validation of time-to-event prognostic models. Int J Epidemiol. (2022) 51:615–25. doi: 10.1093/ije/dyab256

PubMed Abstract | CrossRef Full Text | Google Scholar

5. Saad ED, Zalcberg JR, Péron J, Coart E, Burzykowski T, Buyse M. Understanding and communicating measures of treatment effect on survival: can we do better? J Natl Cancer Inst. (2018) 110:232–40. doi: 10.1093/jnci/djx179

PubMed Abstract | CrossRef Full Text | Google Scholar

6. Blagoev KB, Wilkerson J, Fojo T. Hazard ratios in cancer clinical trials—a primer. Nat Rev Clin Oncol. (2012) 9:178–83. doi: 10.1038/nrclinonc.2011.217

PubMed Abstract | CrossRef Full Text | Google Scholar

7. Zhao L, Tian L, Claggett B, Pfeffer M, Kim DH, Solomon S, et al. Estimating treatment effect with clinical interpretation from a comparative clinical trial with an end point subject to competing risks. JAMA Cardiol. (2018) 3:357–8. doi: 10.1001/jamacardio.2018.0127

PubMed Abstract | CrossRef Full Text | Google Scholar

8. Wu H, Yuan H, Yang Z, Hou Y, Chen Z. Implementation of an alternative method for assessing competing risks: restricted mean time lost. Am J Epidemiol. (2022) 191:163–72. doi: 10.1093/aje/kwab235

PubMed Abstract | CrossRef Full Text | Google Scholar

9. Lyu J, Hou Y, Chen Z. Combined tests based on restricted mean time lost for competing risks data. Stat Biopharmaceutical Res. (2023) 15:332–9. doi: 10.1080/19466315.2021.1994456

CrossRef Full Text | Google Scholar

10. Lyu J, Hou Y, Chen Z. The use of restricted mean time lost under competing risks data. BMC Med Res Method. (2020) 20:197. doi: 10.1186/s12874-020-01040-9

CrossRef Full Text | Google Scholar

11. Andersen PK. Decomposition of number of life years lost according to causes of death. Stat Med. (2013) 32:5278–85. doi: 10.1002/sim.5903

PubMed Abstract | CrossRef Full Text | Google Scholar

12. Conner SC, Trinquart L. Estimation and modeling of the restricted mean time lost in the presence of competing risks. Stat Med. (2021) 40:2177–96. doi: 10.1002/sim.8896

PubMed Abstract | CrossRef Full Text | Google Scholar

13. Zhang M, Peng P, Gu K, Cai H, Qin G, Shu XO, et al. Time-varying effects of prognostic factors associated with long-term survival in breast cancer. Endocr Relat Cancer. (2018) 25:509–21. doi: 10.1530/ERC-17-0502

PubMed Abstract | CrossRef Full Text | Google Scholar

14. Baulies S, Belin L, Mallon P, Senechal C, Pierga J-Y, Cottu P, et al. Time-varying effect and long-term survival analysis in breast cancer patients treated with neoadjuvant chemotherapy. Br J Cancer. (2015) 113:30–6. doi: 10.1038/bjc.2015.174

PubMed Abstract | CrossRef Full Text | Google Scholar

15. Zhang C, Huang B, Wu H, Yuan H, Hou Y, Chen Z. Restricted mean survival time regression model with time-dependent covariates. Stat Med. (2022) 41:4081–90. doi: 10.1002/sim.9495

PubMed Abstract | CrossRef Full Text | Google Scholar

16. Zhong Y, Schaubel DE. Restricted mean survival time as a function of restriction time. Biometrics. (2022) 78:192–201. doi: 10.1111/biom.13414

PubMed Abstract | CrossRef Full Text | Google Scholar

17. Fontein DBY, Klinten Grand M, Nortier JWR, Seynaeve C, Meershoek-Klein Kranenbarg E, Dirix LY, et al. Dynamic prediction in breast cancer: proving feasibility in clinical practice using the TEAM trial. Ann Oncol. (2015) 26:1254–62. doi: 10.1093/annonc/mdv146

PubMed Abstract | CrossRef Full Text | Google Scholar

18. van Geloven N, Giardiello D, Bonneville EF, Teece L, Ramspek CL, van Smeden M, et al. Validation of prediction models in the presence of competing risks: a guide through modern methods. BMJ. (2022) 377:e069249. doi: 10.1136/bmj-2021-069249

PubMed Abstract | CrossRef Full Text | Google Scholar

19. Dignam JJ, Zhang Q, Kocherginsky M. The use and interpretation of competing risks regression models. Clin Cancer Res. (2012) 18:2301–8. doi: 10.1158/1078-0432.CCR-11-2097

PubMed Abstract | CrossRef Full Text | Google Scholar

20. Yao N, Shi W, Liu T, Siyin ST, Wang W, Duan N, et al. Clinicopathologic characteristics and prognosis for male breast cancer compared to female breast cancer. Sci Rep. (2022) 12:220. doi: 10.1038/s41598-021-04342-0

PubMed Abstract | CrossRef Full Text | Google Scholar

21. Xu F, Yang J, Han D, Huang Q, Li C, Zheng S, et al. Nomograms for estimating cause-specific death rates of patients with inflammatory breast cancer: A competing-risks analysis. Technol Cancer Res Treat. (2021) 20:1–12. doi: 10.1177/15330338211016371

CrossRef Full Text | Google Scholar

22. Natarajan L, Pu M, Parker BA, Thomson CA, Caan BJ, Flatt SW, et al. Time-varying effects of prognostic factors associated with disease-free survival in breast cancer. Am J Epidemiol. (2009) 169:1463–70. doi: 10.1093/aje/kwp077

PubMed Abstract | CrossRef Full Text | Google Scholar

23. Ren Y, Black DM, Mittendorf EA, Liu P, Li X, Du XL, et al. Crossover effects of estrogen receptor status on breast cancer-specific hazard rates by age and race. PloS One. (2014) 9:e110281. doi: 10.1371/journal.pone.0110281

PubMed Abstract | CrossRef Full Text | Google Scholar

24. De la Cruz Ku G, Karamchandani M, Chambergo-Michilot D, Narvaez-Rojas AR, Jonczyk M, Príncipe-Meneses FS, et al. Does Breast-Conserving Surgery with Radiotherapy have a Better Survival than Mastectomy? A Meta-Analysis of More than 1,500,000 Patients. Ann Surg Oncol. (2022) 29:6163–88. doi: 10.1245/s10434-022-12133-8

PubMed Abstract | CrossRef Full Text | Google Scholar

25. Kim H, Lee SB, Nam S-J, Lee ES, Park B-W, Park HY, et al. Survival of breast-conserving surgery plus radiotherapy versus total mastectomy in early breast cancer. Ann Surg Oncol. (2021) 28:5039–47. doi: 10.1245/s10434-021-09591-x

PubMed Abstract | CrossRef Full Text | Google Scholar

26. Rakovitch E, Sutradhar R, Hallett M, Thompson AM, Gu S, Dumeaux V, et al. The time-varying effect of radiotherapy after breast-conserving surgery for DCIS. Breast Cancer Res Treat. (2019) 178:221–30. doi: 10.1007/s10549-019-05377-8

PubMed Abstract | CrossRef Full Text | Google Scholar

27. Schumacher M, Hieke S, Ihorst G, Engelhardt M. Dynamic prediction: A challenge for biostatisticians, but greatly needed by patients, physicians and the public. Biometrical J. (2020) 62:822–35. doi: 10.1002/bimj.201800248

CrossRef Full Text | Google Scholar

28. Cui Y, Peng L. Assessing dynamic covariate effects with survival data. Lifetime Data Anal. (2022) 28:675–99. doi: 10.1007/s10985-022-09571-7

PubMed Abstract | CrossRef Full Text | Google Scholar

29. Brenner H, Hakulinen T. Are patients diagnosed with breast cancer before age 50 years ever cured? J Clin Oncol. (2004) 22:432–8. doi: 10.1200/JCO.2004.04.067

PubMed Abstract | CrossRef Full Text | Google Scholar