import numpy as np, pandas as pd
from scipy import stats
d = pd.read_csv("marker.csv")
pos, neg = d.marker[d.disease == 1].to_numpy(), d.marker[d.disease == 0].to_numpy()
# AUC = probability a random diseased mouse has a higher value than a random healthy one (ties count half)
psi = (pos[:, None] > neg[None, :]) + 0.5 * (pos[:, None] == neg[None, :])
auc = psi.mean()
# DeLong variance (DeLong et al. 1988)
v10, v01 = psi.mean(axis=1), psi.mean(axis=0)
var = v10.var(ddof=1) / len(pos) + v01.var(ddof=1) / len(neg)
z = stats.norm.ppf(0.975)
print(f"n = {len(pos)} diseased, {len(neg)} healthy; AUC = {auc:.3f} (95% CI {auc - z * np.sqrt(var):.3f} to {auc + z * np.sqrt(var):.3f}, DeLong)")
# Youden cut-off chosen on these data, and its in-sample sensitivity and specificity
cuts = np.sort(np.unique(d.marker))
best = max(cuts, key=lambda c: (pos >= c).mean() + (neg < c).mean() - 1)
print(f"Youden cut-off (chosen on these data): {best:.2f} ng/mL; sensitivity {(pos >= best).mean():.2f}, specificity {(neg < best).mean():.2f}")