Why Event Order and Gap Time Matter
A recurrent-event analysis must specify both who is at risk for each event and what time scale describes their risk. The Prentice-Williams-Peterson (PWP) gap-time model makes two deliberate choices: a subject enters the risk set for event order \(k\) only after experiencing \(k-1\) events, and the clock for that interval is reset to zero at the preceding event. It therefore models the waiting time to the next event conditional on reaching that event order.
This setup is useful when the first, second, and later events have different meanings or when time since the most recent event is clinically relevant. For example, the hazard of a subsequent hospitalization might depend on how long it has been since the previous hospitalization. A PWP fit does not treat all event records as interchangeable: it conditions on event history through the event-order-specific risk sets.
Model and Assumptions
A common specification is a stratified Cox model. For subject \(i\)'s \(k\)th event, let \(h_{ik}(u\mid\mathcal H_{i,k-1},x_i)\) denote the hazard at gap time \(u\), conditional on the observed history through event \(k-1\), \(\mathcal H_{i,k-1}\), and covariates \(x_i\). The model is
The baseline hazard \(h_{0k}(u)\) is allowed to differ by event order. The coefficient \(\beta_k\) may also differ by event order, as written here. If scientific and statistical considerations support a common covariate effect, one may instead constrain \(\beta_k=\beta\) across orders while retaining separate baseline hazards. These are distinct modeling choices: stratifying by order alone does not make covariate effects order-specific.
Interpretation requires a correctly defined event history and risk intervals, proportional hazards for the modeled covariate effects within each event-order stratum, and censoring that is non-informative for the event process conditional on the modeled covariates and observed history. Subjects must not enter an order-\(k\) risk set before actually reaching that order. If event times are tied, the software’s tie method (such as Efron’s approximation) should be stated; with few ties, the choice is often less influential, but it does not repair incorrect risk sets.
Constructing the Gap-Time Risk Sets
Each event-order interval has a start at gap time zero and an endpoint at the next event or censoring. For the first event, follow-up begins at the study origin. For the second event, follow-up begins at the first event; for the third, it begins at the second event. Censoring after a subject’s first event still supplies information about that subject’s second-event hazard up to the censoring time, but it does not make the subject eligible for the third-event risk set.
Worked Example: Who Is at Risk for the Second and Third Events?
Suppose five hypothetical participants have the following event histories, with times measured in months from study entry. A dash means no further event before censoring.
| Participant | First event | Second event | Third event | Censoring time |
|---|---|---|---|---|
| A | 2 | 5 | — | 8 |
| B | 3 | — | — | 6 |
| C | 1 | 4 | 7 | 9 |
| D | 4 | — | — | 9 |
| E | 2 | 6 | — | 7 |
For event order 2, subtract the first-event time from the second-event or censoring time. Participant A has a second-event gap of \(5-2=3\) months; B is censored at gap \(6-3=3\); C has a second-event gap of \(4-1=3\); D is censored at gap \(9-4=5\); and E has a second-event gap of \(6-2=4\). All five have experienced a first event, so all can enter the second-event risk set. At gap time \(u=2\), all five remain under observation and at risk for event 2.
For event order 3, only A, C, and E have experienced two events. Their third-event intervals start at the second event: A is censored after \(8-5=3\) months; C has a third event after \(7-4=3\) months; and E is censored after \(7-6=1\) month. B and D are not eligible for this risk set because neither experienced a second event. At \(u=1.5\), A and C are at risk for event 3; E has already been censored.
This calculation checks both defining features of the model: the event-order condition determines eligibility, and the elapsed time is measured from the preceding event rather than from study entry.
Fitting and Interpreting the Model
In a long-format data set, create one row for each subject and event-order interval that is observed. Include an event-order variable, the gap-time endpoint, and an indicator equal to 1 if the event of that order occurred and 0 if the interval ended in censoring. A Cox fit then uses gap time as its analysis time and event order as the stratum. If covariate effects may vary by order, include interactions between the covariates and event order, or fit separate models by order when that is appropriate.
For a binary treatment \(x\), the hazard ratio for event order \(k\) is \(\exp(\beta_k)\). It compares the instantaneous hazard of the \(k\)th event at the same gap time among people who have already experienced \(k-1\) events and remain at risk, for otherwise equal modeled covariates. It is not a comparison of the probability of ever having that event, nor of the cumulative number of events over follow-up.
Worked Example: Interpreting an Event-Order-Specific Estimate
Consider a hypothetical PWP gap-time fit with treatment coded 1 for the intervention and 0 for control. For the second-event stratum, suppose the fitted treatment coefficient is \(\hat\beta_2=-0.223\) with standard error \(SE=0.100\). Assume the gap-time risk sets and censoring definition are correct, and that the proportional-hazards assumption is reasonable for treatment in this stratum.
The estimated treatment hazard ratio for the second event is
Thus, among participants who have had a first event and are still at risk for a second, the fitted instantaneous second-event hazard under treatment is about 0.80 times the control hazard at the same gap time, conditional on the modeled covariates. This is an estimated 20.0% lower hazard, not evidence that treatment reduces the probability of a second event by 20.0%.
A Wald statistic for \(H_0:\beta_2=0\) is
At a two-sided significance level of 0.05, this result is inconsistent with \(H_0\) under the model. A 95% Wald interval for the log hazard ratio is \(-0.223\pm1.96(0.100)=[-0.419,-0.027]\). Exponentiating the endpoints gives a 95% confidence interval for the hazard ratio of \([\exp(-0.419),\exp(-0.027)]=[0.6576,0.9734]\). In repeated applications under the model, this procedure has 95% coverage; the interval is not a 95% probability statement about the fixed parameter.
Worked Example: Allowing Treatment Effects to Vary by Event Order
Suppose a hypothetical model includes treatment-by-order interactions and reports treatment coefficients \(\hat\beta_1=-0.105\), \(\hat\beta_2=-0.223\), and \(\hat\beta_3=0.182\). Treat these as event-order-specific estimates, not as one common treatment coefficient.
The corresponding hazard ratios are \(\exp(-0.105)=0.9003\) for the first event, \(\exp(-0.223)=0.8001\) for the second, and \(\exp(0.182)=1.200\) for the third. The fitted treatment association is therefore lower hazard for the first and second events and higher hazard for the third, each comparison conditional on being eligible and at risk in that order’s risk set. These point estimates alone do not establish that effects differ across orders: a formal contrast or interaction test, together with uncertainty estimates, is needed for that conclusion.
Why the Partial Likelihood Uses These Risk Sets
Within event order \(k\), suppose an event occurs at gap time \(u\), and let \(R_k(u)\) be the subjects who have reached order \(k\) and remain at risk immediately before \(u\). Under the Cox model, the instantaneous event rate for subject \(i\) is proportional to \(\exp(\beta_k^\mathsf{T}x_i)\). Conditional on one event occurring at that time in this stratum, the probability that subject \(i\) is the one who experiences it is its relative hazard within \(R_k(u)\):
Multiplying these contributions across event times gives the event-order-specific partial likelihood. The baseline hazard cancels from each conditional probability, which is why the covariate coefficients can be estimated without specifying \(h_{0k}(u)\). Importantly, the denominator must use the conditional PWP risk set, not everyone still enrolled in the study.
Order each subject’s recurrent events and confirm that each event has the correct preceding event and time origin.
For each observed order, calculate the gap from the previous event to the next event or censoring, and mark whether that endpoint is an event.
Stratify the Cox baseline hazard by event order; add covariate-by-order interactions only when event-order-specific effects are intended.
Inspect the number of events and risk-set membership by order, assess proportional hazards within strata, and report the time scale and conditional interpretation.
Diagnostics, Alternatives, and Limitations
A practical diagnostic begins with counts by event order: how many subjects become eligible, how many events occur, and how many intervals are censored. Sparse high-order strata can produce unstable estimates or wide intervals. Check that calculated gap times are nonnegative, that no subject enters order \(k\) without \(k-1\) prior events, and that an event is not accidentally coded as censoring. Plot or test scaled Schoenfeld residuals within strata to assess proportional hazards; a systematic time pattern suggests that a constant hazard ratio may be inadequate.
A common mistake is to use time since study entry for every event while claiming to fit a gap-time model. That gives later-event risk intervals an origin at entry, not at the preceding event, and changes the estimand and risk-set comparisons. Another mistake is to include a subject in the second-event risk set before their first event. The correction is to calculate the interval from the immediately preceding event and include only subjects who have reached the relevant order.
The Andersen-Gill model uses a common event-process time scale and a risk set that can include subjects according to their observation status, without the same event-order conditioning. PWP is preferable when inference should condition on having reached each event order; Andersen-Gill can be more suitable when a common recurrent-event intensity is the target. Neither approach is automatically best. PWP can lose precision at high event orders, and results depend on the chosen time scale, event-order structure, and censoring assumptions. A PWP coefficient is a conditional hazard association, not a causal effect without additional design and causal assumptions.
Repeated intervals from one person are not independent observations in the ordinary sampling sense. Variance estimation should reflect the subject-level structure appropriate to the analysis; the fitted conditional model and the variance estimator are separate choices. Also assess whether treatment or other covariates change over time and whether terminal events prevent future recurrence: these features may require a more specific model or estimand rather than being handled by mechanically adding rows.
Check Your Understanding
Answer each question, then compare your reasoning with the solutions below.
- A subject has events at months 2 and 9, then is censored at month 12. What are the observed gap time and status for the second-event interval?
- Can a participant with one observed event contribute to the third-event risk set? Explain.
- A fitted order-3 treatment coefficient is \(\log(1.25)\). State the hazard-ratio interpretation in context.
- Why does stratifying the baseline hazard by event order not, by itself, give different covariate effects by order?
Solutions. (1) The second-event gap is \(9-2=7\) months, and its status is event. Censoring at month 12 occurs after that event and does not change this completed interval. (2) No. The third-event risk set requires two prior events, so a participant with only one event has not reached order 3. (3) The order-3 treatment hazard ratio is \(\exp\{\log(1.25)\}=1.25\): among subjects who have experienced two prior events and remain at risk for a third, the treated group has an estimated 1.25 times the instantaneous third-event hazard of the control group at the same gap time, conditional on modeled covariates. (4) Stratification permits a different baseline hazard \(h_{0k}(u)\) for each order, but a shared coefficient \(\beta\) still imposes a common covariate hazard ratio; order-specific effects require separate coefficients, such as treatment-by-order interactions.