Capstone - Complete Analytics System

# Python: Player Projection Framework Overview
import pandas as pd
import numpy as np
from scipy import stats

# Define projection system components
projection_components = {
    "historical_performance": "Weighted average of past seasons",
    "aging_curve": "Expected change based on age",
    "regression": "Regression to the mean based on sample size",
    "context": "League/team quality adjustments",
    "uncertainty": "Confidence intervals for projections"
}

# Sample projection workflow
def project_player(player_history, metric="npxg_90"):
    """Project future performance for a player."""
    # 1. Weight recent seasons more heavily
    weighted_avg = calculate_weighted_average(player_history, metric)

    # 2. Apply aging curve adjustment
    age_adjusted = apply_aging_curve(weighted_avg, player_history["age"])

    # 3. Regress to the mean based on minutes played
    regressed = regress_to_mean(age_adjusted, player_history["minutes"])

    # 4. Calculate projection uncertainty
    uncertainty = calculate_uncertainty(player_history, metric)

    return {
        "projection": regressed,
        "lower_95": regressed - 1.96 * uncertainty,
        "upper_95": regressed + 1.96 * uncertainty
    }

print("Projection Framework Components:")
for name, description in projection_components.items():
    print(f"- {name}: {description}")

# Python: Building Aging Curves with Delta Method
import pandas as pd
import numpy as np
from scipy.interpolate import UnivariateSpline

# Load multi-season player data
player_seasons = pd.read_csv("player_season_data.csv")

def build_aging_curve_delta(data, metric):
    """Build aging curve using the delta method."""
    # Sort by player and season
    data = data.sort_values(["player_id", "season"])

    # Create player-season pairs
    paired = []
    for player_id, group in data.groupby("player_id"):
        group = group.reset_index(drop=True)
        for i in range(len(group) - 1):
            current = group.iloc[i]
            next_row = group.iloc[i + 1]
            delta = next_row[metric] - current[metric]
            paired.append({
                "player_id": player_id,
                "age": current["age"],
                "current_value": current[metric],
                "next_value": next_row[metric],
                "delta": delta
            })

    paired_df = pd.DataFrame(paired)

    # Calculate average delta at each age
    age_deltas = paired_df.groupby("age").agg(
        mean_delta=("delta", "mean"),
        se_delta=("delta", lambda x: x.std() / np.sqrt(len(x))),
        n=("delta", "count")
    ).reset_index()

    return age_deltas

# Build curve for npxG/90
npxg_aging = build_aging_curve_delta(player_seasons, "npxg_90")

def fit_smooth_aging(age_deltas):
    """Fit smooth aging curve with spline."""
    # Fit weighted spline
    spline = UnivariateSpline(
        age_deltas["age"],
        age_deltas["mean_delta"],
        w=np.sqrt(age_deltas["n"]),
        s=0.5
    )

    # Predict across age range
    ages = np.arange(18, 39)
    predicted_delta = spline(ages)

    # Cumulate to get aging curve
    cumulative = np.cumsum(predicted_delta)
    peak_value = cumulative.max()
    relative = cumulative - peak_value

    return pd.DataFrame({
        "age": ages,
        "predicted_delta": predicted_delta,
        "cumulative": cumulative,
        "relative": relative
    })

npxg_curve = fit_smooth_aging(npxg_aging)
print(npxg_curve)

# Python: Position-Specific Aging Analysis
# Different positions age differently

def build_position_curves(data):
    """Build aging curves for each position and metric."""
    positions = ["Forward", "Midfielder", "Defender", "Goalkeeper"]
    metrics = ["npxg_90", "xa_90", "progressive_carries", "aerial_wins"]

    curves = {}

    for pos in positions:
        pos_data = data[data["position"] == pos]
        curves[pos] = {}

        for metric in metrics:
            if metric in pos_data.columns:
                age_deltas = build_aging_curve_delta(pos_data, metric)
                if len(age_deltas) > 5:  # Need sufficient data
                    curves[pos][metric] = fit_smooth_aging(age_deltas)

    return curves

def find_peak_ages(curves):
    """Find peak age for each position-metric combination."""
    results = []

    for pos, metrics in curves.items():
        for metric, curve in metrics.items():
            peak_idx = curve["cumulative"].idxmax()
            peak_age = curve.loc[peak_idx, "age"]

            results.append({
                "position": pos,
                "metric": metric,
                "peak_age": peak_age
            })

    return pd.DataFrame(results)

position_curves = build_position_curves(player_seasons)
peak_ages = find_peak_ages(position_curves)
print(peak_ages)

# Python: Regression to the Mean Framework
import pandas as pd
import numpy as np

def calculate_regression_factor(data, metric, min_minutes=900):
    """Calculate year-to-year correlation and regression factor."""
    # Filter for sufficient playing time
    qualified = data[data["minutes"] >= min_minutes].copy()

    # Create consecutive season pairs
    qualified = qualified.sort_values(["player_id", "season"])
    qualified["next_value"] = qualified.groupby("player_id")[metric].shift(-1)
    paired = qualified.dropna(subset=["next_value"])

    # Calculate year-to-year correlation
    correlation = paired[metric].corr(paired["next_value"])

    # Regression factor = 1 - correlation
    regression_factor = 1 - correlation

    return {
        "correlation": correlation,
        "regression_factor": regression_factor,
        "n_pairs": len(paired)
    }

# Calculate for different metrics
metrics = ["goals_90", "npxg_90", "xa_90", "pass_completion", "tackles_90"]

regression_results = []
for m in metrics:
    result = calculate_regression_factor(player_seasons, m)
    regression_results.append({
        "metric": m,
        "correlation": result["correlation"],
        "regression_factor": result["regression_factor"],
        "n": result["n_pairs"]
    })

regression_factors = pd.DataFrame(regression_results)
print(regression_factors)

def regress_to_mean(observed, league_mean, reliability, sample_size,
                    typical_sample=2000):
    """Apply regression to the mean for projection."""
    # Reliability increases with sample size
    adjusted_reliability = reliability * min(1, sample_size / typical_sample)

    # Regressed estimate
    projection = (adjusted_reliability * observed) + \
                 ((1 - adjusted_reliability) * league_mean)

    return projection

# Example: Player with 0.50 npxG/90 in 1500 minutes
regressed_npxg = regress_to_mean(
    observed=0.50,
    league_mean=0.30,
    reliability=0.70,
    sample_size=1500
)

print(f"\nObserved: 0.50 npxG/90")
print(f"Regressed projection: {regressed_npxg:.3f} npxG/90")

# Python: Marcel-Style Football Projection System
import pandas as pd
import numpy as np

# Marcel weights: 5/4/3 for most recent three seasons
MARCEL_WEIGHTS = [5, 4, 3]

def calculate_marcel_projection(player_history, metric, league_avg,
                                 aging_curve, reliability):
    """Calculate Marcel-style projection for a player."""
    # Sort by recency (most recent first)
    history = player_history.sort_values("season", ascending=False).head(3)
    n_seasons = len(history)

    if n_seasons == 0:
        return {"projection": league_avg, "uncertainty": np.nan}

    # Apply Marcel weights
    weights = np.array(MARCEL_WEIGHTS[:n_seasons])
    weights = weights / weights.sum()  # Normalize

    # Weighted average of past performance
    weighted_avg = (history[metric].values * weights).sum()

    # Total weighted minutes for regression
    total_minutes = (history["minutes"].values * weights).sum()

    # Current age and aging adjustment
    current_age = history["age"].max() + 1  # Project for next season
    age_row = aging_curve[aging_curve["age"] == current_age]
    age_adjustment = age_row["relative"].values[0] if len(age_row) > 0 else 0

    # Age-adjusted weighted average
    age_adjusted = weighted_avg + age_adjustment

    # Regression to the mean
    regression_weight = min(1, total_minutes / 3000) * reliability
    projection = (regression_weight * age_adjusted) + \
                 ((1 - regression_weight) * league_avg)

    # Uncertainty based on sample size and metric stability
    uncertainty = history[metric].std() / np.sqrt(n_seasons) * \
                  (1 + 0.5 * (1 - regression_weight))

    return {
        "projection": projection,
        "uncertainty": uncertainty,
        "weighted_avg": weighted_avg,
        "age_adjustment": age_adjustment,
        "regression_weight": regression_weight
    }

# Example projection
player_history = pd.DataFrame({
    "season": [2024, 2023, 2022],
    "age": [26, 25, 24],
    "npxg_90": [0.45, 0.38, 0.32],
    "minutes": [2800, 2400, 1800]
})

projection = calculate_marcel_projection(
    player_history=player_history,
    metric="npxg_90",
    league_avg=0.30,
    aging_curve=npxg_curve,
    reliability=0.71
)

print("Marcel Projection Results:")
print(f"Weighted average: {projection['weighted_avg']:.3f}")
print(f"Age adjustment: {projection['age_adjustment']:.3f}")
print(f"Final projection: {projection['projection']:.3f} npxG/90")
lower = projection["projection"] - 1.96 * projection["uncertainty"]
upper = projection["projection"] + 1.96 * projection["uncertainty"]
print(f"95% CI: [{lower:.3f}, {upper:.3f}]")

# Python: Project Entire Squad
def project_squad(squad_data, metrics, aging_curves,
                  league_averages, reliabilities):
    """Generate projections for entire squad."""
    projections = []

    for player_id in squad_data["player_id"].unique():
        player_history = squad_data[squad_data["player_id"] == player_id]
        player_name = player_history["player_name"].iloc[0]

        player_proj = {
            "player_id": player_id,
            "player_name": player_name
        }

        for metric in metrics:
            proj = calculate_marcel_projection(
                player_history=player_history,
                metric=metric,
                league_avg=league_averages[metric],
                aging_curve=aging_curves.get(metric, npxg_curve),
                reliability=reliabilities[metric]
            )

            player_proj[f"{metric}_proj"] = proj["projection"]
            player_proj[f"{metric}_lower"] = proj["projection"] - 1.96 * proj["uncertainty"]
            player_proj[f"{metric}_upper"] = proj["projection"] + 1.96 * proj["uncertainty"]

        projections.append(player_proj)

    return pd.DataFrame(projections)

# Define parameters
metrics = ["npxg_90", "xa_90", "progressive_passes_90"]
league_averages = {"npxg_90": 0.30, "xa_90": 0.15, "progressive_passes_90": 3.5}
reliabilities = {"npxg_90": 0.71, "xa_90": 0.68, "progressive_passes_90": 0.75}

# Generate squad projections
squad_projections = project_squad(
    squad_data=squad_historical_data,
    metrics=metrics,
    aging_curves=position_curves.get("Forward", {}),
    league_averages=league_averages,
    reliabilities=reliabilities
)

print(squad_projections.head())

# Python: ML-Based Player Projections
import pandas as pd
import numpy as np
from sklearn.model_selection import train_test_split
from xgboost import XGBRegressor
from sklearn.metrics import mean_squared_error, r2_score

def prepare_projection_features(player_data):
    """Prepare features for ML projection model."""
    data = player_data.sort_values(["player_id", "season"]).copy()

    # Create lag and rolling features
    for player_id, group in data.groupby("player_id"):
        idx = group.index

        # Lag features
        data.loc[idx, "metric_lag1"] = group["npxg_90"].shift(1)
        data.loc[idx, "metric_lag2"] = group["npxg_90"].shift(2)
        data.loc[idx, "metric_lag3"] = group["npxg_90"].shift(3)

        # Rolling averages
        data.loc[idx, "metric_ma2"] = group["npxg_90"].shift(1).rolling(2).mean()
        data.loc[idx, "metric_ma3"] = group["npxg_90"].shift(1).rolling(3).mean()

        # Trend
        data.loc[idx, "metric_trend"] = group["npxg_90"].shift(1) - group["npxg_90"].shift(2)

        # Minutes features
        data.loc[idx, "minutes_lag1"] = group["minutes"].shift(1)
        data.loc[idx, "minutes_ma2"] = group["minutes"].shift(1).rolling(2).mean()

        # Target
        data.loc[idx, "target"] = group["npxg_90"].shift(-1)

    # Age features
    data["age_squared"] = data["age"] ** 2

    # Filter complete cases
    data = data.dropna(subset=["target", "metric_lag1"])

    return data

def train_projection_model(train_data):
    """Train XGBoost projection model."""
    feature_cols = ["metric_lag1", "metric_lag2", "metric_ma2", "metric_trend",
                    "minutes_lag1", "minutes_ma2", "age", "age_squared"]

    # Remove rows with missing features
    complete_data = train_data.dropna(subset=feature_cols + ["target"])

    X = complete_data[feature_cols]
    y = complete_data["target"]

    # Train-validation split
    X_train, X_val, y_train, y_val = train_test_split(
        X, y, test_size=0.2, random_state=42
    )

    model = XGBRegressor(
        max_depth=4,
        learning_rate=0.05,
        n_estimators=500,
        min_child_weight=10,
        subsample=0.8,
        colsample_bytree=0.8,
        early_stopping_rounds=50,
        random_state=42
    )

    model.fit(
        X_train, y_train,
        eval_set=[(X_val, y_val)],
        verbose=False
    )

    # Evaluate
    val_pred = model.predict(X_val)
    print(f"Validation RMSE: {np.sqrt(mean_squared_error(y_val, val_pred)):.4f}")
    print(f"Validation R²: {r2_score(y_val, val_pred):.4f}")

    return model

# Train model
projection_features = prepare_projection_features(player_seasons)
ml_model = train_projection_model(projection_features)

# Feature importance
importance = pd.DataFrame({
    "feature": feature_cols,
    "importance": ml_model.feature_importances_
}).sort_values("importance", ascending=False)
print("\nFeature Importance:")
print(importance)

# Python: Ensemble Projection System
def create_ensemble_projection(player_history, ml_model,
                                league_avg, aging_curve, reliability):
    """Create ensemble projection from multiple methods."""
    # Method 1: Marcel projection
    marcel = calculate_marcel_projection(
        player_history, "npxg_90", league_avg, aging_curve, reliability
    )

    # Method 2: ML projection
    features = prepare_projection_features(player_history)
    latest = features.sort_values("season", ascending=False).iloc[[0]]
    feature_cols = ["metric_lag1", "metric_lag2", "metric_ma2", "metric_trend",
                    "minutes_lag1", "minutes_ma2", "age", "age_squared"]
    ml_pred = ml_model.predict(latest[feature_cols])[0]

    # Method 3: Simple weighted average
    recent = player_history.sort_values("season", ascending=False).head(2)
    simple_avg = np.average(recent["npxg_90"], weights=[0.7, 0.3])

    # Combine methods (weighted by historical accuracy)
    weights = {"marcel": 0.35, "ml": 0.40, "simple": 0.25}

    ensemble = (weights["marcel"] * marcel["projection"] +
                weights["ml"] * ml_pred +
                weights["simple"] * simple_avg)

    # Uncertainty from method disagreement
    methods = [marcel["projection"], ml_pred, simple_avg]
    uncertainty = np.std(methods) + marcel["uncertainty"]

    return {
        "ensemble": ensemble,
        "marcel": marcel["projection"],
        "ml": ml_pred,
        "simple": simple_avg,
        "uncertainty": uncertainty,
        "lower_95": ensemble - 1.96 * uncertainty,
        "upper_95": ensemble + 1.96 * uncertainty
    }

# Generate ensemble projection
ensemble_result = create_ensemble_projection(
    player_history=player_history,
    ml_model=ml_model,
    league_avg=0.30,
    aging_curve=npxg_curve,
    reliability=0.71
)

print("Ensemble Projection:")
print(f"Marcel: {ensemble_result['marcel']:.3f}")
print(f"ML: {ensemble_result['ml']:.3f}")
print(f"Simple: {ensemble_result['simple']:.3f}")
print(f"Ensemble: {ensemble_result['ensemble']:.3f} "
      f"(95% CI: [{ensemble_result['lower_95']:.3f}, {ensemble_result['upper_95']:.3f}])")

# Python: In-Season Total Projections
def calculate_season_projection(preseason_proj, current_stats,
                                 games_played, total_games=38):
    """Project final season totals from partial data."""
    games_remaining = total_games - games_played
    season_progress = games_played / total_games

    # Calculate pace
    current_pace = current_stats * (total_games / games_played)

    # Dynamic weighting based on games played
    # More games = more weight on observed performance
    reliability_by_games = min(0.95, 0.2 + 0.8 * (games_played / 30))

    # Blend preseason projection with current pace
    blended_projection = (reliability_by_games * current_pace +
                         (1 - reliability_by_games) * preseason_proj)

    # Uncertainty decreases as season progresses
    pace_variance = current_stats * 0.3  # Rough variance estimate
    uncertainty = np.sqrt(games_remaining / total_games) * pace_variance

    return {
        "projection": blended_projection,
        "current_pace": current_pace,
        "reliability": reliability_by_games,
        "uncertainty": uncertainty,
        "lower_90": blended_projection - 1.645 * uncertainty,
        "upper_90": blended_projection + 1.645 * uncertainty
    }

# Example: Player projections at different points in season
preseason_xg = 15.0  # Preseason projection: 15 xG
current_xg = 6.5     # Current: 6.5 xG in 14 games

season_proj = calculate_season_projection(
    preseason_proj=preseason_xg,
    current_stats=current_xg,
    games_played=14
)

print("In-Season Projection (14 games):")
print(f"Preseason projection: {preseason_xg:.1f} xG")
print(f"Current pace: {season_proj['current_pace']:.1f} xG")
print(f"Blended projection: {season_proj['projection']:.1f} xG "
      f"(90% CI: [{season_proj['lower_90']:.1f}, {season_proj['upper_90']:.1f}])")
print(f"Current performance weight: {season_proj['reliability'] * 100:.0f}%")

# Python: Track Projection Updates Through Season
import matplotlib.pyplot as plt

def track_season_projections(player_game_log, preseason_proj,
                              metric="xg"):
    """Track how projections evolve through the season."""
    updates = []
    total_games = 38

    cumulative = 0
    for i, (_, game) in enumerate(player_game_log.iterrows(), 1):
        cumulative += game[metric]

        proj = calculate_season_projection(
            preseason_proj=preseason_proj,
            current_stats=cumulative,
            games_played=i,
            total_games=total_games
        )

        updates.append({
            "matchweek": i,
            "cumulative": cumulative,
            "pace": proj["current_pace"],
            "projection": proj["projection"],
            "lower": proj["lower_90"],
            "upper": proj["upper_90"],
            "reliability": proj["reliability"]
        })

    return pd.DataFrame(updates)

def plot_projection_evolution(updates):
    """Visualize how projections change through season."""
    fig, ax = plt.subplots(figsize=(12, 6))

    # Confidence interval
    ax.fill_between(updates["matchweek"],
                    updates["lower"], updates["upper"],
                    alpha=0.3, color="blue", label="90% CI")

    # Projection line
    ax.plot(updates["matchweek"], updates["projection"],
            color="blue", linewidth=2, label="Projection")

    # Current pace
    ax.plot(updates["matchweek"], updates["pace"],
            color="gray", linestyle="--", label="Current Pace")

    # Actual cumulative
    ax.scatter(updates["matchweek"], updates["cumulative"],
               color="red", s=30, label="Actual")

    ax.set_xlabel("Matchweek")
    ax.set_ylabel("Expected Goals (xG)")
    ax.set_title("Season Projection Evolution")
    ax.legend()
    ax.grid(True, alpha=0.3)

    plt.tight_layout()
    plt.savefig("projection_evolution.png", dpi=150)
    plt.show()

# Generate and plot
projection_updates = track_season_projections(player_game_log, 15.0, "xg")
plot_projection_evolution(projection_updates)

# Python: Career Trajectory Projections
def project_career(player_history, aging_curve, years_ahead=5):
    """Project performance over multiple future seasons."""
    # Get current state
    latest = player_history.sort_values("season", ascending=False).iloc[0]

    current_age = latest["age"]
    current_level = latest["npxg_90"]

    # Project each future season
    projections = []

    for year in range(1, years_ahead + 1):
        future_age = current_age + year

        # Get aging adjustment
        age_row = aging_curve[aging_curve["age"] == future_age]
        if len(age_row) > 0:
            age_adj = age_row["relative"].values[0]
        else:
            age_adj = -0.05 * year  # Default decline

        # Calculate projection with increasing uncertainty
        base_proj = current_level + age_adj
        uncertainty = 0.05 * np.sqrt(year)  # Uncertainty grows with time

        projections.append({
            "season": int(latest["season"]) + year,
            "age": future_age,
            "projection": base_proj,
            "lower_80": base_proj - 1.28 * uncertainty,
            "upper_80": base_proj + 1.28 * uncertainty,
            "lower_95": base_proj - 1.96 * uncertainty,
            "upper_95": base_proj + 1.96 * uncertainty
        })

    return pd.DataFrame(projections)

# Project career for player
career_proj = project_career(
    player_history=player_history,
    aging_curve=npxg_curve,
    years_ahead=8
)

print(career_proj)

# Python: Project Contract Value
def calculate_contract_value(career_proj, wage_per_npxg=50000,
                              discount_rate=0.05):
    """Calculate present value of projected future performance."""
    proj = career_proj.copy()

    # Calculate expected value of each season
    proj["season_value"] = proj["projection"] * wage_per_npxg * 38  # 38 games

    # Discount future seasons
    proj["year_num"] = range(len(proj))
    proj["discount_factor"] = 1 / (1 + discount_rate) ** proj["year_num"]

    # Present value
    proj["present_value"] = proj["season_value"] * proj["discount_factor"]

    # Uncertainty in value
    proj["value_uncertainty"] = ((proj["upper_95"] - proj["lower_95"]) / 2) * \
                                 wage_per_npxg * 38 * proj["discount_factor"]

    # Total contract value
    total_value = proj["present_value"].sum()
    total_uncertainty = np.sqrt((proj["value_uncertainty"] ** 2).sum())

    return {
        "yearly_values": proj[["season", "age", "projection",
                               "season_value", "present_value"]],
        "total_value": total_value,
        "lower_bound": total_value - 1.96 * total_uncertainty,
        "upper_bound": total_value + 1.96 * total_uncertainty
    }

# Calculate value for 5-year contract
contract_value = calculate_contract_value(
    career_proj.head(5),
    wage_per_npxg=50000,  # £50k per npxG per season
    discount_rate=0.05
)

print("5-Year Contract Valuation:")
print(f"Total Present Value: £{contract_value['total_value']/1e6:.2f}M")
print(f"95% CI: [£{contract_value['lower_bound']/1e6:.2f}M, "
      f"£{contract_value['upper_bound']/1e6:.2f}M]")
print(contract_value["yearly_values"])

# Python: Bootstrapped Confidence Intervals
import numpy as np

def bootstrap_projection(player_history, metric, n_boot=1000,
                          league_avg=0.30, reliability=0.70):
    """Generate bootstrapped confidence intervals."""
    boot_projections = []

    for _ in range(n_boot):
        # Bootstrap resample of seasons
        boot_sample = player_history.sample(n=len(player_history), replace=True)

        # Calculate projection on resampled data
        proj = calculate_marcel_projection(
            player_history=boot_sample,
            metric=metric,
            league_avg=league_avg,
            aging_curve=npxg_curve,
            reliability=reliability
        )

        boot_projections.append(proj["projection"])

    boot_projections = np.array(boot_projections)

    # Calculate percentiles
    percentiles = np.percentile(boot_projections,
                                 [5, 10, 25, 50, 75, 90, 95])

    return {
        "mean": np.mean(boot_projections),
        "sd": np.std(boot_projections),
        "percentiles": dict(zip([5, 10, 25, 50, 75, 90, 95], percentiles)),
        "distribution": boot_projections
    }

# Run bootstrap
boot_result = bootstrap_projection(
    player_history=player_history,
    metric="npxg_90",
    n_boot=1000
)

print("Bootstrap Projection Distribution:")
print(f"Mean: {boot_result['mean']:.3f}")
print(f"Std Dev: {boot_result['sd']:.3f}")
print("\nPercentiles:")
for p, v in boot_result["percentiles"].items():
    print(f"  {p}th: {v:.3f}")

# Python: Complete Projection Dashboard
import matplotlib.pyplot as plt
from matplotlib.gridspec import GridSpec

def create_projection_dashboard(player_data, player_id):
    """Create comprehensive projection dashboard for a player."""
    player = player_data[player_data["player_id"] == player_id]
    player_name = player["player_name"].iloc[0]

    # Generate all projections
    marcel_proj = calculate_marcel_projection(player, "npxg_90", 0.30, npxg_curve, 0.71)
    career_proj = project_career(player, npxg_curve, 5)
    boot_proj = bootstrap_projection(player, "npxg_90")

    # Create figure with subplots
    fig = plt.figure(figsize=(14, 10))
    gs = GridSpec(2, 2, figure=fig)

    # 1. Historical + Projection chart
    ax1 = fig.add_subplot(gs[0, :])
    ax1.plot(player["season"], player["npxg_90"], "bo-", label="Historical")
    ax1.fill_between(career_proj["season"],
                     career_proj["lower_95"], career_proj["upper_95"],
                     alpha=0.3, color="blue", label="95% CI")
    ax1.plot(career_proj["season"], career_proj["projection"],
             "b--", label="Projection")
    ax1.set_xlabel("Season")
    ax1.set_ylabel("npxG/90")
    ax1.set_title(f"{player_name} - npxG/90 Trajectory")
    ax1.legend()
    ax1.grid(True, alpha=0.3)

    # 2. Projection distribution
    ax2 = fig.add_subplot(gs[1, 0])
    ax2.hist(boot_proj["distribution"], bins=30, color="steelblue",
             alpha=0.7, edgecolor="white")
    ax2.axvline(marcel_proj["projection"], color="red",
                linestyle="--", linewidth=2, label="Point Estimate")
    ax2.set_xlabel("npxG/90 Projection")
    ax2.set_ylabel("Frequency")
    ax2.set_title("Projection Distribution")
    ax2.legend()

    # 3. Percentile comparison
    ax3 = fig.add_subplot(gs[1, 1])
    percentiles = [5, 25, 50, 75, 95]
    values = [boot_proj["percentiles"][p] for p in percentiles]
    colors = ["#d73027", "#fc8d59", "#fee08b", "#91cf60", "#1a9850"]
    ax3.barh(range(len(percentiles)), values, color=colors)
    ax3.set_yticks(range(len(percentiles)))
    ax3.set_yticklabels([f"{p}th" for p in percentiles])
    ax3.set_xlabel("npxG/90")
    ax3.set_title("Projection Percentiles")
    ax3.axvline(marcel_proj["projection"], color="black",
                linestyle="--", label="Point Estimate")

    plt.tight_layout()
    plt.savefig(f"projection_dashboard_{player_id}.png", dpi=150)
    plt.show()

    # Summary table
    comparison = pd.DataFrame({
        "Metric": ["Historical Avg", "Marcel Projection", "10th Percentile",
                   "Median", "90th Percentile"],
        "Value": [
            player["npxg_90"].mean(),
            marcel_proj["projection"],
            boot_proj["percentiles"][10],
            boot_proj["percentiles"][50],
            boot_proj["percentiles"][90]
        ]
    })

    return {
        "comparison": comparison,
        "career_proj": career_proj,
        "boot_dist": boot_proj["distribution"]
    }

# Generate dashboard
dashboard = create_projection_dashboard(squad_data, "player_123")
print(dashboard["comparison"])