DeLong Test in Python: Comparing Two Models' AUC Honestly

If you search for a DeLong test in Python, you mostly find GitHub repositories and package pages. They give you a function but don't say when to trust its p-value. As an AI researcher at the University of Alberta, much of my work is comparing clinical prediction models, and the evaluation side of that work leans on two tools: DeLong testing and patient-level bootstrap confidence intervals. This essay gives you both as runnable code, shows how to check that the DeLong implementation is right, and walks through a simulation of the mistake that makes a weak result look significant.

Everything below runs on synthetic data. None of it comes from a clinical dataset, and none of the numbers are research results.

DeLong test or bootstrap: which question are you asking?

Both tools answer some version of "is model A better than model B?" They make different assumptions, and the assumptions matter more than the code.

My rule of thumb: if every patient contributes exactly one row, DeLong is fast and fine. If any patient can contribute more than one row, which is the norm when a patient has several hospital episodes, report the patient-level bootstrap. Use DeLong as a quick check at most, and expect it to be optimistic.

The paired DeLong test in Python

The DeLong method estimates the variance of an AUC without resampling. Each positive case gets a value describing how it ranks against the negatives, and each negative case gets the mirror-image value. The spread of those per-case values gives the variance. Because both models are scored on the same cases, the method also estimates the covariance between the two AUCs, which is what makes it a paired test. Without that covariance term, two correlated models look less distinguishable than they are.

The core fits in about thirty lines. It needs NumPy and SciPy, and it handles tied scores because rankdata assigns ties their average rank.

# auc_compare.py
import numpy as np
from scipy import stats


def delong_paired(y, score_a, score_b):
    """Paired DeLong test for two models scored on the same cases.

    Returns (auc_a, auc_b, z, two-sided p-value).
    """
    y = np.asarray(y).astype(bool)
    pos, neg = y, ~y
    m, n = pos.sum(), neg.sum()
    aucs, v10, v01 = [], [], []
    # rankdata gives tied scores their average rank (midranks), which DeLong needs.
    for s in (np.asarray(score_a), np.asarray(score_b)):
        r_all = stats.rankdata(s)
        r_pos = stats.rankdata(s[pos])
        r_neg = stats.rankdata(s[neg])
        aucs.append((r_all[pos].sum() - m * (m + 1) / 2) / (m * n))
        v10.append((r_all[pos] - r_pos) / n)          # one value per positive case
        v01.append(1 - (r_all[neg] - r_neg) / m)      # one value per negative case
    s10 = np.cov(np.vstack(v10))
    s01 = np.cov(np.vstack(v01))
    cov = s10 / m + s01 / n
    diff = aucs[0] - aucs[1]
    var = cov[0, 0] + cov[1, 1] - 2 * cov[0, 1]
    z = diff / np.sqrt(var)
    p = 2 * stats.norm.sf(abs(z))
    return aucs[0], aucs[1], z, p


def fast_auc(y, s):
    """AUC from ranks: same value as roc_auc_score, without the input checks."""
    r = stats.rankdata(s)
    m = y.sum()
    n = len(y) - m
    return (r[y == 1].sum() - m * (m + 1) / 2) / (m * n)


def patient_bootstrap(y, score_a, score_b, patient_id, n_boot=2000, seed=0):
    """95% CI for AUC(a) - AUC(b), resampling patients rather than rows."""
    rng = np.random.default_rng(seed)
    patients = np.unique(patient_id)
    rows_of = {p: np.flatnonzero(patient_id == p) for p in patients}
    diffs = []
    for _ in range(n_boot):
        picked = rng.choice(patients, size=len(patients), replace=True)
        idx = np.concatenate([rows_of[p] for p in picked])
        if y[idx].min() == y[idx].max():
            continue  # a resample with one class has no AUC
        diffs.append(fast_auc(y[idx], score_a[idx]) - fast_auc(y[idx], score_b[idx]))
    return np.percentile(diffs, [2.5, 97.5])

If you would rather not maintain your own, the roc_comparison repository is an existing DeLong implementation that cites Sun and Xu's 2014 fast algorithm. Whichever you use, run the next check first.

Check the implementation before you trust it

A significance test that is quietly wrong is worse than none, because people will believe it. Two cheap checks catch most bugs. First, the AUCs must match sklearn.metrics.roc_auc_score exactly. Second, when rows really are independent, the DeLong interval should agree with a bootstrap interval. Here the data has one row per patient, and the scores are rounded on purpose so that ties are exercised:

# validate.py
import numpy as np
from auc_compare import delong_paired, patient_bootstrap

# One row per patient, so every row really is independent.
rng = np.random.default_rng(3)
n = 600
risk = rng.normal(0, 1.5, n)
y = (risk + rng.normal(0, 1, n) > 1.5).astype(int)
a = np.round(risk + rng.normal(0, 1.0, n), 1)   # rounding creates ties on purpose
b = np.round(risk + rng.normal(0, 1.4, n), 1)
pid = np.arange(n)

auc_a, auc_b, z, p = delong_paired(y, a, b)
se = (auc_a - auc_b) / z
lo, hi = patient_bootstrap(y, a, b, pid, n_boot=4000)
print(f"DeLong 95% CI:    [{auc_a - auc_b - 1.96 * se:.3f}, {auc_a - auc_b + 1.96 * se:.3f}]")
print(f"Bootstrap 95% CI: [{lo:.3f}, {hi:.3f}]")

DeLong 95% CI:    [0.008, 0.100]
Bootstrap 95% CI: [0.007, 0.100]

The two intervals agree to within a thousandth. If yours don't, the bug is usually in the tie handling or in the sign of one of the per-case terms.

A patient-level bootstrap for AUC confidence intervals

Now the realistic case. In clinical data, one patient can have several hospital episodes, each with its own row and label. The synthetic cohort below has 400 patients with one to eight episodes each. Episodes of the same patient share a patient-level effect, and each model makes a patient-level error too: a model that misjudges a patient tends to misjudge every episode of theirs.

# demo.py
import numpy as np
from sklearn.metrics import roc_auc_score
from auc_compare import delong_paired, patient_bootstrap

rng = np.random.default_rng(42)

# Synthetic cohort: 400 patients, 1 to 8 episodes each.
# Episodes of the same patient share a patient-level effect, so they are correlated.
n_patients = 400
episodes = rng.integers(1, 9, size=n_patients)
patient_id = np.repeat(np.arange(n_patients), episodes)
patient_effect = rng.normal(0, 1.5, size=n_patients)[patient_id]

risk = patient_effect + rng.normal(0, 1, size=len(patient_id))
y = (risk + rng.normal(0, 1, size=len(risk)) > 1.5).astype(int)

# Two models that see the same risk. Each also makes a patient-level error:
# a model that misjudges a patient tends to misjudge every episode of theirs.
err_a = rng.normal(0, 0.8, size=n_patients)[patient_id]
err_b = rng.normal(0, 0.8, size=n_patients)[patient_id]
score_a = risk + err_a + rng.normal(0, 0.8, size=len(risk))
score_b = risk + err_b + rng.normal(0, 1.0, size=len(risk))

print(f"rows: {len(y)}, patients: {n_patients}, positives: {y.sum()}")
auc_a, auc_b, z, p = delong_paired(y, score_a, score_b)
print(f"AUC A = {auc_a:.3f}  (sklearn: {roc_auc_score(y, score_a):.3f})")
print(f"AUC B = {auc_b:.3f}  (sklearn: {roc_auc_score(y, score_b):.3f})")
print(f"DeLong: z = {z:.2f}, p = {p:.4f}")

lo, hi = patient_bootstrap(y, score_a, score_b, patient_id)
print(f"Patient-level bootstrap 95% CI for A - B: [{lo:.3f}, {hi:.3f}]")

rows: 1789, patients: 400, positives: 391
AUC A = 0.879  (sklearn: 0.879)
AUC B = 0.848  (sklearn: 0.848)
DeLong: z = 3.09, p = 0.0020
Patient-level bootstrap 95% CI for A - B: [0.011, 0.050]

In this run, both methods agree that A is better, because the true gap is large enough. The bootstrap resamples whole patients: if a patient is drawn twice, all of their episodes come along twice. The resample keeps the dependence structure of the real data, so the interval reflects how much the result would move with a different set of patients, not a different set of rows. That is the question a reviewer is asking.

Report the interval, not just the p-value. "A beats B by 0.011 to 0.050 AUC" tells a reader how much the improvement is worth. "p = 0.002" doesn't.

The mistake that inflates significance: rows are not patients

The DeLong variance formula above treats every row as an independent case. With one row per patient, that holds. With several correlated rows per patient, the method thinks it has more evidence than it does, so its p-values come out too small.

You can measure this directly. The simulation below builds 200 cohorts in which the two models are equally good by construction, so any "significant" difference is a false positive. A correct test at the 5% level should flag about 5% of them.

# null_sim.py
import numpy as np
from auc_compare import delong_paired, patient_bootstrap


def make_cohort(rng, n_patients=300):
    """Two equally good models: the true AUC difference is zero."""
    episodes = rng.integers(1, 9, size=n_patients)
    pid = np.repeat(np.arange(n_patients), episodes)
    risk = rng.normal(0, 1.5, n_patients)[pid] + rng.normal(0, 1, len(pid))
    y = (risk + rng.normal(0, 1, len(pid)) > 1.5).astype(int)
    a = risk + rng.normal(0, 1.2, n_patients)[pid] + rng.normal(0, 0.3, len(pid))
    b = risk + rng.normal(0, 1.2, n_patients)[pid] + rng.normal(0, 0.3, len(pid))
    return y, a, b, pid


rng = np.random.default_rng(7)
n_sims = 200
delong_hits = boot_hits = 0
for i in range(n_sims):
    y, a, b, pid = make_cohort(rng)
    if delong_paired(y, a, b)[3] < 0.05:
        delong_hits += 1
    lo, hi = patient_bootstrap(y, a, b, pid, n_boot=500, seed=i)
    if lo > 0 or hi < 0:
        boot_hits += 1

print(f"False positives at the 5% level, {n_sims} cohorts with no real difference:")
print(f"  DeLong, rows treated as independent: {delong_hits / n_sims:.1%}")
print(f"  Patient-level bootstrap:             {boot_hits / n_sims:.1%}")

False positives at the 5% level, 200 cohorts with no real difference:
  DeLong, rows treated as independent: 11.0%
  Patient-level bootstrap:             4.0%

It takes about a minute to run. In this setup, the row-level DeLong test called a nonexistent difference significant roughly twice as often as it should. With 200 simulations these rates are rough, give or take a couple of points, but the direction is the point. How large the inflation gets depends on how many rows each patient has and how strongly a model's errors cluster within a patient. Real clinical data can be worse than this toy.

A test that counts episodes as patients is measuring your data's duplication, not your model.

The same reasoning applies to the split itself. If one patient's episodes land in both training and test, the test set is no longer independent of training. A correct bootstrap can't repair that. That is why the evaluation in my research includes patient-disjoint evaluation and leakage audits, and why the early vs. late fusion essay treats patient-disjoint splits as non-negotiable.

Seed noise: the variance neither test sees

Both tests answer "would this result hold on different patients?" Neither answers "would this result hold if I retrained with a different random seed?" For deep models, that second source of variance can be the same size as the gap you are trying to claim.

# seeds.py
import numpy as np
from sklearn.datasets import make_classification
from sklearn.model_selection import train_test_split
from sklearn.neural_network import MLPClassifier
from sklearn.metrics import roc_auc_score

X, y = make_classification(n_samples=3000, n_features=30, n_informative=8,
                           flip_y=0.1, random_state=0)
X_tr, X_te, y_tr, y_te = train_test_split(X, y, test_size=0.3, random_state=0)

aucs = []
for seed in range(5):
    model = MLPClassifier(hidden_layer_sizes=(64,), max_iter=2000, random_state=seed)
    model.fit(X_tr, y_tr)
    aucs.append(roc_auc_score(y_te, model.predict_proba(X_te)[:, 1]))

print("AUC by seed:", np.round(aucs, 3))
print(f"spread (max - min): {max(aucs) - min(aucs):.3f}")

AUC by seed: [0.911 0.918 0.915 0.916 0.914]
spread (max - min): 0.007

Nothing changed between these five runs except the seed, and the AUC moved by 0.007. If architecture A beats architecture B by 0.005 from a single run each, you haven't shown anything yet. My fusion work involved 8 architectures, and seed-noise analysis is part of its evaluation alongside the patient-level intervals. A gap smaller than the seed spread is not a finding.

A checklist before you report "model A beats model B"

This is the order I would work in. Each step is cheap compared with retracting a claim.

1. Split by patient. No patient appears in both training and test. Check it with a set intersection on patient IDs, not by trusting the split function.
2. Audit for leakage. Look for any feature, note or timestamp that encodes the outcome or comes from after the prediction time.
3. Score both models on the same test cases. A paired comparison needs paired data. Drop or impute consistently for both models.
4. Validate your DeLong code. Match AUCs to roc_auc_score, and match the DeLong interval to a bootstrap on one-row-per-patient data.
5. Count rows per patient. If any patient has more than one row, use the patient-level bootstrap as the primary result.
6. Retrain over seeds. If the gap between models is inside the seed spread, say so.
7. Report the interval and the size of the gap. Write "0.011 to 0.050", not "significant".

None of this makes a weak model look better. It stops a weak comparison from looking stronger than it is. That is the evaluation discipline behind the research listed under my current experience.

The same habit carries over to product work. When an AI feature is about to ship on the strength of one benchmark number, the question to ask is how wide the interval is and what unit was resampled. That is part of what I mean by an AI product having a harder definition of done.