mixle.analysis.epidemiology module

Epidemiological cohort attribution: Cox PH + Aalen-Johansen (work-plan Sec.7-K, K9).

Given a cohort’s covariates (one column marking the exposure of interest), event times, and censoring, cohort_attribution() answers the two questions a health-risk sign-off actually needs: how much does the exposure multiply the hazard (the hazard ratio, HR), and what fraction of the hazard among the exposed is attributable to it (AF = (HR - 1) / HR). Both come with an honest interval, not a bare point estimate:

  • the hazard ratio and its Wald confidence interval fall straight out of the already-frozen mixle.inference.survival.cox_ph() partial-likelihood fit (CoxResult.se) – no new survival math here, this module is a thin, health-domain-shaped wrapper around it (K9’s non-goal list).

  • the attributable-fraction interval is a nonparametric bootstrap over cohort resamples, refitting cox_ph each draw; the bootstrap AF distribution is also handed back in provenance as an object satisfying the IC-1 DerivedQuantity shape (samples + credible_interval + prior_dominated), so it slots into the same decision-quantity plumbing (IC-8) the rest of the codebase uses for “never emit a number without its honesty flag.”

  • competing=True additionally runs mixle.inference.survival.aalen_johansen() on the raw (multi-cause) event labels for cause-specific cumulative incidence – the right absolute-risk curve when a competing cause (e.g. all-other-cause mortality) can remove a subject from the risk set before the outcome of interest can occur, which plain 1 - KM overstates.

  • latency > 0 encodes an exposure-to-effect lag as left-truncation: subjects contribute to the risk set (and to the bootstrap) only once they have survived past the latency period, via cox_ph’s counting-process start argument. This is the standard occupational/environmental-epi device for “the exposure cannot have caused an event that happened before the biological lag had elapsed.”

Every number in the returned CohortAttribution traces back to one cox_ph fit (plus, when competing, one aalen_johansen run) and one RNG’s full reproducible state – both recorded in provenance. provenance["rng_state"] is the bit-generator state (not just its construction seed): a seed alone cannot replay the draws of a generator the caller had already advanced before passing it in, and is unavailable at all for a bit generator built without a SeedSequence.

class CohortAttribution(hazard_ratio, hr_ci, attributable_fraction, af_ci, cif, provenance)[source]

Bases: object

Cox-PH hazard ratio + attributable fraction for one exposure, with intervals and provenance.

Parameters:
hazard_ratio

exp(coef[exposure_col]) from the fitted Cox model.

Type:

float

hr_ci

Wald (lo, hi) confidence interval on hazard_ratio (log scale, from CoxResult.se).

Type:

tuple[float, float]

attributable_fraction

fraction of hazard among the exposed attributable to exposure, (hazard_ratio - 1) / hazard_ratio.

Type:

float

af_ci

bootstrap (lo, hi) confidence interval on attributable_fraction.

Type:

tuple[float, float]

cif

cause-specific Aalen-Johansen cumulative incidence, {cause: array} (empty unless competing=True).

Type:

dict[int, numpy.ndarray]

provenance

fit diagnostics, the RNG’s full reproducible state (rng_state), and af_distribution – an _AFDistribution (IC-1 DerivedQuantity-shaped) when the bootstrap produced adequate evidence, or an _InsufficientBootstrapEvidence when it did not (see _MIN_BOOTSTRAP_DRAWS).

Type:

dict

The record is frozen and owns its cumulative-incidence curves (read-only copies), so the point estimates, the intervals and the provenance that documents how they were fitted cannot drift apart after the fact: rebinding hazard_ratio without hr_ci – or editing a cif array in place through the caller’s own handle – left a record whose parts described different fits.

Construction now also enforces the record’s own domain (MXR-080-1900). It previously enforced none: a hazard_ratio of -2.0 (impossible – it is exp(coef)), an attributable_fraction of 4.2 (impossible – (HR - 1) / HR approaches but never exceeds 1), and an inverted hr_ci of (9.0, 1.0) all constructed silently and were indistinguishable downstream from a real fit. The invariants, each one the producer already guarantees:

  • hazard_ratio finite and strictly positive.

  • hr_ci finite, strictly positive, and ordered lo <= hi.

  • attributable_fraction finite and <= 1. It is NOT bounded below by 0: a protective exposure has HR < 1 and a legitimately negative attributable fraction.

  • af_ci either two finite values ordered lo <= hi, or (nan, nan) – the documented marker for the insufficient-bootstrap-evidence case (see _InsufficientBootstrapEvidence), which is a real state the producer emits. A half-NaN interval is not.

Deliberately NOT checked: that attributable_fraction equals (hazard_ratio - 1) / hazard_ratio exactly, or that either interval brackets its point estimate. Both are true of cohort_attribution()’s own output, but a bootstrap percentile interval is not required to contain the point estimate in general, and pinning the AF identity would turn one documented parameterization into the only constructible one.

provenance is copied DEEPLY: the previous dict(...) detached only the top level, so nested entries – coef/se lists, the rng_state mapping that exists precisely so a fit can be replayed – stayed aliased to the caller’s own objects and could be rewritten after the attribution was reported.

cohort_attribution(covariates, time, event, *, exposure_col=0, competing=False, latency=0.0, n_boot=1000, rng=None)[source]

Attribute a cohort’s hazard to one exposure via Cox PH, with an Aalen-Johansen competing-risks CIF.

Parameters:

  • covariates (ndarray) – (n, p) covariates (finite); column exposure_col is the exposure of interest.

  • time (ndarray) – (n,) event/censoring times (finite, non-negative).

  • event (ndarray) – (n,) event indicator. Binary (0/1) when competing=False; an integer cause label (0 = censored, 1 = the outcome of interest, 2..K = competing causes) when competing=True. Either way, the Cox fit is cause-specific for cause 1: competing-cause events are censored at their time for the hazard-ratio fit (only the CIF, when requested, reports the competing causes’ own cumulative incidence). Values must be exact integers – fractional codes are rejected outright, not truncated (a truncated 1.5 would silently relabel a subject from “the event of interest” to “censored”).

  • exposure_col (int) – column of covariates treated as the exposure (must be a valid column index).

  • competing (bool) – if True, also run aalen_johansen for cause-specific cumulative incidence.

  • latency (float) – exposure-to-effect lag (left-truncation, via cox_ph’s start); 0 disables it. Must be non-negative – negative latency is rejected, not silently treated as 0.

  • n_boot (int) – cohort-resample bootstrap draws for the attributable-fraction interval. If fewer than _MIN_BOOTSTRAP_DRAWS converge (or too small a fraction of n_boot, see _MIN_BOOTSTRAP_FRACTION), af_ci is (nan, nan) and provenance["af_distribution"] is an _InsufficientBootstrapEvidence rather than a real interval.

  • rng – seed, numpy.random.Generator, or None.

Returns:

A CohortAttribution.

Raises:

ValueError – if the cohort is malformed – see _validate_cohort() for the exact checks (mismatched array shapes, non-finite values, negative time, an out-of-vocabulary or fractional event code, negative latency, or an out-of-bounds exposure_col).

Return type:

CohortAttribution