Multi-State Models With Omics Covariates

Why Model Transitions Rather Than Only Time to One Event?

Clinical progression often unfolds through distinct states. A patient may begin in a stable state, develop disease progression, and later die; some patients may die without first having documented progression. A single time-to-event outcome can conceal these different paths. A multi-state model represents the process as movement among defined states, allowing an omics measurement to be associated differently with different transitions.

The key distinction is between a transition rate and a state-occupation probability. A transition rate describes the instantaneous rate of moving from one particular state to another among people currently in the origin state. An occupation probability describes the probability of being in a particular state at a specified time, accounting for all possible paths into and out of it. These quantities answer related but different clinical questions.

Definition: A multi-state model specifies a set of states and the allowed transitions between them. For a subject in state $g$ at time $t$, the transition intensity from state $g$ to state $h$ is $$ \lambda_{gh}(t\mid \mathcal{H}_t, X)=\lim_{\Delta t\downarrow 0}\frac{P\{Z(t+\Delta t)=h\mid Z(t)=g,\mathcal{H}_t,X\}}{\Delta t}, $$ where $Z(t)$ is the state occupied at time $t$, $\mathcal{H}_t$ is the observed history up to time $t$, and $X$ denotes covariates, including omics measurements. An occupation probability is $P\{Z(t)=h\mid X\}$, given a specified starting-state distribution.

The transition intensity has units of inverse time, such as per month. It is not itself a probability: for a short interval $\Delta t$, the transition probability is approximately $\lambda_{gh}(t\mid\cdot)\Delta t$, when the interval is sufficiently small. For a given origin state, the total rate of leaving that state is the sum of the intensities for its allowed outgoing transitions.

Transition-Specific Hazard Models

Let $X$ be a vector of baseline omics measurements and any adjustment covariates. A common specification is a separate proportional-hazards model for each allowed transition:

\lambda_{gh}(t\mid X)=\lambda_{0,gh}(t)\exp\{\beta_{gh}^{\mathsf T}X\}.

Here $\lambda_{0,gh}(t)$ is the baseline intensity for transition $g\to h$, and $\beta_{gh}$ is its transition-specific coefficient vector. For a one-unit increase in covariate $X_j$, holding the other included covariates fixed, $\exp(\beta_{gh,j})$ is the multiplicative change in the instantaneous rate of transition $g\to h$, among subjects currently in state $g$. It is a hazard ratio for that transition, not a risk ratio for being in state $h$ at a later time.

The main assumptions should match the model and data collection. The states and allowed transitions must represent the clinical process adequately, and event times and states must be observed or censored in a way the analysis can handle. In a baseline-covariate proportional-hazards model, the covariates are treated as fixed over follow-up and each transition has a proportional covariate effect over time. Censoring must be independent of future transition experience conditional on the modeled information, or a suitable censoring adjustment is needed. If the intensity is modeled as depending only on the current state and covariates, rather than on earlier history, the process is assumed to have the Markov property: given the current state and covariates, the past provides no additional information about future transitions. This can be too strong when time already spent in the current state matters.

Worked Example: Interpreting Rates From Transition Data

Suppose a hypothetical cohort starts in state 0, stable disease. The allowed transitions are progression from state 0 to state 1 and death from either state 0 or state 1; state 2 is death and is absorbing. During follow-up, the data include 20 observed $0\to1$ transitions and 10 observed $0\to2$ transitions over 2,000 person-months at risk in state 0. There are 15 observed $1\to2$ transitions over 500 person-months at risk in state 1. For illustration, assume constant transition intensities and independent censoring.

Under a constant-rate model, the likelihood contribution for a transition with $d$ observed events over $Y$ units of state-specific person-time is proportional to $\lambda^d\exp(-\lambda Y)$. Maximizing it gives $\hat{\lambda}=d/Y$. Thus:

\hat{\lambda}_{01}=\frac{20}{2000}=0.0100,\qquad \hat{\lambda}_{02}=\frac{10}{2000}=0.00500,\qquad \hat{\lambda}_{12}=\frac{15}{500}=0.0300

The estimated rates are, respectively, 0.0100 progression transitions per person-month among those in state 0, 0.00500 deaths per person-month among those in state 0, and 0.0300 deaths per person-month among those in state 1. The total rate of leaving state 0 is $0.0100+0.00500=0.0150$ per month. These rates describe the observed hypothetical data under the stated constant-rate model; they are not probabilities of an event over an unspecified follow-up period.

From Transition Intensities to State-Occupation Probabilities

A transition-specific model does not by itself give an occupation probability by simply exponentiating a coefficient. To obtain the probability of occupying each state, the model must account for all allowed paths and the time spent in each state. In a finite-state model, the transition intensities determine a transition-intensity matrix $Q(t)$: its off-diagonal entry $q_{gh}(t)$ is $\lambda_{gh}(t)$ for an allowed transition and zero otherwise; the diagonal entry $q_{gg}(t)$ is the negative sum of the outgoing intensities. Occupation probabilities evolve according to the forward equations $p'(t)=p(t)Q(t)$, where $p(t)$ is the row vector of state probabilities.

For a simple illness-death model with states 0 (stable), 1 (progressed), and 2 (death), suppose the only transitions are $0\to1$, $0\to2$, and $1\to2$. Let the corresponding constant rates be $a_{01}$, $a_{02}$, and $a_{12}$. Starting everyone in state 0, define $a=a_{01}+a_{02}$, the total rate of leaving state 0. Then $p_0(t)$ and $p_1(t)$ satisfy

p_0'(t)=-a p_0(t),\qquad p_1'(t)=a_{01}p_0(t)-a_{12}p_1(t), \qquad p_0(0)=1,\quad p_1(0)=0.

The first equation gives $p_0(t)=\exp(-at)$. To solve the second, multiply it by the integrating factor $\exp(a_{12}t)$. The product rule yields

\frac{d}{dt}\{p_1(t)\exp(a_{12}t)\} =a_{01}\exp\{-(a-a_{12})t\}.

Integrating from 0 to $t$, when $a\ne a_{12}$, gives the occupation probability for state 1:

p_1(t)=a_{01}\frac{\exp(-at)-\exp(-a_{12}t)}{a_{12}-a}.

If $a=a_{12}$, the limiting expression is $p_1(t)=a_{01}t\exp(-at)$. Because death is absorbing and these are the only three states, $p_2(t)=1-p_0(t)-p_1(t)$. This is a model-based probability of occupying death state 2 by time $t$; it is not a transition-specific hazard.

Worked Example: Occupation Probabilities at Six Months

Consider a hypothetical illness-death model with constant rates $a_{01}=0.10$, $a_{02}=0.05$, and $a_{12}=0.20$ per month. Everyone begins in state 0, and assume the state definitions, constant rates, and Markov structure are appropriate. The total rate leaving state 0 is $a=0.10+0.05=0.15$ per month.

At six months, the probability of remaining in state 0 is

p_0(6)=\exp(-0.15\times 6)=\exp(-0.9)\approx0.4066.

The probability of occupying state 1 is

p_1(6)=\frac{0.10\{\exp(-0.15\times6)-\exp(-0.20\times6)\}}{0.20-0.15} =\frac{0.10(0.4065697-0.3011942)}{0.05} \approx0.2108.

The absorbing-state probability is the remainder:

p_2(6)=1-0.4065697-0.2107509\approx0.3827.

Thus the hypothetical model assigns about 40.66% probability to state 0, 21.08% to state 1, and 38.27% to state 2 at six months. The unrounded probabilities sum to 1. The state 1 probability reflects both entry into progression and subsequent exit to death; it is not simply the probability of a $0\to1$ transition.

How Omics Covariates Change the Transition Pattern

A covariate can be associated with different transitions in different directions. This is why transition-specific coefficients matter. An omics feature associated with a higher rate of progression need not be associated with a higher rate of death from the stable state or from the progressed state. A clinical interpretation should name the particular transition and condition on occupying its origin state.

Worked Example: Comparing Two Omics Profiles

Continue the hypothetical three-state model above and let $X$ be a standardized omics score, with $X=0$ the reference profile and $X=1$ a profile one standard deviation higher. Suppose the transition hazard ratios are 2.0 for $0\to1$, 0.75 for $0\to2$, and 1.5 for $1\to2$. Assume each association is proportional over time and that the Markov illness-death model is appropriate.

For $X=1$, the transition rates are the reference rates multiplied by their respective hazard ratios:

a_{01}=0.10(2.0)=0.20,\qquad a_{02}=0.05(0.75)=0.0375,\qquad a_{12}=0.20(1.5)=0.30

The total rate leaving state 0 is $a=0.20+0.0375=0.2375$ per month. At six months,

p_0(6)=\exp(-0.2375\times6)=\exp(-1.425)\approx0.2405.

Using the occupation-probability formula with $a_{12}-a=0.30-0.2375=0.0625$,

p_1(6)=\frac{0.20\{\exp(-1.425)-\exp(-0.30\times6)\}}{0.0625} =\frac{0.20(0.240510-0.165299)}{0.0625} \approx0.2407.

Therefore,

p_2(6)=1-0.240510-0.240675\approx0.5188.

For this hypothetical model, the higher score corresponds to a higher modeled probability of occupying state 2 at six months: about 51.88%, compared with 38.27% for the reference profile. This contrast follows from all three transition rates, not just the $0\to2$ hazard ratio. It is a model-based association between profile and occupation probability; by itself it does not establish that changing the omics score would change the clinical outcome.

Diagnostics, Alternatives, and Limitations

Interpretation depends on whether the state structure and intensity assumptions are credible. A useful audit begins with the clinical definitions: can two states be distinguished consistently, are transitions recorded with adequate timing, and are important pathways omitted? Then inspect transition counts and state-specific time at risk. Sparse transitions can make a transition-specific coefficient unstable, even when the total sample size is large. Follow-up and censoring should also be examined within origin state and across relevant covariate profiles.

For a proportional-hazards specification, assess whether the covariate association with each transition is reasonably constant over time; the assumption needs checking transition by transition. Compare predicted transitions or occupation probabilities with observed patterns where possible. Poor agreement, systematic residual patterns, or markedly different estimates under reasonable state definitions can signal misspecification. If risk depends on time already spent in the current state, a model based only on the current state may be inadequate; include duration or history in the intensity specification, or use a semi-Markov formulation. If omics measurements change over follow-up, treating them as fixed baseline covariates may also be misleading.

A nonparametric Aalen-Johansen estimator can estimate transition or occupation probabilities without imposing the proportional-hazards regression form, under its censoring and observation conditions. It is useful for descriptive estimates, but does not by itself provide adjusted transition-specific omics associations. Competing-risk methods focus on event types occurring as first events from a common starting state; a multi-state model generalizes this perspective to later transitions and intermediate states. Transition-specific hazard regression is useful for modeling covariate associations, but occupation probabilities require combining the transitions and are sensitive to the adequacy of the full model.

A common mistake is to read a transition hazard ratio of 2 as meaning that the probability of occupying the destination state at a chosen time doubles. That is incorrect: occupation depends on the rate of entering the state, the rate of leaving it, and all other competing paths. For the same reason, a lower hazard on one transition does not necessarily imply a lower probability of death by a fixed time. Report the transition-specific hazard association and, when clinically relevant, the model-based occupation probability at a stated time as distinct results.

Key takeaway: Define the states and allowed transitions first; model each transition rate with its own covariate association; then derive occupation probabilities from the complete transition system. A coefficient describes a transition-specific rate, while an occupation probability summarizes the resulting paths through the states.

Check Your Understanding

Answer each question using the definitions and formulas above. Solutions follow each question.

In the hypothetical model with rates $a_{01}=0.10$, $a_{02}=0.05$, and $a_{12}=0.20$ per month, what is the total rate of leaving state 0, and what is the probability of remaining in state 0 after six months?
Solution: The total rate is $0.10+0.05=0.15$ per month. Thus $p_0(6)=\exp(-0.15\times6)=\exp(-0.9)\approx0.4066$.
If an omics coefficient for the $1\to2$ transition is $\log(1.5)$, interpret it.
Solution: Holding other modeled covariates fixed, a one-unit increase in the omics feature is associated with a 1.5-fold instantaneous rate of transition from state 1 to state 2 among subjects currently in state 1, under the model assumptions. It does not mean a 50% increase in the probability of occupying state 2 by a specified time.
Why can a higher $0\to1$ rate coexist with a lower probability of occupying state 1 at a later time?
Solution: Occupation of state 1 depends on both entry and exit. If the rate of leaving state 1 is also high, subjects may enter state 1 more often but spend less time there. The occupation probability combines these transition rates and all modeled paths.
For an illness-death model with $a=a_{12}$, give the limiting expression for $p_1(t)$.
Solution: Starting in state 0, $p_1(t)=a_{01}t\exp(-at)$. This is the limit of $a_{01}\{\exp(-at)-\exp(-a_{12}t)\}/(a_{12}-a)$ as $a_{12}$ approaches $a$.
In the one-unit-higher profile example, state 0, state 1, and death probabilities at six months are approximately 0.2405, 0.2407, and 0.5188. What check should these values pass?
Solution: They must be nonnegative and sum to 1, up to rounding. Here $0.2405+0.2407+0.5188=1.0000$, consistent with mutually exclusive, exhaustive states at six months.