brier_score()

IPCW (Graf) Brier score of predicted survival probabilities at specified times.

Usage

Source

brier_score(
    surv,
    survival_prob,
    times,
)

Measures calibration and accuracy of predicted survival probabilities at fixed time points using inverse-probability-of-censoring-weighted (IPCW) averaging. The Brier score is the mean squared difference between predicted and observed outcomes, weighted so censored subjects contribute honestly without bias.

Interpretation:

  • Ranges from 0 (perfect predictions) to 1 (worst possible).
  • Lower is better. A Brier score of 0.25 means, on average, predictions are off by \(\pm 0.5\) in terms of squared deviation.
  • Null model baseline: A model predicting 50% survival at every time has Brier score \(\approx 0.25\). Compare your model to this baseline to assess practical improvement.
  • Score typically increases with time (harder to predict farther into the future).

Practical use: After fitting a survival model (Cox, parametric, flexible), evaluate calibration at important clinical horizons (e.g., 1-year, 5-year survival). Compute Brier scores at multiple times, then use integrated_brier_score() for a single summary.

Parameters

surv: Surv

A right-censored Surv response (time-to-event data).

survival_prob: Any

Predicted survival probabilities, shape (n_subjects, n_times). Each entry is a predicted probability of surviving beyond the corresponding time. Must be between 0 and 1. Example: columns from CoxPH.predict(type="survival", times=...) (excluding the time column), transposed so rows are times and columns are subjects.

times: Any
Evaluation times where Brier scores are computed. 1-D array-like. Must have length equal to the second dimension of survival_prob.

Returns

ndarray
Brier score at each time, shape (len(times),). Lower is better.

Details

Graf (IPCW) Brier score: The unbiased Brier score under censoring is

\[ BS(t) = E\left[(S(t) - \hat{S}(t))^2 \cdot \text{weights}\right], \]

where:

  • \(S(t)\) is the true survival status at time \(t\) (1 if alive, 0 if dead)
  • \(\hat{S}(t)\) is the predicted survival probability
  • \(\text{weights}\) are inverse-probability-of-censoring: inverse of the censoring survival function

Mathematically:

\[ BS(t) = \frac{1}{n} \sum_{\text{dead at } t} \frac{(\hat{S}_i(t))^2}{G(t_i^-)} + \frac{1}{n} \sum_{\text{alive at } t} \frac{(1 - \hat{S}_i(t))^2}{G(t)} \]

where \(G(u)\) is the Kaplan-Meier estimate of the censoring distribution (probability of not being censored).

Advantages over MSE: IPCW weighting makes the Brier score unbiased under censoring, unlike naive MSE which would be biased (censored subjects look artificially “correct”).

Time dependence: Brier scores typically increase with time in chronic-disease settings (longer prediction horizons are harder); they may decrease in acute-illness settings (events cluster early, predictions stabilize).

Examples

Fit a Cox model on the lung dataset and evaluate its calibration (how accurately does it predict survival?) at a few horizons.

import greenwood as gw
import numpy as np

lung = gw.load_dataset("lung", backend="polars")
y = gw.Surv.right(lung["time"], event=(lung["status"] == 2))
cox = gw.CoxPH().fit(y, lung[["age", "sex"]])

times = [180, 365, 540]
surv_pred = cox.predict(lung[["age", "sex"]], type="survival", times=times, format="pandas")
probs = surv_pred.iloc[:, 1:].to_numpy().T
brier = gw.brier_score(y, probs, times)
brier
array([0.19160509, 0.23614214, 0.18644928])

Brier scores at three time points. Scores typically increase over time. Compare to a null model (all subjects at 50% survival) to assess improvement:

null_probs = np.full_like(probs, 0.5)
null_brier = gw.brier_score(y, null_probs, times)
print(f"Null model Brier: {null_brier}")
print(f"Cox model Brier: {brier}")
print(f"Improvement: {null_brier - brier}")
Null model Brier: [0.25 0.25 0.25]
Cox model Brier: [0.19160509 0.23614214 0.18644928]
Improvement: [0.05839491 0.01385786 0.06355072]

Summarize Brier scores across times into a single summary via time-averaged Brier score:

ibs = gw.integrated_brier_score(y, probs, times)
print(f"Integrated Brier Score: {ibs:.3f}")
Integrated Brier Score: 0.213