Analysis code and aggregate result tables for the Lumbar Fusion Digital Twin (LF-DT): an electronic health record–based framework that couples calibrated multi-outcome prediction with g-computation re-simulation of alternative lumbar-fusion operative decisions.
The clinical question LF-DT targets is counterfactual — how would this patient's outcomes differ under a different operative plan? — rather than prognostic. Existing spine risk models estimate baseline vulnerability under the observed plan; LF-DT re-estimates every outcome channel after changing the plan, and then asks whether those changes agree with an independently constructed causal benchmark.
Scope. This is a proof of concept. LF-DT produced calibrated baseline risks and cohort-level counterfactual contrasts that preserved the direction and ranking of effects estimated by a doubly robust method. It does not validate patient-specific treatment effects, and it is not validated for clinical deployment, patient counseling, or operative planning. Discrimination is modest (composite AUROC 0.687; reoperation C-index 0.580). External validation and prospective evaluation are required before any clinical use.
Included — the analysis code and the aggregate tables that back every number reported in the manuscript.
Not included, by design:
- patient records or any patient-level identifiers
- temporal-holdout assignments
- per-operation predicted risks or individual counterfactual contrasts
- trained model objects (
*.joblib) - the analytic feature matrix or any TriNetX-derived data
Aggregate cell counts below 11 are masked as <11 (see Small-cell suppression).
The data analyzed in this study were obtained from the TriNetX federated health research network, which aggregates de-identified electronic health record data contributed by participating healthcare organizations. TriNetX data are available to investigators at subscribing institutions under a license agreement and are not publicly downloadable. The authors are not permitted to redistribute patient-level records, and the identities of contributing healthcare organizations and their individual contributions are not disclosed, consistent with the network's data-governance terms.
Investigators at licensed institutions may request access through https://trinetx.com.
Running this pipeline requires an analytic matrix that this repository does not and cannot ship. Two things are needed beyond TriNetX access:
1. A cohort_features.parquet with the following contract:
| Element | Specification |
|---|---|
| Key | (patient_id, surgery_date) — one row per operation, not per patient |
| Cohort | Adults ≥18 undergoing primary lumbar fusion with ≥730 days follow-up or recorded death; trauma, malignancy, spinal infection, and deformity excluded |
| Predictors | 190 baseline and perioperative columns — demographics, 37 Elixhauser comorbidity indicators, laboratory values, vital signs, history variables, preoperative opioid/steroid/anticoagulant exposure, and surgical descriptors (approach, fusion extent, instrumentation, interbody, decompression, graft) |
| Survival columns | time_to_event_days, event_indicator |
| Outcome columns | The 30 channel columns enumerated in CHANNEL_SPEC in src/config.py |
src/config.py is the machine-readable specification: it defines the channel spec, the reliability-tiering
rule, the pre-specified levers, and every path the pipeline resolves. Read it first.
2. Two upstream helper modules that are not vendored here. src/config.py resolves paths into a
reoperation survival study and a binary-outcome validation toolkit, loaded through the opt-in
reop_src_on_path() and revision_utils_on_path() helpers. These supply the reused survival, calibration,
bootstrap, and propensity-score components described in the manuscript's Scope and reporting subsection.
Substituting equivalent implementations is straightforward — the seams are the PanelPredictor protocol and
the data contracts in src/interfaces.py — but they are not drop-in.
In short: this repository specifies the analysis completely and executes it given the matrix above. It is not a turnkey reproduction from raw TriNetX exports.
src/
├── config.py # channel spec, reliability rule, levers, paths — the single source of truth
├── interfaces.py # Twin / Trajectory / PanelPredictor / Lever / WhatIfResult contracts
├── twin_state.py # assemble the 190-feature pre-operative state vector
├── twin.py # Twin object + instantiate()
├── trajectory.py # reoperation trajectory across 30d / 90d / 1y / 2y
├── calibration.py # Platt scaling, calibration slope, CITL, ECE
├── panel/
│ ├── channels.py # channel construction + data-driven reliability tiering
│ ├── train_panel.py # per-channel logistic + histogram gradient boosting, 5-fold OOF, calibrated
│ └── predict.py # locked application of the refit panel to the temporal holdout
├── whatif.py # g-computation counterfactual engine
├── levers.py # the five pre-specified decision levers
├── validation/
│ ├── fidelity.py # temporal-holdout discrimination + calibration
│ ├── target_trial.py # doubly robust AIPW benchmark (emulated target trial)
│ ├── concordance.py # directional, Spearman, Kendall, Lin CCC agreement
│ ├── negative_controls.py # placebo levers + permuted-treatment AIPW
│ ├── face_validity.py # observed-vs-predicted deciles and subgroups
│ └── baseline_comparison.py # panel vs three-feature parsimonious model
├── figures.py # Figures 2–5
└── figure_schematic.py # Figure 1
outputs/tables/ # the 11 aggregate tables below
mask_small_cells.py # the suppression applied before publication
pip install -r requirements.txt
# 1. Instantiate a twin from the analytic matrix
python -m src.twin --example
# 2. Train the outcome panel (per-channel calibrated models)
python -m src.panel.train_panel
# 3. Re-simulate a decision on one operation
python -m src.whatif --example --lever interbody --to yes
# 4. Validation backbone — writes every table in outputs/tables/
python -m src.validation.fidelity # temporal-holdout performance
python -m src.validation.target_trial # AIPW benchmark across 15 contrasts
python -m src.validation.concordance # agreement statistics
python -m src.validation.negative_controls # placebo levers + permuted treatment
python -m src.validation.face_validity # deciles and subgroups
python -m src.validation.baseline_comparison # panel vs bedside model
# ...or run the whole validation suite
python -m src.validation| Table | Contents |
|---|---|
channel_reliability.csv |
All 30 outcome channels: development event counts, prevalence, and the data-driven reliability tier (17 modeled-and-validated, 13 descriptive) |
fidelity_binary.csv |
Temporal-holdout AUROC with bootstrap CI, ECE, Brier, calibration slope and intercept for the 16 modeled binary channels |
fidelity_survival.csv |
Reoperation Harrell C-index on the temporal holdout |
reoperation_strata.csv |
Predicted-risk quartiles vs observed 2-year reoperation rates |
panel_performance.csv |
Per-channel performance across all 30 channels with tier labels |
target_trial_agreement.csv |
The 15 lever × outcome contrasts: LF-DT g-computation Δ and AIPW Δ with intervals, overlap diagnostics, and max absolute SMD |
concordance.csv |
Spearman ρ, Kendall τ, Lin CCC with bootstrap CI, Pearson r across the 15 contrasts |
negative_control.csv |
Four negative controls: random off-model placebo, ethnicity and marital-status placebo levers, permuted-treatment AIPW |
face_validity_deciles.csv |
Observed vs mean-predicted composite risk by predicted-risk decile |
face_validity_subgroups.csv |
Observed vs predicted by age, sex, BMI, prior surgery, and fusion levels, with calibration gap |
baseline_comparison.csv |
190-feature panel vs the three-feature bedside model on the identical holdout |
Headline numbers, all from the temporal holdout of 7,680 operations (3,224 composite events):
- Composite 30-day complication: AUROC 0.687 (95% CI 0.675–0.699), ECE 0.012, calibration slope 0.889, calibration-in-the-large 0.015, Brier 0.218
- 16 binary channels: median AUROC 0.756 (range 0.601–0.902), median calibration slope 0.95
- Reoperation: C-index 0.580; observed 2-year rates rise 8.3% → 10.1% → 12.9% → 15.6% across predicted quartiles
- Counterfactual concordance vs AIPW across 15 contrasts: directional agreement 14/15, Spearman ρ = 0.743 (P = 0.0015), Kendall τ = 0.619, Lin CCC = 0.594 (95% CI 0.169–0.720)
- All four negative controls null
- Three-feature bedside comparator: AUROC 0.637 (95% CI 0.624–0.649)
Before publication, every table is passed through mask_small_cells.py, which masks
patient and event counts below 11 as <11, masks the derived prevalence on those rows as <0.001 (an exact
prevalence plus a known denominator would recover the suppressed count), and blanks the performance metrics on
those rows. The affected rows are all descriptive_grouped channels that the manuscript states were not
modeled-and-validated; a discrimination metric computed on a single-digit event count is uninterpretable, and
reporting one would contradict the reliability tiering. Four rows across two tables are affected. No masked
value appears in the manuscript. The script is idempotent, so the masking is verifiable by re-running it.
If you use this code, please cite both the software and the manuscript. See
CITATION.cff for machine-readable metadata.
MIT.