From Transition Models to State Probabilities
A transition-specific Cox model describes the instantaneous rate of one permitted move between states. A prediction about a patient’s status at a future time requires more: it must combine all the relevant transitions, including pathways that pass through intermediate states. The result is a vector of state-occupation probabilities, one for each state, for a specified starting state and covariate profile.
This tutorial focuses on turning fitted transition-specific coefficients and baseline hazards into those absolute probabilities. The calculations are conditional on the chosen profile; they are not automatically averages over a population. The earlier material on predicting transition probabilities with probtrans establishes the transition-probability framework. Here the emphasis is on how covariates change the transition hazards before those hazards are combined.
The relative hazard multiplier for transition \(r\to s\) is \(\exp\{\beta_{rs}^{T}x\}\). If \(\Lambda_{rs,0}(t)\) is the baseline cumulative hazard, then the profile-specific cumulative hazard is $$ \Lambda_{rs}(t\mid x)=\Lambda_{rs,0}(t)\exp\{\beta_{rs}^{T}x\}. $$ This scaling applies separately to every transition. A covariate may increase one transition hazard, leave another unchanged, and decrease a third. The resulting state probabilities cannot generally be found by scaling a single baseline probability: transitions compete and can lead through intermediate states.
Assumptions and the Calculation
A covariate-specific prediction is meaningful only relative to a clearly defined model and prediction setting. The transition graph must specify which moves are possible, and each fitted model must correspond to the appropriate transition. The Cox proportional-hazards form assumes that, for a fixed covariate profile, each transition hazard is its baseline hazard multiplied by a time-constant relative hazard. The estimated baseline hazards and coefficients must be appropriate for the population and time scale used for prediction.
For the usual state-occupation calculation, the process is treated as Markov on the chosen time scale: given the current state and covariates, the future transition rates do not depend on earlier state history. Covariates are held fixed at the profile values. If hazards depend on time since entering the current state, for example, a prediction based only on time since study origin and current state may not represent the intended process. The time scale and history used by the fitted models therefore need to agree with the prediction question.
To combine transitions, form a matrix of cumulative transition hazards. At a small time increment, its off-diagonal entry for \(r\to s\) is the corresponding cumulative-hazard increment, and its diagonal entry in row \(r\) is minus the sum of the outgoing increments from state \(r\). The transition-probability matrix over time is obtained by multiplying the incremental transition matrices in time order (the product integral). The row for the specified starting state gives the state-occupation probabilities.
Writing \(I\) for the identity matrix, the transition-probability matrix from time 0 to time \(t\) is $$ P(0,t\mid x)=\prod_{0<u\le t}\{I+dA(u\mid x)\}. $$ In a fitted Cox model, baseline cumulative hazards are commonly step functions, so the calculation uses their estimated increments at event times. Software implements the corresponding product-integral calculation; hand calculations below use constant rates over an interval to make the arithmetic transparent.
Worked Examples
Worked Example: Scaling Each Transition Hazard for a Profile
Suppose a hypothetical three-state illness-death model has transitions \(1\to2\), \(1\to3\), and \(2\to3\). A binary covariate \(x\) has fitted coefficients \(\beta_{12}=\log(2)\), \(\beta_{13}=\log(1.5)\), and \(\beta_{23}=\log(0.5)\). For illustration, suppose the fitted baseline hazards are constant over the first two years at 0.20, 0.10, and 0.15 per year, respectively. Calculate the cumulative hazards for profile \(x=1\).
For transition \(1\to2\), the baseline cumulative hazard by year 2 is \(0.20\times2=0.40\), and the covariate multiplier is \(\exp\{\log(2)\}=2\). Thus $$ \Lambda_{12}(2\mid x=1)=0.40\times2=0.80. $$
For transition \(1\to3\), the baseline cumulative hazard is \(0.10\times2=0.20\), and the multiplier is \(\exp\{\log(1.5)\}=1.5\). Therefore \(\Lambda_{13}(2\mid x=1)=0.20\times1.5=0.30\). For transition \(2\to3\), it is \(0.15\times2=0.30\), multiplied by \(\exp\{\log(0.5)\}=0.5\), giving \(\Lambda_{23}(2\mid x=1)=0.30\times0.5=0.15\).
These are transition-specific cumulative hazards, not probabilities. In particular, \(0.80\) is not an 80% probability of moving from state 1 to state 2 by year 2. The probability depends on all exits from state 1 and, for state occupation, on transitions after leaving state 1.
Worked Example: State Occupancy at Two Years for Two Profiles
Use the same hypothetical model and constant baseline rates as in the previous calculation. Assume the person starts in state 1. For profile \(x=0\), the transition rates are \(a=\lambda_{12}=0.20\), \(b=\lambda_{13}=0.10\), and \(c=\lambda_{23}=0.15\) per year. For profile \(x=1\), the rates are \(a=0.40\), \(b=0.15\), and \(c=0.075\) per year. State 3 is absorbing. Find the probabilities of occupying states 1, 2, and 3 at year 2.
For constant rates, let \(q=a+b\) be the total rate of leaving state 1. The probability of remaining in state 1 at time \(t\) is \(e^{-qt}\). The probability of occupying state 2 at time \(t\) is $$ P_{12}(t)=\int_0^t a e^{-qu}e^{-c(t-u)}\,du =\frac{a\{e^{-ct}-e^{-qt}\}}{q-c}, \qquad q\ne c. $$ The integral counts paths that move from state 1 to state 2 at time \(u\), remain in state 1 until \(u\), and then avoid the transition \(2\to3\) until \(t\). Since the three states exhaust the possibilities, \(P_{13}(t)=1-P_{11}(t)-P_{12}(t)\).
For \(x=0\), \(q=0.20+0.10=0.30\) per year. At \(t=2\), $$ P_{11}(2)=e^{-0.30(2)}=e^{-0.60}=0.5488, $$ and $$ P_{12}(2)=\frac{0.20\{e^{-0.15(2)}-e^{-0.30(2)}\}}{0.30-0.15} =\frac{0.20(0.7408-0.5488)}{0.15}=0.2560. $$ Thus \(P_{13}(2)=1-0.5488-0.2560=0.1952\). The predicted occupation probabilities are 0.5488, 0.2560, and 0.1952 for states 1, 2, and 3, respectively.
For \(x=1\), \(q=0.40+0.15=0.55\) per year. Then $$ P_{11}(2)=e^{-0.55(2)}=e^{-1.10}=0.3329, $$ and $$ P_{12}(2)=\frac{0.40\{e^{-0.075(2)}-e^{-0.55(2)}\}}{0.55-0.075} =\frac{0.40(0.8607-0.3329)}{0.475}=0.4445. $$ Therefore \(P_{13}(2)=1-0.3329-0.4445=0.2226\). Compared with \(x=0\), the \(x=1\) profile has a lower probability of remaining in state 1 and a higher probability of occupying state 2 at year 2. Its probability of being in state 3 is also higher in this example, but that conclusion reflects the combined transition rates, not just the direction of one coefficient.
Worked Example: Changing the Starting State
For the \(x=0\) profile in the same model, suppose instead that the person starts in state 2. State 2 has only one possible exit, to absorbing state 3, at rate \(c=0.15\) per year. At year 2, $$ P_{22}(2)=e^{-0.15(2)}=e^{-0.30}=0.7408, \qquad P_{23}(2)=1-0.7408=0.2592. $$
The year-2 prediction is therefore 0 for state 1, 0.7408 for state 2, and 0.2592 for state 3. This is a different row of the transition-probability matrix from the one used in the preceding example. A prediction must specify its starting state as well as its profile and time.
Practical Workflow, Checks, and Alternatives
For a fitted model, the calculation can be organized as a sequence. For each profile, supply its covariate values to the relevant transition-specific models, scale each baseline-hazard increment by the corresponding relative hazard, and combine the increments in time order. Request or extract the transition-probability row for the intended starting state, then report all state probabilities at the specified follow-up time. In software workflows based on mstate, the profile information supplied to the transition-model prediction step must align with the transition definitions and the fitted covariate structure.
Name the starting state, covariate profile, time origin, follow-up time, and state ordering.
For each allowed transition, combine its estimated baseline hazard with the profile-specific exponential coefficient multiplier.
Use the product integral of the profile-specific transition intensity matrices, respecting the permitted transition graph.
Give the complete occupation-probability vector, check its total and values, and state the profile and time units.
Several checks catch common errors. At each requested time, probabilities across all states should sum to 1, up to numerical rounding, and should lie between 0 and 1. A state that cannot be reached from the starting state should have probability zero. An absorbing state cannot be left, although the probability of occupying it can change over time. Confirm that baseline hazard time units match the follow-up time units, and that profile values are in the same coding and reference categories used to fit the models. The transition-specific proportional-hazards assumption should also be considered; earlier material on Schoenfeld residuals covers one diagnostic for that assumption.
A common mistake is to multiply a baseline state probability by a hazard ratio and call the result the profile-specific state probability. A hazard ratio scales an instantaneous transition hazard; it does not scale the final probability. The correction is to scale each relevant transition hazard separately and recompute the entire multi-state probability calculation. Another mistake is to report a profile-specific prediction as the population’s average risk. An average requires averaging predictions over an explicitly chosen covariate distribution; a prediction at one profile does not do this.
This approach is useful when the question concerns adjusted predictions for meaningful covariate profiles. An empirical Aalen–Johansen estimate instead summarizes observed state transitions without conditioning on a specified regression profile. Competing-risk prediction methods can answer questions about particular event probabilities, but a multi-state model also tracks occupancy across intermediate and absorbing states. These approaches have different targets; choosing among them depends on the requested estimand, not simply on which gives the most convenient number.
Predictions inherit the limitations of the fitted models. Sparse transitions can make baseline hazards unstable, unusual covariate profiles can require extrapolation, and misspecified transition structures or proportional hazards can distort the resulting occupancy probabilities. The product-integral calculation also relies on a coherent set of transition hazards; a coding error in transition ordering can produce probabilities that look numerically plausible but correspond to the wrong states. These are model-based point predictions. Their uncertainty requires accounting for estimation of both the coefficients and baseline hazards.
Check Your Understanding
Try each question, then compare your reasoning with the solutions.
- A transition has baseline cumulative hazard 0.25 by year 3 and a binary covariate coefficient \(\log(2)\). What is the profile-specific cumulative hazard when the covariate equals 1?
- Why is a transition hazard ratio of 2 not generally a state-occupation probability ratio of 2?
- In the three-state example, what is the year-2 probability of state 1 for \(x=0\), and what is the total probability of states 2 and 3?
- For the \(x=0\) profile starting in state 2, what is the probability of being in state 3 at year 2?
- Name two checks to perform before reporting a profile-specific prediction.
Solutions. (1) \(0.25\times\exp\{\log(2)\}=0.25\times2=0.50\). (2) A hazard ratio scales one instantaneous transition hazard; state occupancy depends on competing exits, later transitions, time, and the starting state. (3) \(P_{11}(2)=0.5488\); the probability of states 2 or 3 is \(1-0.5488=0.4512\), also \(0.2560+0.1952=0.4512\). (4) \(1-e^{-0.15(2)}=1-0.7408=0.2592\). (5) For example, verify that the state probabilities sum to 1 and lie in \([0,1]\), and confirm the starting state, time units, covariate coding, and permitted transitions.