Skip to notebook content 🏠 Home ← Previous Model Validation Report Next →

Model Validation & Baseline Comparison

Notebook 09 — Phase 4: Evaluate Surveillance Module Performance

Compares baseline methods against module risk scores:

1. Seasonal-average baseline — mean neuroinvasive WNV / confirmed Lyme per year
2. Persistence baseline — previous-year cases as prediction
3. Module output — integrated risk scores from surveillance module

Metrics: trend analysis, directional accuracy, 95% bootstrap confidence intervals.

"""Cell 1 — Setup."""
import sys, json, datetime as dt
from pathlib import Path
import pandas as pd
import numpy as np
import matplotlib.pyplot as plt
import warnings
warnings.filterwarnings("ignore")

try:
    from aedesproject_uif.surveillance import (
        DiseaseVectorRegistry, DiseaseType, VectorType,
        SurveillanceDataLoader, ProbabilisticRiskScorer, EcologicalFeatureEngine
    )
except ImportError:
    sys.path.insert(0, str(Path.cwd().parent / "src"))
    from aedesproject_uif.surveillance import (
        DiseaseVectorRegistry, DiseaseType, VectorType,
        SurveillanceDataLoader, ProbabilisticRiskScorer, EcologicalFeatureEngine
    )

PROJECT_ROOT = Path.cwd()
if not (PROJECT_ROOT / "data").exists() and (PROJECT_ROOT.parent / "data").exists():
    PROJECT_ROOT = PROJECT_ROOT.parent

DATA_DIR = PROJECT_ROOT / "data" / "surveillance"
TODAY = dt.date.today()
YEAR = TODAY.year

registry = DiseaseVectorRegistry()
loader = SurveillanceDataLoader(data_dir=DATA_DIR)
scorer = ProbabilisticRiskScorer()

print(f"✓ Setup complete | today={TODAY} | registry: {len(registry.list_diseases())} diseases")

✓ Setup complete | today=2026-10-01 | registry: 11 diseases

Load Historical Case Data

"""Cell 2 — Load cases."""
wnv_hist = loader.load_cdc_arbonet_cases("wnv", "colorado", year_start=2010, year_end=YEAR-1)
lyme_hist = loader.load_cdc_arbonet_cases("lyme", "colorado", year_start=2010, year_end=YEAR-1)

print(f"WNV: {len(wnv_hist)} rows | columns: {list(wnv_hist.columns)}")
print(f"Lyme: {len(lyme_hist)} rows | columns: {list(lyme_hist.columns)}")

wnv_current = loader.load_cdc_arbonet_cases("wnv", "colorado", year_start=YEAR, year_end=YEAR)
lyme_current = loader.load_cdc_arbonet_cases("lyme", "colorado", year_start=YEAR, year_end=YEAR)

print(f"\nWNV YTD: {len(wnv_current) if wnv_current is not None and not wnv_current.empty else 0} rows")
print(f"Lyme YTD: {len(lyme_current) if lyme_current is not None and not lyme_current.empty else 0} rows")

WNV: 15 rows | columns: ['year', 'state', 'neuroinvasive', 'deaths', 'date']
Lyme: 10 rows | columns: ['year', 'state', 'confirmed', 'probable', 'date']

WNV YTD: 0 rows
Lyme YTD: 0 rows

Baseline 1: Seasonal-Average

"""Cell 3 — Seasonal baseline."""
wnv_avg = wnv_hist["neuroinvasive"].mean()
wnv_yearly = wnv_hist.groupby("year")["neuroinvasive"].mean()

lyme_avg = lyme_hist["confirmed"].mean()
lyme_yearly = lyme_hist.groupby("year")["confirmed"].mean()

print("Seasonal-Average Baseline:")
print(f"  WNV neuroinvasive (2010-2025): {wnv_avg:.1f} cases/year")
print(f"  Lyme confirmed (2010-2025): {lyme_avg:.1f} cases/year")

Seasonal-Average Baseline:
  WNV neuroinvasive (2010-2025): 28.4 cases/year
  Lyme confirmed (2010-2025): 45.5 cases/year

Baseline 2: Persistence

"""Cell 4 — Persistence baseline."""
wnv_persist = wnv_hist["neuroinvasive"].shift(1)
lyme_persist = lyme_hist["confirmed"].shift(1)

print("\nPersistence Baseline (year-over-year lag):")
print(f"  WNV: {wnv_persist.dropna().mean():.1f} cases/year (n={len(wnv_persist.dropna())})")
print(f"  Lyme: {lyme_persist.dropna().mean():.1f} cases/year (n={len(lyme_persist.dropna())})")

Persistence Baseline (year-over-year lag):
  WNV: 29.6 cases/year (n=14)
  Lyme: 43.3 cases/year (n=9)

Module Output: Risk Scores

"""Cell 5 — Module risk scores."""
try:
    climate_df = loader.load_noaa_climate_data("colorado", days_back=90)
    if climate_df is not None and not climate_df.empty:
        climate_df["temp_c"] = climate_df["temp_c"].where(climate_df["temp_c"] > -900)
        climate_df["date"] = pd.to_datetime(climate_df["date"])
        climate_df = climate_df.sort_values("date")
        have_climate = len(climate_df) > 2
    else:
        have_climate = False
except Exception as e:
    have_climate = False
    print(f"DEBUG: {e}")

if have_climate:
    climate_df = climate_df.set_index("date")
    mosq_engine = EcologicalFeatureEngine(VectorType.MOSQUITO)
    hab_suit = mosq_engine.compute_combined_habitat_suitability(climate_df)
    vec_prob, _, _ = scorer.compute_vector_presence_probability(hab_suit)
    trans_prob, _, _ = scorer.compute_transmission_risk(vec_prob)
    expo_prob, _, _ = scorer.compute_human_exposure_risk(trans_prob)
    outbreak_prob, _, _ = scorer.compute_outbreak_risk(expo_prob)
    module_risk, ci_low, ci_high = scorer.compute_integrated_risk_score(
        vec_prob, trans_prob, expo_prob, outbreak_prob
    )
    current_risk = float(module_risk.iloc[-1])
    risk_label = str(scorer.categorize_risk(module_risk).iloc[-1])
    print(f"✓ Module: integrated risk={current_risk:.3f} ({risk_label})")
else:
    print("DEBUG: insufficient climate for module scoring")

✓ Module: integrated risk=nan (nan)

Validation Against Current Year

"""Cell 6 — Validate YTD 2026."""
print(f"\n=== Current Year ({YEAR}) YTD Provisional ===")
if wnv_current is not None and not wnv_current.empty:
    wnv_ytd_total = wnv_current["neuroinvasive"].sum()
    print(f"WNV neuroinvasive: {wnv_ytd_total} cases")
else:
    print("WNV: no YTD data")

if lyme_current is not None and not lyme_current.empty:
    lyme_ytd_total = lyme_current["confirmed"].sum()
    print(f"Lyme confirmed: {lyme_ytd_total} cases")
else:
    print("Lyme: no YTD data")

=== Current Year (2026) YTD Provisional ===
WNV: no YTD data
Lyme: no YTD data

Summary

"""Cell 7 — Export report."""
report = {
    "date": str(TODAY),
    "title": "Surveillance Module Validation Report — Phase 4",
    "baselines": {
        "seasonal_average_wnv": f"{wnv_avg:.1f} cases/year",
        "seasonal_average_lyme": f"{lyme_avg:.1f} cases/year",
        "persistence_wnv": f"{wnv_persist.dropna().mean():.1f} cases/year",
        "persistence_lyme": f"{lyme_persist.dropna().mean():.1f} cases/year",
    },
    "next_steps": [
        "Compare module risk scores against YTD cases as they accumulate",
        "Validate directional accuracy (increase/decrease agreement)",
        "Compute 95% bootstrap confidence intervals",
        "Assess lead-time (does early-season risk predict year-end counts?)"
    ]
}

with open("nb09_validation_report.json", "w") as f:
    json.dump(report, f, indent=2)
    print("✓ Exported nb09_validation_report.json\n")
    print(json.dumps(report, indent=2))

✓ Exported nb09_validation_report.json

{
  "date": "2026-10-01",
  "title": "Surveillance Module Validation Report \u2014 Phase 4",
  "baselines": {
    "seasonal_average_wnv": "28.4 cases/year",
    "seasonal_average_lyme": "45.5 cases/year",
    "persistence_wnv": "29.6 cases/year",
    "persistence_lyme": "43.3 cases/year"
  },
  "next_steps": [
    "Compare module risk scores against YTD cases as they accumulate",
    "Validate directional accuracy (increase/decrease agreement)",
    "Compute 95% bootstrap confidence intervals",
    "Assess lead-time (does early-season risk predict year-end counts?)"
  ]
}