A Review of the Connections Among Logistic Regression, the Cox Model, Competing Risk Models, and Multi-State Models for Estimating Hazards and Survival Risks

1. Logistic Regression, Pooled Logistic Regression, and Conditional Logistic Regression

Although this review focuses primarily on survival analysis, we begin with logistic regression because of its close connection to survival methods. Specifically, a conditional logistic regression stratified by risk set is mathematically equivalent to the Cox model, and pooled logistic regression has been used to estimate discrete-time hazards and survival risks.

1.1. Logistic Regression

Logistic regression is used to estimate the association between an exposure and a binary outcome (Model 1). It models the probability of a dichotomous outcome and does not account for the timing of the event. The model is specified as

$$logit\left(P\right)={\beta }_0+{\beta }_{1}X —\mathbf{Model} \mathbf{1}$$

where P is the probability of the outcome. The coefficient ${\beta }_1$ is interpreted as the change in the log odds of the outcome associated with a one-unit increase in X, and $exp\left({\beta }_1\right)$ represents the odds ratio.

1.2. Pooled Logistic Regression

Although pooled logistic regression is a form of logistic regression,[4] it is commonly used to analyze survival data. The follow-up time (person-time) for each participant is divided equally into small intervals (e.g., days, months, or years), creating one observation for each person per interval. There is a fundamental assumption underlying the use of these person-time data in survival analysis, including pooled logistic regression (Assumption 1).

Assumption 1.

Given the covariates in each interval, observations divided into small time intervals are independent, even when they originate from the same individual.

Under this assumption, the data structure can be changed from wide to long format (or counting process data), and all observations are treated as independent (Figure 1, a-c). By changing to the long format data, we are modeling risk per time interval (Model 2). In particular, the intervals are created with equal lengths so that a logistic regression model can be fitted to these person-time data to estimate the incidence rate (as each time interval has the same duration, there is no need to account for the duration of time interval using Poisson model).

$$logit\left(P_t\right)={\beta }_0+{\beta }_{1}X+{\beta }_{2}t —\mathbf{Model} \mathbf{2}$$

where $P_t$ is the probability of the outcome in time interval t, or incidence rate. The time interval t is modeled as a counting variable (t=1, 2, 3, …) to allow for changes in the baseline incidence rate over follow-up. $exp\left({\beta }_1\right)$ is interpreted as the incidence rate ratio of the outcome associated with a one-unit increase in X. When the interval is small, the event probability within each interval is low, and the incidence rate approximates the hazard, and $exp\left({\beta }_1\right)$ approximates hazard ratio.

1.3. Conditional Logistic Regression

Conditional logistic regression is widely used in risk-set sampling within nested case-control studies.[5] Conditional logistic regression is an application of the Cox model to matched case-control studies.[6,7] Because matched case-control studies are conceptually easier to understand, we first introduce conditional logistic regression, and hope that its connection to the Cox model will help readers better understand the concept of partial likelihood.

For risk-set sampling, a risk set is created each time a case occurs, and each case is matched to one or more controls from the same risk set. As shown in Figure 1d, risk sets of the underlying population are created when an event occurs. Based on Assumption 1, the risk sets are also independent of each other. This is also why a control can be selected multiple times in risk set sampling: we are sampling person-time observations rather than individuals, and these observations are assumed to be independent of one another, even when they come from the same individual.

Conditional logistic regression can be used to examine the association between exposure and outcome by conditioning on the matched sets (Model 3). Compared to pooled logistic regression adjusting for time as a covariate, Model 3 controls for time as confounding in a more efficient way, as stratification better controls for confounding, and a conditional likelihood cancels the intercept and makes the estimation more efficient.

$$logit\left(P_t\right)={\beta }_{0t}+{\beta }_{1}X, —\mathbf{Model} \mathbf{3}$$

where $P_t$ is the conditional probability of developing the outcome within the risk set created at time t, $R\left(t\right)$, and $P_t$ is essentially a hazard. ${\beta }_{0t}$ is the baseline log odds over time t (i.e., $X$=0), which is canceled out in conditional likelihood estimation. The coefficient ${\beta }_1$ is interpreted as the log odds ratio of developing the outcome per unit increase in X, conditional on the risk set, and ${\beta }_1$ is essentially log hazard ratio.

The conditional likelihood is based on the probability of observing a case, given the number of participants within a risk set. In both scenarios, the intercept term, ${\beta }_{0t}$ is canceled out through conditioning, and only the regression coefficient, ${\beta }_1$, is estimated. The conditional likelihood for the risk set $R\left(t\right)$ is shown below for one case matched with multiple controls and multiple cases matched with multiple controls. See Appendix 1 for the detailed derivation.

Conditional likelihood for one case and $l$ participants in risk set $R\left(t\right)$:

$${\displaystyle \frac{exp\left({\beta }_1{X}_1\right)}{{\sum }_{l}exp\left({\beta }_1{X}_l\right)}}—\mathbf{Formula} \mathbf{1}$$

Conditional likelihood for $d$ cases and $l$ participants in risk set $R\left(t\right)$:

$${\displaystyle \frac{exp\left({\beta }_1{X}_1\right)exp\left({\beta }_1{X}_2\right)\times ...\times exp\left({\beta }_1{X}_d\right)}{{\sum }_{\frac{l!}{d!(l-d)!}}({\prod }_{a\in A}exp\left({\beta }_1{X}_a\right))}},—\mathbf{Formula} \mathbf{2}$$

where A is the set of participants randomly selected from the risk set $R\left(t\right)$, with the number of participants in A equal to the number of cases.

2. Single Endpoint, the Cox Model, and Pooled Logistic Regression

Survival analysis deals with time-to-event data. In this paper, we focus on approaches that model the hazard, specifically the Cox model and pooled logistic regression. The comparison of the two models is summarized in Table 1.

2.1. Hazard and Survival Risks

The hazard is the instantaneous risk of experiencing an event at a given time among individuals who are still at risk. Survival risk is the probability of experiencing the event over a period. The survival risk can be estimated nonparametrically using two methods.

The Kaplan-Meier estimator directly estimates the survival function,⁸ which is given by

$S\left(t\right)={\prod }_1^t(1-h(t\left)\right)$, where $h\left(t\right)$ is the hazard at time t.

This is intuitive because the probability of surviving beyond a given time depends on having survived all preceding time intervals. This approach is commonly used to estimate survival risks from pooled logistic regression.

Another approach is to estimate the cumulative hazard using the Nelson-Aalen estimator,[9,10] which is the cumulative hazard over time: $H\left(t\right)={\sum }_1^{t}h\left(t\right)$. The survival risk can be estimated as $S\left(t\right)=e^{-H\left(t\right)}$. Of note, the cumulative hazard is different from the cumulative incidence of an event, 1-$S\left(t\right)$, as cumulative hazard is not a probability and thus can exceed 1. This approach has been used to estimate survival risks using the Cox model.

2.2. The Cox Model, Cumulative Hazard, and Continuous-Time Survival Risks

The Cox model models the hazard (Model 4).[11] The events are ranked by the time of occurrence, and a risk set is constructed when an event occurs, as illustrated in Figure 1d. Under the Assumption 1, all risk sets are independent of each other. The Cox model estimates the coefficients using the partial likelihood, which estimates the conditional probability of observing the event given the corresponding risk set.

$$log\left(h\right(t\left)\right)=log\left(h_0\right(t\left)\right)+{\beta }_{1}X —\mathbf{Model} \mathbf{4}$$

The partial likelihood is obtained by multiplying the conditional probabilities across all risk sets over time. The partial likelihood for a risk set $R\left(t\right)$ with l participants in the risk set and one case is shown as:

$${\displaystyle \frac{exp\left({\beta }_1{X}_1\right)}{{\sum }_{l}exp\left({\beta }_1{X}_l\right)}} —\mathbf{Formula} \mathbf{3}$$

The Cox model is connected to conditional logistic regression in three aspects. First, instead of randomly sampling controls for each case in the nested case control study with risk set sampling, the Cox model includes all individuals at risk at each event time in the risk set. It is evident that the conditional likelihood of the conditional logistic regression is equivalent to the partial likelihood of the Cox model (Formula 1 and 3 are equivalent). Second, in the presence of tied events, the exact partial likelihood of the Cox model[12] is equivalent to the conditional likelihood of conditional logistic regression with multiple cases matched to multiple controls (Formula 2), where the denominator is a sum over all possible combinations of selecting the observed number of cases from the risk set. Third, both models are semi-parametric. In the Cox model, the baseline hazard $h_0$(t) is non-parametric. Because it is canceled out in the partial likelihood, we do not need to make any assumptions about $h_0$(t). The covariates are modeled parametrically under the proportional hazards assumption, whereby $h\left(t\right)$ is proportional to $h_0$(t) on the log scale. The conditional logistic regression model also shares these properties as discussed in Section 1.3.

Survival risks can be estimated from the Cox model using the Breslow estimator.[13] Although the baseline hazard function is canceled out during partial likelihood estimation, it can be estimated from the observed data after the coefficients have been obtained. Details of estimating ${\widehat{h}}_0\left(\mathrm{t}\right)$ are provided in Appendix 2. Then the Breslow estimator is used to estimate the cumulative baseline hazard, where ${\widehat{H}}_0\left(\mathrm{t}\right)$=${\sum }_t{\widehat{h}}_0\left(\mathrm{t}\right)$, which is then used to obtain the survival function ${\widehat{\mathrm{S}}}_0\left(\mathrm{t}\right)={\mathrm{e}\mathrm{x}\mathrm{p}(-\widehat{\mathrm{H}}}_0\left(\mathrm{t}\right))$. Individual survival risks can subsequently be estimated by incorporating the estimated coefficients, $\widehat{S}\left(\mathrm{t}\right)={{\widehat{S}}_0\left(\mathrm{t}\right)}^{exp\left({\beta }_{1}X\right)}$. This explains why most Cox model software reports the cumulative hazard and survival risks rather than the hazard itself, as hazard is not directly estimated at time t.

2.3. Pooled Logistic Regression, Discrete-Time Hazard, and Survival Risks

As discussed in Section 1.2, pooled logistic regression can be used to estimate discrete-time incidence rates, which approximate the hazard when the event is rare given the short time intervals. Pooled logistic regression includes time as a covariate in the model, and the hazard is modeled parametrically as:

$$P_t={\displaystyle \frac{e^{{\beta }_0+{\beta }_{1}X+{\beta }_{2}t}}{1+e^{{\beta }_0+{\beta }_{1}X+{\beta }_{2}t}}}$$

Survival risks can then be estimated by Kaplan-Meier estimator. Because follow-up time is divided into equally small intervals rather than risk sets being defined at each event time, the survival risks estimated from pooled logistic regression are discrete-time estimates. The 95% CI for the estimated survival risk can be obtained using bootstrap.

Although pooled logistic regression approach is more computationally intensive than the Cox model in estimating survival risks, it offers two advantages. First, pooled logistic regression estimates the hazard, which can be flexibly used to create complex statistics. For example, in causal inference, methods such as the g-formula[14] and marginal structural models[15,16] construct counterfactual survival risks based on discrete-time hazards estimated using pooled logistic regression.[17] Second, the long-format data used to fit pooled logistic regression provides a data structure for integrating machine learning into survival analysis. In addition to pooled logistic regression, machine learning methods can be used to estimate the risks at each time interval, which approximates the discrete-time hazard. This approach has been used to integrate neural networks and random forests to survival analysis for estimating discrete-time survival risks.[18,19]

A special case: Replacing pooled logistic regression with a Cox model in the long-format data to estimate hazards.

Previously, we discussed fitting a pooled logistic regression to long-format data to estimate discrete-time hazards, where the hazards within each interval are constant. In addition to pooled logistic regression, a Cox model can be applied to the long-format data.[20]

$$log\left(h\right(\tau \left)\right)=log\left(h_0\right(\tau \left)\right)+{\beta }_{1}X+f\left(t\right) —\mathbf{Model} \mathbf{5}$$

where $\tau$ is time to event within the time interval of long data format, $h_0\left(\tau \right)$ is the baseline hazard within a time interval, $f\left(t\right)$ is a function of time interval (t=1, 2, 3..).

The survival risks can be estimated within each time interval that approximate hazards. Thus, instead of estimating a constant hazard within a small interval using pooled logistic regression, we can use survival probabilities from the Cox model at a specific time point within the interval (such as the midpoint or end of the interval) to approximate the hazard (Figure 2).

3. Competing risks, Cause-Specific Cox Models, and Multinominal Pooled Logistic Regression

3.1. Competing Risks and the Cumulative Incidence Function (CIF)

Competing risks are mutually exclusive events in which the occurrence of one event prevents the occurrence of other events. Once a competing event occurs, the individual is no longer at risk of experiencing other events. For example, if the outcome is the first atherosclerotic cardiovascular disease (ASCVD) event, competing risks include stroke, coronary artery disease, and non-ASCVD mortality. The first occurrence of either stroke or CAD constitutes the event of interest. Non-ASCVD mortality is also a competing risk, as if a person dies due to causes other than ASCVD, the individual can no longer experience a first ASCVD event.

Survival risks for competing risks can be estimated non-parametrically using CIF function,[21] which estimates the marginal probability of a specific event occurring over time in the presence of competing risks.

$${CIF}_k\left(t\right)={\sum }_t{h}_k\left(t\right)\times S\left(t\right),$$

where $h_k\left(t\right)$ is the hazard due to cause $k$, and ${CIF}_k\left(t\right)$ is the cumulative incidence for cause $k$.

For competing risks, survival risks should be estimated using the CIF rather than the Kaplan-Meier estimator applied separately to each cause. The reason is that the CIF imposes constrains on the cumulative incidence of all competing risks, which ensures that the total probability across all causes does not exceed 1.

We illustrate this using an example in Table 2. Suppose individuals enter the study in an at-risk state (cause 0), with two competing risks, indicated as cause 1 and cause 2. The cause-specific hazards at times $t_1$-$t_3$ are assumed to be known. If the incidence for each cause at each time point is estimated using the Kaplan-Meier estimator while ignoring the occurrence of the other cause (i.e., treating individuals who experience the competing event as still at risk for cause 0), the resulting cumulative incidence of developing either cause 1 or cause 2 may exceed 1 at $t_2$ and $t_3$ . This results in negative values for the survival risks at $t_2$ and $t_3$ (i.e., the probability of remaining at risk), which is not meaningful. In contrast, the CIF estimates the incidence of each cause while accounting for competing events by removing individuals from the risk set once any event occurs, which imposes a constraint ensuring that the cumulative incidence of all competing risks does not exceed 1 at all time points. Because the CIF accounts for individuals who experience competing events by treating them as no longer at risk, it can be viewed as an extension of the Kaplan-Meier estimator to the competing risks setting (Table 3).

3.2. Cause-Specific Cox Models

For competing risks with causes 1 and 2, a common but incorrect approach is to fit two separate Cox models, one for each outcome, and then estimate survival risk for each outcome using the Kaplan-Meier estimator. As discussed above, the problem is that the estimated survival risks for competing events may exceed 1, as no constraints are imposed to ensure the cumulative incidence of developing either competing risk may exceed 1.

Cause-specific Cox (CSC) models can be applied to model competing risks.[22] It fits separate Cox models to model the hazard for each cause.

$$log\left(h_k\right(t\left)\right)=log\left(h_{k,0}\right(t\left)\right)+{\beta }_{k}X —\mathbf{Model} \mathbf{6}$$

where $h_k\left(t\right)$ is the hazard due to cause $k$, $h_{k,0}\left(t\right)$ is the baseline hazard due to cause $k$, and ${\beta }_k$ is the coefficient of X with the log hazard of developing cause $k$.

The survival risks are obtained using CIF rather than the Kaplan-Meier estimator, which avoids the problem of cumulative incidence of developing either competing risk may exceed 1. Like the Cox model, the baseline hazard for each cause, ${\widehat{h}}_{k,0}\left(t\right)$, can be estimated based on the data, the cumulative hazard for each cause can be estimated using Breslow estimator as ${\widehat{H}}_{k,0}\left(\mathrm{t}\right)$= ${\sum }_t{\widehat{h}}_{k,0}\left(\mathrm{t}\right)$, and the survival risk can be estimated as $\widehat{\mathrm{S}}\left(\mathrm{t}\right)=\mathrm{e}\mathrm{x}\mathrm{p}(-{\sum }_{\mathrm{K}}{\widehat{\mathrm{H}}}_{\mathrm{k},0}\left(\mathrm{t}\right)\exp ({\mathsf{\beta }}_{\mathrm{k}}\mathrm{X}))$. The CIF for cause $k$ can be estimated as

$${\widehat{\mathrm{F}}}_k\left(t\right)={\sum }_t{\widehat{h}}_k\left(t\right)\times \widehat{\mathrm{S}}\left(t\right).$$

Another approach for analyzing competing risks is the Fine-Gray model.[23] In brief, for example, when our goal is to focus on cause 1, rather than treating individuals who experience cause 2 as censored and removing them from the risk set, the Fine-Gray model considers individuals who develop cause 2 remain in the risk set and applies inverse probability weighting to impute the hazard of developing cause 1. The subdistribution hazard includes the observed hazard of developing cause 1 among individuals who truly develop cause 1 and the imputed hazard of developing cause 1 among those who develop cause 2. The Fine-Gray model estimates the cumulative incidence function of cause 1 based on the subdistribution hazard. Figure 3 conceptually compares approaches to handling competing risks using separate Cox models for each cause, cause-specific Cox models, and Fine-Gray models.

3.3. Pooled Multinomial Logistic Regression

Like section 2, the hazards can be estimated using multinormal pooled logistic models by considering competing risks in the long data format as a multinormal distribution.[24,25] The advantages are the same as the pooled logistic regression, as summarized in Table 1. Given these small age intervals, machine learning can be flexibly and directly applied by considering competing risks as a multi-class classification problem, where the estimated probability of each class approximates discrete-time cause-specific hazards. By using the long-format survival data, Bayesian regression trees, random forests, and deep neural network have been integrated into competing risk analysis to model discrete-time cause-specific hazards.[26,27,28]

4. Multi-State Process, Markov Cox Modeling, and a Discrete-Time Split-State Framework

Multi-state models are an emerging area of research because of the accumulation of large cohort data and their broad applicability in modeling the dynamics of chronic diseases progression. Figure 4 illustrates a four-state model with states S₀–S₄. A forward-transition process is assumed, in which an individual in a state may transition to subsequent states.

4.2. Transition Rate, Survival Risk, and Aalen-Johansen (A-J) Estimator

Transition rate is the instantaneous risk of moving from one state to another and can be regarded as an extension of the hazard to a multi-state setting. The transition rates can be constructed into a transition rate matrix (a 4$\times$4 matrix is shown as an example), where each off-diagonal value denotes the transition rate from the state corresponding to its row to the state corresponding to its column. As we assume a forward transition, the rates for reverse transitions (e.g., from S₁ to S₀) are 0. The diagonal values are defined as the negative sum of the outgoing transition rates from each state. Consequently, each row of the transition rate matrix sums to 0.

The survival risks, or state occupation probabilities, are the probability of an individual being in each state at a given time point. These probabilities can be represented as a vector, in which each value denotes the probability of being in a state, and the elements of the vector sum to 1. Suppose all participants start from S₀ state, the survival risks can be estimated nonparametrically using the Aalen-Johansen (A-J) estimator.[29]

$$\mathit{P}\left(\mathrm{t}\right)=\mathit{P}(\mathrm{t}=0)\times {\prod }_t(I+\mathit{Q}(\mathrm{t}))$$

Where $\mathit{Q}\left(\mathrm{t}\right)$ is the transition rate matrix, I is the identity matrix, and $\mathit{P}\left(\mathrm{t}\right)$ is the survival risk vector.

Appendix 3 illustrates the estimation of the survival risks for a four-state process using the A-J estimator. The estimated survival risk across disease states are equivalent to estimating survival risks using the Kaplan-Meier estimator by allowing individuals to transition from the starting state to multiple endpoints and by allowing individuals to transition dynamically in and out intermediate states. It is therefore evident that the Aalen-Johansen estimator is a natural extension of the Kaplan-Meier estimator to multi-state settings (Table 3).

4.3. Markov Cox Model

Under the Markov assumption (Assumption 2), the Cox model models hazards between state transition.[30] Based on this assumption, if we divide a person’s follow-up time according to the disease states they occupy, the resulting separated follow-up intervals are independent. Assumption 2 is essentially an extension of Assumption 1, in which an individual’s person-time observations are independent of one another, with the conditioning variables including not only covariates but also disease states.

Assumption 2 (Markov assumption).  By conditioning on the present state, the future transition rates are independent of the past history.

Under Assumption 2, follow-up time can be divided into data subsets according to the disease state occupied. Figure 5 shows how the follow-up data is divided into state-specific subsets. Subset 1 (blue) begins in S₀ and ends in S₀, S₁, S₂, or S₃. Subset 2 (red) begins in S₁ and ends in S₁, S₂, or S₃. Subset 3 (green) begins in S₂ and ends in S₂ or S₃. Within each subset, Cox models are applied to model hazards of transition between states. Specifically, three Cox models are fitted to data subset 1, two Cox models to subset 2, and one Cox model to subset 3.

$$log\left(h_{rs}\right(t\left)\right)=log\left(h_{rs,0}\right(t\left)\right)+{\beta }_{rs}X—\mathbf{Model} \mathbf{7}$$

where $h_{rs}\left(t\right)$ is the hazard of transition from state r to state s, $h_{rs,0}\left(t\right)$ is the baseline hazard of transition from state r to state s, and ${\beta }_{rs}$ is the coefficient of X with the log hazard of developing state $s$ from state r.

For survival risk estimation using the Markov Cox model, the cumulative transition hazard is estimated from the corresponding Cox model using the Breslow estimator (or the Nelson-Aalen estimator). Rather than applying the cumulative transition hazard to estimate survival risks like the Cox model or the cause-specific Cox model, the Markov Cox model estimates transition hazards at each time point as the change in cumulative transition hazard at a time point.[31,32] These estimated transition hazards are constructed into a transition rate matrix, $\widehat{\mathbf{Q}}\left(\mathrm{t}\right)$, and the Aalen-Johansen estimator is used to estimate the survival risks across disease states, $\widehat{\mathbf{P}}\left(\mathrm{t}\right)$

$$\widehat{\mathbf{P}}\left(\mathrm{t}\right)={\prod }_t(I+\widehat{\mathbf{Q}}(\mathrm{t}))$$

4.4. A Discrete-Time Split-State Framework

A discrete-time split-state framework has been proposed.[20,33,34] The motivation for this framework is that the Markov assumption is strong for chronic diseases. Take cardiovascular disease as an example, the diagnosis of one endpoint often persists across subsequent states, and there is complex interplays between metabolic risk factors in accelerating disease progression, both of which violate the Markov assumption.[35] Moreover, Markov models are memoryless and do not allow for multimorbidity, a common feature of chronic diseases.[36,37] As shown in Figure 6, the framework splits state into substates by conditioning on past states and reconstruct a process of disease substates.

This split-state approach makes a disease process both memoryless and memorable. First, split-state makes a disease process memoryless. As aforementioned, the Markov Cox models require the Markov assumption. This memoryless property restricts their ability to capture multimorbidity and the interactions among disease states in influencing disease progression. For the split-state framework, by constructing substates conditioning on past states, the substates are independent of past states and satisfy the Markov assumption even when the underlying process is non-Markovian.[20] Therefore, this framework accommodates multimorbidity and is well suitable for modeling transition of chronic diseases. Second, split-state makes a disease process memorable. As the substates are constructed by conditioning on past states, they inform past history, or temporality of multimorbidity development.[34] Thus, the framework can synthesize survival risks across all substates over time into summary measures: disease trajectories for prediction use, and multimorbidity-adjusted life year (MALY) for comparing treatment strategies.³³

For model estimation, it is a natural extension of pooled logistic regression in single endpoint to a multi-state setting. Within each data subset shown in Figure 5, it considers the transition from one to multiple endpoints as competing risks. It cuts time into small intervals and changes the data structure from wide to long format and fits cause-specific Cox models (or multinormal pooled logistic regression) to estimate survival risks within each interval to approximate hazard (see the special case in Section 1.3 for using Cox models to estimate survival risks in a small time interval to approximate hazard). We include past state or their interaction with the current state as covariates into the models to allow past states to affect the transition rates of future states and relax the Markov assumption. The discrete-time transition rates are estimated from the models, and the A-J estimator is used to estimate the discrete-time survival risks of being in each substate, with 95% CI estimated using bootstrap.

This discrete-time split-state framework has three advantages.

First, the framework can be used to predict trajectories,[34] as the substates capture the temporal sequence of its development. Survival risks of being in substates, $S_0$, $S_1$, $S_{2|0}$_, $S_{2|0\_1},S_{3|0}$, $S_{3|0\_1}$, $S_{3|0\_2}$, and $S_{3|0\_1\_2}$, correspond to trajectories, “remaining in $S_0$“, “remaining in $S_1$“, “$S_0$→$S_2$“, “$S_0$→$S_1$→ $S_2$“, “$S_0$→$S_3$“, “$S_0$→$S_1$→$S_3$“, “$S_0$→$S_2$→$S_3$“, and “$S_0$→$S_1$→$S_2$→$S_3$“, respectively.

Second, the discrete-time survival risks estimated from the framework allows to flexibly derive complex metrics. Specifically, the counterfactual multimorbidity-adjusted life year (MALY) under a hypothetical intervention can be estimated, where MALY is a weighted sum of survival risks in each substates over time.[33]

$${MALY}_i^A\left({\mathit{X}}_{\mathit{i}}\right)={\sum }_{\mathit{j}}{\mathrm{D}\mathrm{W}}_{\mathrm{j}}\times ({\sum }_t{\mathit{P}}_i^A\left(\mathrm{t}\right|{\mathit{X}}_{\mathit{i}}\left)\left[j\right]\right),$$

where ${\mathrm{D}\mathrm{W}}_{\mathrm{j}}$is the disability weight for substate $\mathrm{j}$, A is a hypothetical treatment regime across disease substates, ${\mathit{P}}_i^A\left(\mathrm{t}\right|{\mathit{X}}_{\mathit{i}})$ is the counterfactual survival risk at t had the person $i$ received $A$.

A multi-state causal approach[33] has been proposed to compare effectiveness of treatment regime using MALY, where the optimal intervention is the one associated with the highest MALY.

${\tau }_i\left({\mathit{L}}_{\mathit{i}}\right)={MALY}_i^{A_2}\left({\mathit{X}}_{\mathit{i}}\right)$ - ${MALY}_i^{A_1}\left({\mathit{X}}_{\mathit{i}}\right)$ — ${\tau }_i\left({\mathit{L}}_{\mathit{i}}\right)$ is the treatment effect based on Cox models

Third, the long-format data allows flexible integration of machine learning into the framework. In addition to cause-specific Cox models (or multinominal logistic regression), machine learning methods can be applied by considering competing risks as a multi-class classification problem, where the estimated probability of each class approximates discrete-time cause-specific hazards. As mentioned in Section 3.3, this approach has been used to integrate machine learning into competing risk estimation.[26,27,28] As competing risks can be considered as a special case of a multi-state process (or a multi-state process can be viewed as a dynamic combination of multiple competing risks), it is of high interest to integrate machine learning into multi-state modeling to model discrete-time transition rates in the near future.