diff --git a/parktrack_ml/config.py b/parktrack_ml/config.py index 7b46bca..629df2a 100644 --- a/parktrack_ml/config.py +++ b/parktrack_ml/config.py @@ -82,9 +82,10 @@ "verbose": -1, } -TRAIN_DAYS_BACK = int(os.environ.get("TRAIN_DAYS_BACK", "200")) -ML_MODEL_TYPE = "lightgbm" -ML_MODEL_VERSION = "3.0" +TRAIN_DAYS_BACK = int(os.environ.get("TRAIN_DAYS_BACK", "200")) +FORECAST_HOURS = int(os.environ.get("FORECAST_HOURS", "48")) +ML_MODEL_TYPE = "lightgbm" +ML_MODEL_VERSION = "3.0" # --------------------------------------------------------------------------- # Scheduling diff --git a/parktrack_ml/evaluate.py b/parktrack_ml/evaluate.py new file mode 100644 index 0000000..cb693f8 --- /dev/null +++ b/parktrack_ml/evaluate.py @@ -0,0 +1,233 @@ +""" +Backtest validator — walk-forward evaluation of the saved model. + +For each historical hourly bucket, the model predicts using only data +strictly before that timestamp (no lookahead), matching real inference. +Compares predicted class and predicted_occupied against actual values. + +Usage: + python -m parktrack_ml.evaluate + python -m parktrack_ml.evaluate --days 60 + python -m parktrack_ml.evaluate --days 30 --zone-id 3 +""" +from __future__ import annotations + +import argparse +import json +import logging +import sys +from datetime import datetime, timedelta, timezone +from typing import Optional + +import numpy as np +import pandas as pd + +from .config import ( + MODEL_FILE, ZONE_META_FILE, + FEATURE_NAMES, CLASS_CENTER_RATES, + THRESHOLD_LOW, THRESHOLD_MEDIUM, TRAIN_DAYS_BACK, +) +from .data_loader import load_observations, load_zone_meta, aggregate_hourly +from .features import build_prediction_vector, label_occupancy +from .model import LGBMWrapper + +logging.basicConfig( + level=logging.INFO, + format="%(asctime)s %(levelname)-8s %(message)s", + handlers=[logging.StreamHandler(sys.stdout)], +) +log = logging.getLogger("evaluate") + +CLASS_NAMES = {0: "Low", 1: "Medium", 2: "High"} + + +# --------------------------------------------------------------------------- +# Core evaluation +# --------------------------------------------------------------------------- + +def _predict_row( + zone_id: int, + hour: pd.Timestamp, + history: pd.DataFrame, + zone_meta: dict, + model: LGBMWrapper, + capacity: int, +) -> dict: + """Predict for a single (zone, hour) using history strictly before that hour.""" + past = history[history["hour"] < hour] + vec = build_prediction_vector(zone_id, hour.to_pydatetime(), past, zone_meta) + + X = vec.reshape(1, -1) + class_code = int(model.predict(X)[0]) + probs = model.predict_proba(X)[0] + expected_rate = sum(float(probs[c]) * CLASS_CENTER_RATES[c] for c in range(3)) + predicted_occ = max(0, min(capacity, round(expected_rate * capacity))) + + return { + "class_code": class_code, + "confidence": float(probs[class_code]), + "predicted_occupied": predicted_occ, + } + + +def evaluate( + days: int = 30, + zone_ids: Optional[list[int]] = None, +) -> pd.DataFrame: + """ + Run walk-forward backtest. + Returns a DataFrame with one row per (zone, hour) containing + actual and predicted values. + """ + log.info("Loading model from %s", MODEL_FILE) + model = LGBMWrapper.load(MODEL_FILE) + + with open(ZONE_META_FILE) as f: + zone_meta_all: dict[int, dict] = {int(k): v for k, v in json.load(f).items()} + + to_dt = datetime.now(tz=timezone.utc) + from_dt = to_dt - timedelta(days=days) + + if zone_ids is None: + zone_ids = sorted(zone_meta_all.keys()) + if not zone_ids: + zone_meta_df = load_zone_meta() + zone_ids = zone_meta_df["zone_id"].tolist() + + log.info("Fetching %d days of occupancy for zones %s ...", days, zone_ids) + raw = load_observations(zone_ids=zone_ids, from_dt=from_dt, to_dt=to_dt) + hourly = aggregate_hourly(raw) + + if hourly.empty: + log.error("No occupancy data for the requested period. Try --days with a larger value.") + return pd.DataFrame() + + log.info("Evaluating %d hourly records across %d zones...", + len(hourly), hourly["zone_id"].nunique()) + + results = [] + for zone_id, grp in hourly.groupby("zone_id"): + zone_id = int(zone_id) + meta = zone_meta_all.get(zone_id, {"capacity": 10, "zone_type_standard": 1}) + capacity = int(meta.get("capacity", 10)) + zone_history = grp.sort_values("hour").reset_index(drop=True) + + for _, row in zone_history.iterrows(): + hour = row["hour"] + actual_rate = float(row["occupancy_rate"]) + actual_occ = int(round(actual_rate * capacity)) + actual_label = label_occupancy(actual_rate) + + pred = _predict_row(zone_id, hour, zone_history, meta, model, capacity) + + results.append({ + "zone_id": zone_id, + "hour": hour, + "actual_class": actual_label, + "predicted_class": pred["class_code"], + "correct": int(pred["class_code"] == actual_label), + "actual_occupied": actual_occ, + "predicted_occupied": pred["predicted_occupied"], + "abs_error": abs(pred["predicted_occupied"] - actual_occ), + "confidence": pred["confidence"], + "actual_rate": actual_rate, + }) + + return pd.DataFrame(results) + + +# --------------------------------------------------------------------------- +# Report +# --------------------------------------------------------------------------- + +def _bar(ratio: float, width: int = 20) -> str: + filled = round(ratio * width) + return "█" * filled + "░" * (width - filled) + + +def print_report(df: pd.DataFrame) -> None: + if df.empty: + print("No data to report.") + return + + n = len(df) + acc = df["correct"].mean() + mae = df["abs_error"].mean() + + print() + print("=" * 60) + print(" BACKTEST REPORT") + print("=" * 60) + print(f" Samples : {n:,}") + print(f" Accuracy : {acc:.1%} {_bar(acc)}") + print(f" MAE (spots): {mae:.2f}") + print() + + # Per-class metrics + print(" ── Per-class ─────────────────────────────────────") + print(f" {'Class':<8} {'Prec':>6} {'Recall':>6} {'F1':>6} {'Support':>7}") + for cls, name in CLASS_NAMES.items(): + y_true = (df["actual_class"] == cls) + y_pred = (df["predicted_class"] == cls) + tp = int((y_pred & y_true).sum()) + fp = int((y_pred & ~y_true).sum()) + fn = int((~y_pred & y_true).sum()) + p = tp / (tp + fp) if (tp + fp) else 0.0 + r = tp / (tp + fn) if (tp + fn) else 0.0 + f1 = 2 * p * r / (p + r) if (p + r) else 0.0 + print(f" {name:<8} {p:>6.1%} {r:>6.1%} {f1:>6.1%} {y_true.sum():>7,}") + print() + + # Per-zone + print(" ── Per-zone ──────────────────────────────────────") + print(f" {'Zone':>5} {'Accuracy':>9} {'MAE':>6} {'Samples':>7}") + for zid, g in df.groupby("zone_id"): + z_acc = g["correct"].mean() + z_mae = g["abs_error"].mean() + print(f" {zid:>5} {z_acc:>9.1%} {z_mae:>6.2f} {len(g):>7,}") + print() + + # Per-hour-of-day accuracy (key check for static-forecast bug) + print(" ── Accuracy by hour of day ───────────────────────") + df["hod"] = pd.to_datetime(df["hour"], utc=True).dt.hour + hod = df.groupby("hod")["correct"].mean() + for h, a in hod.items(): + bar = _bar(a, 15) + print(f" {h:02d}:00 {a:5.1%} {bar}") + print() + + # Confusion matrix + print(" ── Confusion matrix (actual → predicted) ─────────") + labels = [0, 1, 2] + header = " " + "".join(f" {CLASS_NAMES[p]:>8}" for p in labels) + print(header) + for a in labels: + row_mask = df["actual_class"] == a + row = " ".join( + f"{int((df.loc[row_mask, 'predicted_class'] == p).sum()):>8}" + for p in labels + ) + print(f" {CLASS_NAMES[a]:<7} {row}") + print() + print("=" * 60) + + +# --------------------------------------------------------------------------- +# Entry point +# --------------------------------------------------------------------------- + +def main() -> None: + parser = argparse.ArgumentParser(description="Backtest ML occupancy model") + parser.add_argument("--days", type=int, default=30, + help="Days of history to evaluate (default: 30)") + parser.add_argument("--zone-id", type=int, default=None, + help="Evaluate a single zone only") + args = parser.parse_args() + + zone_ids = [args.zone_id] if args.zone_id else None + df = evaluate(days=args.days, zone_ids=zone_ids) + print_report(df) + + +if __name__ == "__main__": + main() diff --git a/parktrack_ml/features.py b/parktrack_ml/features.py index 3744fee..9029e18 100644 --- a/parktrack_ml/features.py +++ b/parktrack_ml/features.py @@ -178,6 +178,7 @@ def _hist_avg(hour_of_day: int) -> float: feats[f'occupancy_ma_{w}h'] = float(np.mean([_hist_avg(dt.hour - i) for i in range(1, w + 1)])) else: history = recent_hourly.set_index('hour')['occupancy_rate'].sort_index() + history = history[~history.index.duplicated(keep='last')] pred_hour = dt.floor('h') history = history[history.index < pred_hour] @@ -185,14 +186,25 @@ def get_lag(lag_h: int) -> float: t = pred_hour - pd.Timedelta(hours=lag_h) if t in history.index: return float(history[t]) - # Sparse data: use historical average for that hour rather than - # propagating the last known value across all lag positions. + # No direct match (future slot) — try same hour yesterday. + # Yesterday's 22:00 is far more predictive than the multi-day + # average, because it captures the real daily on/off pattern. + t_yesterday = t - pd.Timedelta(hours=24) + if t_yesterday in history.index: + return float(history[t_yesterday]) return _hist_avg(t.hour) def get_ma(w: int) -> float: window = history[history.index >= pred_hour - pd.Timedelta(hours=w)] if not window.empty: return float(window.mean()) + # Try same window yesterday + window_yesterday = history[ + (history.index >= pred_hour - pd.Timedelta(hours=w + 24)) & + (history.index < pred_hour - pd.Timedelta(hours=24)) + ] + if not window_yesterday.empty: + return float(window_yesterday.mean()) hours = [(dt.hour - i) % 24 for i in range(1, w + 1)] return float(np.mean([_hist_avg(h) for h in hours])) diff --git a/parktrack_ml/forecaster.py b/parktrack_ml/forecaster.py index 6564797..7872f65 100644 --- a/parktrack_ml/forecaster.py +++ b/parktrack_ml/forecaster.py @@ -9,7 +9,7 @@ from concurrent.futures import ThreadPoolExecutor, as_completed from datetime import datetime, timedelta, timezone -from .config import API_URL, API_TOKEN, ML_MODEL_TYPE, ML_MODEL_VERSION +from .config import API_URL, API_TOKEN, ML_MODEL_TYPE, ML_MODEL_VERSION, FORECAST_HOURS from .api_client import ParkTrackClient from .interfaces import predict from .data_loader import load_recent_observations, aggregate_hourly @@ -76,9 +76,9 @@ def run() -> None: base = now.replace(second=0, microsecond=0) slots = [ base.replace(minute=m) + timedelta(hours=h) - for h in range(25) for m in (0, 30) + for h in range(FORECAST_HOURS + 1) for m in (0, 30) if base.replace(minute=m) + timedelta(hours=h) > now - ][:48] + ][:FORECAST_HOURS * 2] # Pre-fetch recent observations ONCE per zone, not once per (zone, slot). # Without this we'd make zones×slots = potentially 720+ API calls per run. diff --git a/parktrack_ml/visualize.py b/parktrack_ml/visualize.py new file mode 100644 index 0000000..834feee --- /dev/null +++ b/parktrack_ml/visualize.py @@ -0,0 +1,279 @@ +""" +Visualize model performance — generates PNG charts from backtest results. + +Usage: + python -m parktrack_ml.visualize + python -m parktrack_ml.visualize --days 60 --out ./charts + python -m parktrack_ml.visualize --days 30 --zone-id 3 +""" +from __future__ import annotations + +import argparse +import json +import logging +import os +import sys +from datetime import datetime + +import numpy as np +import pandas as pd + +logging.basicConfig(level=logging.INFO, format="%(asctime)s %(message)s") +log = logging.getLogger("visualize") + +CLASS_NAMES = ["Low", "Medium", "High"] +COLORS = {"Low": "#4CAF50", "Medium": "#FF9800", "High": "#F44336"} +PALETTE = [COLORS[c] for c in CLASS_NAMES] + + +def _ensure_mpl(): + try: + import matplotlib + matplotlib.use("Agg") + import matplotlib.pyplot as plt + return plt + except ImportError: + log.error("matplotlib not installed. Run: pip install matplotlib") + sys.exit(1) + + +# --------------------------------------------------------------------------- +# Individual charts +# --------------------------------------------------------------------------- + +def plot_confusion_matrix(df: pd.DataFrame, out_dir: str, plt) -> str: + import matplotlib.ticker as ticker + + fig, ax = plt.subplots(figsize=(6, 5)) + cm = np.zeros((3, 3), dtype=int) + for a in range(3): + for p in range(3): + cm[a, p] = int(((df["actual_class"] == a) & (df["predicted_class"] == p)).sum()) + + im = ax.imshow(cm, cmap="Blues") + fig.colorbar(im, ax=ax) + ax.set_xticks(range(3)); ax.set_xticklabels(CLASS_NAMES) + ax.set_yticks(range(3)); ax.set_yticklabels(CLASS_NAMES) + ax.set_xlabel("Predicted"); ax.set_ylabel("Actual") + ax.set_title("Confusion Matrix") + + total = cm.sum() + for i in range(3): + for j in range(3): + v = cm[i, j] + color = "white" if v > total * 0.15 else "black" + ax.text(j, i, f"{v}", ha="center", va="center", color=color, fontsize=12) + + fig.tight_layout() + path = os.path.join(out_dir, "confusion_matrix.png") + fig.savefig(path, dpi=150) + plt.close(fig) + return path + + +def plot_accuracy_by_hour(df: pd.DataFrame, out_dir: str, plt) -> str: + df = df.copy() + df["hod"] = pd.to_datetime(df["hour"], utc=True).dt.hour + hod = df.groupby("hod")["correct"].mean().reindex(range(24), fill_value=np.nan) + + fig, ax = plt.subplots(figsize=(10, 4)) + bars = ax.bar(hod.index, hod.values * 100, color="#2196F3", alpha=0.8, width=0.7) + ax.axhline(df["correct"].mean() * 100, color="red", linestyle="--", + linewidth=1.5, label=f"Overall {df['correct'].mean():.1%}") + ax.set_xlabel("Hour of day") + ax.set_ylabel("Accuracy (%)") + ax.set_title("Accuracy by Hour of Day") + ax.set_xticks(range(24)) + ax.set_ylim(0, 105) + ax.legend() + ax.grid(axis="y", alpha=0.3) + + fig.tight_layout() + path = os.path.join(out_dir, "accuracy_by_hour.png") + fig.savefig(path, dpi=150) + plt.close(fig) + return path + + +def plot_per_zone_metrics(df: pd.DataFrame, out_dir: str, plt) -> str: + zones = sorted(df["zone_id"].unique()) + accs = [df[df["zone_id"] == z]["correct"].mean() * 100 for z in zones] + maes = [df[df["zone_id"] == z]["abs_error"].mean() for z in zones] + + fig, (ax1, ax2) = plt.subplots(1, 2, figsize=(12, 4)) + + ax1.barh([f"Zone {z}" for z in zones], accs, color="#4CAF50", alpha=0.8) + ax1.axvline(df["correct"].mean() * 100, color="red", linestyle="--", linewidth=1.5, + label=f"Overall {df['correct'].mean():.1%}") + ax1.set_xlabel("Accuracy (%)") + ax1.set_title("Accuracy per Zone") + ax1.set_xlim(0, 105) + ax1.legend() + ax1.grid(axis="x", alpha=0.3) + + ax2.barh([f"Zone {z}" for z in zones], maes, color="#FF9800", alpha=0.8) + ax2.axvline(df["abs_error"].mean(), color="red", linestyle="--", linewidth=1.5, + label=f"Overall MAE {df['abs_error'].mean():.2f}") + ax2.set_xlabel("MAE (spots)") + ax2.set_title("MAE per Zone") + ax2.legend() + ax2.grid(axis="x", alpha=0.3) + + fig.tight_layout() + path = os.path.join(out_dir, "per_zone_metrics.png") + fig.savefig(path, dpi=150) + plt.close(fig) + return path + + +def plot_class_distribution(df: pd.DataFrame, out_dir: str, plt) -> str: + fig, (ax1, ax2) = plt.subplots(1, 2, figsize=(10, 4)) + + actual_counts = [int((df["actual_class"] == i).sum()) for i in range(3)] + pred_counts = [int((df["predicted_class"] == i).sum()) for i in range(3)] + + x = np.arange(3) + w = 0.35 + ax1.bar(x - w/2, actual_counts, w, label="Actual", color=PALETTE, alpha=0.9) + ax1.bar(x + w/2, pred_counts, w, label="Predicted", color=PALETTE, alpha=0.5, + edgecolor="black", linewidth=0.8) + ax1.set_xticks(x); ax1.set_xticklabels(CLASS_NAMES) + ax1.set_ylabel("Count") + ax1.set_title("Class Distribution: Actual vs Predicted") + ax1.legend() + ax1.grid(axis="y", alpha=0.3) + + # Per-class F1 + f1s = [] + for cls in range(3): + y_true = (df["actual_class"] == cls) + y_pred = (df["predicted_class"] == cls) + tp = int((y_pred & y_true).sum()) + fp = int((y_pred & ~y_true).sum()) + fn = int((~y_pred & y_true).sum()) + p = tp / (tp + fp) if (tp + fp) else 0.0 + r = tp / (tp + fn) if (tp + fn) else 0.0 + f1s.append(2 * p * r / (p + r) if (p + r) else 0.0) + + bars = ax2.bar(CLASS_NAMES, [v * 100 for v in f1s], color=PALETTE, alpha=0.85) + ax2.set_ylabel("F1 Score (%)") + ax2.set_title("F1 Score per Class") + ax2.set_ylim(0, 105) + ax2.grid(axis="y", alpha=0.3) + for bar, val in zip(bars, f1s): + ax2.text(bar.get_x() + bar.get_width()/2, bar.get_height() + 1, + f"{val:.1%}", ha="center", va="bottom", fontsize=11) + + fig.tight_layout() + path = os.path.join(out_dir, "class_metrics.png") + fig.savefig(path, dpi=150) + plt.close(fig) + return path + + +def plot_feature_importance(out_dir: str, plt) -> str | None: + from .config import MODEL_FILE + from .model import LGBMWrapper + + try: + model = LGBMWrapper.load(MODEL_FILE) + importance = model.feature_importance() # dict sorted by gain desc + except Exception as exc: + log.warning("Could not load model for feature importance: %s", exc) + return None + + names = list(importance.keys())[:15] + vals = list(importance.values())[:15] + + fig, ax = plt.subplots(figsize=(8, 6)) + colors = ["#2196F3"] * len(names) + ax.barh(names[::-1], vals[::-1], color=colors[::-1], alpha=0.85) + ax.set_xlabel("Gain (feature importance)") + ax.set_title("Top-15 Features by Importance (Gain)") + ax.grid(axis="x", alpha=0.3) + fig.tight_layout() + + path = os.path.join(out_dir, "feature_importance.png") + fig.savefig(path, dpi=150) + plt.close(fig) + return path + + +def plot_predicted_vs_actual(df: pd.DataFrame, out_dir: str, plt) -> str: + # Sample up to 2000 points for readability + sample = df.sample(min(2000, len(df)), random_state=42) + + fig, ax = plt.subplots(figsize=(6, 6)) + c_colors = [PALETTE[int(c)] for c in sample["actual_class"]] + ax.scatter(sample["actual_occupied"], sample["predicted_occupied"], + c=c_colors, alpha=0.35, s=15) + lim = max(sample["actual_occupied"].max(), sample["predicted_occupied"].max()) + 1 + ax.plot([0, lim], [0, lim], "k--", linewidth=1, label="Perfect prediction") + ax.set_xlabel("Actual occupied spots") + ax.set_ylabel("Predicted occupied spots") + ax.set_title("Predicted vs Actual (spots)") + ax.set_xlim(0, lim); ax.set_ylim(0, lim) + ax.legend() + + # Legend patches + from matplotlib.patches import Patch + legend_els = [Patch(facecolor=PALETTE[i], label=CLASS_NAMES[i]) for i in range(3)] + ax.legend(handles=legend_els + [ax.get_lines()[0]], loc="upper left") + + fig.tight_layout() + path = os.path.join(out_dir, "predicted_vs_actual.png") + fig.savefig(path, dpi=150) + plt.close(fig) + return path + + +# --------------------------------------------------------------------------- +# Entry point +# --------------------------------------------------------------------------- + +def main() -> None: + parser = argparse.ArgumentParser(description="Generate ML performance charts") + parser.add_argument("--days", type=int, default=30) + parser.add_argument("--zone-id", type=int, default=None) + parser.add_argument("--out", type=str, default="./charts", + help="Output directory for PNG files") + args = parser.parse_args() + + plt = _ensure_mpl() + os.makedirs(args.out, exist_ok=True) + + log.info("Running backtest (%d days)...", args.days) + from .evaluate import evaluate + zone_ids = [args.zone_id] if args.zone_id else None + df = evaluate(days=args.days, zone_ids=zone_ids) + + if df.empty: + log.error("No data — cannot generate charts.") + sys.exit(1) + + generated = [] + + log.info("Generating charts → %s", args.out) + generated.append(plot_confusion_matrix(df, args.out, plt)) + generated.append(plot_accuracy_by_hour(df, args.out, plt)) + generated.append(plot_per_zone_metrics(df, args.out, plt)) + generated.append(plot_class_distribution(df, args.out, plt)) + generated.append(plot_predicted_vs_actual(df, args.out, plt)) + + fi_path = plot_feature_importance(args.out, plt) + if fi_path: + generated.append(fi_path) + + print() + print("Charts saved:") + for p in generated: + if p: + print(f" {p}") + print() + print(f"Overall accuracy : {df['correct'].mean():.1%}") + print(f"MAE (spots) : {df['abs_error'].mean():.2f}") + print(f"Samples : {len(df):,}") + + +if __name__ == "__main__": + main()