Applying homology-based morphometry to longitudinal spinal cord MRI will show that topological invariants (persistent homology Betti numbers) predict MS lesion merging and splitting events over time with a higher area under the ROC curve than conventional semantic segmentation volume-change models.

In patients with multiple sclerosis undergoing serial spinal cord MRI (≥3 timepoints, ≥6-month total follow-up), a persistent-homology feature vector computed from binarized/probabilistic lesion segmentation masks — specifically Betti-0 (connected components), Betti-1 (loops/handles), and their persistence-diagram summary statistics (birth-death lifetimes, total persistence, Wasserstein/bottleneck distances between consecutive timepoints) — will predict the occurrence of a lesion merge or split event at the subsequent timepoint (binary label, defined via voxelwise component-correspondence tracking) with a higher cross-validated AUROC than a baseline model using only volumetric/semantic-segmentation change features (total lesion volume, per-lesion volume delta, count of new/resolved lesions). The claim is falsifiable: it specifies the comparator model, the outcome definition, the feature class, and requires a pre-registered statistical test (DeLong test or bootstrap CI comparison) on held-out patients.

Minimum viable test (Phase A, n≈60–80 patients, retrospective):

1. Acquire longitudinal spinal cord MRI + existing lesion segmentations from one public/consortium dataset.
2. Compute ground-truth merge/split labels via component-tracking across timepoints (connected-component relabeling + IoU-based correspondence, threshold IoU<0.3 triggers manual adjudication).
3. Extract (a) volumetric features (11-dim: total volume, lesion count, per-lesion volume vector stats, volume-change rate) and (b) persistent homology features (cubical complex filtration on signed-distance-transformed masks; Betti-0/1 curves sampled at 50 filtration steps; persistence image vectorization, 20×20 grid → 400-dim, PCA-reduced to 20 components).
4. Train/evaluate both feature sets with identical classifier family (logistic regression + gradient-boosted trees) under nested 5-fold patient-level cross-validation (no patient leakage across folds).
5. Compare AUROC via DeLong test; report 95% CI via 1000-iteration bootstrap.
6. Pre-register hypothesis, feature definitions, and primary endpoint before unblinding test folds.

Full validation (Phase B): repeat on independent external cohort (different scanner/site), n≥150, with prospective-style held-out lock-box evaluation.

# Phase A pipeline (pseudocode)

for patient in cohort:
    scans = load_longitudinal_MRI(patient)  # T2/PSIR, ≥3 timepoints
    scans = harmonize(scans, N4_bias_correct, resample_iso_0.8mm, SCT_crop_cord)

    masks = []
    for t in scans:
        mask_t = existing_expert_mask(t) or nnUNet_infer(t)
        qc_check(mask_t, dice_threshold=0.7)
        masks.append(mask_t)

    components = track_components(masks)   # IoU-based correspondence across t
    events = label_merge_split(components) # binary per lesion-pair per transition

    vol_features = []
    topo_features = []
    for t, mask_t in zip(scans, masks):
        vol_features.append(extract_volumetric(mask_t))
        sdt = signed_distance_transform(mask_t)
        pd = cubical_persistence(sdt, lib="GUDHI")   # persistence diagram
        pi = persistence_image(pd, grid=20x20)
        betti_curve = betti_numbers_over_filtration(pd, steps=50)
        topo_features.append(concat(pi.flatten(), betti_curve, pd_summary_stats(pd)))

    store(patient, vol_features, topo_features, events)

# Modeling
X_vol, X_topo, y = assemble_dataset(all_patients, horizon="next_timepoint")
split = PatientLevelStratifiedKFold(n_splits=5, stratify_on=y)

for model_class in [LogisticRegression, GradientBoostedTrees]:
    auc_vol  = nested_cv_eval(model_class, X_vol, y, split)
    auc_topo = nested_cv_eval(model_class, X_topo, y, split)
    auc_comb = nested_cv_eval(model_class, concat(X_vol, X_topo), y, split)

delong_test(auc_topo_preds, auc_vol_preds, y)
bootstrap_CI(auc_topo - auc_vol, n_iter=1000)
likelihood_ratio_test(auc_comb_model, auc_vol_model)  # incremental value of topology

# Stability/sensitivity
for perturbation in [dilate(1vx), erode(1vx), remove_small_components(thresh)]:
    masks_perturbed = apply(perturbation, masks)
    topo_features_perturbed = recompute_topology(masks_perturbed)
    cv_stability = coefficient_of_variation(topo_features, topo_features_perturbed)