Section 36.3: Survival analysis and hazard models

"A regression tells you when the machine will die. A hazard model admits it might outlive us both, and quantifies exactly how surprised it would be if it did."
A Statistically Humble AI Agent

The Big Picture

Section 36.2 framed remaining useful life (RUL) as a point prediction: a number of cycles, hours, or kilometres until failure. That framing quietly assumes every machine in your training set actually failed while you were watching. Real fleets do not cooperate. A pump is swapped out during scheduled maintenance before it breaks. A wind turbine gearbox is still spinning happily when the study ends. A vehicle is sold and leaves your telemetry. In each case you know the unit survived at least so long, but not when it would have failed. That is censoring, and a naive regression that treats "last seen at 40,000 hours, healthy" as "failed at 40,000 hours" will systematically underestimate life. Survival analysis is the branch of statistics built precisely to learn from partially observed lifetimes, and its central object, the hazard function, turns out to be the most natural language for prognostics: it answers "given that this bearing has lasted this long, how likely is it to fail in the next hour?" This section gives you the survival, hazard, and cumulative-hazard functions; the Kaplan-Meier estimator; the Cox proportional-hazards model driven by sensor features; and the deep survival models that connect all of this to the neural encoders in the rest of Part VIII.

This section assumes you are comfortable with probability distributions, conditional probability, and maximum likelihood at the level of Chapter 4, and that you have the failure-mode and RUL vocabulary from Sections 36.1 and 36.2. The sensor features we feed these models come from the feature-engineering pipeline of Chapter 8, and the leakage discipline we insist on throughout comes from Chapter 5.

Lifetimes, censoring, and three equivalent functions

Let $T \ge 0$ be the random failure time of a unit. Survival analysis describes $T$ through three interchangeable functions, each answering a different operational question. The survival function is the probability the unit is still alive past time $t$:

$$S(t) = P(T > t) = 1 - F(t),$$

where $F$ is the cumulative distribution of $T$. The hazard function is the instantaneous failure rate given survival so far:

$$h(t) = \lim_{\Delta t \to 0} \frac{P(t \le T < t + \Delta t \mid T \ge t)}{\Delta t} = \frac{f(t)}{S(t)},$$

with $f$ the density of $T$. The hazard is the star of the show because it is conditional: it is exactly the question a maintenance planner asks about a machine that has already survived to today. Integrating the hazard gives the cumulative hazard $H(t) = \int_0^t h(u)\,du$, and the three functions are locked together by the identity

$$S(t) = \exp\!\big(-H(t)\big) = \exp\!\Big(-\!\int_0^t h(u)\,du\Big).$$

Know one, know all three. A rising hazard means an aging, wear-out failure mode (the bathtub curve's right wall from Section 36.1); a flat hazard means memoryless, random failure; a falling hazard means infant mortality is being weeded out.

Censoring is what makes this more than a distribution-fitting exercise. For unit $i$ we observe a time $t_i$ and an event indicator $\delta_i$: $\delta_i = 1$ means we saw the actual failure at $t_i$, while $\delta_i = 0$ means the unit was right-censored, still healthy when we last observed it. The likelihood contribution respects this: a failure contributes its density $f(t_i)$, a censored observation contributes only its survival $S(t_i)$, because all we learned is that $T_i > t_i$. Writing $f = h \cdot S$, the log-likelihood becomes the clean and famous form

$$\ell = \sum_i \big[\, \delta_i \log h(t_i) - H(t_i) \,\big].$$

Key Insight

The reason you cannot just drop censored units and run a regression is that censoring is almost never random with respect to lifetime. Units are removed because they were healthy, or replaced on a preventive schedule that correlates with the very wear you are trying to model. Discarding them biases your estimated life downward; treating their last-seen time as a failure biases it downward even harder. The survival likelihood is the only estimator here that stays unbiased under non-informative censoring, and that single property is why prognostics teams reach for hazard models instead of ordinary regression.

Non-parametric baselines: Kaplan-Meier

Before fitting any sensor-driven model, you want a picture of the raw fleet survival curve with no assumptions about its shape. The Kaplan-Meier estimator gives it. Order the distinct failure times $t_{(1)} < t_{(2)} < \dots$; at each, let $d_j$ be the number of failures and $n_j$ the number of units still at risk (alive and uncensored just before $t_{(j)}$). The estimator is the running product of conditional survivals:

$$\hat{S}(t) = \prod_{t_{(j)} \le t} \left(1 - \frac{d_j}{n_j}\right).$$

Censored units stay in $n_j$ up until the moment they leave, contributing their information without ever forcing a downward step. The result is a staircase curve with a confidence band (via Greenwood's formula) that you can stratify by machine type, operating regime, or site. Kaplan-Meier answers "what fraction of this population lasts past 20,000 hours?" and it is your honest baseline: any fancier model that cannot beat a well-stratified KM curve on held-out units is not earning its complexity.

Cox proportional hazards: sensors enter the model

Kaplan-Meier describes a population; prognostics needs to condition on this unit's sensor state. The Cox proportional-hazards model is the workhorse. It factors the hazard into a shared, unspecified baseline $h_0(t)$ and a covariate-dependent multiplier:

$$h(t \mid \mathbf{x}) = h_0(t)\, \exp\!\big(\boldsymbol{\beta}^\top \mathbf{x}\big),$$

where $\mathbf{x}$ is the sensor feature vector for the unit (vibration RMS, bearing temperature, oil-debris count, spectral kurtosis from Chapter 7) and $\boldsymbol{\beta}$ are learned coefficients. The name "proportional hazards" is literal: two units' hazards stay in constant ratio $\exp(\boldsymbol{\beta}^\top(\mathbf{x}_i - \mathbf{x}_j))$ for all time, regardless of the baseline's shape. The magic is that $\boldsymbol{\beta}$ is estimated by maximizing the partial likelihood, which never requires you to specify $h_0(t)$ at all:

$$L(\boldsymbol{\beta}) = \prod_{i:\,\delta_i = 1} \frac{\exp(\boldsymbol{\beta}^\top \mathbf{x}_i)}{\sum_{j \in R(t_i)} \exp(\boldsymbol{\beta}^\top \mathbf{x}_j)},$$

where $R(t_i)$ is the risk set of units still alive just before failure $i$. Each factor asks: among everyone who could have failed at this instant, how much of the total hazard belonged to the one who actually did? A coefficient $\beta_k = 0.7$ means a one-unit rise in feature $k$ multiplies the failure rate by $e^{0.7} \approx 2$, an interpretable hazard ratio that reliability engineers can defend in a design review. When sensor features drift over a unit's life, the time-varying Cox variant lets $\mathbf{x}(t)$ update, turning the model into a streaming hazard estimator that pairs naturally with the online feature extraction of Chapter 60.

import numpy as np
from lifelines import CoxPHFitter, KaplanMeierFitter

# One row per unit: last-observed time, event flag, and two sensor features.
# event=1 means a real failure; event=0 means right-censored (still healthy).
rng = np.random.default_rng(0)
n = 300
vib = rng.gamma(2.0, 1.0, n)          # vibration RMS feature
temp = rng.normal(0.0, 1.0, n)         # standardized bearing temperature
# Higher vibration -> higher hazard -> shorter life (Weibull-ish true model).
scale = 50.0 * np.exp(-0.6 * vib - 0.3 * temp)
true_life = rng.weibull(1.8, n) * scale
censor_at = rng.uniform(10, 80, n)     # study ends / preventive swap
time = np.minimum(true_life, censor_at)
event = (true_life <= censor_at).astype(int)

# Honest fleet baseline first (no covariates).
km = KaplanMeierFitter().fit(time, event)
print("Median survival (KM):", km.median_survival_time_)

# Sensor-driven hazard model.
import pandas as pd
df = pd.DataFrame({"time": time, "event": event, "vib": vib, "temp": temp})
cox = CoxPHFitter().fit(df, duration_col="time", event_col="event")
print(cox.summary[["coef", "exp(coef)", "p"]])   # hazard ratios per feature

Fitting a Kaplan-Meier fleet baseline and a Cox proportional-hazards model on censored lifetimes with two sensor features. The event column is what separates survival analysis from regression: censored rows (event=0) contribute their survival probability, not a fabricated failure time. The printed exp(coef) values are directly readable hazard ratios.

The code above computes the KM baseline and Cox fit referenced throughout this subsection. Notice that a censored unit never lies about when it failed; it simply asserts it lasted at least as long as observed, and the likelihood does the rest.

Library Shortcut

Writing the Cox partial-likelihood optimizer, the risk-set bookkeeping, tie handling (Efron's approximation), the Breslow baseline-hazard estimator, and Greenwood confidence bands from scratch is roughly 250 to 400 lines of careful, easy-to-get-wrong NumPy. lifelines collapses all of it to two .fit() calls, about 6 lines including the dataframe, and throws in the concordance index, proportional-hazards diagnostic tests, and partial-effect plots for free. Reach for the from-scratch version only when you need a custom likelihood a library will not express; otherwise the library is both shorter and more correct.

Deep survival models: hazards meet sensor encoders

Cox is linear in $\mathbf{x}$ and assumes proportional hazards, two assumptions that raw vibration spectra and multi-sensor histories routinely violate. Deep survival models keep the survival likelihood but replace the linear predictor with a neural encoder. DeepSurv swaps $\boldsymbol{\beta}^\top \mathbf{x}$ for a network $g_\theta(\mathbf{x})$ trained on the Cox partial likelihood, so it stays proportional but becomes nonlinear. To escape proportionality, DeepHit discretizes time into bins and directly predicts the probability mass of failure in each bin with a softmax, handling competing risks (multiple failure modes from Section 36.1) in one head. Continuous-time neural models such as neural ODEs for the cumulative hazard round out the family. Crucially, the encoder $g_\theta$ can be any of the sensor-stream backbones from Part IV: a temporal convolutional network from Chapter 14, or a transformer from Chapter 15, ingesting the full degradation window rather than a hand-picked feature snapshot. Section 36.4 builds exactly these sequence encoders; this section supplies the loss they should be trained against when lifetimes are censored.

Research Frontier

The current frontier fuses survival objectives with the foundation-model embeddings of Section 36.5. Instead of training a hazard head on raw sensors, teams freeze a pretrained wearable or time-series encoder (MOMENT, Moirai, or a domain-specific vibration model) and fit a lightweight Cox or DeepHit head on its embeddings, getting calibrated hazards from a few dozen labeled failures. Open problems: proving proportional-hazards or discrete-time calibration holds on embedding spaces the survival model never saw during pretraining, and extending conformal survival prediction (distribution-free survival bands, tied to Chapter 18) to censored, non-exchangeable fleet data. Evaluation is unsettled too: the time-dependent concordance index and the integrated Brier score are the reported standards, but both are sensitive to the censoring distribution.

Practical Example: Wind-Turbine Gearbox Fleet

An operator runs 400 turbines across three sites. Over four years, 61 gearboxes failed with logged failure dates; the other 339 are still turning or were preventively rebuilt on a maintenance schedule. Treating the 339 as "failed on their last SCADA timestamp" and running a gradient-boosted RUL regressor predicted a fleet median life of 31,000 hours, and the maintenance budget was sized to it. A Cox model on oil-temperature trend, gearbox vibration RMS, and power-curve residual, fit on the same data with the 339 correctly censored, put the median at 44,000 hours with a hazard ratio of 2.3 per standard deviation of vibration RMS. The regression had been pulling healthy units into its loss as if they had died, shortening every prediction. Re-planning against the survival curve deferred 90 gearbox rebuilds by an average of nine months, and the vibration hazard ratio told the reliability team which units to prioritize, not just when the fleet as a whole would wear out.

Exercise

Take the synthetic dataset from the code block. (a) Refit the Cox model but first corrupt the labels by setting every event = 1, that is, pretend nothing was censored. Compare the estimated median survival and the vibration hazard ratio to the correctly censored fit, and quantify the bias. (b) Split the units into "high vibration" and "low vibration" by the median of vib, fit a separate Kaplan-Meier curve to each, and check visually whether the two curves stay proportional (a rough proportional-hazards sanity check). (c) Add a third feature that is pure noise, refit Cox, and confirm its coefficient's p-value is not significant. What does that tell you about using p-values as a feature-selection guardrail on sensor data?

Self-Check

You have $S(t)$ for a unit. Write down its cumulative hazard $H(t)$ and its hazard $h(t)$ in terms of $S$. Which of the three is the quantity a maintenance planner most directly acts on, and why?
A colleague drops all right-censored units and fits an ordinary RUL regression on the survivors "to keep it simple." In which direction is the resulting life estimate biased, and what property of the censoring would have to hold for the shortcut to be even approximately valid?
Cox assumes proportional hazards. Give one concrete sensor scenario in a rotating machine where that assumption fails, and name the model from this section you would switch to instead.

What's Next

In Section 36.4, we build the sequence encoders that feed these hazard heads: recurrent, convolutional, and attention-based models that read a unit's full degradation trajectory rather than a snapshot of features. The survival likelihood you just learned becomes their training objective whenever the fleet's lifetimes are censored, which, on any real fleet, is always.