From abf37e6376c75e2bdf321319395dc946d0ffa9f5 Mon Sep 17 00:00:00 2001 From: Emi Griffith Date: Thu, 23 Jul 2026 06:33:54 -0700 Subject: [PATCH] Add value-drift check comparing NASA vs Open-Meteo (ERA5) MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit backend/drift_check.py quantifies the migration's data impact by fetching both live history sources (NASA POWER primary + the Open-Meteo/ERA5 fallback) for a sample of cells and reporting, per variable: mean-absolute-difference / bias / max (raw drift) and grade-band divergence (the share of days whose graded band actually changes). Since Open-Meteo is ERA5, this doubles as a fidelity check for the ERA5 seed. Because grades are percentiles within each source's own distribution, a uniform bias moves MAD but not grades — so the two numbers are read together. Open-Meteo is kept as a dormant fallback (not deleted), so both sources stay fetchable for this comparison. The comparison logic is unit-tested; the per-cell fetch runs on a networked box (python drift_check.py [--limit N] [--days N]). --- backend/drift_check.py | 136 ++++++++++++++++++++++++++++++ backend/tests/test_drift_check.py | 55 ++++++++++++ 2 files changed, 191 insertions(+) create mode 100644 backend/drift_check.py create mode 100644 backend/tests/test_drift_check.py diff --git a/backend/drift_check.py b/backend/drift_check.py new file mode 100644 index 0000000..b91df04 --- /dev/null +++ b/backend/drift_check.py @@ -0,0 +1,136 @@ +"""Value-drift check: quantify how much the migration moved the data, per variable +and by grade band, by comparing the two live history sources head-to-head. + +Open-Meteo is kept as a dormant fallback, so both sources are still fetchable: this +pulls NASA POWER (the new primary) and the Open-Meteo archive (ERA5) for a sample of +cells and reports, per variable: + - mean-absolute-difference / bias / max — the RAW value drift, and + - grade-band divergence — the share of (day, metric) whose graded band actually + changes between sources. +Since Open-Meteo IS ERA5, this doubles as a fidelity check for the ERA5 seed. + + python drift_check.py [--limit N] [--days N] # on a networked box + +Note a real property of the grading model worth reading the two numbers together: +grades are percentiles within each source's OWN distribution, so a uniform bias +between sources shifts MAD but leaves grades unchanged. High MAD with low grade +divergence means "different absolute values, same story" — usually fine. High grade +divergence is the signal that actually matters. + +The comparison logic (compare_values / grade_band_divergence) is unit-tested; the +fetch is not. +""" +import datetime +import sys + +import polars as pl + +from data import cities +from data import climate +from data import grading +from data import grid + +# Variables compared for raw drift (feels is derived; precip/wind/gust/humid optional). +_METRICS = ("tmax", "tmin", "precip", "wind", "gust", "humid", "feels") + + +def compare_values(a: pl.DataFrame, b: pl.DataFrame, metrics=_METRICS) -> dict: + """Per-variable raw-drift stats over the dates the two frames share: n aligned + days, mean-absolute-difference, mean bias (a - b), and max absolute difference. + Days where either source is null for a metric are dropped from that metric.""" + j = a.join(b, on="date", how="inner", suffix="_b") + out = {} + for m in metrics: + if m not in a.columns or m not in b.columns: + continue + d = j.select((pl.col(m) - pl.col(f"{m}_b")).alias("d")).drop_nulls() + if d.height == 0: + out[m] = {"n": 0, "mad": None, "bias": None, "max_abs": None} + continue + diff = d["d"] + out[m] = { + "n": d.height, + "mad": round(float(diff.abs().mean()), 3), + "bias": round(float(diff.mean()), 3), + "max_abs": round(float(diff.abs().max()), 3), + } + return out + + +def grade_band_divergence(a: pl.DataFrame, b: pl.DataFrame, start, end) -> dict: + """Per-metric share of days whose graded band differs between the two sources over + [start, end]. Each source's days are graded against its OWN ±7-day climatology + (grading.grade_range), so this measures whether the migration changes the story a + day tells, not just its raw value.""" + ga = {d["date"]: d for d in grading.grade_range(a, start, end)} + gb = {d["date"]: d for d in grading.grade_range(b, start, end)} + per: dict[str, list[int]] = {} + for date in ga.keys() & gb.keys(): + da, db = ga[date], gb[date] + for m in grading.CLIMO_METRICS: + va, vb = da.get(m), db.get(m) + if va and vb: + slot = per.setdefault(m, [0, 0]) # [compared, differ] + slot[0] += 1 + slot[1] += int(va["g"] != vb["g"]) + return {m: {"compared": c, "differ": d, "pct": round(100.0 * d / c, 1) if c else None} + for m, (c, d) in sorted(per.items())} + + +def check_cell(cell: dict, days: int) -> dict: + """Fetch both sources for a cell and compare them over the most recent `days`.""" + nasa = climate._fetch_history_nasa(cell) # new primary (+ Meteostat gusts) + om = climate._fetch_history(cell) # Open-Meteo archive (ERA5) + end = min(nasa["date"].max(), om["date"].max()) + start = end - datetime.timedelta(days=days) + return {"values": compare_values(nasa, om), + "grades": grade_band_divergence(nasa, om, start, end)} + + +def _mean(xs: list) -> "float | None": + xs = [x for x in xs if x is not None] + return round(sum(xs) / len(xs), 3) if xs else None + + +def main(limit: "int | None" = None, days: int = 365) -> None: + todo = cities.all_cities() + if limit: + todo = todo[:limit] + mad_by_metric: dict[str, list] = {} + grade_by_metric: dict[str, list] = {} + checked = failed = 0 + for i, c in enumerate(todo, 1): + cell = grid.snap(c["lat"], c["lon"]) + try: + rep = check_cell(cell, days) + except Exception as e: # noqa: BLE001 - one cell's outage shouldn't abort the sweep + failed += 1 + print(f"[{i}/{len(todo)}] FAILED {c['slug']}: {e}") + continue + checked += 1 + vals = ", ".join(f"{m} mad={s['mad']}" for m, s in rep["values"].items() + if s.get("mad") is not None) + grds = ", ".join(f"{m} {s['pct']}%" for m, s in rep["grades"].items()) + print(f"[{i}/{len(todo)}] {c['slug']} ({cell['id']})") + print(f" values: {vals}") + print(f" grade drift: {grds}") + for m, s in rep["values"].items(): + mad_by_metric.setdefault(m, []).append(s.get("mad")) + for m, s in rep["grades"].items(): + grade_by_metric.setdefault(m, []).append(s.get("pct")) + + print(f"\n=== summary over {checked} cells (failed={failed}) ===") + print("mean MAD (raw drift):") + for m in _METRICS: + if m in mad_by_metric: + print(f" {m:7s} {_mean(mad_by_metric[m])}") + print("mean grade-band divergence (share of days a grade changes):") + for m in sorted(grade_by_metric): + print(f" {m:7s} {_mean(grade_by_metric[m])}%") + + +if __name__ == "__main__": + argv = sys.argv[1:] + lim = int(argv[argv.index("--limit") + 1]) if "--limit" in argv else None + dys = int(argv[argv.index("--days") + 1]) if "--days" in argv else 365 + main(limit=lim, days=dys) diff --git a/backend/tests/test_drift_check.py b/backend/tests/test_drift_check.py new file mode 100644 index 0000000..698dc95 --- /dev/null +++ b/backend/tests/test_drift_check.py @@ -0,0 +1,55 @@ +"""Unit tests for the drift-check comparison logic (hermetic — no network). The +source fetch (check_cell) is not tested here; it runs on a networked box.""" +import datetime +import math + +import polars as pl +import pytest + +import drift_check +from data import climate + +_D = [datetime.date(2024, 1, d) for d in (1, 2, 3)] + + +def test_compare_values_stats(): + a = pl.DataFrame({"date": _D, "tmax": [70.0, 80.0, 90.0], "wind": [5.0, 6.0, 7.0]}) + b = pl.DataFrame({"date": _D, "tmax": [72.0, 78.0, 90.0], "wind": [5.0, 6.0, 7.0]}) + out = drift_check.compare_values(a, b, metrics=("tmax", "wind")) + assert out["tmax"] == {"n": 3, "mad": pytest.approx(4 / 3, abs=1e-3), + "bias": pytest.approx(0.0), "max_abs": pytest.approx(2.0)} + assert out["wind"] == {"n": 3, "mad": 0.0, "bias": 0.0, "max_abs": 0.0} + + +def test_compare_values_skips_nulls_and_missing_metric(): + a = pl.DataFrame({"date": _D, "tmax": [70.0, 80.0, 90.0]}) + b = pl.DataFrame({"date": _D, "tmax": [71.0, None, 90.0]}) + out = drift_check.compare_values(a, b, metrics=("tmax", "gust")) + assert out["tmax"]["n"] == 2 # the null day is dropped + assert out["tmax"]["mad"] == pytest.approx(0.5) # |70-71| and |90-90| + assert "gust" not in out # metric absent from both frames + + +def _gradeable(n=800): + """A finalized, gradeable multi-year frame (seasonal temp swing so windows have a + real distribution).""" + start = datetime.date(2019, 1, 1) + dates = [start + datetime.timedelta(days=i) for i in range(n)] + tmax = [70.0 + 15.0 * math.sin(i / 58.0) for i in range(n)] + tmin = [t - 15.0 for t in tmax] + return climate._finalize_frame(pl.DataFrame({ + "date": dates, "tmax": tmax, "tmin": tmin, "precip": [0.0] * n, + "wind": [5.0] * n, "gust": [8.0] * n, "humid": [60.0] * n, + "fmax": tmax, "fmin": tmin, + })) + + +def test_grade_band_divergence_identical_is_zero(): + f = _gradeable() + end = f["date"].max() + start = end - datetime.timedelta(days=60) + div = drift_check.grade_band_divergence(f, f, start, end) + assert div # some metrics were compared + assert any(s["compared"] > 0 for s in div.values()) + assert all(s["differ"] == 0 for s in div.values()) # identical frames never diverge + assert all(s["pct"] == 0.0 for s in div.values())