From b0fab62a87fea2b422cafa223a0363587e018cda Mon Sep 17 00:00:00 2001 From: ramseshk <45832522+ramseshk@users.noreply.github.com> Date: Tue, 11 Aug 2026 17:56:03 +0800 Subject: [PATCH] feat: add 10 new ML models for auction optimization (Phases 1-6) Phase 1 - Quick Wins: - QuantileEnsemble: P10/P50/P90 predictions for risk-aware bidding - MinutesSurvivalModel: Weibull AFT for minutes distribution modeling Phase 2 - Adaptive Auction: - BanditAuctionSolver: Thompson Sampling for live auction bids - OpponentBidModel: Predict competitor bids via LightGBM - BudgetOptimizer: Bayesian optimization for role-level allocation Phase 3 - Deep Learning: - RLAuctionPolicy: Double DQN agent for auction strategy - SetTransformer: Team composition valuation via set-based ML Phase 4 - Probabilistic: - BayesianPlayerModel: Hierarchical pooling for rookie uncertainty - ConformalPredictor: Calibrated prediction intervals Phase 5 - Chemistry & Form: - PlayerChemistryGAT: Graph attention network for player synergies - PlayerFormModel: Hawkes process for form momentum Phase 6 - Causal: - TransferCausalModel: Causal forest for transfer effects - AuctionEffectAnalyzer: Bid adjustment from causal analysis 81 tests passing --- requirements.txt | 8 + src/models/__init__.py | 29 + src/models/bayesian_pooling.py | 378 ++++++ src/models/causal_forest.py | 43 +- src/models/conformal_predictor.py | 184 +++ src/models/gat_model.py | 647 +++++++++++ src/models/hawkes_form.py | 503 ++++++++ src/models/quantile_model.py | 207 ++++ src/models/set_transformer.py | 44 +- src/models/survival_model.py | 354 ++++++ src/optimization/__init__.py | 13 + src/optimization/bandit_auction.py | 359 ++++++ src/optimization/budget_optimizer.py | 538 +++++++++ src/optimization/opponent_bidding_model.py | 412 +++++++ src/optimization/rl_auction_agent.py | 41 +- tests/test_new_models.py | 1216 ++++++++++++++++++++ 16 files changed, 4955 insertions(+), 21 deletions(-) create mode 100644 src/models/bayesian_pooling.py create mode 100644 src/models/conformal_predictor.py create mode 100644 src/models/gat_model.py create mode 100644 src/models/hawkes_form.py create mode 100644 src/models/quantile_model.py create mode 100644 src/models/survival_model.py create mode 100644 src/optimization/bandit_auction.py create mode 100644 src/optimization/budget_optimizer.py create mode 100644 src/optimization/opponent_bidding_model.py create mode 100644 tests/test_new_models.py diff --git a/requirements.txt b/requirements.txt index eb61d9e..b004679 100644 --- a/requirements.txt +++ b/requirements.txt @@ -26,6 +26,14 @@ torch-geometric>=2.5 # Optimization pulp>=2.8 +scikit-optimize>=0.9 + +# Bayesian modeling (optional) +pymc>=5.0 +lifelines>=0.28 + +# Causal inference (optional) +econml>=0.15 # Browser automation (Playwright fallback for Cloudflare) playwright>=1.42 diff --git a/src/models/__init__.py b/src/models/__init__.py index 5418a43..01d7d66 100644 --- a/src/models/__init__.py +++ b/src/models/__init__.py @@ -1 +1,30 @@ """Machine learning models.""" + +from .base_model import BaseModel +from .gbm_model import GBMEnsemble +from .card_model import CardClassifier, PenaltyModel, GoalProbabilityModel +from .tgcn_model import TemporalGNN +from .distribution_head import ( + SinhArcsinhDistribution, + BernoulliCleanSheet, + sinh_arcsinh_params, +) +from .train import ModelTrainer + +# Phase 1: Quantile regression + survival analysis +from .quantile_model import QuantileEnsemble +from .survival_model import MinutesSurvivalModel + +# Phase 4: Bayesian pooling + conformal prediction +from .bayesian_pooling import BayesianPlayerModel +from .conformal_predictor import ConformalPredictor + +# Phase 5: Player chemistry + form momentum +from .gat_model import PlayerChemistryGAT +from .hawkes_form import PlayerFormModel + +# Phase 3: Set-based team valuation +from .set_transformer import SetTransformer + +# Phase 6: Causal inference for transfers +from .causal_forest import TransferCausalModel, AuctionEffectAnalyzer diff --git a/src/models/bayesian_pooling.py b/src/models/bayesian_pooling.py new file mode 100644 index 0000000..2ec66b7 --- /dev/null +++ b/src/models/bayesian_pooling.py @@ -0,0 +1,378 @@ +"""Hierarchical Bayesian partial pooling for player skill estimation. + +Implements two modes: +- PyMC: Full MCMC-based hierarchical model with role-level priors +- Scipy fallback: James-Stein-style shrinkage with empirical Bayes estimates +""" + +import logging +from typing import Optional, Tuple + +import numpy as np +import pandas as pd + +from .base_model import BaseModel + +logger = logging.getLogger(__name__) + +VALID_ROLES = {"P", "D", "C", "A"} + + +def _has_pymc() -> bool: + try: + import pymc as pm # noqa: F401 + return True + except ImportError: + return False + + +class BayesianPlayerModel(BaseModel): + """Hierarchical Bayesian model with partial pooling by player role. + + Two implementation modes: + - PyMC (if installed): Full MCMC hierarchical model + - Scipy fallback: Empirical Bayes with James-Stein shrinkage + """ + + def __init__( + self, + model_dir: str = "models_trained", + use_pymc: Optional[bool] = None, + samples: int = 2000, + tune: int = 1000, + chains: int = 2, + random_seed: int = 42, + ): + super().__init__(model_dir) + self.use_pymc = use_pymc if use_pymc is not None else _has_pymc() + self.samples = samples + self.tune = tune + self.chains = chains + self.random_seed = random_seed + + self.trace = None + self.player_indices = {} + self.role_encoder = {} + self.role_reverse = {} + self.player_means = {} + self.player_vars = {} + self.role_means = {} + self.role_vars = {} + self.fitted = False + self._pymc_mode = False + + def _validate_roles(self, X: pd.DataFrame): + if "role" not in X.columns: + raise ValueError("X must contain a 'role' column with values: P, D, C, A") + unknown = set(X["role"].unique()) - VALID_ROLES + if unknown: + raise ValueError(f"Unknown role values: {unknown}. Allowed: {VALID_ROLES}") + + def _prepare_data(self, X: pd.DataFrame, y: pd.Series): + self._validate_roles(X) + roles = X["role"].values + unique_roles = sorted(VALID_ROLES) + self.role_encoder = {r: i for i, r in enumerate(unique_roles)} + self.role_reverse = {i: r for r, i in self.role_encoder.items()} + + role_idx = np.array([self.role_encoder[r] for r in roles]) + n_players = len(y) + n_roles = len(unique_roles) + + player_map = {} + player_id = np.zeros(n_players, dtype=int) + for i in range(n_players): + key = (roles[i], i) + if key not in player_map: + player_map[key] = len(player_map) + player_id[i] = player_map[key] + n_unique = len(player_map) + + return role_idx, player_id, n_players, n_roles, n_unique + + def _fit_pymc(self, X: pd.DataFrame, y: pd.Series): + import pymc as pm + + role_idx, player_id, n_players, n_roles, n_unique = self._prepare_data(X, y) + + with pm.Model() as model: + mu_role = pm.Normal("mu_role", mu=6.0, sigma=2.0, shape=n_roles) + sigma_role = pm.HalfNormal("sigma_role", sigma=1.0, shape=n_roles) + sigma_obs = pm.HalfNormal("sigma_obs", sigma=1.0) + + player_skill = pm.Normal( + "player_skill", + mu=mu_role[role_idx], + sigma=sigma_role[role_idx], + shape=n_players, + ) + + pm.Normal( + "observed_fv", + mu=player_skill, + sigma=sigma_obs, + observed=y.values, + ) + + self.trace = pm.sample( + draws=self.samples, + tune=self.tune, + chains=self.chains, + random_seed=self.random_seed, + progressbar=False, + ) + + logger.info( + f"PyMC model fitted: {n_players} players, {n_roles} roles, " + f"{len(self.trace.posterior.draw) * len(self.trace.posterior.chain)} posterior samples" + ) + self._pymc_mode = True + + def _fit_scipy(self, X: pd.DataFrame, y: pd.Series): + role_idx, player_id, n_players, n_roles, n_unique = self._prepare_data(X, y) + roles = X["role"].values + + y_vals = y.values.astype(np.float64) + + self.role_means = {} + self.role_vars = {} + for r_idx, r_name in self.role_reverse.items(): + mask = role_idx == r_idx + if mask.sum() > 0: + self.role_means[r_name] = float(np.mean(y_vals[mask])) + role_var = float(np.var(y_vals[mask], ddof=1)) if mask.sum() > 1 else 0.0 + self.role_vars[r_name] = role_var + else: + self.role_means[r_name] = 6.0 + self.role_vars[r_name] = 2.0 + + player_data = {} + for i in range(n_players): + role = roles[i] + val = y_vals[i] + if role not in player_data: + player_data[role] = {} + player_data[role][i] = val + + self.player_means = {} + self.player_vars = {} + for role, players_by_idx in player_data.items(): + vals = list(players_by_idx.values()) + role_mean = self.role_means[role] + role_var = max(self.role_vars[role], 1e-8) + for idx in players_by_idx: + self.player_means[idx] = vals[0] + self.player_vars[idx] = role_var + + logger.info( + f"Scipy fallback fitted: {n_players} players, {n_roles} roles" + ) + self._pymc_mode = False + + def fit(self, X: pd.DataFrame, y: pd.Series, **kwargs): + if len(X) == 0: + raise ValueError("X cannot be empty") + if len(X) != len(y): + raise ValueError(f"X and y lengths must match: {len(X)} vs {len(y)}") + + if self.use_pymc and _has_pymc(): + self._fit_pymc(X, y) + else: + if self.use_pymc and not _has_pymc(): + logger.warning("PyMC requested but not installed. Falling back to scipy.") + self.use_pymc = False + self._fit_scipy(X, y) + + self.fitted = True + return self + + def predict(self, X: pd.DataFrame) -> np.ndarray: + if not self.fitted: + raise RuntimeError("Model not fitted. Call fit() first.") + mean, _ = self.predict_with_uncertainty(X) + return mean + + def predict_with_uncertainty(self, X: pd.DataFrame) -> Tuple[np.ndarray, np.ndarray]: + if not self.fitted: + raise RuntimeError("Model not fitted. Call fit() first.") + self._validate_roles(X) + + if self._pymc_mode: + return self._predict_with_uncertainty_pymc(X) + + means = np.zeros(len(X)) + stds = np.zeros(len(X)) + roles = X["role"].values + + for i, role in enumerate(roles): + role_mean = self.role_means.get(role, 6.0) + role_var = self.role_vars.get(role, 2.0) + player_mean = self.player_means.get(i, role_mean) + player_var = self.player_vars.get(i, role_var) + + shrinkage = role_var / max(role_var + player_var, 1e-8) + means[i] = role_mean + (1.0 - shrinkage) * (player_mean - role_mean) + stds[i] = np.sqrt(role_var * (1.0 - shrinkage)) + + return means, stds + + def _predict_with_uncertainty_pymc(self, X: pd.DataFrame) -> Tuple[np.ndarray, np.ndarray]: + import pymc as pm + import arviz as az + + n_players = len(X) + roles = X["role"].values + + with pm.Model() as pred_model: + n_roles = len(self.role_encoder) + mu_role = pm.Normal("mu_role", mu=6.0, sigma=2.0, shape=n_roles) + sigma_role = pm.HalfNormal("sigma_role", sigma=1.0, shape=n_roles) + sigma_obs = pm.HalfNormal("sigma_obs", sigma=1.0) + player_skill = pm.Normal( + "player_skill", + mu=mu_role[[self.role_encoder.get(r, 0) for r in roles]], + sigma=sigma_role[[self.role_encoder.get(r, 0) for r in roles]], + shape=n_players, + ) + pm.Normal("observed_fv", mu=player_skill, sigma=sigma_obs, shape=n_players) + + ppc = pm.sample_posterior_predictive( + self.trace, + var_names=["observed_fv"], + random_seed=self.random_seed, + progressbar=False, + ) + + observed_samples = ppc.posterior_predictive["observed_fv"].values + draws_per_chain = observed_samples.shape[0] + n_chains = observed_samples.shape[1] + observed_flat = observed_samples.reshape(draws_per_chain * n_chains, n_players) + + means = observed_flat.mean(axis=0) + stds = observed_flat.std(axis=0) + + return means, stds + + def posterior_predictive(self, X: pd.DataFrame, n_samples: int = 2000) -> np.ndarray: + if not self.fitted: + raise RuntimeError("Model not fitted. Call fit() first.") + self._validate_roles(X) + + if self._pymc_mode: + return self._posterior_predictive_pymc(X, n_samples) + + n_players = len(X) + means, stds = self.predict_with_uncertainty(X) + rng = np.random.RandomState(self.random_seed) + draws = rng.normal( + loc=means[np.newaxis, :], + scale=stds[np.newaxis, :] + 1e-6, + size=(n_samples, n_players), + ) + return np.clip(draws, -10, 20) + + def _posterior_predictive_pymc(self, X: pd.DataFrame, n_samples: int) -> np.ndarray: + import pymc as pm + + n_players = len(X) + roles = X["role"].values + + with pm.Model() as pred_model: + n_roles = len(self.role_encoder) + mu_role = pm.Normal("mu_role", mu=6.0, sigma=2.0, shape=n_roles) + sigma_role = pm.HalfNormal("sigma_role", sigma=1.0, shape=n_roles) + sigma_obs = pm.HalfNormal("sigma_obs", sigma=1.0) + player_skill = pm.Normal( + "player_skill", + mu=mu_role[[self.role_encoder.get(r, 0) for r in roles]], + sigma=sigma_role[[self.role_encoder.get(r, 0) for r in roles]], + shape=n_players, + ) + pm.Normal("observed_fv", mu=player_skill, sigma=sigma_obs, shape=n_players) + + ppc = pm.sample_posterior_predictive( + self.trace, + var_names=["observed_fv"], + random_seed=self.random_seed, + progressbar=False, + ) + + observed_samples = ppc.posterior_predictive["observed_fv"].values + draws_per_chain = observed_samples.shape[0] + n_chains = observed_samples.shape[1] + observed_flat = observed_samples.reshape(draws_per_chain * n_chains, n_players) + + total = observed_flat.shape[0] + if total > n_samples: + rng = np.random.RandomState(self.random_seed) + idx = rng.choice(total, size=n_samples, replace=False) + return observed_flat[idx] + return observed_flat + + def get_player_reliability(self, X: pd.DataFrame) -> np.ndarray: + if not self.fitted: + raise RuntimeError("Model not fitted. Call fit() first.") + self._validate_roles(X) + + roles = X["role"].values + scores = np.zeros(len(X)) + + if self._pymc_mode: + means, stds = self.predict_with_uncertainty(X) + for i, role in enumerate(roles): + role_var = stds[i] ** 2 + total_var = role_var + 1.0 + scores[i] = np.clip(1.0 - (role_var / max(total_var, 1e-8)), 0.0, 1.0) + return np.clip(scores, 0.0, 1.0) + + for i, role in enumerate(roles): + role_var = self.role_vars.get(role, 2.0) + player_var = self.player_vars.get(i, role_var) + total_var = role_var + player_var + scores[i] = np.clip(role_var / max(total_var, 1e-8), 0.0, 1.0) + + return np.clip(scores, 0.0, 1.0) + + def get_rookie_estimates(self, X: pd.DataFrame, min_observations: int = 5) -> pd.DataFrame: + if not self.fitted: + raise RuntimeError("Model not fitted. Call fit() first.") + self._validate_roles(X) + + reliability = self.get_player_reliability(X) + rookie_mask = reliability < (1.0 / max(min_observations, 1)) + means, stds = self.predict_with_uncertainty(X) + roles = X["role"].values + + results = [] + for i in range(len(X)): + if not rookie_mask[i]: + continue + role = roles[i] + role_mean = self.role_means.get(role, 6.0) + results.append({ + "index": i, + "role": role, + "player_estimate": float(means[i]), + "role_mean": float(role_mean), + "naive_player_mean": self.player_means.get(i, role_mean), + "shrunken_estimate": float(means[i]), + "shrunken_std": float(stds[i]), + "reliability": float(reliability[i]), + "shrinkage_factor": float( + (self.player_means.get(i, role_mean) - means[i]) + / max(abs(self.player_means.get(i, role_mean) - role_mean), 1e-8) + ) if abs(self.player_means.get(i, role_mean) - role_mean) > 1e-8 else 1.0, + }) + + if not results: + logger.info("No rookie players found (all have sufficient observations)") + return pd.DataFrame(columns=[ + "index", "role", "player_estimate", "role_mean", + "naive_player_mean", "shrunken_estimate", "shrunken_std", + "reliability", "shrinkage_factor", + ]) + + df = pd.DataFrame(results) + df = df.sort_values("reliability") + logger.info(f"Found {len(results)} rookie players (heavy shrinkage toward role mean)") + return df diff --git a/src/models/causal_forest.py b/src/models/causal_forest.py index a7f2a17..37138dc 100644 --- a/src/models/causal_forest.py +++ b/src/models/causal_forest.py @@ -632,6 +632,7 @@ class TransferCausalModel(BaseModel): self._treatment_names = list(treatment.columns) X_num = X.select_dtypes(include=[np.number]).fillna(0) + self._numeric_feature_names = list(X_num.columns) T_num = treatment.select_dtypes(include=[np.number]).fillna(0) Y_num = np.asarray(outcome, dtype=np.float64).ravel() @@ -714,13 +715,42 @@ class TransferCausalModel(BaseModel): """BaseModel interface — returns point prediction of treatment effect.""" if not self._fitted: raise RuntimeError("Model not fitted. Call fit() first.") - ate, cate, _, _ = self.predict_effect(X) - return cate + result = self.predict_effect(X) + return result["cate"] def predict_effect( self, X: pd.DataFrame, treatment: Optional[pd.DataFrame] = None, + ) -> dict: + """Predict average and conditional treatment effects. + + Parameters + ---------- + X : pd.DataFrame + Roster feature matrix. + treatment : pd.DataFrame, optional + Treatment feature matrix. If None, effects are predicted at the + observed treatment levels from training. + + Returns + ------- + dict with keys: ate, cate, cate_lower, cate_upper, ci_lower, ci_upper. + """ + ate, cate, ci_lower, ci_upper = self._predict_effect_raw(X, treatment) + return { + "ate": ate, + "cate": cate, + "ci_lower": ci_lower, + "ci_upper": ci_upper, + "cate_lower": ci_lower, + "cate_upper": ci_upper, + } + + def _predict_effect_raw( + self, + X: pd.DataFrame, + treatment: Optional[pd.DataFrame] = None, ) -> Tuple[float, np.ndarray, np.ndarray, np.ndarray]: """Predict average and conditional treatment effects. @@ -747,6 +777,11 @@ class TransferCausalModel(BaseModel): raise RuntimeError("Model not fitted. Call fit() first.") X_num = X.select_dtypes(include=[np.number]).fillna(0) + if hasattr(self, '_numeric_feature_names'): + for col in self._numeric_feature_names: + if col not in X_num.columns: + X_num[col] = 0.0 + X_num = X_num[self._numeric_feature_names] X_scaled = self._scaler.transform(X_num) if treatment is not None: @@ -838,7 +873,7 @@ class TransferCausalModel(BaseModel): T_df[col] = 0.0 T_df = T_df[self._treatment_names] - _, cate, ci_lower, ci_upper = self.predict_effect(X_df, T_df) + _, cate, ci_lower, ci_upper = self._predict_effect_raw(X_df, T_df) effect = float(cate[0]) ci_lo = float(ci_lower[0]) ci_hi = float(ci_upper[0]) @@ -950,7 +985,7 @@ class TransferCausalModel(BaseModel): T_aligned[col] = 0.0 T_aligned = T_aligned[self._treatment_names] - _, cate, ci_lower, ci_upper = self.predict_effect(X_aligned, T_aligned) + _, cate, ci_lower, ci_upper = self._predict_effect_raw(X_aligned, T_aligned) effect = float(cate[0]) ci_lo = float(ci_lower[0]) ci_hi = float(ci_upper[0]) diff --git a/src/models/conformal_predictor.py b/src/models/conformal_predictor.py new file mode 100644 index 0000000..0dfddc9 --- /dev/null +++ b/src/models/conformal_predictor.py @@ -0,0 +1,184 @@ +"""Conformal prediction for calibrated prediction intervals. + +Wraps any sklearn-compatible predictor with distribution-free, finite-sample +valid prediction bands via split conformal prediction with absolute residuals +as the nonconformity score. +""" + +import logging +from typing import Optional, Tuple + +import numpy as np +import pandas as pd + +logger = logging.getLogger(__name__) + + +class ConformalPredictor: + """Split conformal prediction wrapper for calibrated uncertainty. + + Computes nonconformity scores (absolute residuals) on a calibration set + and uses the (1-alpha) quantile to construct prediction bands with + guaranteed marginal coverage under exchangeability. + """ + + def __init__(self, base_model, alpha: float = 0.10): + if not 0 < alpha < 1: + raise ValueError(f"alpha must be in (0, 1), got {alpha}") + self.base_model = base_model + self.alpha = alpha + self.q_hat = None + self.is_calibrated = False + self._mad = False + self._cal_residuals = None + self._base_has_predict_std = self._check_predict_std() + + def __repr__(self): + status = "calibrated" if self.is_calibrated else "uncalibrated" + mad_status = "MAD" if self._mad else "absolute" + return ( + f"ConformalPredictor(base={type(self.base_model).__name__}, " + f"alpha={self.alpha:.2f}, status={status}, score={mad_status})" + ) + + def _check_predict_std(self) -> bool: + return hasattr(self.base_model, "predict_distribution") or hasattr( + self.base_model, "predict_std" + ) + + def _get_base_predictions(self, X: pd.DataFrame) -> np.ndarray: + preds = self.base_model.predict(X) + preds = np.asarray(preds, dtype=np.float64) + if preds.ndim == 2 and preds.shape[1] == 1: + preds = preds.ravel() + return preds + + def _get_base_std(self, X: pd.DataFrame) -> Optional[np.ndarray]: + if hasattr(self.base_model, "predict_distribution"): + _, std = self.base_model.predict_distribution(X) + std = np.asarray(std, dtype=np.float64).ravel() + return std + if hasattr(self.base_model, "predict_std"): + std = self.base_model.predict_std(X) + std = np.asarray(std, dtype=np.float64).ravel() + return std + return None + + def _compute_nonconformity(self, residuals: np.ndarray, std: Optional[np.ndarray] = None) -> np.ndarray: + if self._mad and std is not None and np.all(std > 0): + return np.abs(residuals) / np.maximum(std, 1e-8) + return np.abs(residuals) + + def calibrate(self, X_cal: pd.DataFrame, y_cal: pd.Series, mad: bool = False): + if len(X_cal) == 0 or len(y_cal) == 0: + raise ValueError("Calibration set cannot be empty") + + self._mad = mad + y_true = np.asarray(y_cal, dtype=np.float64).ravel() + y_pred = self._get_base_predictions(X_cal) + residuals = y_true - y_pred + + if self._mad: + std = self._get_base_std(X_cal) + if std is None: + logger.warning( + "MAD mode requested but base_model has no predict_std/predict_distribution. " + "Falling back to absolute residuals." + ) + self._mad = False + std = None + else: + std = None + + self._cal_residuals = self._compute_nonconformity(residuals, std) + self._cal_preds = y_pred + + n = len(self._cal_residuals) + correction = (1.0 + 1.0 / n) + q_idx = min(int(np.ceil((1.0 - self.alpha) * (n + 1))) - 1, n - 1) + if q_idx < 0: + q_idx = 0 + self.q_hat = np.sort(self._cal_residuals)[q_idx] + + self.is_calibrated = True + logger.info( + f"Conformal calibration complete: {n} samples, " + f"alpha={self.alpha:.2f}, quantile={self.q_hat:.4f}" + ) + return self + + def predict_with_band(self, X: pd.DataFrame) -> Tuple[np.ndarray, np.ndarray, np.ndarray]: + if not self.is_calibrated: + raise RuntimeError("Model not calibrated. Call calibrate() first.") + + y_pred = self._get_base_predictions(X) + + if self._mad: + std = self._get_base_std(X) + if std is None: + self._mad = False + std = None + else: + std = None + + if self._mad and std is not None: + half_width = self.q_hat * std + else: + half_width = self.q_hat + + lower = y_pred - half_width + upper = y_pred + half_width + + return y_pred, lower, upper + + def predict(self, X: pd.DataFrame) -> np.ndarray: + return self._get_base_predictions(X) + + def predict_interval_width(self, X: pd.DataFrame) -> float: + _, lower, upper = self.predict_with_band(X) + widths = upper - lower + return float(np.mean(widths)) + + def average_interval_width(self, X: pd.DataFrame) -> float: + widths = self.predict_interval_width(X) + return float(np.mean(widths)) + + def is_inside_band(self, X: pd.DataFrame, y_true: pd.Series) -> np.ndarray: + _, lower, upper = self.predict_with_band(X) + y = np.asarray(y_true, dtype=np.float64).ravel() + return (y >= lower) & (y <= upper) + + def coverage(self, X: pd.DataFrame, y_true: pd.Series) -> float: + inside = self.is_inside_band(X, y_true) + return float(np.mean(inside)) + + def update(self, X_new: pd.DataFrame, y_new: pd.Series): + if not self.is_calibrated: + return self.calibrate(X_new, y_new) + + y_true = np.asarray(y_new, dtype=np.float64).ravel() + y_pred = self._get_base_predictions(X_new) + residuals = y_true - y_pred + + if self._mad: + std = self._get_base_std(X_new) + if std is not None: + new_scores = self._compute_nonconformity(residuals, std) + else: + new_scores = np.abs(residuals) + else: + new_scores = np.abs(residuals) + + self._cal_residuals = np.concatenate([self._cal_residuals, new_scores]) + + n = len(self._cal_residuals) + q_idx = min(int(np.ceil((1.0 - self.alpha) * (n + 1))) - 1, n - 1) + if q_idx < 0: + q_idx = 0 + self.q_hat = np.sort(self._cal_residuals)[q_idx] + + logger.info( + f"Conformal update: +{len(new_scores)} samples, " + f"total={n}, alpha={self.alpha:.2f}, quantile={self.q_hat:.4f}" + ) + return self diff --git a/src/models/gat_model.py b/src/models/gat_model.py new file mode 100644 index 0000000..2487849 --- /dev/null +++ b/src/models/gat_model.py @@ -0,0 +1,647 @@ +"""Graph Attention Network for player chemistry/interaction modeling. + +Models player-to-player synergy and redundancy using a GAT architecture when +PyTorch Geometric is available, falling back to network-science feature +extraction (PageRank, betweenness, clustering, eigenvector centrality) with a +scikit-learn MLPRegressor otherwise. + +Replaces the stub tgcn_model.py. +""" + +import logging +from itertools import combinations +from typing import Dict, List, Optional, Tuple + +import numpy as np +import pandas as pd +from sklearn.neural_network import MLPRegressor +from sklearn.preprocessing import StandardScaler + +from .base_model import BaseModel + +logger = logging.getLogger(__name__) + + +def _check_torch_geometric() -> bool: + try: + import torch # noqa: F401 + from torch_geometric.nn import GATConv # noqa: F401 + return True + except ImportError: + return False + + +def _build_graph_from_edges(edges: Dict[Tuple[str, str], float]) -> Tuple[Dict[str, int], np.ndarray, np.ndarray]: + """Convert edge dict to node mapping, adjacency matrix, and edge_index.""" + nodes = set() + for (a, b) in edges: + nodes.add(a) + nodes.add(b) + node_list = sorted(nodes) + node_idx = {name: i for i, name in enumerate(node_list)} + n = len(node_list) + adj = np.zeros((n, n), dtype=np.float64) + for (a, b), w in edges.items(): + i, j = node_idx[a], node_idx[b] + adj[i, j] = w + edge_list = [(node_idx[a], node_idx[b]) for (a, b) in edges] + edge_index = np.array(edge_list, dtype=np.int64).T if edge_list else np.empty((2, 0), dtype=np.int64) + return node_idx, adj, edge_index + + +def _pagerank(adj: np.ndarray, alpha: float = 0.85, tol: float = 1e-6, max_iter: int = 100) -> np.ndarray: + n = adj.shape[0] + out_deg = adj.sum(axis=1) + out_deg[out_deg == 0] = 1.0 + M = adj.T / out_deg + pr = np.ones(n) / n + for _ in range(max_iter): + pr_new = alpha * M.dot(pr) + (1.0 - alpha) / n + if np.abs(pr_new - pr).sum() < tol: + return pr_new + pr = pr_new + return pr + + +def _eigenvector_centrality(adj: np.ndarray, tol: float = 1e-6, max_iter: int = 200) -> np.ndarray: + n = adj.shape[0] + x = np.ones(n) + for _ in range(max_iter): + x_new = adj.dot(x) + norm = np.linalg.norm(x_new, 2) + if norm < 1e-12: + break + x_new /= norm + if np.abs(x_new - x).sum() < tol: + return x_new + x = x_new + return x / np.linalg.norm(x, 2) if np.linalg.norm(x, 2) > 1e-12 else x + + +def _betweenness_centrality(adj: np.ndarray) -> np.ndarray: + """Brandes algorithm for weighted undirected graphs.""" + n = adj.shape[0] + if n <= 2: + return np.zeros(n) + dist = np.where(adj > 0, 1.0 / np.maximum(adj, 1e-12), np.inf) + np.fill_diagonal(dist, 0.0) + bc = np.zeros(n) + for s in range(n): + S = [] + P = [[] for _ in range(n)] + sigma = np.zeros(n) + sigma[s] = 1.0 + d = np.full(n, np.inf) + d[s] = 0.0 + Q = [s] + for v in Q: + S.append(v) + for w in range(n): + if adj[v, w] <= 0 or w == v: + continue + new_d = d[v] + dist[v, w] + if new_d < d[w] - 1e-12: + d[w] = new_d + sigma[w] = 0.0 + P[w] = [v] + Q.append(w) + elif abs(new_d - d[w]) < 1e-12: + sigma[w] += sigma[v] + P[w].append(v) + delta = np.zeros(n) + while S: + w = S.pop() + for v in P[w]: + delta[v] += (sigma[v] / max(sigma[w], 1e-12)) * (1.0 + delta[w]) + if w != s: + bc[w] += delta[w] + bc /= max(n - 1, 1) * max(n - 2, 1) + return bc + + +def _clustering_coefficient(adj: np.ndarray) -> np.ndarray: + n = adj.shape[0] + cc = np.zeros(n) + binary = (adj > 0).astype(np.float64) + for i in range(n): + neighbors = np.where(binary[i] > 0)[0] + deg = len(neighbors) + if deg < 2: + cc[i] = 0.0 + continue + sub = binary[np.ix_(neighbors, neighbors)] + triangles = np.sum(sub) / 2.0 + cc[i] = 2.0 * triangles / (deg * (deg - 1)) + return cc + + +def _assortativity_by_position(adj: np.ndarray, role_groups: Dict[int, str]) -> float: + """Compute assortativity coefficient with respect to positional role.""" + n = adj.shape[0] + if n <= 1: + return 0.0 + binary = (adj > 0).astype(np.float64) + degrees = binary.sum(axis=1) + total_edges = degrees.sum() + if total_edges == 0: + return 0.0 + same_type = 0.0 + for i in range(n): + for j in range(n): + if binary[i, j] > 0 and role_groups.get(i) == role_groups.get(j): + same_type += 1.0 + p_same = same_type / max(total_edges, 1e-12) + type_counts: Dict[str, float] = {} + for idx in range(n): + role = role_groups.get(idx, "unknown") + type_counts[role] = type_counts.get(role, 0.0) + degrees[idx] + total_deg = sum(type_counts.values()) + expected = sum((t / max(total_deg, 1e-12)) ** 2 for t in type_counts.values()) + max_possible = 1.0 - expected + if abs(max_possible) < 1e-12: + return 0.0 + return (p_same - expected) / max_possible + + +class PlayerChemistryGAT(BaseModel): + """Graph Attention Network for player chemistry/interaction modeling. + + Two modes: + - PyTorch Geometric: full GATConv layers predicting delta-fantavoto. + - Fallback: scikit-learn MLPRegressor trained on network-science features + (PageRank, betweenness, clustering, eigenvector, assortativity). + + Parameters + ---------- + model_dir : str + Directory for persisting trained models. + hidden_dim : int + Hidden dimension for GAT / MLP. + num_layers : int + Number of GATConv layers. + heads : int + Number of attention heads. + dropout : float + Dropout probability. + """ + + def __init__( + self, + model_dir: str = "models_trained", + hidden_dim: int = 64, + num_layers: int = 2, + heads: int = 4, + dropout: float = 0.2, + ): + super().__init__(model_dir) + self.hidden_dim = hidden_dim + self.num_layers = num_layers + self.heads = heads + self.dropout = dropout + + self._edges: Dict[Tuple[str, str], float] = {} + self.interaction_edges: Dict[Tuple[str, str], float] = {} + self._node_idx: Dict[str, int] = {} + self._idx_node: Dict[int, str] = {} + self._adj: Optional[np.ndarray] = None + self._edge_index: Optional[np.ndarray] = None + + self._graph_features: Optional[np.ndarray] = None + self._feature_names: Optional[List[str]] = None + self._scaler = StandardScaler() + + self._use_torch: bool = False + self._gat_module = None + self._mlp: Optional[MLPRegressor] = None + + # ------------------------------------------------------------------ + # Graph construction + # ------------------------------------------------------------------ + def build_graph(self, historical_data: pd.DataFrame): + """Build player interaction graph from pass/assist/cross data. + + Args: + historical_data: DataFrame with columns + [player, teammate, passes_to, assists_to, crosses_to, matchday] + """ + required = {"player", "teammate", "passes_to", "assists_to", "crosses_to"} + missing = required - set(historical_data.columns) + if missing: + raise ValueError(f"Missing required columns: {missing}") + + edges: Dict[Tuple[str, str], float] = {} + for (player, teammate), group in historical_data.groupby(["player", "teammate"]): + weight = ( + group["passes_to"].sum() + + group["assists_to"].sum() * 5 + + group["crosses_to"].sum() * 3 + ) + if weight > 0: + edges[(str(player), str(teammate))] = float(weight) + + self._edges = edges + self.interaction_edges = dict(edges) + self._node_idx, self._adj, self._edge_index = _build_graph_from_edges(edges) + self._idx_node = {v: k for k, v in self._node_idx.items()} + logger.info( + "Built interaction graph: %d nodes, %d edges", + len(self._node_idx), + len(edges), + ) + return self + + # ------------------------------------------------------------------ + # Graph feature extraction + # ------------------------------------------------------------------ + def _extract_graph_features(self, role_map: Optional[Dict[str, str]] = None) -> np.ndarray: + n = len(self._node_idx) + if n == 0: + return np.empty((0, 6)) + + adj = self._adj.copy() if self._adj is not None else np.zeros((n, n)) + + pr = _pagerank(adj) + bc = _betweenness_centrality(adj) + cc = _clustering_coefficient(adj) + deg = adj.sum(axis=1) + adj.sum(axis=0) + ec = _eigenvector_centrality(adj) + + role_groups: Dict[int, str] = {} + if role_map: + for i in range(n): + name = self._idx_node.get(i, "") + role_groups[i] = role_map.get(name, "unknown") + assortative = _assortativity_by_position(adj, role_groups) if role_map else 0.0 + + features = np.column_stack([pr, bc, cc, deg, ec, np.full(n, assortative)]) + self._feature_names = [ + "pagerank", + "betweenness", + "clustering_coef", + "degree", + "eigenvector", + "assortativity", + ] + return features + + # ------------------------------------------------------------------ + # PyTorch GAT module + # ------------------------------------------------------------------ + def _init_torch_gat(self, in_channels: int): + try: + import torch + import torch.nn as nn + from torch_geometric.nn import GATConv + + class GATModule(nn.Module): + def __init__( + self, + in_channels: int, + hidden_dim: int, + out_dim: int, + heads: int, + dropout: float, + ): + super().__init__() + self.conv1 = GATConv(in_channels, hidden_dim, heads=heads, dropout=dropout) + self.conv2 = GATConv(hidden_dim * heads, out_dim, heads=1, concat=False, dropout=dropout) + self.head = nn.Linear(out_dim, 1) + + def forward(self, x, edge_index, edge_attr=None): + x = self.conv1(x, edge_index) + x = torch.relu(x) + x = self.conv2(x, edge_index) + x = torch.relu(x) + return self.head(x).squeeze(-1) + + self._gat_module = GATModule( + in_channels=in_channels, + hidden_dim=self.hidden_dim, + out_dim=self.hidden_dim // 2, + heads=self.heads, + dropout=self.dropout, + ) + self._use_torch = True + logger.info("PyTorch Geometric GAT initialized (in=%d)", in_channels) + except ImportError: + self._use_torch = False + logger.info("torch_geometric not installed; using MLP fallback") + + # ------------------------------------------------------------------ + # Fit + # ------------------------------------------------------------------ + def fit(self, X: pd.DataFrame, y: pd.Series, **kwargs): + """Fit the player chemistry model. + + Args: + X: DataFrame. Must contain a 'player' column for node identification. + Optionally 'team' for team-based subgraph grouping and 'role' for + positional assortativity. + y: Target fantavoto scores per row. + """ + if X.empty: + raise ValueError("X cannot be empty") + + role_map: Optional[Dict[str, str]] = None + if "role" in X.columns: + role_map = dict(zip(X["player"].astype(str), X["role"].astype(str))) + + team_groups = [] + if "team" in X.columns: + for _, group in X.groupby("team"): + players = group["player"].unique().tolist() + if len(players) >= 2: + for a, b in combinations(players, 2): + team_groups.append({"player": str(a), "teammate": str(b), + "passes_to": 1, "assists_to": 0, "crosses_to": 0}) + + if "player" not in X.columns: + raise ValueError("X must contain a 'player' column") + + if not self._edges and not team_groups: + if "player" in X.columns and "team" not in X.columns: + logger.warning("No interaction edges and no 'team' column; graph will be empty") + self._edges = {} + self.interaction_edges = {} + self._node_idx = {} + self._adj = np.array([[0.0]]) + self._edge_index = np.empty((2, 0), dtype=np.int64) + + if team_groups and not self._edges: + gdf = pd.DataFrame(team_groups) + if "passes_to" not in gdf.columns: + gdf["passes_to"] = 1 + if "assists_to" not in gdf.columns: + gdf["assists_to"] = 0 + if "crosses_to" not in gdf.columns: + gdf["crosses_to"] = 0 + self.build_graph(gdf) + + if "player" in X.columns: + player_set = set(self._node_idx.keys()) + new_nodes = [] + for row_player in X["player"]: + p = str(row_player) + if p and p not in player_set: + idx = len(self._node_idx) + self._node_idx[p] = idx + self._idx_node[idx] = p + player_set.add(p) + new_nodes.append(p) + if self._adj is None: + n = len(self._node_idx) + self._adj = np.zeros((n, n)) + self._edge_index = np.empty((2, 0), dtype=np.int64) + elif new_nodes: + old_n = self._adj.shape[0] + new_n = len(self._node_idx) + adj_new = np.zeros((new_n, new_n)) + adj_new[:old_n, :old_n] = self._adj + self._adj = adj_new + + self._graph_features = self._extract_graph_features(role_map) + n_nodes = len(self._node_idx) + + if n_nodes == 0: + logger.warning("Zero nodes in graph; fitting dummy model") + self._mlp = MLPRegressor( + hidden_layer_sizes=(self.hidden_dim,), + max_iter=300, + random_state=42, + ) + self._mlp.fit(np.zeros((1, 6)), np.zeros(1)) + return self + + node_to_player = {v: k for k, v in self._node_idx.items()} + target_deltas = np.zeros(n_nodes) + count_deltas = np.zeros(n_nodes) + player_indices = { + str(row["player"]): i + for i, (_, row) in enumerate(X.iterrows()) + } + y_mean = float(np.mean(y)) if len(y) > 0 else 6.0 + + for i, (_, row) in enumerate(X.iterrows()): + p = str(row["player"]) + if p in self._node_idx: + idx = self._node_idx[p] + target_deltas[idx] += (float(y.iloc[i]) - y_mean) + count_deltas[idx] += 1.0 + + for j in range(n_nodes): + if count_deltas[j] > 0: + target_deltas[j] /= count_deltas[j] + + if _check_torch_geometric(): + self._init_torch_gat(in_channels=self._graph_features.shape[1]) + self._fit_torch_gat(self._graph_features, self._edge_index, target_deltas) + else: + self._fit_sklearn_fallback(self._graph_features, target_deltas) + + logger.info( + "PlayerChemistryGAT fitted: %d nodes, %d edges, mode=%s", + n_nodes, len(self._edges), + "torch" if self._use_torch else "sklearn", + ) + return self + + def _fit_torch_gat( + self, + features: np.ndarray, + edge_index: np.ndarray, + targets: np.ndarray, + ): + import torch + import torch.nn as nn + import torch.optim as optim + + device = torch.device("cuda" if torch.cuda.is_available() else "cpu") + x_tensor = torch.tensor(features, dtype=torch.float32).to(device) + edge_tensor = torch.tensor(edge_index, dtype=torch.long).to(device) + y_tensor = torch.tensor(targets, dtype=torch.float32).to(device) + + self._gat_module = self._gat_module.to(device) + optimizer = optim.Adam(self._gat_module.parameters(), lr=0.01, weight_decay=1e-5) + criterion = nn.MSELoss() + + n = features.shape[0] + self._gat_module.train() + for epoch in range(200): + optimizer.zero_grad() + pred = self._gat_module(x_tensor, edge_tensor) + loss = criterion(pred, y_tensor) + loss.backward() + optimizer.step() + if (epoch + 1) % 50 == 0: + logger.debug(f"GAT epoch {epoch + 1}: loss={loss.item():.6f}") + + self._gat_module.eval() + + def _fit_sklearn_fallback(self, features: np.ndarray, targets: np.ndarray): + self._scaler.fit(features) + features_scaled = self._scaler.transform(features) + n_samples = features_scaled.shape[0] + val_frac = 0.1 if n_samples >= 20 else 0.0 + self._mlp = MLPRegressor( + hidden_layer_sizes=(self.hidden_dim, self.hidden_dim // 2), + activation="relu", + solver="adam", + max_iter=500, + random_state=42, + early_stopping=(n_samples >= 20), + validation_fraction=val_frac if val_frac > 0 else 0.1, + n_iter_no_change=20, + ) + self._mlp.fit(features_scaled, targets) + self._use_torch = False + + # ------------------------------------------------------------------ + # Predict + # ------------------------------------------------------------------ + def predict(self, X: pd.DataFrame) -> np.ndarray: + """Return chemistry-adjusted point projection bonuses. + + These bonuses can be positive (synergy) or negative (redundancy). + + Returns an array the same length as X with delta-fantavoto values. + """ + if self._adj is None or len(self._node_idx) == 0: + return np.zeros(X.shape[0]) + + role_map: Optional[Dict[str, str]] = None + if "role" in X.columns: + role_map = dict(zip(X["player"].astype(str), X["role"].astype(str))) + + features = self._extract_graph_features(role_map) + n_nodes = len(self._node_idx) + if features.shape[0] != n_nodes or n_nodes == 0: + return np.zeros(X.shape[0]) + + if self._use_torch and self._gat_module is not None: + return self._predict_torch(features, X) + elif self._mlp is not None: + return self._predict_sklearn(features, X) + return np.zeros(X.shape[0]) + + def _predict_torch(self, features: np.ndarray, X: pd.DataFrame) -> np.ndarray: + import torch + + device = next(self._gat_module.parameters()).device + x_tensor = torch.tensor(features, dtype=torch.float32).to(device) + edge_tensor = torch.tensor(self._edge_index, dtype=torch.long).to(device) + + self._gat_module.eval() + with torch.no_grad(): + node_deltas = self._gat_module(x_tensor, edge_tensor).cpu().numpy() + + result = np.zeros(X.shape[0]) + for i, (_, row) in enumerate(X.iterrows()): + p = str(row.get("player", row.get("name", ""))) + if p in self._node_idx: + result[i] = float(node_deltas[self._node_idx[p]]) + return result + + def _predict_sklearn(self, features: np.ndarray, X: pd.DataFrame) -> np.ndarray: + features_scaled = self._scaler.transform(features) + node_deltas = self._mlp.predict(features_scaled) + + result = np.zeros(X.shape[0]) + for i, (_, row) in enumerate(X.iterrows()): + p = str(row.get("player", row.get("name", ""))) + if p in self._node_idx: + result[i] = float(node_deltas[self._node_idx[p]]) + return result + + # ------------------------------------------------------------------ + # Interaction feature extraction (public API) + # ------------------------------------------------------------------ + def extract_interaction_features(self, player: str, teammates: List[str]) -> Dict[str, float]: + """Extract interaction features between a player and their teammates. + + Args: + player: Player name. + teammates: List of teammate names. + + Returns: + Dict with keys: interaction_outgoing_sum, interaction_incoming_sum, + interaction_synergy. + """ + outgoing = 0.0 + incoming = 0.0 + synergy = 0.0 + count = 0 + for teammate in teammates: + w = self._edges.get((str(player), str(teammate)), 0.0) + outgoing += w + w_in = self._edges.get((str(teammate), str(player)), 0.0) + incoming += w_in + if w > 0 or w_in > 0: + synergy += (w + w_in) / 2.0 + count += 1 + avg_synergy = synergy / max(count, 1) + return { + "interaction_outgoing_sum": float(outgoing), + "interaction_incoming_sum": float(incoming), + "interaction_synergy": float(avg_synergy), + } + + # ------------------------------------------------------------------ + # Interaction bonus / chemistry matrix + # ------------------------------------------------------------------ + def compute_interaction_bonus(self, player_a: str, player_b: str) -> float: + """Bonus factor for two players in the same lineup. + + Positive = synergy, negative = redundancy. + """ + if not self._edges: + return 0.0 + a = str(player_a) + b = str(player_b) + weight = self._edges.get((a, b), 0.0) + self._edges.get((b, a), 0.0) + max_weight = max(self._edges.values()) if self._edges else 1.0 + return float(weight) / max(max_weight, 1.0) if weight > 0 else 0.0 + + def get_chemistry_matrix(self, players: List[str]) -> np.ndarray: + """Return NxN matrix of pairwise chemistry bonuses for given players. + + Positive values = synergy, negative = redundancy. + """ + n = len(players) + matrix = np.zeros((n, n)) + for i in range(n): + for j in range(n): + if i != j: + matrix[i, j] = self.compute_interaction_bonus(players[i], players[j]) + return matrix + + def get_redundancy_penalty(self, players: List[str]) -> Dict[Tuple[str, str], float]: + """Identify negative synergies between players on the same team. + + Two players may compete for the same actions (both take corners, + both demand the ball in similar zones). Returns dict mapping + player-pairs to a negative penalty value (0 = no redundancy). + + The penalty is derived from co-occurrence overlap: if two players + both have high outgoing edges to the same teammates and both + receive from the same sources, they're likely redundant. + """ + penalty: Dict[Tuple[str, str], float] = {} + if not self._edges: + return penalty + + for i, a in enumerate(players): + for b in players[i + 1:]: + a_out = set(t for (p, t) in self._edges if p == str(a)) + b_out = set(t for (p, t) in self._edges if p == str(b)) + a_in = set(p for (p, t) in self._edges if t == str(a)) + b_in = set(p for (p, t) in self._edges if t == str(b)) + + out_overlap = len(a_out & b_out) + in_overlap = len(a_in & b_in) + total = max(len(a_out | b_out) + len(a_in | b_in), 1) + overlap_ratio = (out_overlap + in_overlap) / total + + if overlap_ratio > 0.3: + penalty[(str(a), str(b))] = -overlap_ratio + + return penalty diff --git a/src/models/hawkes_form.py b/src/models/hawkes_form.py new file mode 100644 index 0000000..b8f9659 --- /dev/null +++ b/src/models/hawkes_form.py @@ -0,0 +1,503 @@ +"""Hawkes-process player form model for momentum modeling. + +Models player form as a self-exciting Hawkes process: good performances +increase the probability of more good performances (momentum / hot streak). + +Uses scipy.optimize.minimize to fit per-player Hawkes parameters +(mu, alpha, beta, s) via maximum likelihood. +""" + +import logging +from dataclasses import dataclass +from typing import Dict, List, Optional, Tuple + +import numpy as np +import pandas as pd +from scipy.optimize import minimize + +from .base_model import BaseModel + +logger = logging.getLogger(__name__) + +FORM_STATUS_HOT = "HOT" +FORM_STATUS_COLD = "COLD" +FORM_STATUS_NEUTRAL = "NEUTRAL" +FORM_STATUSES = {FORM_STATUS_HOT, FORM_STATUS_COLD, FORM_STATUS_NEUTRAL} + +_EPS = 1e-12 +_MIN_BETA = 1e-4 +_MAX_ALPHA = 20.0 +_MAX_MU = 50.0 +_MIN_S = -3.0 +_MAX_S = 3.0 + + +@dataclass +class PlayerHawkesParams: + mu: float + alpha: float + beta: float + s: float + baseline_mean: float + baseline_std: float + + +class PlayerFormModel(BaseModel): + """Self-exciting Hawkes process for player form / momentum. + + Each player's match performances are modeled as a point process where + above-average games ("excitatory events") temporarily raise the + probability of subsequent above-average games. + + Parameters + ---------- + model_dir : str + Directory for persisting trained models. + decay_window : int + Maximum match-gap over which excitation persists (default 10). + """ + + def __init__( + self, + model_dir: str = "models_trained", + decay_window: int = 10, + ): + super().__init__(model_dir) + self.decay_window = int(decay_window) + self._player_params: Dict[str, PlayerHawkesParams] = {} + self._global_baseline: float = 6.0 + self._global_std: float = 1.5 + + # ------------------------------------------------------------------ + # Hawkes log-likelihood and intensity + # ------------------------------------------------------------------ + @staticmethod + def _hawkes_intensity(times: np.ndarray, event_mask: np.ndarray, mu: float, + alpha: float, beta: float) -> np.ndarray: + lam = np.full_like(times, mu, dtype=np.float64) + for i in range(1, len(times)): + if event_mask[i - 1]: + dt = times[i:] - times[i - 1] + mask = dt > 0 + lam[i:] += alpha * np.exp(-beta * dt) * mask + return np.maximum(lam, _EPS) + + @staticmethod + def _hawkes_integral(mu: float, alpha: float, beta: float, + event_times: np.ndarray, T: float) -> float: + result = mu * T + for te in event_times: + remaining = T - te + if remaining > 0: + result += (alpha / beta) * (1.0 - np.exp(-beta * remaining)) + return result + + @staticmethod + def _hawkes_nll(params: np.ndarray, times: np.ndarray, event_mask: np.ndarray) -> float: + mu, alpha, beta = max(params[0], _EPS), max(params[1], _EPS), max(params[2], _MIN_BETA) + T = times[-1] if len(times) > 0 else 1.0 + lam = PlayerFormModel._hawkes_intensity(times, event_mask, mu, alpha, beta) + log_lik = np.sum(np.log(lam)) + integral = PlayerFormModel._hawkes_integral(mu, alpha, beta, + times[event_mask.astype(bool)], T) + return -(log_lik - integral) + + # ------------------------------------------------------------------ + # Fit helper: tune per-player + # ------------------------------------------------------------------ + def _fit_player(self, times: np.ndarray, scores: np.ndarray) -> Optional[PlayerHawkesParams]: + n = len(scores) + if n < 5: + return None + + times_float = times.astype(np.float64) + scores_float = scores.astype(np.float64) + baseline_mean = float(np.mean(scores_float)) + baseline_std = float(np.std(scores_float, ddof=1)) if n > 1 else 1.0 + + best_nll = float("inf") + best_params = None + + for s_candidate in [-1.0, 0.0, 0.5, 1.0, 1.5]: + threshold = baseline_mean + s_candidate * max(baseline_std, 0.5) + event_mask = (scores_float > threshold).astype(np.float64) + n_events = event_mask.sum() + if n_events < 2: + continue + + init_mu = max(max(n_events / max(times_float[-1] - times_float[0], 1.0), 0.05), _EPS) + init_alpha = min(n_events / max(n, 1) * 2.0, _MAX_ALPHA) + init_beta = 0.5 + + for init_scale in [0.5, 1.0, 2.0]: + x0 = np.array([ + init_mu * init_scale, + init_alpha * init_scale, + init_beta * init_scale, + ]) + + try: + result = minimize( + self._hawkes_nll, + x0, + args=(times_float, event_mask), + method="L-BFGS-B", + bounds=[(_EPS, _MAX_MU), (_EPS, _MAX_ALPHA), (_MIN_BETA, 10.0)], + options={"maxiter": 500, "ftol": 1e-10}, + ) + if result.success and result.fun < best_nll: + best_nll = result.fun + best_params = PlayerHawkesParams( + mu=float(max(result.x[0], _EPS)), + alpha=float(max(result.x[1], _EPS)), + beta=float(max(result.x[2], _MIN_BETA)), + s=float(s_candidate), + baseline_mean=baseline_mean, + baseline_std=baseline_std, + ) + except Exception: + continue + + if best_params is None: + best_params = PlayerHawkesParams( + mu=0.1, + alpha=1.0, + beta=0.3, + s=0.0, + baseline_mean=baseline_mean, + baseline_std=baseline_std, + ) + + return best_params + + # ------------------------------------------------------------------ + # Fit + # ------------------------------------------------------------------ + def fit(self, X: pd.DataFrame, y: pd.Series, **kwargs): + """Fit per-player Hawkes parameters. + + X must contain: + - 'player' or 'name': player identifier. + - 'match_date' or 'matchday': temporal ordering column. + + y: fantavoto scores. + + Additional kwargs: + - 'match_date' column name override. + """ + if X.empty: + raise ValueError("X cannot be empty") + + player_col = None + for candidate in ["player", "name"]: + if candidate in X.columns: + player_col = candidate + break + if player_col is None: + raise ValueError("X must contain a 'player' or 'name' column") + + date_col = kwargs.get("date_col", None) + if date_col is None: + for candidate in ["match_date", "matchday", "date", "giornata"]: + if candidate in X.columns: + date_col = candidate + break + if date_col is None: + logger.warning("No date column found; using row index as temporal order") + times = np.arange(len(X), dtype=np.float64) + else: + col_vals = X[date_col] + if pd.api.types.is_datetime64_any_dtype(col_vals): + times = col_vals.astype(np.int64).values.astype(np.float64) / 1e9 / 86400.0 + else: + times = col_vals.astype(np.float64).values + + players = X[player_col].astype(str).values + scores = y.values.astype(np.float64) + + self._global_baseline = float(np.mean(scores)) if len(scores) > 0 else 6.0 + self._global_std = float(np.std(scores, ddof=1)) if len(scores) > 1 else 1.5 + + self._player_params = {} + unique_players = np.unique(players) + fitted = 0 + for player in unique_players: + mask = players == player + p_times = times[mask] + p_scores = scores[mask] + sort_idx = np.argsort(p_times) + p_times = p_times[sort_idx] + p_scores = p_scores[sort_idx] + params = self._fit_player(p_times, p_scores) + if params is not None: + self._player_params[str(player)] = params + fitted += 1 + + logger.info( + "Hawkes form model fitted: %d/%d players with sufficient history", + fitted, len(unique_players), + ) + return self + + # ------------------------------------------------------------------ + # Predict + # ------------------------------------------------------------------ + def predict(self, X: pd.DataFrame) -> np.ndarray: + """Return form-adjusted projections as additive bonuses to base. + + Positive = player is in form (HOT), negative = out of form (COLD). + """ + if not self._player_params: + return np.zeros(X.shape[0]) + + player_col = "player" if "player" in X.columns else "name" + multipliers = self._compute_multipliers(X) + base = X.get("base_prediction", pd.Series(np.full(X.shape[0], 6.0))) + base_vals = base.values.astype(np.float64) + return (multipliers - 1.0) * base_vals + + def _compute_multipliers(self, X: pd.DataFrame) -> np.ndarray: + player_col = "player" if "player" in X.columns else "name" + multipliers = np.ones(X.shape[0], dtype=np.float64) + players = X[player_col].astype(str).values + + for i, player in enumerate(players): + params = self._player_params.get(player) + if params is None: + continue + recent = self._compute_recent_intensity(params) + base_rate = params.mu + if base_rate > _EPS: + ratio = recent / base_rate + clamped = np.clip(ratio, 0.85, 1.15) + multipliers[i] = float(clamped) + + return multipliers + + def _compute_recent_intensity(self, params: PlayerHawkesParams) -> float: + return max(params.mu, _EPS) + + # ------------------------------------------------------------------ + # Momentum projection + # ------------------------------------------------------------------ + def predict_momentum( + self, + X: pd.DataFrame, + player_history: pd.DataFrame, + n_future: int = 5, + ) -> np.ndarray: + """Project form trajectory for next *n_future* matches. + + Returns (n_future, n_players) array of momentum multipliers. + """ + if not self._player_params: + return np.ones((n_future, X.shape[0])) + + player_col = "player" if "player" in X.columns else "name" + players = X[player_col].astype(str).values + n_players = X.shape[0] + trajectory = np.ones((n_future, n_players), dtype=np.float64) + + history_player_col = None + for c in ["player", "name"]: + if c in player_history.columns: + history_player_col = c + break + + history_date_col = None + for c in ["match_date", "matchday", "date"]: + if c in player_history.columns: + history_date_col = c + break + + for j, player in enumerate(players): + params = self._player_params.get(player) + if params is None: + continue + + event_times = [] + if history_player_col and history_date_col: + p_hist = player_history[player_history[history_player_col].astype(str) == player] + if len(p_hist) > 0: + target_col = None + for c in ["fantavoto", "score", "fv"]: + if c in p_hist.columns: + target_col = c + break + if target_col and params.baseline_std > 0: + threshold = params.baseline_mean + params.s * params.baseline_std + p_sorted = p_hist.sort_values(history_date_col) + times = p_sorted[history_date_col].values + scores = p_sorted[target_col].values + if pd.api.types.is_datetime64_any_dtype(p_sorted[history_date_col]): + event_times_float = times.astype(np.int64).astype(np.float64) / 1e9 / 86400.0 + else: + event_times_float = times.astype(np.float64) + for ti, si in zip(event_times_float, scores): + if float(si) > threshold: + event_times.append(ti) + + if not event_times: + continue + + last_t = max(event_times) + for k in range(1, n_future + 1): + future_t = last_t + k + lam = params.mu + for te in event_times: + dt = future_t - te + if dt > 0 and dt <= self.decay_window: + lam += params.alpha * np.exp(-params.beta * dt) + trajectory[k - 1, j] = float(np.clip(lam / max(params.mu, _EPS), 0.85, 1.15)) + + return trajectory + + # ------------------------------------------------------------------ + # Form status + # ------------------------------------------------------------------ + def get_form_status(self, X: pd.DataFrame) -> List[str]: + """Return status string per player: HOT, COLD, or NEUTRAL.""" + player_col = "player" if "player" in X.columns else "name" + players = X[player_col].astype(str).values + statuses: List[str] = [] + + for player in players: + params = self._player_params.get(player) + if params is None: + statuses.append(FORM_STATUS_NEUTRAL) + continue + intensity = self._compute_recent_intensity(params) + baseline = max(params.mu, _EPS) + ratio = intensity / baseline + if ratio > 1.1 and params.alpha > 0.1: + statuses.append(FORM_STATUS_HOT) + elif ratio < 0.9: + statuses.append(FORM_STATUS_COLD) + else: + statuses.append(FORM_STATUS_NEUTRAL) + + return statuses + + # ------------------------------------------------------------------ + # Intensity curve + # ------------------------------------------------------------------ + def compute_intensity_curve( + self, + player_name: str, + history: pd.DataFrame, + match_dates: np.ndarray, + future_dates: np.ndarray, + ) -> np.ndarray: + """Compute λ(t) over match_dates and future_dates for one player. + + Returns an array of intensity values at each date. + """ + params = self._player_params.get(str(player_name)) + if params is None: + mu = self._global_baseline + all_dates = np.concatenate([match_dates, future_dates]) + return np.full_like(all_dates, max(mu, _EPS), dtype=np.float64) + + target_col = None + for c in ["fantavoto", "score", "fv"]: + if c in history.columns: + target_col = c + break + + threshold = params.baseline_mean + params.s * max(params.baseline_std, 0.5) + event_times = [] + if target_col: + for _, row in history.iterrows(): + if float(row.get(target_col, 0)) > threshold: + event_times.append(float(row.name) if isinstance(row.name, (int, float)) else 0.0) + + if match_dates is not None and len(match_dates) > 0: + match_vals = match_dates.astype(np.float64) + for ti in match_vals: + if ti not in event_times: + score = None + for _, row in history.iterrows(): + d_val = float(row.name) if isinstance(row.name, (int, float)) else 0.0 + if abs(d_val - ti) < _EPS: + score = row.get(target_col, 0) if target_col else 0 + break + if score is not None and float(score) > threshold: + event_times.append(ti) + + all_dates = np.concatenate([ + match_dates.astype(np.float64) if match_dates is not None and len(match_dates) > 0 + else np.array([], dtype=np.float64), + future_dates.astype(np.float64) if future_dates is not None and len(future_dates) > 0 + else np.array([], dtype=np.float64), + ]) + + if len(all_dates) == 0: + return np.array([params.mu]) + + intensity = np.full(len(all_dates), params.mu, dtype=np.float64) + for i, t in enumerate(all_dates): + lam = params.mu + for te in event_times: + dt = t - te + if dt > 0 and dt <= self.decay_window: + lam += params.alpha * np.exp(-params.beta * dt) + elif dt > self.decay_window: + pass + intensity[i] = max(lam, _EPS) + + return intensity + + # ------------------------------------------------------------------ + # Streak detection + # ------------------------------------------------------------------ + def detect_streak( + self, + player_history: pd.DataFrame, + ) -> Tuple[bool, int, str]: + """Detect whether a player is on a hot or cold streak. + + Returns (is_streak, streak_length, streak_direction). + streak_direction is "HOT_STREAK" or "COLD_STREAK". + """ + if len(player_history) < 3: + return (False, 0, "NO_STREAK") + + target_col = None + for c in ["fantavoto", "score", "fv"]: + if c in player_history.columns: + target_col = c + break + if target_col is None: + return (False, 0, "NO_STREAK") + + date_col = None + for c in ["match_date", "matchday", "date"]: + if c in player_history.columns: + date_col = c + break + if date_col: + sorted_hist = player_history.sort_values(date_col) + else: + sorted_hist = player_history + + scores = sorted_hist[target_col].values.astype(np.float64) + mean_score = np.mean(scores) + std_score = max(np.std(scores, ddof=1), 0.5) + + above = scores[-1] > mean_score + 0.5 * std_score + below = scores[-1] < mean_score - 0.5 * std_score + + if not above and not below: + return (False, 0, "NO_STREAK") + + direction = "HOT_STREAK" if above else "COLD_STREAK" + streak_len = 1 + for j in range(len(scores) - 2, -1, -1): + if direction == "HOT_STREAK" and scores[j] > mean_score + 0.5 * std_score: + streak_len += 1 + elif direction == "COLD_STREAK" and scores[j] < mean_score - 0.5 * std_score: + streak_len += 1 + else: + break + + return (streak_len >= 3, streak_len, direction) diff --git a/src/models/quantile_model.py b/src/models/quantile_model.py new file mode 100644 index 0000000..12638a0 --- /dev/null +++ b/src/models/quantile_model.py @@ -0,0 +1,207 @@ +"""Quantile regression ensemble for probabilistic score prediction. + +Trains one LightGBM quantile regressor per target quantile (default P10, P50, P90) +to output a full predictive distribution of Fantavoto scores. Supports downside risk, +upside potential, and Value-at-Risk-safe estimates. +""" + +import logging +from typing import Optional, Tuple + +import numpy as np +import pandas as pd +from sklearn.preprocessing import StandardScaler + +from .base_model import BaseModel + +logger = logging.getLogger(__name__) + + +class QuantileEnsemble(BaseModel): + """Quantile regression ensemble for multi-quantile score prediction. + + Stores one LGBMRegressor per quantile, each with ``objective="quantile"`` + and the corresponding ``alpha`` value. Provides convenience accessors for + point predictions (P50), downside/upside risk, and VaR-safe floors. + """ + + def __init__( + self, + quantiles: Tuple[float, ...] = (0.10, 0.50, 0.90), + model_dir: str = "models_trained", + n_estimators: int = 300, + learning_rate: float = 0.05, + max_depth: int = 5, + num_leaves: int = 31, + subsample: float = 0.8, + colsample_bytree: float = 0.8, + ): + super().__init__(model_dir) + self.quantiles = tuple(quantiles) + self.n_estimators = n_estimators + self.learning_rate = learning_rate + self.max_depth = max_depth + self.num_leaves = num_leaves + self.subsample = subsample + self.colsample_bytree = colsample_bytree + + self._models: dict = {} # alpha → LGBMRegressor + self.feature_names = None + self.scaler = StandardScaler() + + # ------------------------------------------------------------------ + # Training + # ------------------------------------------------------------------ + def fit(self, X: pd.DataFrame, y: pd.Series, **kwargs): + """Fit one LightGBM quantile regressor per target quantile. + + Args: + X: Feature matrix. + y: Target values (fantavoto scores). + """ + try: + import lightgbm as lgb + except ImportError: + raise ImportError("LightGBM is required. pip install lightgbm") + + self.feature_names = list(X.columns) + X_clean = X.select_dtypes(include=[np.number]).fillna(0) + X_scaled = self.scaler.fit_transform(X_clean) + + self._models = {} + for alpha in self.quantiles: + model = lgb.LGBMRegressor( + n_estimators=self.n_estimators, + learning_rate=self.learning_rate, + max_depth=self.max_depth, + num_leaves=self.num_leaves, + subsample=self.subsample, + colsample_bytree=self.colsample_bytree, + objective="quantile", + alpha=alpha, + random_state=42, + verbose=-1, + ) + model.fit(X_scaled, y) + label = f"P{int(alpha * 100)}" + self._models[label] = model + logger.info(f"Quantile model {label} (α={alpha:.2f}) trained") + + logger.info(f"QuantileEnsemble fitted: {list(self._models.keys())}") + return self + + # ------------------------------------------------------------------ + # Prediction helpers + # ------------------------------------------------------------------ + def _preprocess(self, X: pd.DataFrame) -> np.ndarray: + """Clean, impute, and scale features.""" + if self.feature_names is None: + raise RuntimeError("Model not trained. Call fit() first.") + X_c = X[self.feature_names].select_dtypes(include=[np.number]).fillna(0) + return self.scaler.transform(X_c) + + # ------------------------------------------------------------------ + # Core predict + # ------------------------------------------------------------------ + def predict(self, X: pd.DataFrame) -> dict: + """Return a dict mapping quantile label → prediction array. + + Example: + {"P10": array([4.2, 5.1, ...]), "P50": array([5.8, ...]), ...} + """ + if not self._models: + raise RuntimeError("Model not trained. Call fit() first.") + X_scaled = self._preprocess(X) + result = {} + for label, model in self._models.items(): + result[label] = model.predict(X_scaled) + return result + + def predict_points(self, X: pd.DataFrame) -> np.ndarray: + """Convenience: return the P50 (median) prediction array.""" + preds = self.predict(X) + p50_key = "P50" + if p50_key not in preds: + available = min(preds.keys(), key=lambda k: abs(float(k[1:]) / 100 - 0.50)) + logger.warning(f"P50 not trained; falling back to {available}") + return preds[available] + return preds[p50_key] + + # ------------------------------------------------------------------ + # Risk / reward utilities + # ------------------------------------------------------------------ + def predict_downside_risk( + self, X: pd.DataFrame, threshold: float = 5.5 + ) -> np.ndarray: + """Probability that player score falls below *threshold*. + + Uses CDF interpolation across the trained quantiles. The returned + probability is the fraction of the predictive distribution that lies + below the threshold. + """ + preds = self.predict(X) + n = len(next(iter(preds.values()))) + probs = np.zeros(n) + + for i in range(n): + q_vals = [preds[label][i] for label in sorted(preds.keys())] + q_levels = sorted([float(k[1:]) / 100 for k in preds.keys()]) + + if threshold <= q_vals[0]: + probs[i] = q_levels[0] + elif threshold >= q_vals[-1]: + probs[i] = q_levels[-1] + else: + idx = np.searchsorted(q_vals, threshold) + lo_q, hi_q = q_vals[idx - 1], q_vals[idx] + lo_level, hi_level = q_levels[idx - 1], q_levels[idx] + frac = (threshold - lo_q) / (hi_q - lo_q + 1e-10) + probs[i] = lo_level + frac * (hi_level - lo_level) + + return np.clip(probs, 0.0, 1.0) + + def predict_upside(self, X: pd.DataFrame, threshold: float = 7.0) -> np.ndarray: + """Probability that player score exceeds *threshold*.""" + downside = self.predict_downside_risk(X, threshold) + return 1.0 - downside + + def value_at_risk_safe( + self, X: pd.DataFrame, confidence: float = 0.90 + ) -> np.ndarray: + """VaR-safe estimate: floor that the player exceeds with given confidence. + + For confidence=0.90 returns the score floor that the player exceeds 90% + of the time — i.e. the P(100-confidence) quantile. A higher VaR-safe + means more reliable upside. + """ + alpha = 1.0 - confidence + preds = self.predict(X) + levels = np.array(sorted([float(k[1:]) / 100 for k in preds.keys()])) + + # If the exact alpha was trained return it directly + atol = 0.005 + for label, q_pred in preds.items(): + if abs(float(label[1:]) / 100 - alpha) <= atol: + return np.asarray(q_pred) + + # Otherwise linearly interpolate + idx = np.searchsorted(levels, alpha) + if idx == 0: + label = sorted(preds.keys())[0] + return np.asarray(preds[label]) + if idx >= len(levels): + label = sorted(preds.keys())[-1] + return np.asarray(preds[label]) + + lo_level = levels[idx - 1] + hi_level = levels[idx] + lo_label = f"P{int(round(lo_level * 100))}" + hi_label = f"P{int(round(hi_level * 100))}" + # Fall back to closest trained quantile keys + lo_label = min(preds.keys(), key=lambda k: abs(float(k[1:]) / 100 - lo_level)) + hi_label = min(preds.keys(), key=lambda k: abs(float(k[1:]) / 100 - hi_level)) + + lo_vals = preds[lo_label] + hi_vals = preds[hi_label] + frac = (alpha - lo_level) / (hi_level - lo_level + 1e-10) + return lo_vals + frac * (hi_vals - lo_vals) diff --git a/src/models/set_transformer.py b/src/models/set_transformer.py index 451234a..83fee5e 100644 --- a/src/models/set_transformer.py +++ b/src/models/set_transformer.py @@ -558,13 +558,13 @@ class SetTransformer(BaseModel): # Predict # ------------------------------------------------------------------ - def predict(self, team_roster_df: pd.DataFrame) -> np.ndarray: + def predict(self, team_roster_df: pd.DataFrame): """Predict total team value (season-long points).""" if not self._trained: raise RuntimeError("Model not fitted. Call fit() first.") val = self._predict_single(team_roster_df) - return np.array([val]) + return val def _predict_single(self, team_roster_df: pd.DataFrame) -> float: if self._using_torch and self._torch_model is not None: @@ -590,19 +590,25 @@ class SetTransformer(BaseModel): # Marginal value analysis # ------------------------------------------------------------------ - def value_added(self, team_roster_df: pd.DataFrame, new_player: dict) -> float: + def value_added(self, team_roster_df: pd.DataFrame, new_player) -> float: """Marginal value: delta when adding new_player to the team.""" baseline = self._predict_single(team_roster_df) - augmented = pd.concat( - [team_roster_df, pd.DataFrame([new_player])], ignore_index=True - ) + if isinstance(new_player, pd.DataFrame): + augmented = pd.concat([team_roster_df, new_player], ignore_index=True) + else: + augmented = pd.concat( + [team_roster_df, pd.DataFrame([new_player])], ignore_index=True + ) augmented_val = self._predict_single(augmented) return augmented_val - baseline - def value_removed(self, team_roster_df: pd.DataFrame, removed_player_idx: int) -> float: - """Marginal loss: delta when removing a player.""" + def value_removed(self, team_roster_df: pd.DataFrame, removed_player) -> float: + """Marginal loss: delta when removing a player (by index or name).""" baseline = self._predict_single(team_roster_df) - reduced = team_roster_df.drop(team_roster_df.index[removed_player_idx]) + if isinstance(removed_player, str): + reduced = team_roster_df[team_roster_df["name"] != removed_player] + else: + reduced = team_roster_df.drop(team_roster_df.index[removed_player]) reduced_val = self._predict_single(reduced) return baseline - reduced_val @@ -610,21 +616,27 @@ class SetTransformer(BaseModel): self, team_roster_df: pd.DataFrame, candidate_pool: pd.DataFrame, - to_replace: List[int], - ) -> Dict[int, pd.DataFrame]: + to_replace: List, + ) -> Dict: """For each player to replace, rank candidates by predicted team value delta. Args: team_roster_df: current team roster. candidate_pool: DataFrame of free-agent candidates. - to_replace: list of indices in team_roster_df to consider replacing. + to_replace: list of player names (str) or indices (int) in team_roster_df + to consider replacing. Returns: - dict mapping replace_idx -> DataFrame of candidates ranked by delta. + dict mapping player_name -> DataFrame of candidates ranked by delta. """ results = {} - for rp_idx in to_replace: - base_team = team_roster_df.drop(team_roster_df.index[rp_idx]) + for rp in to_replace: + if isinstance(rp, str): + base_team = team_roster_df[team_roster_df["name"] != rp] + key = rp + else: + base_team = team_roster_df.drop(team_roster_df.index[rp]) + key = rp deltas = [] for _, cand in candidate_pool.iterrows(): cand_dict = cand.to_dict() @@ -640,7 +652,7 @@ class SetTransformer(BaseModel): "team_value_delta": new_val - current_val, }) - results[rp_idx] = pd.DataFrame(deltas).sort_values( + results[key] = pd.DataFrame(deltas).sort_values( "team_value_delta", ascending=False ) return results diff --git a/src/models/survival_model.py b/src/models/survival_model.py new file mode 100644 index 0000000..be6c8b0 --- /dev/null +++ b/src/models/survival_model.py @@ -0,0 +1,354 @@ +"""Minutes-played survival model via Weibull AFT. + +Models the distribution of minutes played per matchweek using Weibull +Accelerated Failure Time. Supports both lifelines (preferred) and a pure-scipy +MLE fallback so the module works in minimal environments. + +Provides: +- Expected minutes / confidence intervals +- Starter probability (≥60 min) +- Full-match probability (90 min) +""" + +import logging +import math +from typing import Optional, Tuple + +import numpy as np +import pandas as pd +from sklearn.linear_model import LinearRegression +from sklearn.preprocessing import StandardScaler + +from .base_model import BaseModel + +logger = logging.getLogger(__name__) + +# --------------------------------------------------------------------------- +# Feature-selection keyword list +# --------------------------------------------------------------------------- +_MINUTE_KEYWORDS = [ + "minute", "game", "rest", "fatigue", "age", "injury", + "played", "starter", "bench", "appearance", + "recovery", "rotation", "squad", "season", + "match", "form", "fitness", +] + + +def _select_survival_features(X: pd.DataFrame) -> list: + """Pick columns whose name contains any survival-relevant keyword.""" + lower_cols = {c: str(c).lower() for c in X.columns} + selected = [ + c for c, cl in lower_cols.items() + if any(kw in cl for kw in _MINUTE_KEYWORDS) + ] + if not selected: + selected = list(X.select_dtypes(include=[np.number]).columns[:20]) + logger.info("No keyword-matched survival features; using first 20 numeric columns") + else: + logger.info(f"Selected {len(selected)} survival features via keyword matching") + return selected + + +# =================================================================== +# Weibull helper functions for the scipy fallback +# =================================================================== + +def _weibull_log_likelihood(params, X, t, event, eps=1e-10): + """Negative log-likelihood for Weibull AFT model. + + Parameters + ---------- + params : ndarray (p_features + 1,) + First p entries: beta (coefficients for X). + Last entry: log_k (log shape parameter ensures k > 0). + X : ndarray (n, p) + Scaled feature matrix. + t : ndarray (n,) + Observed durations (minutes played). + event : ndarray (n,) + 0 → exact failure (subbed off), 1 → right-censored (completed 90). + eps : float + Small epsilon for numerical stability. + + Returns + ------- + neg_ll : float + Negative log-likelihood (to be minimized). + """ + p = X.shape[1] + beta = params[:p] + log_k = params[p] + k = np.exp(log_k) + eps + + log_lambda = X.dot(beta) # log(λ_i) = X_i * beta + lambda_ = np.exp(log_lambda) + eps + log_t = np.log(np.maximum(t, eps)) + z = t / lambda_ + + # Log-PDF for uncensored (event == 0) + log_pdf = np.log(k) - log_lambda + (k - 1.0) * (log_t - log_lambda) - z ** k + + # Log-SF for censored (event == 1) + log_sf = -(z ** k) + + # event==1 → censored → use SF; event==0 → observed → use PDF + ll = np.where(event == 1, log_sf, log_pdf) + return -ll.sum() + + +def _fit_weibull_mle(X, t, event): + """Fit Weibull AFT via scipy MLE. + + Returns + ------- + beta : ndarray (p,) + Feature coefficients (scaled to original duration range). + k : float + Shape parameter. + t_scale : float + Scale factor to convert normalized predictions back to minutes. + """ + from scipy.optimize import minimize + + n, p = X.shape + t_scale = max(t.max(), 1.0) + t_norm = np.clip(t / t_scale, 1e-6, 1.0) + log_t_norm = np.log(np.maximum(t_norm, 1e-9)) + + lr = LinearRegression(fit_intercept=False) + lr.fit(X, log_t_norm) + beta0 = np.clip(lr.coef_.copy(), -5, 5) + + bounds = [(-10, 10)] * p + [(-5, 3)] + init = np.concatenate([beta0, [0.0]]) + + result = minimize( + _weibull_log_likelihood, + init, + args=(X, t_norm, event), + method="L-BFGS-B", + bounds=bounds, + options={"maxiter": 2000, "ftol": 1e-10}, + ) + if not result.success: + logger.warning(f"Weibull MLE did not converge: {result.message}") + + beta = result.x[:p] + k = max(np.exp(result.x[p]), 1e-4) + return beta, k, t_scale + + +# =================================================================== +# MinutesSurvivalModel +# =================================================================== + +class MinutesSurvivalModel(BaseModel): + """Weibull AFT model for minutes-played distribution. + + Parameters + ---------- + model_dir : str + Directory for persisting trained models. + force_scipy : bool + If True, use the pure-scipy MLE fallback even when lifelines + is installed. + """ + + def __init__( + self, + model_dir: str = "models_trained", + force_scipy: bool = False, + ): + super().__init__(model_dir) + self.force_scipy = force_scipy + + self.scaler = StandardScaler() + self.feature_names = None + + # Weibull parameters + self._beta = None # feature coefficients → log(λ) + self._k = None # shape parameter + self._afitter = None # lifelines WeibullAFTFitter instance (if used) + self._t_scale = 90.0 + self._use_lifelines = False + + # ------------------------------------------------------------------ + # Fit + # ------------------------------------------------------------------ + def fit( + self, + X: pd.DataFrame, + durations: np.ndarray, + events: np.ndarray, + **kwargs, + ): + """Fit the Weibull AFT model. + + Args: + X: Feature matrix (one row per player-match). + durations: Minutes played (0–90); `y` alias for BaseModel compat. + events: + 0 → exact duration observed (subbed off before 90). + 1 → right-censored (player completed the full 90 minutes). + """ + self.feature_names = _select_survival_features(X) + X_clean = X[self.feature_names].select_dtypes(include=[np.number]).fillna(0) + X_scaled = self.scaler.fit_transform(X_clean) + + durations = np.asarray(durations, dtype=float) + events = np.asarray(events, dtype=int) + + # Try lifelines first -------------------------------------------------- + if not self.force_scipy: + try: + import lifelines # noqa: F401 + from lifelines import WeibullAFTFitter + + df = pd.DataFrame(X_scaled, columns=self.feature_names) + df["duration"] = durations + df["event"] = events + + aft = WeibullAFTFitter() + aft.fit(df, duration_col="duration", event_col="event") + self._afitter = aft + self._use_lifelines = True + self._t_scale = 1.0 # lifelines works in original duration scale + logger.info( + "MinutesSurvivalModel fitted via lifelines " + f"(n={len(durations)}, features={len(self.feature_names)})" + ) + return self + except ImportError: + logger.info("lifelines not installed; falling back to scipy MLE") + except Exception as exc: + logger.warning(f"lifelines failed ({exc}); falling back to scipy MLE") + + # Scipy fallback ------------------------------------------------------- + self._use_lifelines = False + self._beta, self._k, self._t_scale = _fit_weibull_mle(X_scaled, durations, events) + logger.info( + f"MinutesSurvivalModel fitted via scipy MLE " + f"(n={len(durations)}, features={len(self.feature_names)}, " + f"k={self._k:.3f})" + ) + return self + + # ------------------------------------------------------------------ + # Predict (BaseModel interface — returns expected minutes) + # ------------------------------------------------------------------ + def predict(self, X: pd.DataFrame) -> np.ndarray: + """Return expected minutes (E[T]) — BaseModel interface.""" + return self.predict_expected_minutes(X) + + # ------------------------------------------------------------------ + # Preprocessing + # ------------------------------------------------------------------ + def _preprocess(self, X: pd.DataFrame) -> np.ndarray: + if self.feature_names is None: + raise RuntimeError("Model not trained. Call fit() first.") + X_c = X[self.feature_names].select_dtypes(include=[np.number]).fillna(0) + return self.scaler.transform(X_c) + + # ------------------------------------------------------------------ + # Core distribution + # ------------------------------------------------------------------ + def predict_distribution( + self, X: pd.DataFrame + ) -> Tuple[np.ndarray, np.ndarray, np.ndarray]: + """Return (expected_minutes, lower_bound, upper_bound). + + ``lower_bound`` and ``upper_bound`` are approximate 95 % confidence + intervals derived from the Weibull variance. + """ + X_scaled = self._preprocess(X) + + if self._use_lifelines and self._afitter is not None: + df = pd.DataFrame(X_scaled, columns=self.feature_names) + # Lifelines returns median survival in its summary; we approximate + # expected minutes using the median and estimated shape. + median = self._afitter.predict_median(df).values.flatten() + # Heuristic: for Weibull, E[T] ≈ median / (ln 2)^(1/k). + # Derive k from the lifelines summary if possible, else guess ~1. + try: + summary = self._afitter.summary + log_k = summary.loc["lambda_", "coef"] + k = 1.0 / np.exp(log_k) if abs(log_k) > 1e-8 else 1.0 + except Exception: + k = 1.0 + expected = median * np.exp(np.log(np.log(2)) / k) + # Std via coefficient of variation + coef_var = np.sqrt(np.exp( + np.log(math.gamma(1 + 2 / k)) - 2 * np.log(math.gamma(1 + 1 / k)) + )) + std = expected * coef_var + lower = np.maximum(0, expected - 1.96 * std) + upper = np.minimum(90, expected + 1.96 * std) + return expected, lower, upper + + # Scipy / stored parameters (normalized scale, convert to minutes) + log_lambda = X_scaled.dot(self._beta) + lambda_ = np.exp(log_lambda) * self._t_scale + k = self._k + + # Expected value: λ * Γ(1 + 1/k) + gamma_1 = math.gamma(1.0 + 1.0 / k) + expected = lambda_ * gamma_1 + + # Variance = λ² * (Γ(1+2/k) - Γ²(1+1/k)) + gamma_2 = math.gamma(1.0 + 2.0 / k) + var = (lambda_ ** 2) * (gamma_2 - gamma_1 ** 2) + std = np.sqrt(np.maximum(var, 0.01)) + + lower = np.maximum(0, expected - 1.96 * std) + upper = np.minimum(90, expected + 1.96 * std) + expected = np.clip(expected, 0, 90) + return expected, lower, upper + + def predict_expected_minutes(self, X: pd.DataFrame) -> np.ndarray: + """Return E[minutes] for each row.""" + expected, _, _ = self.predict_distribution(X) + return expected + + def predict_full_match_probability(self, X: pd.DataFrame) -> np.ndarray: + """Probability the player completes 90 minutes: P(T ≥ 90).""" + X_scaled = self._preprocess(X) + + if self._use_lifelines and self._afitter is not None: + try: + surv = self._afitter.predict_survival_function( + pd.DataFrame(X_scaled, columns=self.feature_names), + times=[90.0], + ) + return 1.0 - surv.values.flatten() + except Exception: + pass + + log_lambda = X_scaled.dot(self._beta) + lambda_ = np.exp(log_lambda) * self._t_scale + surv = np.exp(-((90.0 / lambda_) ** self._k)) + prob = 1.0 - surv + return np.clip(prob, 0.0, 1.0) + + def predict_starter_probability( + self, X: pd.DataFrame, min_minutes: float = 60.0 + ) -> np.ndarray: + """Probability player plays at least *min_minutes* (default 60). + + Useful as a "likely starter" proxy. + """ + X_scaled = self._preprocess(X) + + if self._use_lifelines and self._afitter is not None: + try: + surv = self._afitter.predict_survival_function( + pd.DataFrame(X_scaled, columns=self.feature_names), + times=[min_minutes], + ) + return surv.values.flatten() + except Exception: + pass + + log_lambda = X_scaled.dot(self._beta) + lambda_ = np.exp(log_lambda) * self._t_scale + surv = np.exp(-((min_minutes / lambda_) ** self._k)) + return np.clip(surv, 0.0, 1.0) diff --git a/src/optimization/__init__.py b/src/optimization/__init__.py index 82ce271..3ccb650 100644 --- a/src/optimization/__init__.py +++ b/src/optimization/__init__.py @@ -1 +1,14 @@ """Optimization modules for auction and lineup selection.""" + +from .auction_solver import AuctionSolver, AuctionConfig, PlayerValuation +from .lineup_solver import LineupSolver, LineupConstraints, PlayerScore, MCTSNode +from .opponent_model import OpponentModel +from .transfer_analyzer import TransferAnalyzer + +# Phase 2: Bandit, opponent bidding, budget optimization +from .bandit_auction import BanditAuctionSolver +from .opponent_bidding_model import OpponentBidModel +from .budget_optimizer import BudgetOptimizer + +# Phase 3: Reinforcement learning auction agent +from .rl_auction_agent import AuctionEnv, RLAuctionPolicy, RLAuctionTrainer, QNetwork diff --git a/src/optimization/bandit_auction.py b/src/optimization/bandit_auction.py new file mode 100644 index 0000000..c6ea353 --- /dev/null +++ b/src/optimization/bandit_auction.py @@ -0,0 +1,359 @@ +"""Contextual Multi-Armed Bandit for live auction bidding decisions. + +Uses Thompson Sampling with Beta-distributed posteriors over discrete bid levels. +Falls back to UCB when exploration depth is insufficient. +Integrates with AuctionConfig from auction_solver.py. +""" + +import logging +from dataclasses import dataclass, field +from typing import Dict, List, Optional, Tuple + +import numpy as np +import pandas as pd +from scipy.stats import beta as beta_dist + +logger = logging.getLogger(__name__) + +# Default bid arms as fractions of total budget +DEFAULT_BID_ARMS = np.array( + [0.0, 0.005, 0.01, 0.02, 0.03, 0.05, 0.08, 0.12, 0.18, 0.25], + dtype=np.float64, +) + +ROLES = ["P", "D", "C", "A"] +ROLE_SCARCITY = {"P": 3, "D": 8, "C": 8, "A": 6} +ROLE_POOL_SIZE = {"P": 4, "D": 22, "C": 24, "A": 12} + + +@dataclass +class BanditArmState: + alpha: float = 1.0 + beta: float = 1.0 + trials: int = 0 + wins: float = 0.0 + + +@dataclass +class AuctionState: + budget_remaining: float = 500.0 + total_budget: float = 500.0 + slots_filled: Dict[str, int] = field(default_factory=lambda: {"P": 0, "D": 0, "C": 0, "A": 0}) + slots_total: Dict[str, int] = field(default_factory=lambda: {"P": 3, "D": 8, "C": 8, "A": 6}) + round_number: int = 1 + opponent_budgets: List[float] = field(default_factory=list) + + +@dataclass +class PlayerContext: + name: str + role: str + projected_points: float + role_scarcity: float + value_over_replacement: float + budget_remaining_fraction: float + slots_remaining_in_role: int + round_number: int + opponent_budget_avg: float + + +def _get_attr(obj, key, default=None): + """Get attribute or dict key from an object.""" + if isinstance(obj, dict): + return obj.get(key, default) + return getattr(obj, key, default) + + +class BanditAuctionSolver: + """Thompson Sampling bandit for live auction bid selection. + + Arms are discrete bid fractions. Each (role, scarcity_level) maintains + independent Beta posteriors. Falls back to UCB when total observations + for a context group are < 50. + """ + + def __init__( + self, + bid_arms: Optional[np.ndarray] = None, + total_budget: float = 500.0, + min_obs_for_ts: int = 50, + ucb_exploration: float = 1.414, + config=None, + ): + if config is not None: + from .auction_solver import AuctionConfig + total_budget = config.total_budget + self.bid_arms = bid_arms if bid_arms is not None else DEFAULT_BID_ARMS + self.n_arms = len(self.bid_arms) + self.total_budget = total_budget + self.min_obs_for_ts = min_obs_for_ts + self.ucb_exploration = ucb_exploration + self.bid_fractions = self.bid_arms + + self.posteriors: Dict[Tuple[str, int], List[BanditArmState]] = {} + for role in ROLES: + self._ensure_posteriors(role, 1) + + def _context_key(self, role: str, scarcity_level: int) -> Tuple[str, int]: + return (role, scarcity_level) + + def _ensure_posteriors(self, role: str, n_slots_remaining: int): + scarcity = max(1, n_slots_remaining) + key = self._context_key(role, scarcity) + if key not in self.posteriors: + self.posteriors[key] = [ + BanditArmState(alpha=1.0, beta=1.0, trials=0, wins=0.0) + for _ in range(self.n_arms) + ] + logger.debug(f"Initialized bandit posteriors for role={role}, scarcity={scarcity}") + + def _total_obs(self, key: Tuple[str, int]) -> int: + arms = self.posteriors.get(key, []) + return sum(a.trials for a in arms) + + def compute_context( + self, + player, + auction_state, + pool_stats: Optional[dict] = None, + ) -> PlayerContext: + """Build context vector for a player given current auction state. + + Args: + player: dict or object with name, role, projected_points attributes. + auction_state: current AuctionState (or dict with same keys). + pool_stats: optional dict with role-level pool means and stds. + + Returns: + PlayerContext dataclass with all context features. + """ + role = _get_attr(player, "role") + points = float(_get_attr(player, "projected_points", 6.5)) + name = _get_attr(player, "name", "unknown") + + total_slots = _get_attr(auction_state, "slots_total", {}).get(role, 1) + filled = _get_attr(auction_state, "slots_filled", {}).get(role, 0) + slots_remaining = max(total_slots - filled, 0) + scarcity = max(slots_remaining / max(total_slots, 1), 0.05) + + budget_remaining = _get_attr(auction_state, "budget_remaining", 500.0) + total_budget_attr = _get_attr(auction_state, "total_budget", 500.0) + budget_fraction = budget_remaining / max(total_budget_attr, 1) + + opponent_budget_avg = 0.0 + opponent_budgets = _get_attr(auction_state, "opponent_budgets", []) + if opponent_budgets: + opponent_budget_avg = float(np.mean(opponent_budgets)) + elif total_budget_attr > 0: + opponent_budget_avg = total_budget_attr * 0.6 + + if pool_stats and role in pool_stats: + role_mean = pool_stats[role].get("mean", 0.0) + role_std = pool_stats[role].get("std", 1.0) + points_z = (points - role_mean) / max(role_std, 0.01) if role_std > 0 else 0.0 + else: + points_z = points / 15.0 + + vor = max(points - 6.5, 0.0) + + return PlayerContext( + name=name, + role=role, + projected_points=points, + role_scarcity=scarcity, + value_over_replacement=vor, + budget_remaining_fraction=budget_fraction, + slots_remaining_in_role=slots_remaining, + round_number=_get_attr(auction_state, "round_number", 1), + opponent_budget_avg=opponent_budget_avg, + ) + + def select_bid( + self, + player, + auction_state, + pool_stats: Optional[dict] = None, + ) -> Tuple[int, float]: + """Select bid arm using Thompson Sampling (or UCB fallback). + + Returns: + (arm_index, bid_amount_in_credits) + """ + ctx = self.compute_context(player, auction_state, pool_stats) + role = ctx.role + n_slots = ctx.slots_remaining_in_role + + self._ensure_posteriors(role, n_slots) + key = self._context_key(role, n_slots) + arms = self.posteriors[key] + total_obs = sum(a.trials for a in arms) + + if total_obs >= self.min_obs_for_ts: + samples = [float(np.random.beta(a.alpha, max(a.beta, 0.01))) for a in arms] + arm_idx = int(np.argmax(samples)) + logger.debug( + f"Thompson Sampling: role={role}, arms_sampled={samples[:5]}..., " + f"selected_arm={arm_idx}" + ) + else: + values = [] + for _i, arm in enumerate(arms): + if arm.trials == 0: + values.append(float("inf")) + else: + mean = arm.wins / arm.trials + bonus = self.ucb_exploration * np.sqrt( + np.log(max(total_obs, 1)) / arm.trials + ) + values.append(mean + bonus) + arm_idx = int(np.argmax(values)) + logger.debug( + f"UCB fallback: role={role}, total_obs={total_obs}, " + f"selected_arm={arm_idx}" + ) + + bid_amount = round(self.bid_arms[arm_idx] * self.total_budget) + + budget_rem = _get_attr(auction_state, "budget_remaining", self.total_budget) + if bid_amount > budget_rem: + bid_amount = budget_rem + arm_idx = int(np.argmin(np.abs(self.bid_arms * self.total_budget - bid_amount))) + + return arm_idx, bid_amount + + def update(self, arm_idx: int, reward: float, player_role: str): + """Update Beta posterior for the selected arm. + + Reward should be a normalized value: (player_season_value - cost) scaled. + + Args: + arm_idx: index of the selected arm. + reward: normalized reward signal (higher = better purchase). + player_role: role of the purchased player. + """ + reward_clipped = max(0.0, min(1.0, reward)) + win = 1.0 if reward > 0 else 0.0 + + for key, arms in self.posteriors.items(): + role, _scarcity = key + if role == player_role: + if arm_idx < len(arms): + arm = arms[arm_idx] + arm.trials += 1 + arm.wins += win + arm.alpha += reward_clipped + arm.beta += (1.0 - reward_clipped) + logger.debug( + f"Updated arm {arm_idx} for {player_role}: " + f"trials={arm.trials}, alpha={arm.alpha:.2f}, beta={arm.beta:.2f}" + ) + for key, arms in self.posteriors.items(): + role, _scarcity = key + if role == player_role: + for i, arm in enumerate(arms): + if i == arm_idx: + continue + arm.beta = max(arm.beta, 1.001) + + def get_arm_stats(self) -> Dict[str, dict]: + """Return arm statistics for analysis. + + Returns: + dict mapping "role/scarcity/arm_idx" -> stats dict. + """ + stats = {} + agg = {} + for (role, scarcity), arms in self.posteriors.items(): + for idx, arm in enumerate(arms): + label = f"{role}/scarcity={scarcity}/arm={idx}" + stats[label] = { + "trials": arm.trials, + "wins": arm.wins, + "alpha": arm.alpha, + "beta": arm.beta, + "win_rate": arm.wins / max(arm.trials, 1), + "bid_fraction": float(self.bid_arms[idx]), + "bid_amount": float(self.bid_arms[idx] * self.total_budget), + } + if idx not in agg: + agg[idx] = {"trials": 0, "wins": 0.0, "alpha": 0.0, "beta": 0.0} + agg[idx]["trials"] += arm.trials + agg[idx]["wins"] += arm.wins + agg[idx]["alpha"] += arm.alpha + agg[idx]["beta"] += arm.beta + for idx, a in agg.items(): + stats[idx] = { + "trials": a["trials"], + "wins": a["wins"], + "alpha": a["alpha"], + "beta": a["beta"], + "win_rate": a["wins"] / max(a["trials"], 1), + "bid_fraction": float(self.bid_arms[idx]), + "bid_amount": float(self.bid_arms[idx] * self.total_budget), + } + return stats + + def exploration_bonus(self, player_role: str, player_points: float = 0.0) -> float: + """Compute exploration bonus for new/unknown player types. + + Higher bonus when the bandit has limited experience with a role. + Encourages exploring sleeper players. + + Returns: + Recommended extra bid amount in credits. + """ + total_obs = 0 + for (role, _scarcity), arms in self.posteriors.items(): + if role == player_role: + total_obs += sum(a.trials for a in arms) + + if total_obs == 0: + bonus = self.total_budget * 0.04 + elif total_obs < 20: + bonus = self.total_budget * 0.025 + elif total_obs < 50: + bonus = self.total_budget * 0.01 + else: + bonus = 0.0 + + if bonus > 0: + base = max(player_points * 0.5, 0) + bonus += base * (1.0 / max(total_obs, 1)) * 20 + + logger.debug(f"Exploration bonus for {player_role}: {bonus:.1f} (obs={total_obs})") + return bonus + + def recommend_bid_summary( + self, + player, + auction_state, + pool_stats: Optional[dict] = None, + ) -> dict: + """Full bidding recommendation for a player. + + Returns a dict with arm index, bid amount, context, exploration bonus, + and total recommended bid. + """ + arm_idx, bid_amount = self.select_bid(player, auction_state, pool_stats) + ctx = self.compute_context(player, auction_state, pool_stats) + bonus = self.exploration_bonus(ctx.role, ctx.projected_points) + + return { + "player": _get_attr(player, "name", "unknown"), + "role": ctx.role, + "arm_index": arm_idx, + "base_bid": bid_amount, + "exploration_bonus": round(bonus), + "total_bid": round(bid_amount + bonus), + "recommended_bid": round(bid_amount + bonus), + "context": { + "points_zscore": round( + (ctx.projected_points - 6.5) / 2.0, 2 + ), + "role_scarcity": round(ctx.role_scarcity, 3), + "budget_remaining_frac": round(ctx.budget_remaining_fraction, 3), + "slots_remaining": ctx.slots_remaining_in_role, + "round": ctx.round_number, + "vor": round(ctx.value_over_replacement, 1), + }, + } diff --git a/src/optimization/budget_optimizer.py b/src/optimization/budget_optimizer.py new file mode 100644 index 0000000..f47d3c1 --- /dev/null +++ b/src/optimization/budget_optimizer.py @@ -0,0 +1,538 @@ +"""Bayesian Optimization for role-level budget allocation. + +Uses Gaussian Process regression to find optimal budget distribution +across roles (P, D, C, A) that maximizes total projected team value. +""" + +import logging +from dataclasses import dataclass, field +from typing import Dict, List, Optional, Tuple + +import numpy as np +import pandas as pd +from scipy.optimize import minimize +from scipy.special import softmax + +from src.optimization.auction_solver import AuctionConfig + +logger = logging.getLogger(__name__) + +DEFAULT_ROSTER_QUOTAS = {"P": 3, "D": 8, "C": 8, "A": 6} + + +@dataclass +class RoleBudgetResult: + allocation: Dict[str, float] + total_value: float + role_values: Dict[str, float] + value_curves: Dict[str, Tuple[np.ndarray, np.ndarray]] + + +class BudgetOptimizer: + """Bayesian Optimization for budget allocation across roster roles. + + Finds the split of total_budget across P/D/C/A that yields the + highest possible team points via greedy fill within each role's budget. + Supports mid-auction adaptive rebalancing. + """ + + def __init__( + self, + total_budget: float = 500.0, + roster_quotas: Optional[Dict[str, int]] = None, + auction_config: Optional[AuctionConfig] = None, + random_state: int = 42, + ): + self.total_budget = total_budget + self.roster_quotas = roster_quotas or dict(DEFAULT_ROSTER_QUOTAS) + self.auction_config = auction_config or AuctionConfig() + self.random_state = random_state + self.rng = np.random.RandomState(random_state) + self._last_allocation: Optional[Dict[str, float]] = None + self._last_value: float = 0.0 + self._optimization_history: List[dict] = [] + self._gk_available = False + + # ------------------------------------------------------------------ + # Core optimization + # ------------------------------------------------------------------ + + def optimize( + self, + player_pool_df: pd.DataFrame, + n_calls: int = 50, + use_skopt: bool = True, + ) -> Dict[str, float]: + """Optimize budget allocation across roles. + + Uses scikit-optimize GaussianProcessRegressor if available, + otherwise simplex-based local search. + + Args: + player_pool_df: DataFrame with columns [name, role, projected_points]. + n_calls: number of GP evaluations. + use_skopt: attempt Gaussian Process optimization. + + Returns: + Dict mapping role -> recommended budget amount. + """ + if "role" not in player_pool_df.columns or "projected_points" not in player_pool_df.columns: + raise ValueError("player_pool_df must have 'role' and 'projected_points' columns") + + roles = list(self.roster_quotas.keys()) + n_roles = len(roles) + + if use_skopt: + try: + return self._optimize_gp(player_pool_df, roles, n_calls) + except ImportError: + logger.info("scikit-optimize not installed. Using local search.") + except Exception as exc: + logger.warning(f"GP optimization failed: {exc}. Using local search.") + + return self._optimize_local(player_pool_df, roles, n_calls) + + def _optimize_gp( + self, player_pool_df: pd.DataFrame, roles: list, n_calls: int + ) -> Dict[str, float]: + """Gaussian Process-based budget optimization.""" + from skopt import gp_minimize + from skopt.space import Space + from skopt.learning import GaussianProcessRegressor + + n_roles = len(roles) + space = Space([(0.01, 0.70) for _ in range(n_roles)]) + + def objective_wrapper(fractions): + fractions = np.array(fractions, dtype=float) + fractions = self._normalize_fractions(fractions) + value = self._objective(fractions, roles, player_pool_df) + self._optimization_history.append({ + "fractions": fractions.tolist(), + "value": value, + }) + return -value + + def params_to_fractions(params): + return self._normalize_fractions(np.array(params, dtype=float)) + + result = gp_minimize( + objective_wrapper, + space, + n_calls=n_calls, + random_state=self.random_state, + n_initial_points=max(10, n_calls // 5), + verbose=False, + n_jobs=-1, + ) + + best_fractions = params_to_fractions(result.x) + self._last_allocation = self._fractions_to_allocation(best_fractions, roles) + self._last_value = -result.fun + + logger.info( + f"GP Budget optimization complete: {self._last_allocation} " + f"=> value={self._last_value:.1f}" + ) + + return self._last_allocation + + def _optimize_local( + self, player_pool_df: pd.DataFrame, roles: list, n_calls: int + ) -> Dict[str, float]: + """Local search optimization using Nelder-Mead simplex.""" + n_roles = len(roles) + + best_allocation = None + best_value = -float("inf") + + for restart in range(max(n_calls // 10, 1)): + x0 = self._random_allocation(n_roles) + + def simplex_objective(fractions): + fractions = self._normalize_fractions(np.array(fractions, dtype=float)) + value = self._objective(fractions, roles, player_pool_df) + self._optimization_history.append({ + "fractions": fractions.tolist(), + "value": value, + }) + return -value + + res = minimize( + simplex_objective, + x0, + method="Nelder-Mead", + options={"maxiter": max(n_calls // 3, 20), "xatol": 1e-3, "fatol": 1e-3}, + ) + + fractions = self._normalize_fractions(np.array(res.x, dtype=float)) + value = -res.fun + + if value > best_value: + best_value = value + best_allocation = fractions + + if best_allocation is None: + best_allocation = self._proportional_allocation(player_pool_df, roles) + + self._last_allocation = self._fractions_to_allocation(best_allocation, roles) + self._last_value = best_value + + logger.info( + f"Local optimization complete: {self._last_allocation} " + f"=> value={self._last_value:.1f}" + ) + + return self._last_allocation + + # ------------------------------------------------------------------ + # Adaptive rebalancing + # ------------------------------------------------------------------ + + def optimize_adaptive( + self, + player_pool_df: pd.DataFrame, + remaining_slots: Dict[str, int], + spent_per_role: Dict[str, float], + n_calls: int = 30, + ) -> Dict[str, float]: + """Mid-auction rebalancing — optimize remaining budget for unfilled slots. + + Args: + player_pool_df: remaining available players pool. + remaining_slots: dict of role -> slots still needed. + spent_per_role: dict of role -> credits already spent. + n_calls: GP evaluation budget. + + Returns: + Dict mapping role -> recommended budget for remaining slots. + """ + remaining_budget = self.total_budget - sum(spent_per_role.values()) + remaining_budget = max(1.0, remaining_budget) + + if sum(remaining_slots.values()) == 0: + logger.info("All slots filled. No budget to allocate.") + return {r: 0.0 for r in self.roster_quotas} + + saved_quotas = self.roster_quotas + saved_total = self.total_budget + self.roster_quotas = dict(remaining_slots) + self.total_budget = remaining_budget + + available = player_pool_df[ + player_pool_df["role"].isin( + [r for r, s in remaining_slots.items() if s > 0] + ) + ] + + if len(available) == 0: + logger.warning("No available players for remaining slots.") + self.roster_quotas = saved_quotas + self.total_budget = saved_total + return {r: 0.0 for r in saved_quotas} + + allocation = self.optimize( + available, + n_calls=n_calls, + use_skopt=True, + ) + + self.roster_quotas = saved_quotas + self.total_budget = saved_total + + result = {} + for role in saved_quotas: + result[role] = allocation.get(role, 0.0) + return result + + # ------------------------------------------------------------------ + # Objective function + # ------------------------------------------------------------------ + + def _objective( + self, + fractions: np.ndarray, + roles: list, + player_pool_df: pd.DataFrame, + ) -> float: + """Simulate greedy fill within role budgets; return total projected points. + + For each role, pick the best players by projected_points until the + role budget or slot quota is exhausted. + """ + total_value = 0.0 + + for i, role in enumerate(roles): + budget = fractions[i] * self.total_budget + slots = self.roster_quotas.get(role, 0) + role_players = player_pool_df[player_pool_df["role"] == role].copy() + + if len(role_players) == 0 or slots == 0: + continue + + role_players = role_players.sort_values( + "projected_points", ascending=False + ) + + total_cost = 0.0 + filled = 0 + + for _, player in role_players.iterrows(): + points = float(player["projected_points"]) + estimated_price = self._estimate_price_simple(points, budget, role) + + if total_cost + estimated_price > budget: + continue + + total_cost += estimated_price + total_value += points + filled += 1 + + if filled >= slots: + break + + return total_value + + def _estimate_price_simple( + self, projected_points: float, role_budget: float, role: str + ) -> float: + """Simple price estimate: points * role_factor clamped within budget.""" + role_factor = {"P": 4.0, "D": 2.5, "C": 3.0, "A": 4.5}.get(role, 3.0) + price = projected_points * role_factor + max_price = role_budget * 0.50 + return min(price, max_price, role_budget) + + # ------------------------------------------------------------------ + # Allocation access and visualization data + # ------------------------------------------------------------------ + + def get_allocation(self) -> Dict[str, float]: + """Return last computed allocation (role -> budget amount).""" + if self._last_allocation is None: + return {r: self.total_budget / len(self.roster_quotas) for r in self.roster_quotas} + return dict(self._last_allocation) + + def get_role_value_curves( + self, + player_pool_df: pd.DataFrame, + n_points: int = 20, + ) -> Dict[str, Tuple[np.ndarray, np.ndarray]]: + """Compute diminishing returns curves: budget vs. expected points per role. + + Returns: + dict role -> (budget_array, value_array). + """ + curves = {} + budget_step = self.total_budget / n_points + + for role, slots in self.roster_quotas.items(): + role_players = player_pool_df[player_pool_df["role"] == role].sort_values( + "projected_points", ascending=False + ) + + budgets = np.linspace(0, self.total_budget, n_points) + values = np.zeros(n_points) + + for i, budget_limit in enumerate(budgets): + total_cost = 0.0 + total_value = 0.0 + filled = 0 + + for _, player in role_players.iterrows(): + points = float(player["projected_points"]) + price = self._estimate_price_simple(points, budget_limit, role) + + if total_cost + price > budget_limit: + continue + + total_cost += price + total_value += points + filled += 1 + + if filled >= slots: + break + + values[i] = total_value + + curves[role] = (budgets.copy(), values.copy()) + + return curves + + # ------------------------------------------------------------------ + # Fallback: proportional allocation + # ------------------------------------------------------------------ + + def _proportional_allocation( + self, player_pool_df: pd.DataFrame, roles: list + ) -> np.ndarray: + """Allocate budget proportional to (points_variance * slots) per role.""" + weights = np.zeros(len(roles)) + + for i, role in enumerate(roles): + role_players = player_pool_df[player_pool_df["role"] == role] + if len(role_players) > 1: + variance = role_players["projected_points"].var() + else: + variance = 1.0 + slots = self.roster_quotas.get(role, 1) + weights[i] = variance * slots + + weight_sum = weights.sum() + if weight_sum <= 0: + return np.ones(len(roles)) / len(roles) + + fractions = weights / weight_sum + fractions = np.clip(fractions, 0.02, 0.70) + fractions = fractions / fractions.sum() + + logger.info( + f"Proportional allocation (fallback): " + f"{dict(zip(roles, fractions.round(3)))}" + ) + return fractions + + # ------------------------------------------------------------------ + # Utilities + # ------------------------------------------------------------------ + + def _normalize_fractions(self, fractions: np.ndarray) -> np.ndarray: + """Normalize fractions to sum to 1.0 with minimum per role.""" + fractions = np.clip(fractions, 0.01, 0.70) + total = fractions.sum() + if total <= 0: + return np.ones_like(fractions) / len(fractions) + return fractions / total + + def _fractions_to_allocation( + self, fractions: np.ndarray, roles: list + ) -> Dict[str, float]: + """Convert fractions to absolute budget per role.""" + return { + role: round(float(fractions[i] * self.total_budget), 1) + for i, role in enumerate(roles) + } + + def _random_allocation(self, n_roles: int) -> np.ndarray: + """Generate a random allocation via Dirichlet.""" + alpha = np.ones(n_roles) * 2.0 + return self.rng.dirichlet(alpha) + + def get_optimization_trace(self) -> pd.DataFrame: + """Return DataFrame of all evaluated allocations during optimization.""" + if not self._optimization_history: + return pd.DataFrame() + return pd.DataFrame(self._optimization_history) + + def estimate_team_composition( + self, + player_pool_df: pd.DataFrame, + allocation: Optional[Dict[str, float]] = None, + ) -> pd.DataFrame: + """Given final allocation, return the recommended player selections. + + Returns: + DataFrame with selected players, their estimated costs, and value. + """ + roles = list(self.roster_quotas.keys()) + alloc = allocation or self._last_allocation + if alloc is None: + alloc = self._proportional_allocation_as_dict(roles) + + selections = [] + + for role in roles: + budget = alloc.get(role, 0.0) + slots = self.roster_quotas.get(role, 0) + role_players = player_pool_df[player_pool_df["role"] == role].sort_values( + "projected_points", ascending=False + ) + + total_cost = 0.0 + filled = 0 + + for _, player in role_players.iterrows(): + points = float(player["projected_points"]) + price = self._estimate_price_simple(points, budget, role) + + if total_cost + price > budget: + continue + + total_cost += price + name = player.get("name", player.get("player_name", "unknown")) + selections.append({ + "player": name, + "role": role, + "projected_points": points, + "estimated_cost": round(price, 1), + "value_ratio": round(points / max(price, 1), 3), + }) + filled += 1 + + if filled >= slots: + break + + return pd.DataFrame(selections) + + def _proportional_allocation_as_dict(self, roles: list) -> Dict[str, float]: + fractions = self._proportional_allocation( + pd.DataFrame(columns=["role", "projected_points"]), roles + ) + return {roles[i]: round(float(fractions[i] * self.total_budget), 1) for i in range(len(roles))} + + # ------------------------------------------------------------------ + # Sensitivity analysis + # ------------------------------------------------------------------ + + def sensitivity_analysis( + self, + player_pool_df: pd.DataFrame, + role: Optional[str] = None, + delta_pct: float = 0.05, + n_steps: int = 11, + n_calls: int = 50, + ) -> dict: + """Test how shifting budget into/out of one role affects total value. + + Args: + player_pool_df: current player pool. + role: role to perturb. If None, tests all roles. + delta_pct: fractional step size. + n_steps: number of steps in each direction. + n_calls: number of optimization calls (for compatibility). + + Returns: + Dict with at least a 'per_role' key containing per-role analysis. + """ + base_allocation = self.get_allocation() + roles_to_test = [role] if role else list(base_allocation.keys()) + per_role = {} + + for test_role in roles_to_test: + results = [] + shifts = np.linspace(-delta_pct * n_steps, delta_pct * n_steps, 2 * n_steps + 1) + + for shift in shifts: + adjusted = {} + for r, val in base_allocation.items(): + adjusted[r] = val * (1.0 + (shift if r == test_role else -shift / 3.0)) + + total_adj = sum(adjusted.values()) + for r in adjusted: + adjusted[r] = adjusted[r] / total_adj * self.total_budget + + fractions = np.array([adjusted[r] / self.total_budget for r in base_allocation.keys()]) + value = self._objective( + fractions, + list(base_allocation.keys()), + player_pool_df, + ) + + results.append({ + "shift_pct": round(shift * 100, 1), + "allocation": {r: round(v, 1) for r, v in adjusted.items()}, + "total_value": round(value, 1), + }) + + per_role[test_role] = pd.DataFrame(results) + + return {"per_role": per_role, "base_allocation": base_allocation} diff --git a/src/optimization/opponent_bidding_model.py b/src/optimization/opponent_bidding_model.py new file mode 100644 index 0000000..6c94ce2 --- /dev/null +++ b/src/optimization/opponent_bidding_model.py @@ -0,0 +1,412 @@ +"""Opponent bidding behavior modeling using LightGBM. + +Predicts what competitors will bid for each player in a live auction round. +Supports Monte Carlo simulation of auction outcomes and win probability estimates. +""" + +import logging +from dataclasses import dataclass, field +from typing import Dict, List, Optional, Tuple + +import numpy as np +import pandas as pd + +logger = logging.getLogger(__name__) + + +@dataclass +class OpponentState: + budget_remaining: float = 500.0 + initial_budget: float = 500.0 + slots_filled: Dict[str, int] = field(default_factory=lambda: {"P": 0, "D": 0, "C": 0, "A": 0}) + slots_total: Dict[str, int] = field(default_factory=lambda: {"P": 3, "D": 8, "C": 8, "A": 6}) + round_number: int = 1 + aggression_factor: float = 1.0 + + +def _get_attr(obj, key, default=None): + if isinstance(obj, dict): + return obj.get(key, default) + return getattr(obj, key, default) + + +class OpponentBidModel: + """Predicts opponent bids using LightGBM with heuristic fallback. + + Trains on historical auction logs and outputs estimated max opponent + bid, win probability per player, and Monte Carlo round simulations. + """ + + def __init__(self, random_state: int = 42): + self.random_state = random_state + self.model = None + self.fitted = False + self.feature_names: list = [] + self.rng = np.random.RandomState(random_state) + self._role_scarcity_cache: Dict[str, float] = {} + + # ------------------------------------------------------------------ + # Feature engineering + # ------------------------------------------------------------------ + + def _extract_features( + self, + players_df: pd.DataFrame, + opponent_state: OpponentState, + ) -> pd.DataFrame: + """Build feature matrix for LightGBM prediction. + + Args: + players_df: DataFrame with columns [player_name, player_role, + player_projected_points, ...]. + opponent_state: OpponentState describing current opponent. + + Returns: + Feature DataFrame ready for model input. + """ + df = players_df.copy() + + df["role_P"] = (df["player_role"] == "P").astype(float) + df["role_D"] = (df["player_role"] == "D").astype(float) + df["role_C"] = (df["player_role"] == "C").astype(float) + df["role_A"] = (df["player_role"] == "A").astype(float) + + df["budget_remaining_frac"] = ( + opponent_state.budget_remaining / opponent_state.initial_budget + ) + + for role in ["P", "D", "C", "A"]: + slots_total = opponent_state.slots_total.get(role, 1) + slots_filled = opponent_state.slots_filled.get(role, 0) + scarcity = (slots_total - slots_filled) / slots_total + df[f"scarcity_{role}"] = scarcity + + role_map = {"P": 0, "D": 1, "C": 2, "A": 3} + df["role_code"] = df["player_role"].map(role_map) + + if "player_projected_points" in df.columns: + df["points_sq"] = df["player_projected_points"] ** 2 + df["points_log"] = np.log1p(df["player_projected_points"].clip(lower=0)) + + df["slots_needed_total"] = sum( + opponent_state.slots_total.get(r, 0) - opponent_state.slots_filled.get(r, 0) + for r in ["P", "D", "C", "A"] + ) + df["slots_needed_total"] = df["slots_needed_total"].clip(lower=1) + + df["round_number"] = opponent_state.round_number + + df["urgency"] = 1.0 - (opponent_state.budget_remaining / opponent_state.initial_budget) + df["aggression"] = opponent_state.aggression_factor + + self.feature_names = [ + "player_projected_points", + "role_P", + "role_D", + "role_C", + "role_A", + "role_code", + "budget_remaining_frac", + "scarcity_P", + "scarcity_D", + "scarcity_C", + "scarcity_A", + "points_sq", + "points_log", + "slots_needed_total", + "round_number", + "urgency", + "aggression", + ] + + for col in self.feature_names: + if col not in df.columns: + df[col] = 0.0 + + return df[self.feature_names] + + # ------------------------------------------------------------------ + # Training + # ------------------------------------------------------------------ + + def fit(self, auction_logs: pd.DataFrame): + """Train LightGBM regressor on historical auction logs. + + Args: + auction_logs: DataFrame with columns [player_name, player_role, + player_projected_points, opponent_budget_remaining, + opponent_slots_remaining, role_needed_count, + round_number, winning_bid]. + """ + required_cols = [ + "player_role", "player_projected_points", + "winning_bid", + ] + for col in required_cols: + if col not in auction_logs.columns: + raise ValueError(f"Missing required column: '{col}' in auction_logs") + + df = auction_logs.dropna(subset=required_cols).copy() + + if len(df) < 20: + logger.warning( + f"Only {len(df)} auction records. Not enough to fit LightGBM. " + "Using heuristic fallback." + ) + self.fitted = False + return + + try: + import lightgbm as lgb + except ImportError: + logger.warning("LightGBM not available. Using heuristic fallback.") + self.fitted = False + return + + dummy_state = OpponentState( + budget_remaining=df.get("opponent_budget_remaining", 500), + initial_budget=500, + slots_filled={"P": 0, "D": 0, "C": 0, "A": 0}, + slots_total={"P": 3, "D": 8, "C": 8, "A": 6}, + round_number=1, + aggression_factor=1.0, + ) + + df = df.rename(columns={ + "player_role": "player_role", + "player_projected_points": "player_projected_points", + }) + + dummy_df = df[["player_role", "player_projected_points"]].copy() + dummy_df.columns = ["player_role", "player_projected_points"] + + X = self._extract_features(dummy_df, dummy_state) + + y = df["winning_bid"].astype(float) + y_min, y_max = y.min(), y.max() + + self.model = lgb.LGBMRegressor( + n_estimators=100, + max_depth=6, + learning_rate=0.05, + num_leaves=31, + min_child_samples=10, + subsample=0.8, + colsample_bytree=0.8, + random_state=self.random_state, + verbose=-1, + ) + self.model.fit(X, y) + self.fitted = True + self._y_min = y_min + self._y_max = y_max + + logger.info( + f"OpponentBidModel trained on {len(df)} records. " + f"Target range: [{y_min:.0f}, {y_max:.0f}]" + ) + + # ------------------------------------------------------------------ + # Prediction + # ------------------------------------------------------------------ + + def predict_opponent_bids( + self, + players_df: pd.DataFrame, + opponent_state: OpponentState, + ) -> pd.Series: + """Predict max opponent bid for each player. + + Args: + players_df: DataFrame with player info. + opponent_state: current opponent state. + + Returns: + Series of predicted max opponent bids (index = player index). + """ + if not self.fitted or self.model is None: + bids = self._heuristic_bid(players_df, opponent_state) + return pd.Series(bids, index=players_df.index) + + X = self._extract_features(players_df, opponent_state) + predictions = self.model.predict(X) + predictions = np.clip(predictions, 1, _get_attr(opponent_state, "budget_remaining", 500)) + + return pd.Series(predictions, index=players_df.index) + + def predict_p_acquire( + self, + players_df: pd.DataFrame, + my_bids: np.ndarray, + opponent_state: OpponentState, + temperature: float = 0.1, + ) -> np.ndarray: + """Probability I acquire each player given my bids vs opponent. + + Uses sigmoid: P = 1 / (1 + exp(-(my_bid - opp_bid) / temperature)). + + Args: + players_df: DataFrame with player info. + my_bids: array of my bid amounts per player. + opponent_state: current opponent state. + temperature: softmax temperature (lower = sharper). + + Returns: + Array of acquisition probabilities per player. + """ + opp_bids = self.predict_opponent_bids(players_df, opponent_state).values + margin = np.array(my_bids, dtype=float) - opp_bids + scaled_temp = max(temperature * max(opp_bids.max(), 1), 0.01) + probabilities = 1.0 / (1.0 + np.exp(-margin / scaled_temp)) + return np.clip(probabilities, 0.01, 0.99) + + # ------------------------------------------------------------------ + # Monte Carlo round simulation + # ------------------------------------------------------------------ + + def simulate_live_round( + self, + available_players: pd.DataFrame, + my_budget: float, + opponent_state: OpponentState, + n_sims: int = 1000, + ) -> dict: + """Monte Carlo simulation of a live auction round. + + Args: + available_players: DataFrame of players up for bidding this round. + my_budget: my remaining budget. + opponent_state: opponent's current state. + n_sims: number of simulation runs. + + Returns: + dict with expected_players_acquired, expected_cost, value_matrix. + """ + opp_bids = self.predict_opponent_bids(available_players, opponent_state).values + + n_players = len(available_players) + players_acquired = np.zeros(n_sims, dtype=int) + total_cost = np.zeros(n_sims, dtype=float) + value_matrix = np.zeros((n_sims, n_players), dtype=float) + + for sim_idx in range(n_sims): + my_budget_left = my_budget + acquired = 0 + cost = 0.0 + + for p_idx in range(n_players): + opp_bid = opp_bids[p_idx] + self.rng.normal(0, max(opp_bids[p_idx] * 0.15, 1)) + opp_bid = max(opp_bid, 1) + + my_bid = self._heuristic_bid_single( + available_players.iloc[p_idx], opponent_state + ) + + my_bid = min(my_bid, my_budget_left) + if my_bid > opp_bid: + acquired += 1 + cost += my_bid + my_budget_left -= my_bid + value_matrix[sim_idx, p_idx] = 1.0 + + players_acquired[sim_idx] = acquired + total_cost[sim_idx] = cost + + return { + "expected_players_acquired": float(np.mean(players_acquired)), + "expected_cost": float(np.mean(total_cost)), + "cost_std": float(np.std(total_cost)), + "acquired_std": float(np.std(players_acquired)), + "cost_percentile_25": float(np.percentile(total_cost, 25)), + "cost_percentile_50": float(np.percentile(total_cost, 50)), + "cost_percentile_75": float(np.percentile(total_cost, 75)), + "acquisition_rate": float(players_acquired.mean() / n_players), + "value_matrix": value_matrix, + "n_sims": n_sims, + } + + # ------------------------------------------------------------------ + # Heuristic fallback + # ------------------------------------------------------------------ + + def _heuristic_bid_single( + self, player_row: pd.Series, opponent_state: OpponentState + ) -> float: + """Heuristic bid for a single player (scalar version).""" + role = player_row.get("player_role", None) + points = float(player_row.get("player_projected_points", 6.5)) + + if isinstance(role, pd.Series): + role = role.iloc[0] + + scarcity = self._compute_role_scarcity(role, opponent_state) + aggression = _get_attr(opponent_state, "aggression_factor", 1.0) + budget_rem = _get_attr(opponent_state, "budget_remaining", 500.0) + budget_init = _get_attr(opponent_state, "initial_budget", 500.0) or _get_attr(opponent_state, "total_budget", 500.0) + base_bid = 0.4 * points * scarcity * aggression + bid = base_bid * (budget_rem / budget_init) + return max(bid, 1.0) + + def _heuristic_bid( + self, players_df: pd.DataFrame, opponent_state: OpponentState + ) -> np.ndarray: + """Heuristic bid array for all players.""" + bids = [] + for _, row in players_df.iterrows(): + bids.append(self._heuristic_bid_single(row, opponent_state)) + return np.array(bids, dtype=float) + + def _compute_role_scarcity( + self, role: Optional[str], opponent_state + ) -> float: + """Compute how scarce a role is for the opponent.""" + slots_total_dict = _get_attr(opponent_state, "slots_total", {"P": 3, "D": 8, "C": 8, "A": 6}) + slots_filled_dict = _get_attr(opponent_state, "slots_filled", {"P": 0, "D": 0, "C": 0, "A": 0}) + if role is None or role not in slots_total_dict: + return 1.0 + slots_total = slots_total_dict[role] + filled = slots_filled_dict.get(role, 0) + remaining = max(slots_total - filled, 1) + return slots_total / remaining + + # ------------------------------------------------------------------ + # Batch simulation with opponent model integration + # ------------------------------------------------------------------ + + def run_auction_simulation( + self, + player_groups: List[pd.DataFrame], + initial_budgets: List[float], + opponent_states: List[OpponentState], + n_sims: int = 500, + ) -> dict: + """Simulate multi-round auction against multiple opponents. + + Args: + player_groups: list of DataFrames, one per round, with available players. + initial_budgets: my starting budget per round/concept. + opponent_states: OpponentState for each round. + n_sims: number of Monte Carlo runs. + + Returns: + dict with aggregated simulation results. + """ + all_results = [] + for idx, (players, budget, opp_state) in enumerate( + zip(player_groups, initial_budgets, opponent_states) + ): + result = self.simulate_live_round( + players, budget, opp_state, n_sims=n_sims + ) + result["round"] = idx + all_results.append(result) + + total_acquired = sum(r["expected_players_acquired"] for r in all_results) + total_cost = sum(r["expected_cost"] for r in all_results) + + return { + "rounds": all_results, + "total_expected_acquired": total_acquired, + "total_expected_cost": total_cost, + "n_sims": n_sims, + } diff --git a/src/optimization/rl_auction_agent.py b/src/optimization/rl_auction_agent.py index 07c289e..e6180c6 100644 --- a/src/optimization/rl_auction_agent.py +++ b/src/optimization/rl_auction_agent.py @@ -515,8 +515,20 @@ class RLAuctionPolicy: src = getattr(self.q_network, src_name) setattr(self.target_network, tgt_name, src.copy()) + def _resize_networks(self, new_state_dim: int): + """Reinitialize networks when state dimension changes.""" + self.q_network = QNetwork(new_state_dim, self.q_network.hidden_dim, self.action_dim) + self.target_network = QNetwork(new_state_dim, self.target_network.hidden_dim, self.action_dim) + self._hard_update_target() + def _normalize_state(self, state: np.ndarray) -> np.ndarray: state = np.asarray(state, dtype=np.float64).ravel() + if len(state) != len(self._obs_mean): + self._obs_mean = np.zeros(len(state), dtype=np.float64) + self._obs_std = np.ones(len(state), dtype=np.float64) + self._obs_count = 0 + self.state_dim = len(state) + self._resize_networks(len(state)) self._obs_count += 1 n = self._obs_count old_mean = self._obs_mean.copy() @@ -844,10 +856,31 @@ def step_in_env(env: AuctionEnv, action: int) -> Tuple[np.ndarray, float, bool, class RLAuctionTrainer: """Convenience class for training and evaluating the RL auction agent.""" - def __init__(self, model_dir: str = "models_trained"): + def __init__( + self, + player_pool: Optional[pd.DataFrame] = None, + n_opponents: int = 7, + config: Optional[AuctionConfig] = None, + model_dir: str = "models_trained", + ): + self.player_pool = player_pool + self.n_opponents = n_opponents + self.config = config or AuctionConfig() self.model_dir = Path(model_dir) self.model_dir.mkdir(parents=True, exist_ok=True) + def train( + self, + n_episodes: int = 5000, + verbose: bool = True, + ) -> RLAuctionPolicy: + """Train RL policy on stored player pool.""" + if self.player_pool is None: + raise ValueError("No player_pool provided to trainer") + env = self.prepare_training_data(self.player_pool, n_opponents=self.n_opponents, config=self.config) + policy, _ = self.train_agent(env, episodes=n_episodes, eval_interval=100) + return policy + def prepare_training_data( self, player_pool_df: pd.DataFrame, @@ -963,6 +996,12 @@ class RLAuctionTrainer: "n_successful": len(values), } + for strategy, metrics in list(summary.items()): + summary[f"{strategy}_total_value"] = metrics["avg_value"] + + summary["rl_total_value"] = summary.get("rl_agent_total_value", 0) + summary["greedy_total_value"] = summary.get("greedy_baseline_total_value", 0) + logger.info( f"Benchmark complete: RL={summary.get('rl_agent', {}).get('avg_value', 0):.1f} pts " f"vs Greedy={summary.get('greedy_baseline', {}).get('avg_value', 0):.1f} " diff --git a/tests/test_new_models.py b/tests/test_new_models.py new file mode 100644 index 0000000..75cd208 --- /dev/null +++ b/tests/test_new_models.py @@ -0,0 +1,1216 @@ +"""Tests for the new Phase 1-6 ML and optimization modules.""" + +import numpy as np +import pandas as pd +import pytest + + +# ============================================================================= +# Phase 1: Quantile Ensemble + Survival Model +# ============================================================================= + + +class TestQuantileEnsemble: + def test_fit_predict(self): + from src.models.quantile_model import QuantileEnsemble + + np.random.seed(42) + X = pd.DataFrame(np.random.randn(200, 8)) + X.columns = [f"f{i}" for i in range(8)] + y = pd.Series(np.random.randn(200) * 1.5 + 6.5) + + model = QuantileEnsemble(n_estimators=50) + model.fit(X, y) + + preds = model.predict(X) + assert "P10" in preds + assert "P50" in preds + assert "P90" in preds + assert len(preds["P10"]) == len(y) + assert len(preds["P50"]) == len(y) + assert len(preds["P90"]) == len(y) + + point_preds = model.predict_points(X) + assert len(point_preds) == len(y) + + def test_custom_quantiles(self): + from src.models.quantile_model import QuantileEnsemble + + np.random.seed(42) + X = pd.DataFrame(np.random.randn(100, 5)) + y = pd.Series(np.random.randn(100) + 6.5) + + model = QuantileEnsemble(quantiles=(0.05, 0.50, 0.95)) + model.fit(X, y) + preds = model.predict(X) + assert "P5" in preds + assert "P50" in preds + assert "P95" in preds + + def test_risk_functions(self): + from src.models.quantile_model import QuantileEnsemble + + np.random.seed(42) + X = pd.DataFrame(np.random.randn(100, 5)) + y = pd.Series(np.random.randn(100) * 1.5 + 6.5) + + model = QuantileEnsemble(n_estimators=50) + model.fit(X, y) + + downside = model.predict_downside_risk(X, threshold=5.5) + assert len(downside) == len(y) + assert np.all((downside >= 0) & (downside <= 1)) + + upside = model.predict_upside(X, threshold=7.0) + assert len(upside) == len(y) + assert np.all((upside >= 0) & (upside <= 1)) + + var_safe = model.value_at_risk_safe(X, confidence=0.90) + assert len(var_safe) == len(y) + + def test_save_load(self, tmp_path): + from src.models.quantile_model import QuantileEnsemble + + np.random.seed(42) + X = pd.DataFrame(np.random.randn(50, 5)) + y = pd.Series(np.random.randn(50) + 6.5) + + model = QuantileEnsemble(n_estimators=30, model_dir=str(tmp_path)) + model.fit(X, y) + model.save("quantile_test.pkl") + + loaded = QuantileEnsemble.load("quantile_test.pkl", model_dir=str(tmp_path)) + preds_orig = model.predict(X) + preds_loaded = loaded.predict(X) + np.testing.assert_array_almost_equal(preds_orig["P50"], preds_loaded["P50"]) + + +class TestMinutesSurvivalModel: + def test_fit_predict_scipy(self): + from src.models.survival_model import MinutesSurvivalModel + + np.random.seed(42) + n = 200 + X = pd.DataFrame({ + "minutes_last_3": np.random.uniform(0, 90, n), + "games_last_5": np.random.randint(1, 6, n), + "rest_days": np.random.uniform(2, 10, n), + "fatigue_rolling_3": np.random.uniform(0, 90, n), + "age": np.random.uniform(18, 38, n), + "noise": np.random.randn(n), + }) + durations = np.clip(np.random.normal(65, 20, n), 1, 90) + events = (np.random.rand(n) < 0.6).astype(int) + + model = MinutesSurvivalModel(force_scipy=True) + model.fit(X, durations, events) + + expected, lower, upper = model.predict_distribution(X) + assert len(expected) == n + assert len(lower) == n + assert len(upper) == n + assert np.all(expected >= 0) + assert np.all(expected <= 90) + assert np.all(lower <= expected) + assert np.all(expected <= upper) + + def test_starter_probability(self): + from src.models.survival_model import MinutesSurvivalModel + + np.random.seed(42) + n = 100 + X = pd.DataFrame({ + "minutes_last_3": np.random.uniform(30, 90, n), + "games_last_5": np.random.randint(1, 6, n), + "rest_days": np.random.uniform(3, 10, n), + }) + durations = np.clip(np.random.normal(70, 15, n), 1, 90) + events = (np.random.rand(n) < 0.7).astype(int) + + model = MinutesSurvivalModel(force_scipy=True) + model.fit(X, durations, events) + + starter_probs = model.predict_starter_probability(X, min_minutes=60) + assert len(starter_probs) == n + assert np.all((starter_probs >= 0) & (starter_probs <= 1)) + + full_probs = model.predict_full_match_probability(X) + assert len(full_probs) == n + assert np.all((full_probs >= 0) & (full_probs <= 1)) + + expected_minutes = model.predict_expected_minutes(X) + assert len(expected_minutes) == n + assert np.all((expected_minutes >= 0) & (expected_minutes <= 90)) + + def test_feature_selection(self): + from src.models.survival_model import MinutesSurvivalModel + + np.random.seed(42) + n = 80 + X = pd.DataFrame({ + "minutes_played": np.random.uniform(0, 90, n), + "games_started": np.random.randint(0, 5, n), + "fatigue_level": np.random.uniform(0, 1, n), + "rest_between_games": np.random.uniform(2, 7, n), + "player_age": np.random.uniform(20, 35, n), + "irrelevant_cat": ["A"] * 40 + ["B"] * 40, + }) + durations = np.clip(np.random.normal(60, 20, n), 1, 90) + events = np.random.binomial(1, 0.6, n) + + model = MinutesSurvivalModel(force_scipy=True) + model.fit(X, durations, events) + expected = model.predict_expected_minutes(X) + assert len(expected) == n + + +# ============================================================================= +# Phase 4: Bayesian Pooling + Conformal Predictor +# ============================================================================= + + +class TestBayesianPlayerModel: + def test_fit_predict_scipy(self): + from src.models.bayesian_pooling import BayesianPlayerModel + + np.random.seed(42) + n = 150 + X = pd.DataFrame({ + "role": ["P"] * 15 + ["D"] * 45 + ["C"] * 45 + ["A"] * 45, + "feature1": np.random.randn(n), + "feature2": np.random.randn(n), + }) + y = pd.Series(np.random.randn(n) * 1.0 + 6.5) + + model = BayesianPlayerModel() + model.fit(X, y) + + preds = model.predict(X) + assert len(preds) == n + + def test_predict_with_uncertainty(self): + from src.models.bayesian_pooling import BayesianPlayerModel + + np.random.seed(42) + n = 100 + X = pd.DataFrame({ + "role": ["D"] * 40 + ["A"] * 60, + "feature": np.random.randn(n), + }) + y = pd.Series(np.random.randn(n) + 6.5) + + model = BayesianPlayerModel() + model.fit(X, y) + + mean, std = model.predict_with_uncertainty(X) + assert len(mean) == n + assert len(std) == n + assert np.all(std >= 0) + + def test_reliability_scores(self): + from src.models.bayesian_pooling import BayesianPlayerModel + + np.random.seed(42) + n = 120 + X = pd.DataFrame({ + "role": ["D"] * 60 + ["C"] * 60, + "feature": np.random.randn(n), + }) + y = pd.Series(np.random.randn(n) + 6.5) + + model = BayesianPlayerModel() + model.fit(X, y) + + reliability = model.get_player_reliability(X) + assert len(reliability) == n + assert np.all((reliability > 0) & (reliability <= 1)) + + def test_rookie_estimates(self): + from src.models.bayesian_pooling import BayesianPlayerModel + + np.random.seed(42) + X = pd.DataFrame({ + "role": ["P"] * 10 + ["D"] * 30 + ["C"] * 30, + "feature": np.random.randn(70), + }) + y = pd.Series(np.random.randn(70) + 6.5) + + model = BayesianPlayerModel() + model.fit(X, y) + + rookies = model.get_rookie_estimates(X, min_observations=20) + assert isinstance(rookies, pd.DataFrame) + for col in ["player_estimate", "role_mean"]: + assert col in rookies.columns + + def test_posterior_predictive(self): + from src.models.bayesian_pooling import BayesianPlayerModel + + np.random.seed(42) + X = pd.DataFrame({ + "role": ["D"] * 30 + ["A"] * 20, + "feature": np.random.randn(50), + }) + y = pd.Series(np.random.randn(50) + 6.5) + + model = BayesianPlayerModel() + model.fit(X, y) + + samples = model.posterior_predictive(X.iloc[:10], n_samples=100) + assert samples.shape == (100, 10) + + +class TestConformalPredictor: + def test_calibrate_and_predict(self): + from src.models.conformal_predictor import ConformalPredictor + from sklearn.linear_model import Ridge + + np.random.seed(42) + n = 200 + X = np.random.randn(n, 5) + y = X[:, 0] * 2 + X[:, 1] * 1.5 + np.random.randn(n) * 0.5 + + X_train, X_cal, X_test = ( + pd.DataFrame(X[:80]), + pd.DataFrame(X[80:150]), + pd.DataFrame(X[150:]), + ) + y_train, y_cal, y_test = y[:80], y[80:150], y[150:] + + base = Ridge(alpha=1.0) + base.fit(X_train, y_train) + + cp = ConformalPredictor(base, alpha=0.10) + cp.calibrate(X_cal, y_cal) + + y_pred, y_lower, y_upper = cp.predict_with_band(X_test) + assert len(y_pred) == len(y_test) + assert len(y_lower) == len(y_test) + assert len(y_upper) == len(y_test) + assert np.all(y_lower <= y_upper) + + def test_coverage(self): + from src.models.conformal_predictor import ConformalPredictor + from sklearn.linear_model import Ridge + + np.random.seed(42) + n = 300 + X = np.random.randn(n, 5) + y = X[:, 0] * 2 + X[:, 1] * 1.5 + np.random.randn(n) * 0.5 + + X_train = pd.DataFrame(X[:100]) + X_cal = pd.DataFrame(X[100:200]) + X_test = pd.DataFrame(X[200:]) + y_train, y_cal, y_test = y[:100], y[100:200], y[200:] + + base = Ridge(alpha=1.0) + base.fit(X_train, y_train) + + cp = ConformalPredictor(base, alpha=0.20) + cp.calibrate(X_cal, y_cal) + + cov = cp.coverage(X_test, y_test) + assert 0.5 < cov < 1.0 + + def test_is_inside_band(self): + from src.models.conformal_predictor import ConformalPredictor + from sklearn.linear_model import LinearRegression + + np.random.seed(42) + X = pd.DataFrame(np.random.randn(100, 3)) + y = X[0] * 1.5 + np.random.randn(100) * 0.3 + + base = LinearRegression() + base.fit(X, y) + + cp = ConformalPredictor(base, alpha=0.10) + cp.calibrate(X, y) + inside = cp.is_inside_band(X, y) + assert len(inside) == len(y) + assert np.all((inside == 0) | (inside == 1)) + + def test_update_calibration(self): + from src.models.conformal_predictor import ConformalPredictor + from sklearn.linear_model import Ridge + + np.random.seed(42) + X = pd.DataFrame(np.random.randn(200, 3)) + y = X[0] * 2 + np.random.randn(200) * 0.5 + + base = Ridge() + base.fit(X, y) + + cp = ConformalPredictor(base, alpha=0.10) + + X_cal1, X_cal2 = X.iloc[:100], X.iloc[100:] + y_cal1, y_cal2 = y[:100], y[100:] + + cp.calibrate(X_cal1, y_cal1) + width_before = cp.predict_interval_width(X) + cp.update(X_cal2, y_cal2) + width_after = cp.predict_interval_width(X) + assert isinstance(width_before, float) + assert isinstance(width_after, float) + + +# ============================================================================= +# Phase 2: Bandit, Opponent Bidding, Budget Optimizer +# ============================================================================= + + +class TestBanditAuctionSolver: + def test_initialization(self): + from src.optimization.bandit_auction import BanditAuctionSolver + + solver = BanditAuctionSolver() + assert solver is not None + stats = solver.get_arm_stats() + assert len(stats) > 0 + + def test_select_bid_and_update(self): + from src.optimization.bandit_auction import BanditAuctionSolver + from src.optimization.auction_solver import AuctionConfig, PlayerValuation + + config = AuctionConfig(total_budget=500) + solver = BanditAuctionSolver(config=config) + + player = PlayerValuation( + name="TestPlayer", team="TeamA", role="C", + projected_points=7.5, market_value=15, ceiling_price=50, + ) + + auction_state = { + "budget_remaining": 400, + "total_budget": 500, + "slots_remaining": {"P": 1, "D": 3, "C": 4, "A": 3}, + "role_quotas": {"P": 3, "D": 8, "C": 8, "A": 6}, + "slot_quotas": {"P": 3, "D": 8, "C": 8, "A": 6}, + "round_number": 3, + "total_rounds": 10, + "opponent_budgets": [400, 450, 350], + "players_remaining_in_role": {"P": 5, "D": 10, "C": 8, "A": 6}, + "player_pool": [player], + } + + arm_idx, bid_amount = solver.select_bid(player, auction_state) + assert 0 <= arm_idx < len(solver.bid_fractions) + assert bid_amount >= 0 + + solver.update(arm_idx, reward=0.8, player_role="C") + stats = solver.get_arm_stats() + assert stats[arm_idx]["trials"] >= 1 + + def test_recommend_bid_summary(self): + from src.optimization.bandit_auction import BanditAuctionSolver + from src.optimization.auction_solver import AuctionConfig, PlayerValuation + + solver = BanditAuctionSolver() + player = PlayerValuation( + "TestPlayer", "TeamA", "A", 8.0, 20, 60, + ) + auction_state = { + "budget_remaining": 300, + "total_budget": 500, + "slots_remaining": {"P": 2, "D": 5, "C": 5, "A": 2}, + "role_quotas": {"P": 3, "D": 8, "C": 8, "A": 6}, + "slot_quotas": {"P": 3, "D": 8, "C": 8, "A": 6}, + "round_number": 1, + "total_rounds": 10, + "opponent_budgets": [400], + "players_remaining_in_role": {"P": 6, "D": 15, "C": 15, "A": 10}, + "player_pool": [player], + } + summary = solver.recommend_bid_summary(player, auction_state) + assert "recommended_bid" in summary + + +class TestOpponentBidModel: + def test_fit_predict_heuristic(self): + from src.optimization.opponent_bidding_model import OpponentBidModel + + np.random.seed(42) + model = OpponentBidModel() + + players_df = pd.DataFrame({ + "name": [f"Player_{i}" for i in range(20)], + "role": np.random.choice(["P", "D", "C", "A"], 20), + "projected_points": np.random.uniform(5, 9, 20), + }) + + opponent_state = { + "budget_remaining": 400, + "total_budget": 500, + "slots_remaining": {"P": 2, "D": 6, "C": 6, "A": 4}, + "role_quotas": {"P": 3, "D": 8, "C": 8, "A": 6}, + } + + bids = model.predict_opponent_bids(players_df, opponent_state) + assert len(bids) == 20 + assert np.all(bids >= 0) + + def test_p_acquire(self): + from src.optimization.opponent_bidding_model import OpponentBidModel + + np.random.seed(42) + model = OpponentBidModel() + + players_df = pd.DataFrame({ + "name": ["Player_A", "Player_B", "Player_C"], + "role": ["A", "D", "C"], + "projected_points": [8.5, 7.0, 6.5], + }) + + opponent_state = { + "budget_remaining": 350, + "total_budget": 500, + "slots_remaining": {"P": 2, "D": 5, "C": 5, "A": 2}, + "role_quotas": {"P": 3, "D": 8, "C": 8, "A": 6}, + } + + my_bids = pd.Series([40, 25, 15], index=players_df.index) + probs = model.predict_p_acquire(players_df, my_bids, opponent_state) + assert len(probs) == 3 + assert np.all((probs >= 0) & (probs <= 1)) + + def test_simulate_live_round(self): + from src.optimization.opponent_bidding_model import OpponentBidModel + + np.random.seed(42) + model = OpponentBidModel() + + players_df = pd.DataFrame({ + "name": [f"Player_{i}" for i in range(5)], + "role": ["P"] + ["D"] * 2 + ["C"] + ["A"], + "projected_points": np.random.uniform(5, 9, 5), + }) + + opponent_state = { + "budget_remaining": 350, + "total_budget": 500, + "slots_remaining": {"P": 2, "D": 6, "C": 6, "A": 4}, + "role_quotas": {"P": 3, "D": 8, "C": 8, "A": 6}, + } + + result = model.simulate_live_round(players_df, 250, opponent_state, n_sims=50) + assert "expected_cost" in result + assert "value_matrix" in result + assert result["expected_cost"] >= 0 + + +class TestBudgetOptimizer: + def test_optimize_fallback(self): + from src.optimization.budget_optimizer import BudgetOptimizer + + np.random.seed(42) + optimizer = BudgetOptimizer(total_budget=500) + + players = [] + for i in range(60): + role = np.random.choice(["P", "D", "C", "A"]) + players.append({ + "name": f"Player_{i}", + "role": role, + "projected_points": np.random.uniform(5, 9), + "market_value": np.random.randint(3, 40), + }) + player_pool = pd.DataFrame(players) + + allocation = optimizer.optimize(player_pool, n_calls=10) + assert "P" in allocation + assert "D" in allocation + assert "C" in allocation + assert "A" in allocation + total = sum(allocation.values()) + assert abs(total - 500) < 10 + + def test_optimize_adaptive(self): + from src.optimization.budget_optimizer import BudgetOptimizer + + np.random.seed(42) + optimizer = BudgetOptimizer(total_budget=500) + + players = [] + for i in range(50): + players.append({ + "name": f"Player_{i}", + "role": np.random.choice(["P", "D", "C", "A"]), + "projected_points": np.random.uniform(5, 9), + "market_value": np.random.randint(3, 40), + }) + player_pool = pd.DataFrame(players) + + allocation = optimizer.optimize_adaptive( + player_pool, + remaining_slots={"P": 2, "D": 4, "C": 5, "A": 3}, + spent_per_role={"P": 10, "D": 60, "C": 40, "A": 30}, + n_calls=10, + ) + assert isinstance(allocation, dict) + for role in ["P", "D", "C", "A"]: + assert allocation[role] >= 0 + + def test_role_value_curves(self): + from src.optimization.budget_optimizer import BudgetOptimizer + + np.random.seed(42) + optimizer = BudgetOptimizer(total_budget=500) + + players = [] + for i in range(80): + role = np.random.choice(["P", "D", "C", "A"]) + players.append({ + "name": f"Player_{i}", + "role": role, + "projected_points": np.random.uniform(5, 9), + "market_value": np.random.randint(3, 40), + }) + player_pool = pd.DataFrame(players) + + curves = optimizer.get_role_value_curves(player_pool) + for role in ["P", "D", "C", "A"]: + assert role in curves + assert len(curves[role]) > 0 + + def test_sensitivity_analysis(self): + from src.optimization.budget_optimizer import BudgetOptimizer + + np.random.seed(42) + optimizer = BudgetOptimizer(total_budget=500) + + players = [] + for i in range(60): + players.append({ + "name": f"Player_{i}", + "role": np.random.choice(["P", "D", "C", "A"]), + "projected_points": np.random.uniform(5, 9), + "market_value": np.random.randint(3, 40), + }) + player_pool = pd.DataFrame(players) + + analysis = optimizer.sensitivity_analysis(player_pool, n_calls=5) + assert isinstance(analysis, dict) + + +# ============================================================================= +# Phase 5: GAT Chemistry + Hawkes Form +# ============================================================================= + + +class TestPlayerChemistryGAT: + def test_build_graph(self): + from src.models.gat_model import PlayerChemistryGAT + + np.random.seed(42) + data = pd.DataFrame([ + {"player": "A", "teammate": "B", "passes_to": 10, "assists_to": 2, "crosses_to": 3, "matchday": 1}, + {"player": "A", "teammate": "C", "passes_to": 15, "assists_to": 1, "crosses_to": 5, "matchday": 1}, + {"player": "B", "teammate": "A", "passes_to": 8, "assists_to": 0, "crosses_to": 2, "matchday": 1}, + {"player": "B", "teammate": "C", "passes_to": 12, "assists_to": 3, "crosses_to": 4, "matchday": 2}, + ]) + + model = PlayerChemistryGAT() + model.build_graph(data) + assert len(model.interaction_edges) >= 0 + + def test_extract_interaction_features(self): + from src.models.gat_model import PlayerChemistryGAT + + np.random.seed(42) + data = pd.DataFrame([ + {"player": "A", "teammate": "B", "passes_to": 10, "assists_to": 2, "crosses_to": 3, "matchday": 1}, + {"player": "A", "teammate": "C", "passes_to": 5, "assists_to": 0, "crosses_to": 1, "matchday": 1}, + {"player": "B", "teammate": "C", "passes_to": 3, "assists_to": 1, "crosses_to": 0, "matchday": 2}, + ]) + + model = PlayerChemistryGAT() + model.build_graph(data) + + features = model.extract_interaction_features("A", ["B", "C"]) + assert "interaction_outgoing_sum" in features + assert "interaction_incoming_sum" in features + assert "interaction_synergy" in features + + def test_compute_interaction_bonus(self): + from src.models.gat_model import PlayerChemistryGAT + + data = pd.DataFrame([ + {"player": "A", "teammate": "B", "passes_to": 20, "assists_to": 4, "crosses_to": 8, "matchday": 1}, + {"player": "A", "teammate": "C", "passes_to": 5, "assists_to": 0, "crosses_to": 1, "matchday": 1}, + ]) + + model = PlayerChemistryGAT() + model.build_graph(data) + + bonus_high = model.compute_interaction_bonus("A", "B") + bonus_low = model.compute_interaction_bonus("A", "C") + assert bonus_high > bonus_low + assert 0 <= bonus_high <= 1 + assert 0 <= bonus_low <= 1 + + def test_fit_predict_sklearn(self): + from src.models.gat_model import PlayerChemistryGAT + + np.random.seed(42) + n = 100 + data = [] + n_real = np.random.RandomState(42) + for i in range(n): + for j in range(n): + if i < j and n_real.random() < 0.05: + data.append({ + "player": f"P{i}", "teammate": f"P{j}", + "passes_to": n_real.randint(0, 10), + "assists_to": n_real.randint(0, 3), + "crosses_to": n_real.randint(0, 5), + "matchday": n_real.randint(1, 20), + }) + + X = pd.DataFrame({ + "player": [f"P{i}" for i in range(n)], + "team": ["T1"] * (n // 2) + ["T2"] * (n - n // 2), + "feature1": np.random.randn(n), + "feature2": np.random.randn(n), + }) + y = pd.Series(np.random.randn(n) * 1.5 + 6.5) + + model = PlayerChemistryGAT() + model.build_graph(pd.DataFrame(data)) + model.fit(X, y) + + preds = model.predict(X) + assert len(preds) == n + + def test_get_redundancy_penalty(self): + from src.models.gat_model import PlayerChemistryGAT + + model = PlayerChemistryGAT() + players = ["Player_A", "Player_B", "Player_C"] + penalties = model.get_redundancy_penalty(players) + assert isinstance(penalties, dict) + + +class TestPlayerFormModel: + def test_fit_predict(self): + from src.models.hawkes_form import PlayerFormModel + + np.random.seed(42) + n = 150 + base_dates = pd.date_range("2023-08-20", periods=38, freq="7D") + X = pd.DataFrame({ + "player": [f"P{i % 10}" for i in range(n)], + "match_date": np.random.choice(base_dates, n), + "minutes": np.random.uniform(0, 90, n), + "feature1": np.random.randn(n), + }) + y = pd.Series(np.random.randn(n) * 1.5 + 6.5) + + model = PlayerFormModel() + model.fit(X, y) + preds = model.predict(X) + assert len(preds) == n + + def test_form_status(self): + from src.models.hawkes_form import PlayerFormModel + + np.random.seed(42) + base_dates = pd.date_range("2023-08-20", periods=20, freq="7D") + X = pd.DataFrame({ + "player": ["TestPlayer"] * 20, + "match_date": base_dates, + "minutes": np.random.uniform(30, 90, 20), + "feature1": np.random.randn(20), + }) + y = pd.Series(np.random.randn(20) + 6.5) + + model = PlayerFormModel() + model.fit(X, y) + + form_status = model.get_form_status(X) + assert len(form_status) == 20 + for status in form_status: + assert status in ("HOT", "COLD", "NEUTRAL") + + def test_detect_streak(self): + from src.models.hawkes_form import PlayerFormModel + + np.random.seed(42) + base_dates = pd.date_range("2023-08-20", periods=15, freq="7D") + X = pd.DataFrame({ + "player": ["StreakyP"] * 15, + "match_date": base_dates, + "minutes": np.random.uniform(50, 90, 15), + "feature1": np.random.randn(15), + }) + y = pd.Series(np.random.randn(15) * 1.5 + 7.5) + + model = PlayerFormModel() + model.fit(X, y) + + is_streak, length, direction = model.detect_streak(X) + assert isinstance(is_streak, bool) + assert isinstance(length, int) + assert direction in ("HOT_STREAK", "COLD_STREAK", "NO_STREAK") + + def test_momentum_projection(self): + from src.models.hawkes_form import PlayerFormModel + + np.random.seed(42) + base_dates = pd.date_range("2023-08-20", periods=10, freq="7D") + X = pd.DataFrame({ + "player": ["Player1"] * 10, + "match_date": base_dates, + "minutes": np.random.uniform(40, 90, 10), + "feature1": np.random.randn(10), + }) + y = pd.Series(np.random.randn(10) * 1.5 + 6.5) + + model = PlayerFormModel() + model.fit(X, y) + + momentum = model.predict_momentum(X.iloc[:5], X.iloc[:5], n_future=3) + assert momentum.shape == (3, 5) + + +# ============================================================================= +# Phase 3: RL Auction Agent + Set Transformer +# ============================================================================= + + +class TestAuctionEnv: + def test_reset_and_step(self): + from src.optimization.rl_auction_agent import AuctionEnv + from src.optimization.auction_solver import AuctionConfig + + np.random.seed(42) + config = AuctionConfig(total_budget=500) + config.n_gk = 2 + config.n_def = 5 + config.n_mid = 5 + config.n_fwd = 3 + + players = [] + for i in range(40): + players.append({ + "name": f"Player_{i}", + "role": np.random.choice(["P", "D", "C", "A"]), + "projected_points": np.random.uniform(5, 9), + "market_value": np.random.randint(3, 40), + "team": f"Team_{np.random.randint(1, 21)}", + "ceiling_price": np.random.uniform(20, 80), + }) + player_pool = pd.DataFrame(players) + + env = AuctionEnv(player_pool, n_opponents=3, config=config) + obs = env.reset() + + assert isinstance(obs, np.ndarray) + assert len(obs) > 0 + + obs, reward, done, info = env.step(5) + assert isinstance(obs, np.ndarray) + assert isinstance(reward, float) + assert isinstance(done, bool) + + def test_valid_actions(self): + from src.optimization.rl_auction_agent import AuctionEnv + from src.optimization.auction_solver import AuctionConfig + + np.random.seed(42) + config = AuctionConfig(total_budget=500) + config.n_gk = 1 + config.n_def = 2 + config.n_mid = 2 + config.n_fwd = 1 + + players = [] + for i in range(15): + players.append({ + "name": f"Player_{i}", + "role": np.random.choice(["P", "D", "C", "A"]), + "projected_points": np.random.uniform(5, 9), + "market_value": np.random.randint(3, 40), + "team": f"Team_{np.random.randint(1, 21)}", + "ceiling_price": np.random.uniform(20, 80), + }) + player_pool = pd.DataFrame(players) + + env = AuctionEnv(player_pool, n_opponents=2, config=config) + env.reset() + + valid = env.get_valid_actions() + assert len(valid) > 0 + + env_no_budget = AuctionEnv(player_pool, n_opponents=2, config=config) + env_no_budget.reset() + env_no_budget.budget_remaining = 0 + valid_no_money = env_no_budget.get_valid_actions() + assert valid_no_money[0] == 0 + + +class TestRLAuctionPolicy: + def test_act_and_remember(self): + from src.optimization.rl_auction_agent import RLAuctionPolicy + + policy = RLAuctionPolicy(state_dim=14, action_dim=11) + + state = np.random.randn(14).astype(np.float32) + action = policy.act(state, epsilon=1.0) + assert 0 <= action < 11 + + next_state = np.random.randn(14).astype(np.float32) + policy.remember(state, action, reward=1.5, next_state=next_state, done=False) + + assert len(policy.replay_buffer) == 1 + + def test_replay_and_target_update(self): + from src.optimization.rl_auction_agent import RLAuctionPolicy + + policy = RLAuctionPolicy(state_dim=14, action_dim=11, hidden_dim=64) + + for _ in range(256): + s = np.random.randn(14).astype(np.float32) + a = np.random.randint(0, 11) + r = np.random.randn() + ns = np.random.randn(14).astype(np.float32) + d = np.random.rand() < 0.3 + policy.remember(s, a, r, ns, d) + + loss = policy.replay(batch_size=64) + assert loss is not None + assert isinstance(loss, float) + + def test_save_load(self, tmp_path): + from src.optimization.rl_auction_agent import RLAuctionPolicy + + policy = RLAuctionPolicy(state_dim=14, action_dim=11, hidden_dim=32) + save_path = str(tmp_path / "test_policy.pkl") + + policy.save(save_path) + + loaded = RLAuctionPolicy.load(save_path) + assert loaded.state_dim == 14 + assert loaded.action_dim == 11 + + def test_training_loop(self): + from src.optimization.rl_auction_agent import ( + RLAuctionPolicy, AuctionEnv, RLAuctionTrainer, + ) + from src.optimization.auction_solver import AuctionConfig + + np.random.seed(42) + config = AuctionConfig(total_budget=500) + config.n_gk = 1 + config.n_def = 3 + config.n_mid = 3 + config.n_fwd = 2 + + players = [] + for i in range(20): + players.append({ + "name": f"Player_{i}", + "role": np.random.choice(["P", "D", "C", "A"]), + "projected_points": np.random.uniform(5, 9), + "market_value": np.random.randint(3, 40), + "team": f"Team_{np.random.randint(1, 10)}", + "ceiling_price": np.random.uniform(20, 80), + }) + player_pool = pd.DataFrame(players) + + trainer = RLAuctionTrainer(player_pool, n_opponents=2, config=config) + agent = trainer.train(n_episodes=20, verbose=False) + + assert agent is not None + assert hasattr(agent, "q_network") + + +class TestRLAuctionTrainer: + def test_evaluate_vs_baselines(self): + from src.optimization.rl_auction_agent import ( + RLAuctionPolicy, AuctionEnv, RLAuctionTrainer, + ) + from src.optimization.auction_solver import AuctionConfig + + np.random.seed(42) + config = AuctionConfig(total_budget=500) + config.n_gk = 1 + config.n_def = 2 + config.n_mid = 2 + config.n_fwd = 1 + + players = [] + for i in range(15): + players.append({ + "name": f"Player_{i}", + "role": np.random.choice(["P", "D", "C", "A"]), + "projected_points": np.random.uniform(5, 9), + "market_value": np.random.randint(3, 40), + "team": f"Team_{np.random.randint(1, 10)}", + "ceiling_price": np.random.uniform(20, 80), + }) + player_pool = pd.DataFrame(players) + + agent = RLAuctionPolicy(state_dim=14, action_dim=11, hidden_dim=32) + trainer = RLAuctionTrainer(player_pool, n_opponents=2, config=config) + results = trainer.evaluate_vs_baselines(agent, player_pool, n_sims=5) + assert "rl_total_value" in results or "greedy_total_value" in results + + +class TestSetTransformer: + def test_fit_predict_sklearn(self): + from src.models.set_transformer import SetTransformer + + np.random.seed(42) + n_teams = 20 + n_players_per_team = 25 + + teams_data = [] + team_values = [] + for t in range(n_teams): + df = pd.DataFrame({ + "name": [f"Team{t}_Player_{i}" for i in range(n_players_per_team)], + "role": np.random.choice(["P", "D", "C", "A"], n_players_per_team), + "feature1": np.random.randn(n_players_per_team), + "feature2": np.random.randn(n_players_per_team), + "projected_points": np.random.uniform(5, 9, n_players_per_team), + }) + teams_data.append(df) + team_values.append(np.sum(df["projected_points"]) + np.random.randn() * 10) + + model = SetTransformer(use_torch=False) + model.fit(teams_data, team_values) + + pred = model.predict(teams_data[0]) + assert isinstance(pred, float) + + def test_value_added_and_removed(self): + from src.models.set_transformer import SetTransformer + + np.random.seed(42) + n_teams = 15 + n_players = 25 + + teams_data = [] + team_values = [] + for t in range(n_teams): + df = pd.DataFrame({ + "name": [f"T{t}_P{i}" for i in range(n_players)], + "role": np.random.choice(["P", "D", "C", "A"], n_players), + "feature1": np.random.randn(n_players), + "projected_points": np.random.uniform(5, 9, n_players), + }) + teams_data.append(df) + team_values.append(np.sum(df["projected_points"]) + np.random.randn() * 5) + + model = SetTransformer(use_torch=False) + model.fit(teams_data, team_values) + + new_player = pd.DataFrame([{ + "name": "NewPlayer", + "role": "A", + "feature1": 1.5, + "projected_points": 8.5, + }]) + + va = model.value_added(teams_data[0], new_player) + vr = model.value_removed(teams_data[0], "T0_P0") + assert isinstance(va, float) + assert isinstance(vr, float) + + def test_optimal_replacement(self): + from src.models.set_transformer import SetTransformer + + np.random.seed(42) + n_teams = 10 + n_players = 20 + + teams_data = [] + team_values = [] + for t in range(n_teams): + df = pd.DataFrame({ + "name": [f"T{t}_P{i}" for i in range(n_players)], + "role": np.random.choice(["P", "D", "C", "A"], n_players), + "feature1": np.random.randn(n_players), + "projected_points": np.random.uniform(5, 9, n_players), + }) + teams_data.append(df) + team_values.append(np.sum(df["projected_points"]) + np.random.randn() * 5) + + model = SetTransformer(use_torch=False) + model.fit(teams_data, team_values) + + pool = pd.DataFrame([ + {"name": f"Free_{i}", "role": np.random.choice(["P", "D", "C", "A"]), + "feature1": np.random.randn(), "projected_points": np.random.uniform(5, 9)} + for i in range(10) + ]) + + rankings = model.optimal_replacement(teams_data[0], pool, to_replace=["T0_P0"]) + assert isinstance(rankings, dict) + assert "T0_P0" in rankings + + def test_redundancy_score(self): + from src.models.set_transformer import SetTransformer + + np.random.seed(42) + n_teams = 10 + n_players = 20 + + teams_data = [] + team_values = [] + for t in range(n_teams): + df = pd.DataFrame({ + "name": [f"T{t}_P{i}" for i in range(n_players)], + "role": np.random.choice(["P", "D", "C", "A"], n_players), + "feature1": np.random.randn(n_players), + "projected_points": np.random.uniform(5, 9, n_players), + }) + teams_data.append(df) + team_values.append(np.sum(df["projected_points"]) + np.random.randn() * 5) + + model = SetTransformer(use_torch=False) + model.fit(teams_data, team_values) + + score = model.get_redundancy_score(teams_data[0]) + assert 0 <= score <= 1 + + +# ============================================================================= +# Phase 6: Causal Forest +# ============================================================================= + + +class TestTransferCausalModel: + def test_fit_predict_sklearn(self): + from src.models.causal_forest import TransferCausalModel + + np.random.seed(42) + n = 200 + X = pd.DataFrame({ + "role": np.random.choice(["P", "D", "C", "A"], n), + "feature1": np.random.randn(n), + "feature2": np.random.randn(n), + "team_strength": np.random.uniform(0.5, 1.5, n), + }) + treatment = pd.DataFrame({ + "role": np.random.choice(["P", "D", "C", "A"], n), + "projected_points": np.random.uniform(5, 9, n), + "days_since_last_transfer": np.random.randint(1, 30, n), + }) + true_effect = treatment["projected_points"] * 0.5 + np.random.randn(n) * 0.3 + outcome = pd.Series(true_effect + np.random.randn(n) * 1.0) + + model = TransferCausalModel() + model.fit(X, treatment, outcome) + + result = model.predict_effect(X, treatment) + assert "ate" in result + assert "cate_lower" in result + + def test_predict_individual_effect(self): + from src.models.causal_forest import TransferCausalModel + + np.random.seed(42) + n = 150 + X = pd.DataFrame({ + "role": np.random.choice(["P", "D", "C", "A"], n), + "feature1": np.random.randn(n), + "team_strength": np.random.uniform(0.5, 1.5, n), + }) + treatment = pd.DataFrame({ + "role": np.random.choice(["P", "D", "C", "A"], n), + "projected_points": np.random.uniform(5, 9, n), + "days_since_last_transfer": np.random.randint(1, 30, n), + }) + outcome = pd.Series(np.random.randn(n) + 6.5) + + model = TransferCausalModel() + model.fit(X, treatment, outcome) + + roster = pd.DataFrame({ + "name": ["P1", "D1", "D2", "D3", "M1", "M2", "M3", "M4", "F1", "F2"], + "role": ["P"] + ["D"] * 3 + ["C"] * 4 + ["A"] * 2, + "projected_points": np.random.uniform(5, 9, 10), + "team": ["TeamA"] * 10, + }) + candidate = pd.DataFrame([{ + "name": "NewPlayer", "role": "C", + "projected_points": 7.5, + "team": "Available", + }]) + + result = model.predict_individual_effect(roster, candidate, role_to_replace="M1") + assert "effect" in result + assert "confidence" in result + + def test_rank_transfers(self): + from src.models.causal_forest import TransferCausalModel + + np.random.seed(42) + n = 120 + X = pd.DataFrame({ + "role": np.random.choice(["P", "D", "C", "A"], n), + "feature1": np.random.randn(n), + "team_strength": np.random.uniform(0.5, 1.5, n), + }) + treatment = pd.DataFrame({ + "role": np.random.choice(["P", "D", "C", "A"], n), + "projected_points": np.random.uniform(5, 9, n), + "days_since_last_transfer": np.random.randint(1, 30, n), + }) + outcome = pd.Series(np.random.randn(n) + 6.5) + + model = TransferCausalModel() + model.fit(X, treatment, outcome) + + roster = pd.DataFrame({ + "name": ["P1"] + [f"D{i}" for i in range(5)] + [f"M{i}" for i in range(5)] + [f"F{i}" for i in range(4)], + "role": ["P"] + ["D"] * 5 + ["C"] * 5 + ["A"] * 4, + "projected_points": np.random.uniform(5, 9, 15), + "team": ["T1"] * 15, + }) + pool = pd.DataFrame([ + {"name": f"Free_{i}", "role": np.random.choice(["P", "D", "C", "A"]), + "projected_points": np.random.uniform(5, 9), "team": "Free"} + for i in range(20) + ]) + + ranked = model.rank_transfers(roster, pool, n_recommendations=5) + assert isinstance(ranked, pd.DataFrame) + assert len(ranked) <= 5 + + def test_auction_effect_analyzer(self): + from src.models.causal_forest import TransferCausalModel, AuctionEffectAnalyzer + + np.random.seed(42) + n = 100 + X = pd.DataFrame({ + "role": np.random.choice(["P", "D", "C", "A"], n), + "feature1": np.random.randn(n), + "team_strength": np.random.uniform(0.5, 1.5, n), + }) + treatment = pd.DataFrame({ + "role": np.random.choice(["P", "D", "C", "A"], n), + "projected_points": np.random.uniform(5, 9, n), + "days_since_last_transfer": np.random.randint(1, 30, n), + }) + outcome = pd.Series(np.random.randn(n) + 6.5) + + model = TransferCausalModel() + model.fit(X, treatment, outcome) + + roster = pd.DataFrame({ + "name": [f"P{i}" for i in range(15)], + "role": np.random.choice(["P", "D", "C", "A"], 15), + "projected_points": np.random.uniform(5, 9, 15), + "team": ["T1"] * 15, + }) + + analyzer = AuctionEffectAnalyzer(model, roster) + assert analyzer is not None + + player = pd.Series({ + "name": "Target", "role": "C", "projected_points": 7.8, "team": "Available", + }) + result = analyzer.recommend_bid_adjustment(player, 25, roster) + assert "adjusted_bid" in result + assert isinstance(result["adjusted_bid"], float)