Why Examine Covariate Form?
A transition-specific Cox model can represent the effect of a continuous covariate \(x\) with a linear term, so that a one-unit increase in \(x\) multiplies the transition hazard by the same factor at every value of \(x\). This is a strong functional-form assumption. If the effect bends, levels off, or changes more sharply in part of the covariate range, a linear term can obscure that structure and give an inadequate description of the transition hazard.
Martingale residuals compare observed transition events with the number expected under a fitted model. Their patterns against a continuous covariate can help identify a poor functional form. They are diagnostics, not a standalone test: a smooth curve can suggest where to investigate, but it does not establish that a particular nonlinear model is correct.
Definition for a Transition-Specific Cox Model
Let \(r\) denote a transition and \(i\) a subject. Let \(N_{ir}(t)\) count the subject’s observed transitions of type \(r\) by time \(t\), and let \(Y_{ir}(t)\) indicate whether that subject is at risk for transition \(r\) at time \(t\). The transition-specific Cox model specifies the event intensity as
Here \(\mathcal{F}_{t-}\) represents the observed history just before time \(t\), \(\lambda_{0r}(t)\) is the baseline hazard for transition \(r\), and \(\mathbf{Z}_i(t)\) contains the covariates in the model. For a linear-form check, one component of \(\mathbf{Z}_i(t)\) is the continuous covariate \(x_i(t)\), entered linearly. The model assumes that the covariate effects have this specified form and are proportional over time; it also requires risk sets and covariate histories to be correctly defined for the transition being modeled.
In practice, the baseline cumulative hazard is estimated from the transition-specific risk sets, and the regression coefficients are estimated from the transition-specific Cox partial likelihood. Residuals are calculated using those fitted quantities. A residual is not an individual probability of experiencing the transition: it is a difference between an observed count and a model-based expected count.
Why a Residual Plot Can Reveal Curvature
If the model’s event intensity is correctly specified, the observed counting process minus its cumulative intensity is a martingale. In particular, its increments have conditional mean zero given the observed history up to the preceding time. A persistent residual pattern against a fixed covariate is therefore a warning that the fitted intensity may not describe event accumulation equally well across its values.
There are two useful versions of the diagnostic. First, fit a model that includes the other relevant covariates but omits the continuous covariate being examined. Plot those reduced-model martingale residuals against the omitted covariate. A smooth trend suggests that the omitted variable may matter; bending in the trend can suggest that a linear contribution is insufficient. Second, plot residuals from the model that includes the covariate linearly. A remaining smooth pattern can signal that the linear specification has not captured the covariate relationship. Neither plot directly estimates a hazard ratio curve, and the reduced-model plot in particular should not be read as an exact estimate of the omitted effect.
A smooth curve, such as a locally weighted smoother, is commonly added to the scatterplot. Its purpose is to summarize broad structure rather than follow every residual. The zero line is a reference: repeated departures above or below it, or a systematic curved shape, merit investigation. Individual residuals can be noisy, especially when events are few.
Worked Example: Calculating Residuals for One Transition
Hypothetical data. Six subjects are followed for transition \(r\), and each can experience this transition at most once during follow-up. A reduced Cox model includes an adjustment covariate but omits continuous covariate \(x\). Suppose the fitted cumulative intensities at the subjects’ observed follow-up times are \(0.10, 0.16, 0.25, 0.40, 0.55,\) and \(0.80\), respectively. These values already account for each subject’s at-risk time. The observed event indicators and residual calculations are:
| Subject | \(x\) | Observed event \(N_{ir}\) | Fitted cumulative intensity | Residual \(N_{ir}-\widehat{A}_{ir}\) |
|---|---|---|---|---|
| 1 | 1 | 0 | 0.10 | \(-0.10\) |
| 2 | 2 | 0 | 0.16 | \(-0.16\) |
| 3 | 3 | 1 | 0.25 | \(0.75\) |
| 4 | 4 | 0 | 0.40 | \(-0.40\) |
| 5 | 5 | 1 | 0.55 | \(0.45\) |
| 6 | 6 | 1 | 0.80 | \(0.20\) |
For subject 5, the calculation is \(1-0.55=0.45\). For subject 4, it is \(0-0.40=-0.40\). A positive residual means that more events were observed than the fitted cumulative intensity would suggest for that subject; a negative residual means fewer were observed.
If a smooth through these points showed a broad rise and fall across \(x\), it would prompt examination of whether a linear term for \(x\) is adequate. These six observations alone do not establish curvature: the residuals are noisy, and the apparent shape must be judged in relation to the risk sets, the amount of data, and the smooth’s stability.
Using the Diagnostic in Practice
For a focused check, retain the other model terms and examine one continuous covariate at a time. In multi-state data, the risk set is transition-specific: a subject contributes only while eligible for the transition under consideration. Compute and plot residuals for the same transition whose covariate form is in question. Combining residuals from distinct transitions can blur different baseline hazards, covariate effects, and risk processes.
Confirm which subjects are at risk for transition \(r\), when they enter and leave that risk set, and how the continuous covariate is measured or updated.
Retain the other prespecified covariates. To screen for a possible effect, omit the target covariate; to check the adequacy of its linear form, include it linearly and examine remaining residual structure.
For each subject or appropriate observation interval, subtract fitted cumulative intensity from the observed transition count. Plot the resulting residuals against the target covariate and add a restrained smooth.
Check how many events and at-risk observations support each part of the curve. If a stable departure from linearity remains plausible, compare with a flexible representation, such as a spline form, and assess the fitted relationship on the hazard scale.
The justification for subtracting the fitted cumulative intensity is direct. Under the specified intensity model, the increment \(dN_{ir}(t)\) has conditional expectation \(Y_{ir}(t)\lambda_{0r}(t)\exp\{\boldsymbol{\beta}_r^{\mathsf T}\mathbf{Z}_i(t)\}\,dt\) over an infinitesimal interval, in the continuous-time formulation. Subtracting this conditional expected increment yields a mean-zero increment:
Accumulating these increments gives the martingale residual. If the covariate form is wrong, the fitted cumulative intensity can be systematically too large in some covariate ranges and too small in others. The resulting residuals may show structure against \(x\). This reasoning explains the diagnostic, but it does not guarantee that every misspecified form produces a clear visual pattern in a finite sample.
Worked Example: Interpreting a Curved Residual Smooth
Hypothetical diagnostic. For transition \(r\), a linear model includes age in years and other prespecified covariates. A smooth of the fitted-model martingale residuals against age lies near zero at ages 40–55, below zero around 60–65, and above zero around 70–75. Assume the plot represents adequate at-risk support across these ranges and is not driven by one observation.
The curve is a departure from a flat residual pattern: the model’s residuals vary systematically with age. A possible interpretation is that the linear age term does not capture the age relationship for this transition. The plot does not show that the hazard itself decreases and then increases in the same way; it shows a residual pattern that warrants checking whether a more flexible age effect improves the model’s description.
A reasonable next step is to fit a model with a flexible age function, using the same transition risk sets and adjustment terms, and compare the resulting relationship and uncertainty with the linear specification. The conclusion should be based on the combined evidence, not on the smooth alone. If the apparent bends occur only where there are very few events, they may be unstable rather than persuasive evidence of nonlinearity.
Worked Example: Why Diagnostics Must Be Transition-Specific
Hypothetical data. Consider two transitions, \(r=1\) and \(r=2\), and a subject eligible for both at different times. For transition 1, the subject has one observed event and a fitted cumulative intensity of \(0.70\), giving \(\widehat{M}_{i1}=1-0.70=0.30\). For transition 2, the subject has no observed event and a fitted cumulative intensity of \(0.20\), giving \(\widehat{M}_{i2}=0-0.20=-0.20\).
These calculations concern different event processes and their own at-risk histories. Adding the residuals would give \(0.30+(-0.20)=0.10\), but that sum is not a residual for either transition and could conceal the opposing signs. To assess whether \(x\) is adequately modeled for transition 1, plot transition-1 residuals against \(x\); repeat separately for transition 2 if its covariate form is also in question.
If the same covariate has transition-specific effects, a pattern for one transition does not establish a pattern for another. The transition-specific model and residual calculation must match the question being asked.
Diagnostics, Alternatives, and Limitations
A martingale residual plot is most informative when considered alongside checks of the data structure and the fitted model. Before interpreting a smooth, verify that transition status and at-risk intervals were coded correctly, that the covariate values used correspond to the correct times, and that enough events support the covariate range. Sparse support at an endpoint can make a smoother bend sharply even when the underlying relationship is uncertain. Also inspect whether a few observations dominate the pattern.
A residual pattern can have explanations other than a nonlinear covariate effect. An inadequate risk-set definition, omitted relevant covariates, or a time-related departure from the model can create structure. Schoenfeld residual diagnostics examine proportional-hazards behavior over time; they answer a different question from a martingale residual plot against \(x\). A time-related issue should not automatically be “fixed” by adding curvature in \(x\), and apparent curvature should not be treated as evidence about time variation.
| Approach | What it contributes | Main caution |
|---|---|---|
| Martingale residual plot | Visual screen for residual structure against a continuous covariate | Shape is noisy and does not identify a unique correct model |
| Linear term for the covariate | Simple, interpretable constant log-hazard change per unit increase | Can miss bends or plateaus in the covariate relationship |
| Flexible function, such as a spline | Allows the log hazard to change nonlinearly with the covariate | Needs adequate data support; a highly flexible fit can follow noise |
A residual plot is a visual diagnostic, not a formal hypothesis test with a universal cutoff. It is also not a substitute for estimating and interpreting the fitted hazard relationship. Flexible forms add parameters and can increase uncertainty, especially where events are sparse. Consider whether the suggested shape is stable and scientifically interpretable, and report the chosen covariate form and its uncertainty.
Common Mistake
Wrong reasoning: “Residuals are positive for high values of \(x\), so high \(x\) must increase the transition hazard.” A positive residual only says that the observed count exceeded the fitted cumulative intensity for that observation. It does not, by itself, quantify an effect, establish its direction over the whole covariate range, or demonstrate causation.
Correction: Look for a reproducible pattern across the supported covariate range, check that the transition-specific model and risk sets are correct, and compare a flexible covariate form with the linear specification. Interpret the fitted function on the hazard scale and in context; treat the residual plot as evidence to investigate, not as the final estimate.
Check Your Understanding
Answer each question, then compare your reasoning with the short solutions.
- A subject has no observed transition-\(r\) event and fitted cumulative intensity \(0.35\). What is the martingale residual, and what does its sign mean?
- Why should residuals from two different transitions generally be plotted separately when checking a covariate form?
- A smooth of full-model residuals against a continuous covariate bends away from zero only where there are very few events. What should you conclude?
- What is the difference between a residual plot against a covariate and a Schoenfeld residual diagnostic?
- Under the correct intensity model, why does subtracting cumulative intensity from the observed event count produce a martingale?
Solutions. (1) The residual is \(0-0.35=-0.35\). Fewer events were observed than the fitted cumulative intensity would suggest for this subject; this is not a probability. (2) Each transition has its own event process and at-risk history, so combining their residuals can hide transition-specific patterns. (3) Treat the bend as uncertain: sparse event support can make a smooth unstable, so do not infer nonlinearity without further evidence. (4) A martingale residual plot against a covariate screens for structure related to covariate form; Schoenfeld residual diagnostics examine whether effects vary over time. (5) Under the specified intensity, the conditional expected event increment equals the intensity increment. Subtracting it gives a conditional mean-zero increment, and accumulating those increments yields the martingale residual.