mixle.inference.survival module

Time-to-event (survival) estimators and hazard regression.

Survival analysis models the time until an event in the presence of right-censoring (subjects who leave the study before the event – their time is a lower bound, not a missing value). The toolkit here covers the estimators and the regression layer:

  • kaplan_meier() / nelson_aalen() – nonparametric survival and cumulative-hazard curves with Greenwood / Poisson variance and confidence bands.

  • cox_ph() – the Cox proportional-hazards regression: how covariates multiply the hazard, estimated from the partial likelihood (Efron or Breslow tie handling), with stratification and time-varying covariates (counting-process start, stop input), Breslow baseline hazard, and the concordance index.

  • discrete_time_hazard() (+ to_person_period()) – discrete-time hazard models fit as a binary GLM on the person-period array (logit or complementary-log-log), supporting offsets and fixed effects through the design matrix.

  • aalen_johansen() – competing-risks cumulative incidence functions (cause-specific).

  • aalen_additive() – Aalen’s additive-hazards regression (cumulative covariate effects).

  • frailty_cox() – shared gamma-frailty Cox for clustered survival (random effect per group), fit by EM.

Event indicators are 1 for an observed event and 0 for right-censoring (for competing risks, an integer cause label with 0 = censored).

kaplan_meier(time, event=None, *, ci_level=0.95)[source]

Kaplan–Meier product-limit estimate of the survival function S(t).

Parameters:
  • time (ndarray) – (n,) observed times (event or censoring).

  • event (ndarray | None) – (n,) 1 = event, 0 = right-censored (defaults to all events).

  • ci_level (float) – confidence level for the log–log survival band.

Returns:

{'time', 'survival', 'se', 'ci_low', 'ci_high', 'at_risk', 'n_events', 'median'}.

Return type:

dict[str, ndarray]

nelson_aalen(time, event=None)[source]

Nelson–Aalen estimate of the cumulative hazard H(t) = sum d_i / Y_i.

Returns:

{'time', 'cumhaz', 'se'} with the Poisson-type standard error of the cumulative hazard.

Parameters:
Return type:

dict[str, ndarray]

class CoxResult(coef, se, cov, loglik, baseline_time, baseline_cumhaz, concordance, n_iter)[source]

Bases: object

Fitted Cox proportional-hazards model.

Parameters:
coef

(p,) log-hazard-ratio coefficients.

Type:

numpy.ndarray

se

(p,) standard errors (inverse observed information).

Type:

numpy.ndarray

cov

(p, p) covariance.

Type:

numpy.ndarray

loglik

maximised partial log-likelihood.

Type:

float

baseline_time / baseline_cumhaz

Breslow baseline cumulative hazard.

concordance

Harrell’s C-index.

Type:

float

n_iter

Newton iterations.

Type:

int

coef: ndarray
se: ndarray
cov: ndarray
loglik: float
baseline_time: ndarray
baseline_cumhaz: ndarray
concordance: float
n_iter: int
hazard_ratios()[source]
Return type:

ndarray

z_values()[source]
Return type:

ndarray

p_values()[source]
Return type:

ndarray

cox_ph(x, time, event, *, start=None, strata=None, ties='efron', max_iter=50, tol=1e-9)[source]

Cox proportional-hazards regression by Newton–Raphson on the partial likelihood.

The hazard is h(t | x) = h0(t) exp(x' beta); only the ordering of event times enters, so the baseline h0 is left unspecified (semi-parametric). Time-varying covariates are supported through the counting-process form: pass start so each row is an at-risk interval (start, stop] (a subject contributes several rows), and the risk set at an event time is every interval covering it.

Parameters:
  • x (ndarray) – (n, p) covariates (no intercept – it is absorbed into the baseline).

  • time (ndarray) – (n,) event/censoring times (the interval stop times).

  • event (ndarray) – (n,) 1 = event, 0 = censored.

  • start (ndarray | None) – optional (n,) interval start times for time-varying covariates / left truncation.

  • strata (ndarray | None) – optional (n,) labels; each stratum gets its own baseline hazard (coefficients shared).

  • ties (str) – "efron" (default, more accurate) or "breslow" tie handling.

  • max_iter (int) – Newton controls.

  • tol (float) – Newton controls.

Returns:

A CoxResult.

Return type:

CoxResult

to_person_period(time, event, covariates=None)[source]

Expand right-censored durations into a person-period (long) array for discrete-time models.

Each subject contributes one row per discrete period they were at risk; the binary outcome is 1 in the period the event occurred and 0 otherwise. Integer time is the number of periods observed.

Returns:

{'period', 'outcome', 'subject', 'covariates'} (covariates repeated per period if given).

Parameters:
Return type:

dict[str, ndarray]

discrete_time_hazard(x, outcome, *, link='cloglog', offset=None)[source]

Discrete-time hazard model: a binary GLM on the person-period array.

Fit on the long-format data from to_person_period() (the design x typically holds period indicators / a time trend plus covariates). cloglog gives the grouped-proportional-hazards (interval-censored Cox) interpretation; logit gives the proportional-odds hazard.

Returns:

a mixle.inference.glm.GLMResult (binomial family with the chosen link).

Parameters:
aalen_johansen(time, event, *, causes=None)[source]

Aalen–Johansen cumulative incidence functions for competing risks.

With several mutually exclusive event types, the cause-specific CIF F_k(t) is the probability of failing from cause k by time t accounting for the competing causes (it is not 1 - KM on the cause, which overstates incidence).

Parameters:
  • time (ndarray) – (n,) event/censoring times.

  • event (ndarray) – (n,) integer cause label, 0 = censored, 1..K = causes.

  • causes (ndarray | None) – optional explicit list of cause labels; inferred from event if None.

Returns:

{'time', 'cif': {cause: array}, 'overall_survival'}.

Return type:

dict

aalen_additive(x, time, event, *, intercept=True)[source]

Aalen’s additive-hazards regression: cumulative regression functions B(t).

Models h(t | x) = b0(t) + sum_j x_j b_j(t) with time-varying additive effects. At each event time the increment dB is the least-squares solution over the risk set; the cumulative B(t) (returned) has interpretable slopes – a rising B_j means covariate j adds hazard.

Returns:

{'time', 'cum_coef'} where cum_coef is (n_event_times, p[+1]) cumulative coefficients (the first column is the baseline when intercept is True).

Parameters:
Return type:

dict

class FrailtyCoxResult(coef, se, theta, frailties, groups, n_iter=0)[source]

Bases: object

Shared gamma-frailty Cox result.

Parameters:
coef / se

fixed-effect log-hazard-ratios and standard errors.

theta

estimated frailty variance (0 means no clustering signal).

Type:

float

frailties

posterior mean random effect per group.

Type:

numpy.ndarray

groups

group labels aligned to frailties.

Type:

numpy.ndarray

n_iter

EM iterations.

Type:

int

coef: ndarray
se: ndarray
theta: float
frailties: ndarray
groups: ndarray
n_iter: int = 0
frailty_cox(x, time, event, groups, *, max_iter=50, tol=1e-5, ties='breslow')[source]

Shared gamma-frailty Cox model for clustered survival, by EM.

Subjects in the same group share an unobserved frailty w_g ~ Gamma(1/theta, 1/theta) (mean 1, variance theta) that multiplies the hazard, capturing within-group correlation. The E-step takes the posterior-mean frailties; the M-step refits Cox with log w_g as an offset and updates theta. theta -> 0 indicates no detectable clustering.

Returns:

A FrailtyCoxResult.

Parameters:
Return type:

FrailtyCoxResult