Aalen-Johansen estimator of multi-state transition and occupancy probabilities.
Given counting-process intervals (start, stop] each labelled with the state occupied (state) and the state transitioned to at stop (event, or a censoring marker), this forms the Aalen-Johansen product \(P(0, t) = \prod (I + dA(s))\) and reports the state occupancy probabilities over time. Occupancy probabilities are validated to tolerance against R’s survfit multi-state pstate. (Competing risks and Kaplan-Meier are special cases handled by AalenJohansen and KaplanMeier.)
Returns
Fitted estimator-
Call fit() to produce a fitted estimator with cached results (
states_, time_, occupancy_, and internal transition matrices), accessible as tidy DataFrames.
Examples
The mgus2 patients occupy three states in turn: "mgus" at entry, then possibly "pcm" (plasma-cell malignancy), then "death". Reshape the wide dataset into counting-process intervals (start, stop], one interval per state occupied, labelled with the state entered next. A patient who progresses before dying contributes two intervals; everyone else contributes one. Fitting reports the occupancy probability of each state over time.
import greenwood as gw
mg = gw.load_dataset("mgus2", backend="pandas")
start, stop, state, event = [], [], [], []
for i in range(len(mg)):
pt, ft = mg["ptime"][i], mg["futime"][i]
progressed, died = mg["pstat"][i] == 1, mg["death"][i] == 1
if progressed and pt < ft:
start += [0, pt]; stop += [pt, ft]; state += ["mgus", "pcm"]
event += ["pcm", "death" if died else None]
else:
start += [0]; stop += [ft]; state += ["mgus"]
event += ["death" if died else ("pcm" if progressed else None)]
rows = [(a, b, s, e) for a, b, s, e in zip(start, stop, state, event) if b > a]
start, stop, state, event = map(list, zip(*rows))
ms = gw.MultiState().fit(start, stop, state, event, states=("mgus", "pcm", "death"))
ms.to_frame(format="polars")
shape: (276, 4)| time | mgus | pcm | death |
|---|
| f64 | f64 | f64 | f64 |
| 1.0 | 0.969653 | 0.0 | 0.030347 |
| 2.0 | 0.947961 | 0.001446 | 0.050593 |
| 3.0 | 0.937114 | 0.001446 | 0.061439 |
| 4.0 | 0.924822 | 0.002169 | 0.073009 |
| 5.0 | 0.916868 | 0.002892 | 0.080239 |
| … | … | … | … |
| 356.0 | 0.08175 | 0.0 | 0.91825 |
| 373.0 | 0.0545 | 0.02725 | 0.91825 |
| 376.0 | 0.0545 | 0.02725 | 0.91825 |
| 394.0 | 0.0545 | 0.02725 | 0.91825 |
| 424.0 | 0.0 | 0.02725 | 0.97275 |
Methods
|
Name
|
Description
|
|
fit()
|
Fit a multi-state model using counting-process intervals.
|
|
predict()
|
State occupancy probabilities at specified times.
|
|
to_frame()
|
Return occupancy probabilities over time as a DataFrame.
|
fit()
Fit a multi-state model using counting-process intervals.
fit(start, stop, state, event, *, states=None)
Estimates state-occupancy probabilities and transition dynamics over time from counting-process data (multiple overlapping intervals per subject). The model tracks how subjects move between states and computes the probability of being in each state at any given time, accounting for censoring and competing transitions.
This estimator is ideal for:
- Panel data: Subjects observed at discrete times, with state changes recorded between observations.
- Chronic-disease progression: Modeling progression through stages (e.g., MGUS → PCM → death).
- Multi-event data: Non-absorbing or semi-absorbing intermediate states.
- Irregular follow-up: Each subject’s observation times may differ.
The model estimates without distributional assumptions via non-parametric maximum likelihood. Occupancy probabilities are computed as a product of transition matrices evaluated at each event time.
Parameters
start: Any-
Start time of each interval. Can be a 1-D array-like (or Polars/Pandas Series). Intervals are half-open: (start, stop].
stop: Any-
Stop (end) time of each interval. Must have the same length as start. Intervals define subject-time windows.
state: Any-
The state occupied during each interval (the “from” state). Can be string, int, or other hashable label. Must have the same length as start and stop.
event: Any-
The state transitioned to at the stop time. If None, NaN, or 0, the subject was censored (no transition). Otherwise, must be a valid state label. Must have the same length as start and stop.
states: Any = None-
Optional ordered sequence of all state labels (default: auto-detected from data). If provided, must include all unique states in
state and event. Useful for enforcing a specific state ordering (e.g., disease progression order) or including states with no observed transitions.
Returns
MultiState-
The fitted estimator object itself (for method chaining) with cached results (
states_, time_, occupancy_, transition_) accessible via to_frame() (optionally format=). Occupancy probabilities and transition probabilities can be queried at any time via predict().
Details
Data format: Intervals are half-open (start, stop]. Each row represents a subject-interval: the period during which the subject was in state and either remained (censored) or transitioned to event at stop.
State labels: States can be strings, integers, or other hashable types (e.g., tuples). Mixed types are not allowed. Transitions between the same state (self-loops) treated as censoring.
Handling censoring: Censored intervals (event = None/NaN/0) contribute right-censored data. Subjects re-enter their original state after censoring (common in discrete-time or periodic follow-up studies).
Computational method: Non-parametric maximum likelihood. At each distinct event time, transition intensities are estimated from risk sets (subjects at risk to transition from each state), and occupancy is updated by matrix multiplication of transition probabilities.
Examples
Build a multi-state model from counting-process intervals. First, prepare interval data from a chronic-disease cohort:
import greenwood as gw
mg = gw.load_dataset("mgus2", backend="pandas")
start, stop, state, event = [], [], [], []
for i in range(len(mg)):
pt, ft = mg["ptime"][i], mg["futime"][i]
progressed, died = mg["pstat"][i] == 1, mg["death"][i] == 1
if progressed and pt < ft:
start += [0, pt]; stop += [pt, ft]; state += ["mgus", "pcm"]
event += ["pcm", "death" if died else None]
else:
start += [0]; stop += [ft]; state += ["mgus"]
event += ["death" if died else ("pcm" if progressed else None)]
rows = [(a, b, s, e) for a, b, s, e in zip(start, stop, state, event) if b > a]
start, stop, state, event = map(list, zip(*rows))
ms = gw.MultiState().fit(start, stop, state, event, states=("mgus", "pcm", "death"))
ms
MultiState (Aalen-Johansen multi-state model)
states: mgus, pcm, death
times: 276
final occupancy
mgus 0
pcm 0.02725
death 0.9727
Query occupancy probabilities at specific follow-up times (60, 120, and 240 months); pass format= to choose the backend (here, Polars):
ms.predict([60, 120, 240], format="polars")
shape: (3, 4)| time | mgus | pcm | death |
|---|
| f64 | f64 | f64 | f64 |
| 60.0 | 0.645529 | 0.016006 | 0.338464 |
| 120.0 | 0.40446 | 0.012027 | 0.583513 |
| 240.0 | 0.176158 | 0.011489 | 0.812353 |
predict()
State occupancy probabilities at specified times.
predict(times, *, format=None)
Evaluates the state occupancy probabilities (probability of being in each state) at requested times. The occupancy probability for each state is a right-continuous step function, defined at the fitted time points and interpolated or held constant elsewhere.
Parameters
times: Any-
Time points at which to evaluate state occupancy. Can be a scalar or array-like of floats. Values before the first transition time use the initial distribution; values after the last transition time use the final distribution.
format: str | None = None-
Output format:
None (default), "pandas", "polars", or "pyarrow". When None, a backend is auto-detected (Polars, then Pandas, then PyArrow).
Returns
pandas.DataFrame, polars.DataFrame, or pyarrow.Table-
A table with a
time column and one column per state, containing occupancy probabilities (values between 0 and 1) at each query time.
Examples
Evaluate occupancy probabilities at specific times. The result shows how the probability of being in each state changes over time:
import greenwood as gw
mg = gw.load_dataset("mgus2", backend="pandas")
start, stop, state, event = [], [], [], []
for i in range(len(mg)):
pt, ft = mg["ptime"][i], mg["futime"][i]
progressed, died = mg["pstat"][i] == 1, mg["death"][i] == 1
if progressed and pt < ft:
start += [0, pt]; stop += [pt, ft]; state += ["mgus", "pcm"]
event += ["pcm", "death" if died else None]
else:
start += [0]; stop += [ft]; state += ["mgus"]
event += ["death" if died else ("pcm" if progressed else None)]
rows = [(a, b, s, e) for a, b, s, e in zip(start, stop, state, event) if b > a]
start, stop, state, event = map(list, zip(*rows))
ms = gw.MultiState().fit(start, stop, state, event, states=("mgus", "pcm", "death"))
ms.predict([60, 120, 240], format="polars")
shape: (3, 4)| time | mgus | pcm | death |
|---|
| f64 | f64 | f64 | f64 |
| 60.0 | 0.645529 | 0.016006 | 0.338464 |
| 120.0 | 0.40446 | 0.012027 | 0.583513 |
| 240.0 | 0.176158 | 0.011489 | 0.812353 |
to_frame()
Return occupancy probabilities over time as a DataFrame.
Exports one row per distinct time and one column per state, where each state column contains its occupancy probability at that time.
Parameters
format: str | None = None-
Output format:
None (default), "pandas", "polars", or "pyarrow". When None, a backend is auto-detected (Polars, then Pandas, then PyArrow).
Returns
pandas.DataFrame, polars.DataFrame, or pyarrow.Table-
A tidy table with a
time column and one probability column per state.
Raises
ImportError-
If the requested (or, when auto-detecting, any) DataFrame library is not installed.
Examples
Fit a multi-state model on the mgus2 cohort and export the state-occupancy probabilities as a Polars frame:
import greenwood as gw
mg = gw.load_dataset("mgus2", backend="pandas")
start, stop, state, event = [], [], [], []
for i in range(len(mg)):
pt, ft = mg["ptime"][i], mg["futime"][i]
progressed, died = mg["pstat"][i] == 1, mg["death"][i] == 1
if progressed and pt < ft:
start += [0, pt]; stop += [pt, ft]; state += ["mgus", "pcm"]
event += ["pcm", "death" if died else None]
else:
start += [0]; stop += [ft]; state += ["mgus"]
event += ["death" if died else ("pcm" if progressed else None)]
rows = [(a, b, s, e) for a, b, s, e in zip(start, stop, state, event) if b > a]
start, stop, state, event = map(list, zip(*rows))
ms = gw.MultiState().fit(start, stop, state, event, states=("mgus", "pcm", "death"))
ms.to_frame(format="polars")
shape: (276, 4)| time | mgus | pcm | death |
|---|
| f64 | f64 | f64 | f64 |
| 1.0 | 0.969653 | 0.0 | 0.030347 |
| 2.0 | 0.947961 | 0.001446 | 0.050593 |
| 3.0 | 0.937114 | 0.001446 | 0.061439 |
| 4.0 | 0.924822 | 0.002169 | 0.073009 |
| 5.0 | 0.916868 | 0.002892 | 0.080239 |
| … | … | … | … |
| 356.0 | 0.08175 | 0.0 | 0.91825 |
| 373.0 | 0.0545 | 0.02725 | 0.91825 |
| 376.0 | 0.0545 | 0.02725 | 0.91825 |
| 394.0 | 0.0545 | 0.02725 | 0.91825 |
| 424.0 | 0.0 | 0.02725 | 0.97275 |
Request a different backend with format=:
ms.to_frame(format="pandas")
|
time |
mgus |
pcm |
death |
| 0 |
1.0 |
0.969653 |
0.000000 |
0.030347 |
| 1 |
2.0 |
0.947961 |
0.001446 |
0.050593 |
| 2 |
3.0 |
0.937114 |
0.001446 |
0.061439 |
| 3 |
4.0 |
0.924822 |
0.002169 |
0.073009 |
| 4 |
5.0 |
0.916868 |
0.002892 |
0.080239 |
| ... |
... |
... |
... |
... |
| 271 |
356.0 |
0.081750 |
0.000000 |
0.918250 |
| 272 |
373.0 |
0.054500 |
0.027250 |
0.918250 |
| 273 |
376.0 |
0.054500 |
0.027250 |
0.918250 |
| 274 |
394.0 |
0.054500 |
0.027250 |
0.918250 |
| 275 |
424.0 |
0.000000 |
0.027250 |
0.972750 |
276 rows × 4 columns