{ "cells": [ { "cell_type": "markdown", "id": "d8e5c163", "metadata": {}, "source": [ "# Leakage-aware pre-match volleyball outcome modelling — BMC reproducibility notebook (v5)\n", "\n", "This notebook reproduces the primary analysis from the uploaded/public PlusLiga CSV and adds the final sensitivity analyses used in the BMC submission.\n", "\n", "Key safeguards:\n", "- the target match never contributes to its own feature vector;\n", "- preprocessing is fitted inside training pipelines;\n", "- hyperparameters are selected with expanding time-series folds inside the development cohort;\n", "- the source-labelled 2022/23 cohort is excluded from fitting and tuning, while historical states update prospectively after completed evaluation matches;\n", "- 5,000 outcome-stratified bootstrap resamples quantify fixed-prediction uncertainty;\n", "- Elo parameter and rolling-origin sensitivity analyses are reported separately from the primary analysis.\n", "\n", "The source file contains 14 season labels and no 2019/20 rows; it ends on 5 March 2023, so the 204 source-labelled 2022/23 observations are not described as a complete season.\n" ] }, { "cell_type": "code", "execution_count": 1, "id": "3eea965c", "metadata": { "execution": { "iopub.execute_input": "2026-08-31T18:54:37.320676Z", "iopub.status.busy": "2026-08-31T18:54:37.320571Z", "iopub.status.idle": "2026-08-31T18:55:05.603540Z", "shell.execute_reply": "2026-08-31T18:55:05.603140Z" } }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "SUMMARY\n", "{\n", " \"raw_rows\": 2639,\n", " \"raw_first_date\": \"2008-10-10 18:00:00\",\n", " \"raw_last_date\": \"2023-03-05 14:45:00\",\n", " \"season_labels\": 14,\n", " \"has_2019_20\": false,\n", " \"eligible_rows\": 628,\n", " \"development_n\": 424,\n", " \"test_n\": 204,\n", " \"feature_count\": 17,\n", " \"LR_test_AUC\": 0.7941262277545463,\n", " \"LR_sensitivity_AUC_min\": 0.7786638140620442,\n", " \"LR_sensitivity_AUC_max\": 0.8011280754643586\n", "}\n", "\n", "CV\n", " Model CV AUC mean CV AUC SD Best parameters\n", " Gradient Boosting 0.744819 0.074613 {\"model__learning_rate\": 0.02, \"model__max_depth\": 2, \"model__n_estimators\": 100}\n", " Random Forest 0.733182 0.070224 {\"model__max_depth\": 3, \"model__min_samples_leaf\": 3, \"model__n_estimators\": 200}\n", " XGBoost 0.732552 0.068565 {\"model__colsample_bytree\": 0.8, \"model__learning_rate\": 0.03, \"model__max_depth\": 2, \"model__n_estimators\": 100, \"model__subsample\": 1.0}\n", "Logistic Regression 0.732490 0.086029 {\"model__C\": 0.05, \"model__penalty\": \"l2\"}\n", " MLP 0.693307 0.190221 {\"model__alpha\": 1.0, \"model__hidden_layer_sizes\": [16]}\n", "\n", "TEST\n", " Model Accuracy Precision Recall F1 ROC-AUC Brier Log loss ECE (10 bins)\n", "Logistic Regression 0.7206 0.7593 0.7257 0.7421 0.7941 0.1880 0.5614 0.0736\n", " Random Forest 0.7206 0.7545 0.7345 0.7444 0.7811 0.1963 0.5803 0.0864\n", " Elo baseline 0.7059 0.7431 0.7168 0.7297 0.7629 0.2030 0.5938 0.0942\n", " Gradient Boosting 0.7108 0.7368 0.7434 0.7401 0.7479 0.2001 0.5883 0.0557\n", " XGBoost 0.7157 0.7391 0.7522 0.7456 0.7453 0.2014 0.5916 0.0743\n", " MLP 0.6520 0.6419 0.8407 0.7280 0.7217 0.2149 0.6244 0.0739\n", "Prevalence baseline 0.5539 0.5539 1.0000 0.7129 0.5000 0.2495 0.6922 0.0492\n", "\n", "LR BOOTSTRAP 5000\n", " Model Metric CI_low CI_high\n", "Logistic Regression Accuracy 0.6569 0.7794\n", "Logistic Regression F1 0.6792 0.8000\n", "Logistic Regression AUC 0.7296 0.8528\n", "Logistic Regression Brier 0.1652 0.2127\n", "\n", "PAIRED\n", " Comparator Acc_diff Acc_low Acc_high AUC_diff AUC_low AUC_high\n", " Random Forest 0.0000 -0.0392 0.0441 0.0130 -0.0138 0.0403\n", "Gradient Boosting 0.0098 -0.0343 0.0539 0.0462 0.0134 0.0810\n", " XGBoost 0.0049 -0.0294 0.0441 0.0489 0.0187 0.0804\n", " MLP 0.0686 0.0049 0.1324 0.0724 0.0237 0.1219\n", " Elo baseline 0.0147 -0.0294 0.0588 0.0312 0.0047 0.0584\n", "\n", "ELO SENSITIVITY RANGE\n", " LR_AUC Elo_AUC\n", "min 0.778664 0.736944\n", "max 0.801128 0.786930\n", "\n", "ROLLING ORIGIN\n", " Development Evaluation n_dev n_eval LR_C CV_AUC_mean CV_AUC_SD LR_AUC LR_Acc LR_Brier Elo_AUC Elo_Acc Elo_Brier\n", " 2020/21 2021/22 209 215 0.0500 0.6793 0.1111 0.7244 0.6791 0.2139 0.7719 0.6884 0.2033\n", "2020/21+2021/22 2022/23 424 204 0.0500 0.7325 0.0860 0.7941 0.7206 0.1880 0.7629 0.7059 0.2030\n", "\n", "WARM-START\n", " LR_AUC LR_Accuracy LR_Brier Elo_AUC Elo_Accuracy Elo_Brier\n", " 0.7970 0.7304 0.1872 0.7754 0.7255 0.1976\n" ] } ], "source": [ "from __future__ import annotations\n", "import json, math, warnings\n", "from pathlib import Path\n", "from collections import defaultdict, deque\n", "import numpy as np\n", "import pandas as pd\n", "import matplotlib.pyplot as plt\n", "from scipy.stats import rankdata\n", "from sklearn.base import clone\n", "from sklearn.pipeline import Pipeline\n", "from sklearn.impute import SimpleImputer\n", "from sklearn.preprocessing import StandardScaler\n", "from sklearn.model_selection import TimeSeriesSplit, ParameterGrid\n", "from sklearn.linear_model import LogisticRegression\n", "from sklearn.ensemble import RandomForestClassifier, GradientBoostingClassifier\n", "from sklearn.neural_network import MLPClassifier\n", "from sklearn.metrics import (accuracy_score, precision_score, recall_score, f1_score,\n", " roc_auc_score, brier_score_loss, log_loss, roc_curve)\n", "from sklearn.inspection import permutation_importance\n", "from sklearn.calibration import calibration_curve\n", "from xgboost import XGBClassifier\n", "\n", "warnings.filterwarnings('ignore')\n", "SEED=42\n", "HERE=Path.cwd()\n", "CSV_CANDIDATES=[HERE/'Mens-Volleyball-PlusLiga-2008-2023.csv',\n", " HERE/'Mens-Volleyball-PlusLiga-2008-2023(1).csv',\n", " Path('/mnt/data/Mens-Volleyball-PlusLiga-2008-2023.csv')]\n", "CSV=next((p for p in CSV_CANDIDATES if p.exists()),None)\n", "if CSV is None:\n", " raise FileNotFoundError('Place Mens-Volleyball-PlusLiga-2008-2023.csv beside this notebook.')\n", "OUT=HERE/'revised_outputs_v5'\n", "OUT.mkdir(exist_ok=True)\n", "\n", "def parse_pct(x):\n", " if pd.isna(x): return np.nan\n", " s=str(x).strip().replace('%','').replace(',','.')\n", " try: return float(s)/100.0\n", " except: return np.nan\n", "\n", "def parse_num(x):\n", " if pd.isna(x): return np.nan\n", " s=str(x).strip().replace(',','.')\n", " try: return float(s)\n", " except: return np.nan\n", "\n", "def season_label(dt):\n", " y=dt.year\n", " return f'{y}/{str(y+1)[-2:]}' if dt.month>=7 else f'{y-1}/{str(y)[-2:]}'\n", "\n", "def safe_div(a,b):\n", " return float(a/b) if pd.notna(a) and pd.notna(b) and b!=0 else np.nan\n", "\n", "raw=pd.read_csv(CSV)\n", "raw['dt']=pd.to_datetime(raw['Date'],dayfirst=True,errors='coerce')\n", "raw['season']=raw['dt'].map(season_label)\n", "raw['year_start']=raw['season'].str[:4].astype(int)\n", "srv1=pd.to_numeric(raw['T1_Srv_Sum'],errors='coerce'); ace1=pd.to_numeric(raw['T1_Srv_Ace'],errors='coerce')\n", "srv2=pd.to_numeric(raw['T2_Srv_Sum'],errors='coerce'); ace2=pd.to_numeric(raw['T2_Srv_Ace'],errors='coerce')\n", "raw['serve_header_violation']=(ace1>srv1)|(ace2>srv2)\n", "\n", "# score/outcome agreement across the full source file\n", "calc_winner=np.where(raw['T1_Score']>raw['T2_Score'],0,1)\n", "assert np.array_equal(calc_winner,raw['Winner'].to_numpy())\n", "assert not (raw['season']=='2019/20').any()\n", "\n", "season_audit=raw.groupby('season').agg(matches=('Winner','size'),first_date=('dt','min'),last_date=('dt','max'),serve_ace_header_violations=('serve_header_violation','sum')).reset_index()\n", "season_audit.to_csv(OUT/'season_schema_audit.csv',index=False)\n", "\n", "eligible=raw[raw['year_start']>=2020].copy().sort_values(['dt','Team_1','Team_2']).reset_index(drop=True)\n", "assert len(eligible)==628\n", "assert eligible['serve_header_violation'].sum()==0\n", "\n", "# Further elementary constraints in retained rows\n", "for side in ['T1','T2']:\n", " srv=pd.to_numeric(eligible[f'{side}_Srv_Sum'],errors='coerce')\n", " serr=pd.to_numeric(eligible[f'{side}_Srv_Err'],errors='coerce')\n", " ace=pd.to_numeric(eligible[f'{side}_Srv_Ace'],errors='coerce')\n", " rec=pd.to_numeric(eligible[f'{side}_Rec_Sum'],errors='coerce')\n", " rerr=pd.to_numeric(eligible[f'{side}_Rec_Err'],errors='coerce')\n", " att=pd.to_numeric(eligible[f'{side}_Att_Sum'],errors='coerce')\n", " aerr=pd.to_numeric(eligible[f'{side}_Att_Err'],errors='coerce')\n", " kill=pd.to_numeric(eligible[f'{side}_Att_Kill'],errors='coerce')\n", " assert (serr<=srv).all() and (ace<=srv).all() and (rerr<=rec).all() and (aerr<=att).all() and (kill<=att).all()\n", "\n", "pct_cols=['T1_Srv_Eff','T1_Rec_Pos','T1_Rec_Perf','T1_Att_Kill_Perc','T1_Att_Eff','T2_Srv_Eff','T2_Rec_Pos','T2_Rec_Perf','T2_Att_Kill_Perc','T2_Att_Eff']\n", "for c in pct_cols: eligible[c+'_p']=eligible[c].map(parse_pct)\n", "num_cols=['T1_Sum','T1_BP','T1_Srv_Sum','T1_Srv_Err','T1_Srv_Ace','T1_Rec_Sum','T1_Rec_Err','T1_Att_Sum','T1_Att_Err','T1_Att_Kill','T1_Blk_Sum','T2_Sum','T2_BP','T2_Srv_Sum','T2_Srv_Err','T2_Srv_Ace','T2_Rec_Sum','T2_Rec_Err','T2_Att_Sum','T2_Att_Err','T2_Att_Kill','T2_Blk_Sum']\n", "for c in num_cols: eligible[c+'_n']=eligible[c].map(parse_num)\n", "\n", "metric_names=['win','set_diff','point_diff_per_set','bp_per_set','serve_ace_rate','serve_error_rate','serve_eff','rec_pos','rec_perf','rec_error_rate','attack_kill_rate','attack_eff','attack_error_rate','blocks_per_set']\n", "feature_cols=['elo_diff']+[f'l5_{m}_diff' for m in metric_names]+['prior_matches_log_diff','rest_days_diff']\n", "\n", "def build_features(elo_k=20.0,season_retention=0.75,warmstart_old=False):\n", " hist=defaultdict(lambda:deque(maxlen=5)); elos=defaultdict(lambda:1500.0); last_date={}; last_season=None\n", " if warmstart_old:\n", " older=raw[raw['year_start']<2020].copy().sort_values(['dt','Team_1','Team_2'])\n", " old_season=None\n", " for _,r in older.iterrows():\n", " if old_season is None: old_season=r['season']\n", " elif r['season']!=old_season:\n", " for t in list(elos): elos[t]=1500+season_retention*(elos[t]-1500)\n", " old_season=r['season']\n", " t1,t2=r['Team_1'],r['Team_2']; e1=1/(1+10**((elos[t2]-elos[t1])/400)); a1=1.0 if r['Winner']==0 else 0.0\n", " d=elo_k*(a1-e1); elos[t1]+=d; elos[t2]-=d\n", " # additional shrink across absent 2019/20 period\n", " for t in list(elos): elos[t]=1500+season_retention*(elos[t]-1500)\n", " rows=[]\n", " def agg(team):\n", " recs=list(hist[team])\n", " if not recs: return {m:np.nan for m in metric_names}\n", " out={}\n", " for m in metric_names:\n", " vals=np.array([r[m] for r in recs],dtype=float)\n", " out[m]=float(np.nanmean(vals)) if np.any(~np.isnan(vals)) else np.nan\n", " return out\n", " for _,r in eligible.iterrows():\n", " season=r['season']\n", " if last_season is None: last_season=season\n", " elif season!=last_season:\n", " for t in list(elos): elos[t]=1500+season_retention*(elos[t]-1500)\n", " last_season=season\n", " t1,t2=r['Team_1'],r['Team_2']; h1,h2=agg(t1),agg(t2)\n", " feat={'dt':r['dt'],'season':season,'Team_1':t1,'Team_2':t2,'y':1-int(r['Winner'])}\n", " feat['elo_diff']=(elos[t1]-elos[t2])/400.0\n", " for m in metric_names: feat['l5_'+m+'_diff']=h1[m]-h2[m] if pd.notna(h1[m]) and pd.notna(h2[m]) else np.nan\n", " feat['prior_matches_log_diff']=math.log1p(len(hist[t1]))-math.log1p(len(hist[t2]))\n", " d1=(r['dt']-last_date[t1]).days if t1 in last_date else np.nan; d2=(r['dt']-last_date[t2]).days if t2 in last_date else np.nan\n", " d1=min(d1,60) if pd.notna(d1) else np.nan; d2=min(d2,60) if pd.notna(d2) else np.nan\n", " feat['rest_days_diff']=d1-d2 if pd.notna(d1) and pd.notna(d2) else np.nan\n", " rows.append(feat)\n", " sets=float(r['T1_Score']+r['T2_Score'])\n", " rec1={'win':1.0 if r['Winner']==0 else 0.0,'set_diff':float(r['T1_Score']-r['T2_Score']),'point_diff_per_set':safe_div(r['T1_Sum_n']-r['T2_Sum_n'],sets),'bp_per_set':safe_div(r['T1_BP_n'],sets),'serve_ace_rate':safe_div(r['T1_Srv_Ace_n'],r['T1_Srv_Sum_n']),'serve_error_rate':safe_div(r['T1_Srv_Err_n'],r['T1_Srv_Sum_n']),'serve_eff':r['T1_Srv_Eff_p'],'rec_pos':r['T1_Rec_Pos_p'],'rec_perf':r['T1_Rec_Perf_p'],'rec_error_rate':safe_div(r['T1_Rec_Err_n'],r['T1_Rec_Sum_n']),'attack_kill_rate':safe_div(r['T1_Att_Kill_n'],r['T1_Att_Sum_n']),'attack_eff':r['T1_Att_Eff_p'],'attack_error_rate':safe_div(r['T1_Att_Err_n'],r['T1_Att_Sum_n']),'blocks_per_set':safe_div(r['T1_Blk_Sum_n'],sets)}\n", " rec2={'win':1.0 if r['Winner']==1 else 0.0,'set_diff':float(r['T2_Score']-r['T1_Score']),'point_diff_per_set':safe_div(r['T2_Sum_n']-r['T1_Sum_n'],sets),'bp_per_set':safe_div(r['T2_BP_n'],sets),'serve_ace_rate':safe_div(r['T2_Srv_Ace_n'],r['T2_Srv_Sum_n']),'serve_error_rate':safe_div(r['T2_Srv_Err_n'],r['T2_Srv_Sum_n']),'serve_eff':r['T2_Srv_Eff_p'],'rec_pos':r['T2_Rec_Pos_p'],'rec_perf':r['T2_Rec_Perf_p'],'rec_error_rate':safe_div(r['T2_Rec_Err_n'],r['T2_Rec_Sum_n']),'attack_kill_rate':safe_div(r['T2_Att_Kill_n'],r['T2_Att_Sum_n']),'attack_eff':r['T2_Att_Eff_p'],'attack_error_rate':safe_div(r['T2_Att_Err_n'],r['T2_Att_Sum_n']),'blocks_per_set':safe_div(r['T2_Blk_Sum_n'],sets)}\n", " hist[t1].append(rec1); hist[t2].append(rec2)\n", " e1=1/(1+10**((elos[t2]-elos[t1])/400)); a1=1.0 if r['Winner']==0 else 0.0; d=elo_k*(a1-e1)\n", " elos[t1]+=d; elos[t2]-=d; last_date[t1]=r['dt']; last_date[t2]=r['dt']\n", " return pd.DataFrame(rows)\n", "\n", "feat_df=build_features(); feat_df.to_csv(OUT/'leakage_safe_features.csv',index=False)\n", "train=feat_df[feat_df['season'].isin(['2020/21','2021/22'])].reset_index(drop=True)\n", "test=feat_df[feat_df['season']=='2022/23'].reset_index(drop=True)\n", "assert len(train)==424 and len(test)==204 and len(feature_cols)==17\n", "Xtr,ytr=train[feature_cols],train['y'].astype(int); Xte,yte=test[feature_cols],test['y'].astype(int)\n", "\n", "pipelines={\n", "'Logistic Regression':Pipeline([('imp',SimpleImputer(strategy='median',add_indicator=True)),('scale',StandardScaler()),('model',LogisticRegression(max_iter=3000,random_state=SEED))]),\n", "'Random Forest':Pipeline([('imp',SimpleImputer(strategy='median',add_indicator=True)),('model',RandomForestClassifier(random_state=SEED,n_jobs=-1))]),\n", "'Gradient Boosting':Pipeline([('imp',SimpleImputer(strategy='median',add_indicator=True)),('model',GradientBoostingClassifier(random_state=SEED))]),\n", "'MLP':Pipeline([('imp',SimpleImputer(strategy='median',add_indicator=True)),('scale',StandardScaler()),('model',MLPClassifier(max_iter=2000,early_stopping=True,random_state=SEED))]),\n", "'XGBoost':Pipeline([('imp',SimpleImputer(strategy='median',add_indicator=True)),('model',XGBClassifier(random_state=SEED,eval_metric='logloss',n_jobs=-1,tree_method='hist'))])}\n", "param_grids={\n", "'Logistic Regression':{'model__C':[0.05,0.2,1.0,5.0],'model__penalty':['l2']},\n", "'Random Forest':{'model__n_estimators':[200],'model__max_depth':[3,4,6,None],'model__min_samples_leaf':[3,8]},\n", "'Gradient Boosting':{'model__n_estimators':[100,200],'model__learning_rate':[0.02,0.05],'model__max_depth':[1,2]},\n", "'MLP':{'model__hidden_layer_sizes':[(16,),(32,),(32,16)],'model__alpha':[0.1,1.0]},\n", "'XGBoost':{'model__n_estimators':[100],'model__max_depth':[2,3],'model__learning_rate':[0.03,0.05],'model__subsample':[0.8,1.0],'model__colsample_bytree':[0.8]}}\n", "\n", "tscv=TimeSeriesSplit(n_splits=5); best_models={}; cv_rows=[]\n", "for name,pipe in pipelines.items():\n", " best=(-np.inf,None,None)\n", " for params in ParameterGrid(param_grids[name]):\n", " fs=[]\n", " for tri,vai in tscv.split(Xtr):\n", " m=clone(pipe).set_params(**params).fit(Xtr.iloc[tri],ytr.iloc[tri]); p=m.predict_proba(Xtr.iloc[vai])[:,1]\n", " fs.append(roc_auc_score(ytr.iloc[vai],p))\n", " if np.mean(fs)>best[0]: best=(float(np.mean(fs)),params,fs)\n", " score,params,fs=best; final=clone(pipe).set_params(**params).fit(Xtr,ytr); best_models[name]=final\n", " cv_rows.append({'Model':name,'CV AUC mean':score,'CV AUC SD':float(np.std(fs,ddof=1)),'Best parameters':json.dumps(params,sort_keys=True)})\n", "cv_df=pd.DataFrame(cv_rows).sort_values('CV AUC mean',ascending=False); cv_df.to_csv(OUT/'cv_model_selection.csv',index=False)\n", "\n", "probs={name:m.predict_proba(Xte)[:,1] for name,m in best_models.items()}\n", "prev=float(ytr.mean()); probs['Prevalence baseline']=np.full(len(yte),prev); probs['Elo baseline']=1/(1+10**(-Xte['elo_diff'].fillna(0).to_numpy()))\n", "preds={n:(p>=0.5).astype(int) for n,p in probs.items()}\n", "\n", "def ece(y,p,n_bins=10):\n", " bins=np.linspace(0,1,n_bins+1); inds=np.digitize(p,bins[1:-1],right=False); ans=0.0\n", " for b in range(n_bins):\n", " m=inds==b\n", " if m.sum(): ans+=m.mean()*abs(y.to_numpy()[m].mean()-p[m].mean())\n", " return float(ans)\n", "\n", "metrics=[]\n", "for n,p in probs.items():\n", " pr=preds[n]\n", " metrics.append({'Model':n,'Accuracy':accuracy_score(yte,pr),'Precision':precision_score(yte,pr,zero_division=0),'Recall':recall_score(yte,pr,zero_division=0),'F1':f1_score(yte,pr,zero_division=0),'ROC-AUC':roc_auc_score(yte,p) if np.unique(p).size>1 else 0.5,'Brier':brier_score_loss(yte,p),'Log loss':log_loss(yte,np.clip(p,1e-8,1-1e-8)),'ECE (10 bins)':ece(yte,p,10)})\n", "metrics_df=pd.DataFrame(metrics).sort_values('ROC-AUC',ascending=False); metrics_df.to_csv(OUT/'test_metrics.csv',index=False)\n", "\n", "# Fast 5,000-resample outcome-stratified bootstrap on fixed predictions\n", "rng=np.random.default_rng(SEED); B=5000; yarr=yte.to_numpy(); i0=np.where(yarr==0)[0]; i1=np.where(yarr==1)[0]\n", "def fast_auc(y,p):\n", " n1=int(y.sum()); n0=len(y)-n1; r=rankdata(p,method='average'); return (r[y==1].sum()-n1*(n1+1)/2)/(n1*n0)\n", "def fast_f1(y,pr):\n", " tp=np.sum((y==1)&(pr==1)); fp=np.sum((y==0)&(pr==1)); fn=np.sum((y==1)&(pr==0)); den=2*tp+fp+fn; return 2*tp/den if den else 0.0\n", "names=list(probs); boot={n:{k:np.empty(B) for k in ['Accuracy','F1','AUC','Brier']} for n in names}; others=['Random Forest','Gradient Boosting','XGBoost','MLP','Elo baseline']; pairA={o:np.empty(B) for o in others}; pairAcc={o:np.empty(B) for o in others}\n", "for b in range(B):\n", " idx=np.concatenate([rng.choice(i0,len(i0),True),rng.choice(i1,len(i1),True)]); yy=yarr[idx]\n", " for n in names:\n", " pp=probs[n][idx]; pr=preds[n][idx]; boot[n]['Accuracy'][b]=np.mean(yy==pr); boot[n]['F1'][b]=fast_f1(yy,pr); boot[n]['AUC'][b]=fast_auc(yy,pp) if np.unique(pp).size>1 else .5; boot[n]['Brier'][b]=np.mean((yy-pp)**2)\n", " lrA=boot['Logistic Regression']['AUC'][b]; lrAcc=boot['Logistic Regression']['Accuracy'][b]\n", " for o in others:\n", " pairA[o][b]=lrA-fast_auc(yy,probs[o][idx]); pairAcc[o][b]=lrAcc-np.mean(yy==preds[o][idx])\n", "ci_rows=[]\n", "for n in names:\n", " for met,a in boot[n].items():\n", " lo,hi=np.quantile(a,[.025,.975]); ci_rows.append({'Model':n,'Metric':met,'CI_low':lo,'CI_high':hi})\n", "ci_df=pd.DataFrame(ci_rows); ci_df.to_csv(OUT/'bootstrap_ci_5000.csv',index=False)\n", "pair_rows=[]\n", "for o in others:\n", " pair_rows.append({'Comparator':o,'Acc_diff':accuracy_score(yarr,preds['Logistic Regression'])-accuracy_score(yarr,preds[o]),'Acc_low':np.quantile(pairAcc[o],.025),'Acc_high':np.quantile(pairAcc[o],.975),'AUC_diff':fast_auc(yarr,probs['Logistic Regression'])-fast_auc(yarr,probs[o]),'AUC_low':np.quantile(pairA[o],.025),'AUC_high':np.quantile(pairA[o],.975)})\n", "pair_df=pd.DataFrame(pair_rows); pair_df.to_csv(OUT/'paired_comparisons_5000.csv',index=False)\n", "\n", "perm=permutation_importance(best_models['Logistic Regression'],Xte,yte,n_repeats=20,random_state=SEED,scoring='roc_auc',n_jobs=-1)\n", "imp_df=pd.DataFrame({'Feature':feature_cols,'AUC decrease mean':perm.importances_mean,'AUC decrease SD':perm.importances_std}).sort_values('AUC decrease mean',ascending=False); imp_df.to_csv(OUT/'permutation_importance.csv',index=False)\n", "\n", "# Elo parameter sensitivity\n", "sens=[]\n", "for K in [10,20,30,40]:\n", " for ret in [0.0,0.5,0.75,1.0]:\n", " f=build_features(K,ret,False); tr=f[f.season.isin(['2020/21','2021/22'])].reset_index(drop=True); te=f[f.season=='2022/23'].reset_index(drop=True)\n", " lr=clone(pipelines['Logistic Regression']).set_params(model__C=0.05,model__penalty='l2').fit(tr[feature_cols],tr.y.astype(int)); plr=lr.predict_proba(te[feature_cols])[:,1]; pelo=1/(1+10**(-te['elo_diff'].fillna(0).to_numpy()))\n", " sens.append({'K':K,'Season_retention':ret,'LR_AUC':roc_auc_score(te.y,plr),'LR_Acc':accuracy_score(te.y,plr>=.5),'LR_Brier':brier_score_loss(te.y,plr),'Elo_AUC':roc_auc_score(te.y,pelo),'Elo_Acc':accuracy_score(te.y,pelo>=.5),'Elo_Brier':brier_score_loss(te.y,pelo)})\n", "sens_df=pd.DataFrame(sens); sens_df.to_csv(OUT/'elo_parameter_sensitivity.csv',index=False)\n", "\n", "# Proper rolling-origin LR sensitivity: tune C only within each development window\n", "rolling=[]\n", "for trseas,teseas in [(['2020/21'],'2021/22'),(['2020/21','2021/22'],'2022/23')]:\n", " tr=feat_df[feat_df.season.isin(trseas)].reset_index(drop=True); te=feat_df[feat_df.season==teseas].reset_index(drop=True); ts=TimeSeriesSplit(n_splits=5)\n", " best=(-np.inf,None,None)\n", " for C in [0.05,0.2,1.0,5.0]:\n", " fs=[]\n", " for tri,vai in ts.split(tr):\n", " m=clone(pipelines['Logistic Regression']).set_params(model__C=C,model__penalty='l2').fit(tr.loc[tri,feature_cols],tr.loc[tri,'y'].astype(int)); pp=m.predict_proba(tr.loc[vai,feature_cols])[:,1]; fs.append(roc_auc_score(tr.loc[vai,'y'],pp))\n", " if np.mean(fs)>best[0]: best=(float(np.mean(fs)),C,fs)\n", " cvmean,C,fs=best; m=clone(pipelines['Logistic Regression']).set_params(model__C=C,model__penalty='l2').fit(tr[feature_cols],tr.y.astype(int)); pp=m.predict_proba(te[feature_cols])[:,1]; pe=1/(1+10**(-te['elo_diff'].fillna(0).to_numpy()))\n", " rolling.append({'Development':'+'.join(trseas),'Evaluation':teseas,'n_dev':len(tr),'n_eval':len(te),'LR_C':C,'CV_AUC_mean':cvmean,'CV_AUC_SD':np.std(fs,ddof=1),'LR_AUC':roc_auc_score(te.y,pp),'LR_Acc':accuracy_score(te.y,pp>=.5),'LR_Brier':brier_score_loss(te.y,pp),'Elo_AUC':roc_auc_score(te.y,pe),'Elo_Acc':accuracy_score(te.y,pe>=.5),'Elo_Brier':brier_score_loss(te.y,pe)})\n", "rolling_df=pd.DataFrame(rolling); rolling_df.to_csv(OUT/'rolling_origin_sensitivity.csv',index=False)\n", "\n", "# Outcome-only warm-start sensitivity using older reliable score/outcome fields\n", "fw=build_features(20,0.75,True); tr=fw[fw.season.isin(['2020/21','2021/22'])]; te=fw[fw.season=='2022/23']; lr=clone(pipelines['Logistic Regression']).set_params(model__C=.05,model__penalty='l2').fit(tr[feature_cols],tr.y.astype(int)); pp=lr.predict_proba(te[feature_cols])[:,1]; pe=1/(1+10**(-te['elo_diff'].fillna(0).to_numpy()))\n", "warm_df=pd.DataFrame([{'LR_AUC':roc_auc_score(te.y,pp),'LR_Accuracy':accuracy_score(te.y,pp>=.5),'LR_Brier':brier_score_loss(te.y,pp),'Elo_AUC':roc_auc_score(te.y,pe),'Elo_Accuracy':accuracy_score(te.y,pe>=.5),'Elo_Brier':brier_score_loss(te.y,pe)}]); warm_df.to_csv(OUT/'outcome_only_warmstart_sensitivity.csv',index=False)\n", "\n", "pred_df=test[['dt','season','Team_1','Team_2','y']].copy()\n", "for n,p in probs.items(): pred_df[n+'_p']=p\n", "pred_df.to_csv(OUT/'heldout_predictions.csv',index=False)\n", "\n", "summary={'raw_rows':len(raw),'raw_first_date':str(raw.dt.min()),'raw_last_date':str(raw.dt.max()),'season_labels':raw.season.nunique(),'has_2019_20':bool((raw.season=='2019/20').any()),'eligible_rows':len(eligible),'development_n':len(train),'test_n':len(test),'feature_count':len(feature_cols),'LR_test_AUC':float(roc_auc_score(yte,probs['Logistic Regression'])),'LR_sensitivity_AUC_min':float(sens_df.LR_AUC.min()),'LR_sensitivity_AUC_max':float(sens_df.LR_AUC.max())}\n", "with open(OUT/'analysis_summary_v5.json','w') as f: json.dump(summary,f,indent=2)\n", "\n", "print('SUMMARY'); print(json.dumps(summary,indent=2))\n", "print('\\nCV'); print(cv_df.to_string(index=False))\n", "print('\\nTEST'); print(metrics_df.to_string(index=False,float_format=lambda x:f'{x:.4f}'))\n", "print('\\nLR BOOTSTRAP 5000'); print(ci_df[ci_df.Model=='Logistic Regression'].to_string(index=False,float_format=lambda x:f'{x:.4f}'))\n", "print('\\nPAIRED'); print(pair_df.to_string(index=False,float_format=lambda x:f'{x:.4f}'))\n", "print('\\nELO SENSITIVITY RANGE'); print(sens_df[['LR_AUC','Elo_AUC']].agg(['min','max']))\n", "print('\\nROLLING ORIGIN'); print(rolling_df.to_string(index=False,float_format=lambda x:f'{x:.4f}'))\n", "print('\\nWARM-START'); print(warm_df.to_string(index=False,float_format=lambda x:f'{x:.4f}'))\n" ] }, { "cell_type": "markdown", "id": "4ebfef0e", "metadata": {}, "source": [ "## Interpretation guardrails\n", "\n", "The primary held-out estimates apply only to the 204 source-labelled 2022/23 matches available in this file. They are not evidence of live in-game performance or cross-league external validity. The Elo and rolling-origin analyses are sensitivity checks; they do not convert the single-league study into external validation." ] } ], "metadata": { "kernelspec": { "display_name": "Python 3", "language": "python", "name": "python3" }, "language_info": { "codemirror_mode": { "name": "ipython", "version": 3 }, "file_extension": ".py", "mimetype": "text/x-python", "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", "version": "3.13.5" } }, "nbformat": 4, "nbformat_minor": 5 }