From d5da0cec2b9b53220906bb497cb5886fa78a2bfc Mon Sep 17 00:00:00 2001 From: Graham McIntire Date: Tue, 28 Apr 2026 14:38:02 -0500 Subject: [PATCH] feat(ml): add prop.compare task for algorithm-vs-ML-vs-reality calibration MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit `mix prop.compare` runs both scorers against a recent contact sample and measures both against achieved contact distance. Each run appends to priv/calibration/history.jsonl, so trends in alg-ML divergence and algorithm-distance correlation surface over weeks of running. When drift exceeds threshold (alg-ML RMSE > 8 score points, or current correlation > 0.10 below the history median), the task writes a recommendation pointing at the right remediation: mix recalibrate_scorer # algorithm drift → refit weights mix propagation_train # ML drift → retrain on the new algorithm That is the feedback loop: measure → recommend → recalibrate/retrain → loop. Operational rather than automatic, so a bad data window can't silently corrupt the model. Pure analysis lives in Microwaveprop.Propagation.Calibration with its own unit tests; the Mix task only handles data loading, file I/O, and console formatting. --- .gitignore | 3 + lib_ml/calibration.ex | 227 ++++++++ lib_ml/prop_compare.ex | 533 ++++++++++++++++++ .../propagation/calibration_test.exs | 141 +++++ 4 files changed, 904 insertions(+) create mode 100644 lib_ml/calibration.ex create mode 100644 lib_ml/prop_compare.ex create mode 100644 test/microwaveprop/propagation/calibration_test.exs diff --git a/.gitignore b/.gitignore index 87335bc1..cd593324 100644 --- a/.gitignore +++ b/.gitignore @@ -47,6 +47,9 @@ microwaveprop-*.tar # Dev-only propagation score binaries written by PropagationGridWorker /priv/dev_scores/ +# Calibration reports written by `mix prop.compare` — local artifacts only +/priv/calibration/ + # Dialyzer PLT files /priv/plts/ diff --git a/lib_ml/calibration.ex b/lib_ml/calibration.ex new file mode 100644 index 00000000..6a6678ed --- /dev/null +++ b/lib_ml/calibration.ex @@ -0,0 +1,227 @@ +defmodule Microwaveprop.Propagation.Calibration do + @moduledoc """ + Pure analysis functions for the algorithm-vs-ML-vs-reality + comparison reports written by `mix prop.compare`. + + Every function here takes plain maps and lists, so the Mix task can + load data however it wants and these helpers stay easy to test + without a database. + + Sample shape (one map per QSO): + + %{ + algorithm_score: 0..100, # composite score from Scorer + ml_score: 0..100, # ML model prediction (or nil if unavailable) + distance_km: float(), # actual contact distance — the empirical signal + band: integer(), # MHz + # optional, copied through for outlier inspection: + contact_id: term(), conditions: map() + } + """ + + @score_buckets [ + {"0-20", 0, 20}, + {"20-40", 20, 40}, + {"40-60", 40, 60}, + {"60-80", 60, 80}, + {"80-100", 80, 101} + ] + + @doc """ + Spearman rank correlation. Returns `nil` for fewer than 3 pairs. + Tied ranks are averaged. Returns `0` (by convention) when one input + has zero variance. + """ + @spec spearman([number()], [number()]) :: float() | nil + def spearman(xs, ys) when length(xs) != length(ys), do: nil + def spearman(xs, _) when length(xs) < 3, do: nil + + def spearman(xs, ys) do + n = length(xs) + rx = ranks(xs) + ry = ranks(ys) + + if zero_variance?(rx) or zero_variance?(ry) do + 0.0 + else + sum_d_sq = + rx + |> Enum.zip(ry) + |> Enum.reduce(0.0, fn {a, b}, acc -> + d = a - b + acc + d * d + end) + + 1.0 - 6.0 * sum_d_sq / (n * (n * n - 1)) + end + end + + @doc "Root-mean-squared error between two equal-length numeric lists. Nil for empty." + @spec rmse([number()], [number()]) :: float() | nil + def rmse([], []), do: nil + + def rmse(xs, ys) when length(xs) == length(ys) do + n = length(xs) + + sum_sq = + xs + |> Enum.zip(ys) + |> Enum.reduce(0.0, fn {x, y}, acc -> + d = x - y + acc + d * d + end) + + :math.sqrt(sum_sq / n) + end + + @doc """ + Splits samples into the five fixed score buckets ([0-20, 20-40, ..., + 80-100]) and reports `n` + median observed distance for each. Always + returns 5 entries, even when buckets are empty. + + A well-calibrated scorer shows monotonically increasing median + distance across the buckets. + """ + @spec bucket_by_score([map()], atom()) :: [%{bucket: String.t(), n: non_neg_integer(), median_distance_km: float() | nil}] + def bucket_by_score(samples, score_key) do + Enum.map(@score_buckets, fn {label, low, high} -> + in_bucket = + Enum.filter(samples, fn s -> + score = Map.get(s, score_key) + is_number(score) and score >= low and score < high + end) + + %{ + bucket: label, + n: length(in_bucket), + median_distance_km: in_bucket |> Enum.map(& &1.distance_km) |> median() + } + end) + end + + @doc """ + One-band summary: row count, mean scores, algorithm-vs-ML RMSE, and + Spearman correlations of each score against actual contact distance. + """ + @spec band_summary([map()]) :: %{ + n: non_neg_integer(), + mean_algorithm_score: float() | nil, + mean_ml_score: float() | nil, + alg_ml_rmse: float() | nil, + alg_distance_spearman: float() | nil, + ml_distance_spearman: float() | nil + } + def band_summary([]) do + %{ + n: 0, + mean_algorithm_score: nil, + mean_ml_score: nil, + alg_ml_rmse: nil, + alg_distance_spearman: nil, + ml_distance_spearman: nil + } + end + + def band_summary(samples) do + alg = Enum.map(samples, & &1.algorithm_score) + ml = Enum.map(samples, & &1.ml_score) + distances = Enum.map(samples, & &1.distance_km) + + %{ + n: length(samples), + mean_algorithm_score: mean(alg), + mean_ml_score: mean(ml), + alg_ml_rmse: rmse(alg, ml), + alg_distance_spearman: spearman(alg, distances), + ml_distance_spearman: spearman(ml, distances) + } + end + + @doc "Top `n` samples by absolute (algorithm − ML) score difference, descending." + @spec disagreements([map()], non_neg_integer()) :: [map()] + def disagreements(samples, n) do + samples + |> Enum.map(fn s -> Map.put(s, :divergence, abs(s.algorithm_score - s.ml_score)) end) + |> Enum.sort_by(& &1.divergence, :desc) + |> Enum.take(n) + end + + @doc """ + Drift detector. Two modes: + + * `:ml_drift` — passes the current alg-vs-ML RMSE; flags when it + exceeds `:threshold` (default 8.0 score points). + * `:algorithm_drift` — passes the current Spearman correlation + between algorithm score and contact distance plus a `:history` + list of past correlations; flags when current is more than + `:margin` (default 0.10) below the history median. + """ + @spec detect_drift(:ml_drift | :algorithm_drift, number() | nil, keyword()) :: + :ok | {:drift, String.t()} + def detect_drift(:ml_drift, current_rmse, opts) do + threshold = Keyword.get(opts, :threshold, 8.0) + + cond do + is_nil(current_rmse) -> :ok + current_rmse > threshold -> + {:drift, "ML model has drifted: alg-vs-ML RMSE = #{Float.round(current_rmse, 2)} > #{threshold}"} + true -> :ok + end + end + + def detect_drift(:algorithm_drift, current_corr, opts) do + history = Keyword.get(opts, :history, []) + margin = Keyword.get(opts, :margin, 0.10) + + cond do + is_nil(current_corr) -> :ok + history == [] -> :ok + current_corr < median(history) - margin -> + {:drift, + "Algorithm correlation dropped: current #{format(current_corr)} vs history median " <> + "#{format(median(history))} (margin #{margin})"} + + true -> + :ok + end + end + + ## ── helpers ────────────────────────────────────────────────────── + + defp ranks(values) do + values + |> Enum.with_index() + |> Enum.sort_by(&elem(&1, 0)) + |> Enum.chunk_by(&elem(&1, 0)) + |> Enum.flat_map_reduce(1, fn group, start -> + avg_rank = start + (length(group) - 1) / 2.0 + tagged = Enum.map(group, fn {_v, idx} -> {idx, avg_rank} end) + {tagged, start + length(group)} + end) + |> elem(0) + |> Enum.sort_by(&elem(&1, 0)) + |> Enum.map(&elem(&1, 1)) + end + + defp zero_variance?(list), do: list |> Enum.uniq() |> length() == 1 + + defp mean([]), do: nil + defp mean(list), do: Enum.sum(list) / length(list) + + defp median([]), do: nil + + defp median(list) do + sorted = Enum.sort(list) + n = length(sorted) + mid = div(n, 2) + + if rem(n, 2) == 0 do + (Enum.at(sorted, mid - 1) + Enum.at(sorted, mid)) / 2.0 + else + Enum.at(sorted, mid) * 1.0 + end + end + + defp format(nil), do: "n/a" + defp format(n) when is_number(n), do: :erlang.float_to_binary(n * 1.0, decimals: 3) +end diff --git a/lib_ml/prop_compare.ex b/lib_ml/prop_compare.ex new file mode 100644 index 00000000..dad2cd3d --- /dev/null +++ b/lib_ml/prop_compare.ex @@ -0,0 +1,533 @@ +defmodule Mix.Tasks.Prop.Compare do + @shortdoc "Compare algorithm scoring vs ML model vs actual contact distances" + + @moduledoc """ + Runs the static weighted algorithm and the ML model against a recent + contact sample, then measures both against the empirical signal — + achieved contact distance. + + Three things land on disk every run: + + * `priv/calibration/latest.json` — full breakdown of the most + recent run (overwritten each time). + * `priv/calibration/history.jsonl` — single-line summary appended + so trends show up over weeks/months. + * `priv/calibration/recommendation.txt` — human-readable next + action, written only when drift is detected. + + ## Why this is the feedback loop + + The ML model was trained to predict the algorithm's composite score. + The algorithm itself is a static weighted blend that was last + recalibrated on a snapshot of the contact corpus. Both can drift: + + * **ML drift** — atmospheric data shifts past the training + distribution and the ML model starts diverging from the + algorithm. Surfaced by alg-vs-ML RMSE. + * **Algorithm drift** — physical conditions change (new + stations, new bands in heavy use, climatology shifts) and the + Spearman correlation between algorithm score and observed + contact distance starts dropping. + + Each comparison run logs both metrics to `history.jsonl`. When the + rolling history shows degradation, this task emits a recommendation + pointing at the right remediation: + + mix recalibrate_scorer # algorithm drift → refit weights + mix propagation_train # ML drift → retrain model + + Run on a cron (weekly is plenty) and the loop converges as long as + the underlying contact data keeps flowing. + + ## Usage + + mix prop.compare + mix prop.compare --days 30 --samples 5000 + mix prop.compare --output priv/calibration/ + + ## Options + + * `--days N` — only consider contacts from the last N days + (default 365). Contact density varies a lot: tighten this to + a few weeks once the corpus is dense again. + * `--samples N` — cap the sampled pairs to N rows (default 5000). + * `--output DIR` — directory for report files (default `priv/calibration`). + """ + + use Mix.Task + + alias Microwaveprop.Propagation.BandConfig + alias Microwaveprop.Propagation.Calibration + alias Microwaveprop.Propagation.Model + alias Microwaveprop.Propagation.Scorer + alias Microwaveprop.Repo + + require Logger + + @bands [10_000, 24_000, 47_000, 75_000] + @query_timeout 600_000 + + @impl Mix.Task + def run(args) do + {opts, _, _} = + OptionParser.parse(args, + strict: [days: :integer, samples: :integer, output: :string] + ) + + days = Keyword.get(opts, :days, 365) + samples_cap = Keyword.get(opts, :samples, 5_000) + output_dir = Keyword.get(opts, :output, Path.join(["priv", "calibration"])) + + Application.put_env( + :microwaveprop, + Oban, + Keyword.put(Application.get_env(:microwaveprop, Oban, []), :queues, false) + ) + + Mix.Task.run("app.start") + File.mkdir_p!(output_dir) + + print_header(days, samples_cap) + + rows = load_dataset(days, samples_cap) + IO.puts("Loaded #{length(rows)} contact-profile pairs from the last #{days} days.\n") + + if rows == [] do + IO.puts("No samples — nothing to compare.") + System.halt(0) + end + + {predict_fn, params} = load_ml_or_halt() + + samples = Enum.map(rows, &score_sample(&1, predict_fn, params)) + + summaries = summarize(samples) + print_summaries(summaries) + + history = read_history(output_dir) + recommendations = build_recommendations(summaries, history) + + write_latest(output_dir, samples, summaries, recommendations) + append_history(output_dir, summaries) + write_recommendations(output_dir, recommendations) + + print_recommendations(recommendations) + + "-" |> String.duplicate(80) |> IO.puts() + IO.puts("Report files:") + IO.puts(" #{Path.join(output_dir, "latest.json")}") + IO.puts(" #{Path.join(output_dir, "history.jsonl")}") + + if recommendations != [] do + IO.puts(" #{Path.join(output_dir, "recommendation.txt")}") + end + end + + ## ── Data loading ───────────────────────────────────────────────── + + defp load_dataset(days, samples_cap) do + sql = """ + SELECT + q.id, + q.band::int AS band, + q.distance_km::float AS distance_km, + q.qso_timestamp, + EXTRACT(MONTH FROM q.qso_timestamp)::int AS month, + EXTRACT(HOUR FROM q.qso_timestamp)::int AS utc_hour, + ((q.pos1->>'lat')::float + (q.pos2->>'lat')::float) / 2.0 AS lat, + ((q.pos1->>'lon')::float + (q.pos2->>'lon')::float) / 2.0 AS lon, + h.surface_temp_c, + h.surface_dewpoint_c, + h.surface_pressure_mb, + h.min_refractivity_gradient, + h.hpbl_m, + h.pwat_mm, + h.surface_refractivity, + h.ducting_detected + FROM contacts q + INNER JOIN hrrr_profiles h + ON h.lat = ROUND(((q.pos1->>'lat')::numeric + (q.pos2->>'lat')::numeric) / 2 * 8) / 8 + AND h.lon = ROUND(((q.pos1->>'lon')::numeric + (q.pos2->>'lon')::numeric) / 2 * 8) / 8 + AND h.valid_time = date_trunc('hour', q.qso_timestamp) + WHERE q.distance_km > 0 + AND q.distance_km < 3000 + AND q.qso_timestamp >= NOW() - ($1::int || ' days')::interval + AND q.flagged_invalid IS NOT TRUE + AND q.pos1 IS NOT NULL + AND q.pos2 IS NOT NULL + AND q.band::int = ANY($2::int[]) + ORDER BY q.qso_timestamp DESC + LIMIT $3 + """ + + %{rows: rows, columns: cols} = + Repo.query!(sql, [days, @bands, samples_cap], timeout: @query_timeout) + + Enum.map(rows, fn row -> + cols + |> Enum.zip(row) + |> Map.new(fn {k, v} -> {String.to_atom(k), to_native(k, v)} end) + end) + end + + defp to_native("id", <<_::binary-size(16)>> = bin), do: Ecto.UUID.cast!(bin) + defp to_native(_, %Decimal{} = d), do: Decimal.to_float(d) + defp to_native(_, other), do: other + + ## ── Scoring ────────────────────────────────────────────────────── + + defp score_sample(row, predict_fn, params) do + band_config = BandConfig.get(row.band) + conditions = build_scorer_conditions(row) + + %{score: alg_score} = Scorer.composite_score(conditions, band_config) + + ml_score = Model.predict_score(predict_fn, params, ml_features(row)) + + %{ + contact_id: row.id, + band: row.band, + distance_km: row.distance_km, + algorithm_score: alg_score, + ml_score: ml_score, + conditions: %{ + temp_c: row.surface_temp_c, + dewpoint_c: row.surface_dewpoint_c, + pressure_mb: row.surface_pressure_mb, + gradient: row.min_refractivity_gradient, + hpbl_m: row.hpbl_m, + pwat_mm: row.pwat_mm + } + } + end + + defp build_scorer_conditions(row) do + %{ + abs_humidity: Scorer.absolute_humidity(safe_float(row.surface_temp_c), safe_float(row.surface_dewpoint_c)), + temp_f: row.surface_temp_c |> safe_float() |> Scorer.c_to_f(), + dewpoint_f: row.surface_dewpoint_c |> safe_float() |> Scorer.c_to_f(), + pressure_mb: safe_float(row.surface_pressure_mb), + prev_pressure_mb: nil, + min_refractivity_gradient: safe_float(row.min_refractivity_gradient), + bl_depth_m: safe_float(row.hpbl_m), + pwat_mm: safe_float(row.pwat_mm), + sky_cover_pct: nil, + wind_speed_kts: nil, + rain_rate_mmhr: 0.0, + utc_hour: row.utc_hour, + utc_minute: 0, + month: row.month, + longitude: row.lon, + latitude: row.lat, + best_duct_band_ghz: nil + } + end + + defp ml_features(row) do + %{ + surface_temp_c: safe_float(row.surface_temp_c, 20.0), + surface_dewpoint_c: safe_float(row.surface_dewpoint_c, 10.0), + surface_pressure_mb: safe_float(row.surface_pressure_mb, 1013.0), + min_refractivity_gradient: safe_float(row.min_refractivity_gradient, -70.0), + hpbl_m: safe_float(row.hpbl_m, 500.0), + pwat_mm: safe_float(row.pwat_mm, 20.0), + surface_refractivity: safe_float(row.surface_refractivity, 320.0), + ducting_detected: row.ducting_detected || false, + latitude: row.lat, + longitude: row.lon, + utc_hour: row.utc_hour, + month: row.month, + freq_mhz: row.band + } + end + + defp safe_float(nil), do: 0.0 + defp safe_float(%Decimal{} = d), do: Decimal.to_float(d) + defp safe_float(v) when is_number(v), do: v * 1.0 + defp safe_float(_), do: 0.0 + + defp safe_float(nil, default), do: default + defp safe_float(%Decimal{} = d, _default), do: Decimal.to_float(d) + defp safe_float(v, _default) when is_number(v), do: v * 1.0 + defp safe_float(_, default), do: default + + ## ── Summaries ──────────────────────────────────────────────────── + + defp summarize(samples) do + overall = Calibration.band_summary(samples) + overall_buckets = Calibration.bucket_by_score(samples, :algorithm_score) + + by_band = + @bands + |> Map.new(fn band -> + band_samples = Enum.filter(samples, &(&1.band == band)) + {band, Calibration.band_summary(band_samples)} + end) + + top_disagreements = Calibration.disagreements(samples, 10) + + %{ + n_samples: length(samples), + overall: overall, + overall_buckets: overall_buckets, + by_band: by_band, + top_disagreements: top_disagreements + } + end + + ## ── Drift detection + recommendations ──────────────────────────── + + defp build_recommendations(summaries, history) do + ml_drift = + Calibration.detect_drift(:ml_drift, summaries.overall.alg_ml_rmse, threshold: 8.0) + + history_corrs = + Enum.flat_map(history, fn entry -> + case entry["overall"]["alg_distance_spearman"] do + n when is_number(n) -> [n] + _ -> [] + end + end) + + algorithm_drift = + Calibration.detect_drift( + :algorithm_drift, + summaries.overall.alg_distance_spearman, + history: history_corrs, + margin: 0.10 + ) + + [ + to_recommendation(:ml_drift, ml_drift, "mix propagation_train"), + to_recommendation(:algorithm_drift, algorithm_drift, "mix recalibrate_scorer") + ] + |> Enum.reject(&is_nil/1) + end + + defp to_recommendation(_kind, :ok, _action), do: nil + defp to_recommendation(kind, {:drift, msg}, action), do: %{kind: kind, reason: msg, action: action} + + ## ── File output ────────────────────────────────────────────────── + + defp write_latest(dir, samples, summaries, recommendations) do + payload = %{ + generated_at: DateTime.utc_now() |> DateTime.to_iso8601(), + n_samples: length(samples), + overall: jsonable(summaries.overall), + overall_buckets: summaries.overall_buckets, + by_band: + Map.new(summaries.by_band, fn {band, summary} -> {Integer.to_string(band), jsonable(summary)} end), + top_disagreements: + Enum.map(summaries.top_disagreements, fn s -> + %{ + contact_id: s.contact_id, + band: s.band, + distance_km: s.distance_km, + algorithm_score: s.algorithm_score, + ml_score: s.ml_score, + divergence: s.divergence, + conditions: s.conditions + } + end), + recommendations: recommendations + } + + File.write!(Path.join(dir, "latest.json"), Jason.encode_to_iodata!(payload, pretty: true)) + end + + defp append_history(dir, summaries) do + line = + Jason.encode!(%{ + date: DateTime.utc_now() |> DateTime.to_iso8601(), + n_samples: summaries.n_samples, + overall: jsonable(summaries.overall), + by_band: + Map.new(summaries.by_band, fn {band, summary} -> + {Integer.to_string(band), jsonable(summary)} + end) + }) + + File.write!(Path.join(dir, "history.jsonl"), line <> "\n", [:append]) + end + + defp read_history(dir) do + path = Path.join(dir, "history.jsonl") + + case File.read(path) do + {:ok, content} -> + content + |> String.split("\n", trim: true) + |> Enum.flat_map(fn line -> + case Jason.decode(line) do + {:ok, entry} -> [entry] + {:error, _} -> [] + end + end) + + {:error, _} -> + [] + end + end + + defp write_recommendations(_dir, []), do: :ok + + defp write_recommendations(dir, recs) do + body = + [ + "Calibration drift detected at #{DateTime.utc_now() |> DateTime.to_iso8601()}.", + "", + "Recommended actions (run in order):", + "" + | Enum.flat_map(recs, fn rec -> [" $ #{rec.action}", " reason: #{rec.reason}", ""] end) + ] + |> Enum.join("\n") + + File.write!(Path.join(dir, "recommendation.txt"), body) + end + + ## ── Console output ─────────────────────────────────────────────── + + defp print_header(days, samples_cap) do + "=" |> String.duplicate(80) |> IO.puts() + IO.puts("PROPAGATION CALIBRATION — algorithm vs ML vs empirical contact distance") + "=" |> String.duplicate(80) |> IO.puts() + IO.puts("Window: last #{days} days Sample cap: #{samples_cap}") + IO.puts("Started: #{DateTime.to_string(DateTime.utc_now())}\n") + end + + defp print_summaries(summaries) do + IO.puts(String.duplicate("-", 80)) + IO.puts("OVERALL (n=#{summaries.overall.n})") + IO.puts(String.duplicate("-", 80)) + print_summary_row(summaries.overall) + + IO.puts("\nAlgorithm-score buckets vs achieved distance:") + IO.puts(" #{pad("bucket", 10)} #{pad("n", 8)} #{pad("median km", 12)}") + IO.puts(" " <> String.duplicate("-", 30)) + + Enum.each(summaries.overall_buckets, fn b -> + IO.puts( + " #{pad(b.bucket, 10)} #{pad(Integer.to_string(b.n), 8)} " <> + pad(format_num(b.median_distance_km, 1), 12) + ) + end) + + IO.puts("\n" <> String.duplicate("-", 80)) + IO.puts("PER BAND") + IO.puts(String.duplicate("-", 80)) + + IO.puts( + " #{pad("band", 8)} #{pad("n", 6)} #{pad("mean(alg)", 11)} #{pad("mean(ml)", 11)} " <> + "#{pad("alg-ml RMSE", 13)} #{pad("alg↔dist ρ", 13)} #{pad("ml↔dist ρ", 12)}" + ) + + IO.puts(" " <> String.duplicate("-", 78)) + + Enum.each(@bands, fn band -> + summary = summaries.by_band[band] + label = format_band(band) + + IO.puts( + " #{pad(label, 8)} #{pad(Integer.to_string(summary.n), 6)} " <> + "#{pad(format_num(summary.mean_algorithm_score, 1), 11)} " <> + "#{pad(format_num(summary.mean_ml_score, 1), 11)} " <> + "#{pad(format_num(summary.alg_ml_rmse, 2), 13)} " <> + "#{pad(format_num(summary.alg_distance_spearman, 3), 13)} " <> + "#{pad(format_num(summary.ml_distance_spearman, 3), 12)}" + ) + end) + + IO.puts("\n" <> String.duplicate("-", 80)) + IO.puts("TOP DISAGREEMENTS — biggest |algorithm − ML| score gaps") + IO.puts(String.duplicate("-", 80)) + + IO.puts( + " #{pad("band", 7)} #{pad("dist km", 10)} #{pad("alg", 6)} " <> + "#{pad("ml", 6)} #{pad("Δ", 6)} notable conditions" + ) + + IO.puts(" " <> String.duplicate("-", 70)) + + summaries.top_disagreements + |> Enum.take(10) + |> Enum.each(fn s -> + grad = s.conditions[:gradient] + pwat = s.conditions[:pwat_mm] + + IO.puts( + " #{pad(format_band(s.band), 7)} " <> + "#{pad(format_num(s.distance_km, 1), 10)} " <> + "#{pad(Integer.to_string(s.algorithm_score), 6)} " <> + "#{pad(Integer.to_string(s.ml_score), 6)} " <> + "#{pad(Integer.to_string(s.divergence), 6)} " <> + "grad=#{format_num(grad, 0)} pwat=#{format_num(pwat, 1)}" + ) + end) + + IO.puts("") + end + + defp print_summary_row(s) do + IO.puts( + " mean alg=#{format_num(s.mean_algorithm_score, 2)} " <> + "mean ml=#{format_num(s.mean_ml_score, 2)} " <> + "alg-ml RMSE=#{format_num(s.alg_ml_rmse, 2)} " <> + "alg↔dist ρ=#{format_num(s.alg_distance_spearman, 3)} " <> + "ml↔dist ρ=#{format_num(s.ml_distance_spearman, 3)}" + ) + end + + defp print_recommendations([]) do + IO.puts("\n✓ No drift detected — algorithm + ML are in sync with reality.") + IO.puts(" Re-run periodically (weekly is plenty) to keep watching.") + end + + defp print_recommendations(recs) do + IO.puts("\n" <> String.duplicate("=", 80)) + IO.puts("⚠ DRIFT DETECTED — recommended actions:") + IO.puts(String.duplicate("=", 80)) + + Enum.each(recs, fn rec -> + IO.puts("\n $ #{rec.action}") + IO.puts(" " <> rec.reason) + end) + + IO.puts("") + end + + ## ── Misc ───────────────────────────────────────────────────────── + + defp load_ml_or_halt do + case Model.load() do + {:ok, params} -> + {Model.compile_predict(), params} + + :error -> + IO.puts(:stderr, "ERROR: no trained ML model at #{Model.default_path()}.") + IO.puts(:stderr, "Run `mix propagation_train` first.") + System.halt(1) + end + end + + defp jsonable(map) when is_map(map) do + Map.new(map, fn {k, v} -> + {k, + case v do + n when is_float(n) -> Float.round(n, 6) + other -> other + end} + end) + end + + defp pad(s, n), do: String.pad_trailing(to_string(s), n) + + defp format_num(nil, _), do: "n/a" + defp format_num(n, decimals) when is_number(n), do: :erlang.float_to_binary(n * 1.0, decimals: decimals) + defp format_num(_, _), do: "n/a" + + defp format_band(10_000), do: "10G" + defp format_band(24_000), do: "24G" + defp format_band(47_000), do: "47G" + defp format_band(75_000), do: "75G" + defp format_band(b), do: "#{b}MHz" +end diff --git a/test/microwaveprop/propagation/calibration_test.exs b/test/microwaveprop/propagation/calibration_test.exs new file mode 100644 index 00000000..dc5ea87d --- /dev/null +++ b/test/microwaveprop/propagation/calibration_test.exs @@ -0,0 +1,141 @@ +defmodule Microwaveprop.Propagation.CalibrationTest do + use ExUnit.Case, async: true + + alias Microwaveprop.Propagation.Calibration + + describe "spearman/2" do + test "perfect monotonic correlation is 1.0" do + assert_in_delta Calibration.spearman([1, 2, 3, 4, 5], [10, 20, 30, 40, 50]), 1.0, 1.0e-9 + end + + test "perfect inverse correlation is -1.0" do + assert_in_delta Calibration.spearman([1, 2, 3, 4, 5], [50, 40, 30, 20, 10]), -1.0, 1.0e-9 + end + + test "uncorrelated input lands near zero" do + xs = Enum.to_list(1..100) + # constant ys → undefined ranking; we expect 0 by convention. + ys = List.duplicate(42, 100) + result = Calibration.spearman(xs, ys) + assert result == nil or abs(result) < 0.05 + end + + test "returns nil for fewer than 3 pairs" do + assert Calibration.spearman([1], [1]) == nil + assert Calibration.spearman([1, 2], [3, 4]) == nil + end + end + + describe "rmse/2" do + test "zero for identical lists" do + assert Calibration.rmse([1, 2, 3], [1, 2, 3]) == 0.0 + end + + test "matches the textbook formula" do + # mean of [1, 4, 9] = 14/3 ; sqrt = 2.16 + assert_in_delta Calibration.rmse([1, 2, 3], [2, 4, 6]), :math.sqrt(14 / 3), 1.0e-9 + end + + test "returns nil for empty input" do + assert Calibration.rmse([], []) == nil + end + end + + describe "bucket_by_score/2" do + test "splits samples into 5 fixed score buckets and reports n + median distance" do + samples = [ + %{algorithm_score: 5, distance_km: 50}, + %{algorithm_score: 15, distance_km: 60}, + %{algorithm_score: 25, distance_km: 100}, + %{algorithm_score: 45, distance_km: 200}, + %{algorithm_score: 65, distance_km: 300}, + %{algorithm_score: 85, distance_km: 500} + ] + + buckets = Calibration.bucket_by_score(samples, :algorithm_score) + assert length(buckets) == 5 + + # First bucket is 0-20 with two contacts + assert hd(buckets) == %{ + bucket: "0-20", + n: 2, + median_distance_km: 55.0 + } + + # Last bucket is 80-100 with one contact + last = List.last(buckets) + assert last.bucket == "80-100" + assert last.n == 1 + assert last.median_distance_km == 500.0 + end + + test "returns 5 buckets even when some are empty" do + buckets = Calibration.bucket_by_score([%{algorithm_score: 50, distance_km: 100}], :algorithm_score) + assert length(buckets) == 5 + empty = Enum.filter(buckets, &(&1.n == 0)) + assert length(empty) == 4 + end + end + + describe "band_summary/1" do + test "computes n, mean scores, RMSE, and Spearman correlations" do + samples = + for i <- 1..10 do + %{algorithm_score: i * 10, ml_score: i * 10 + 2, distance_km: i * 50} + end + + summary = Calibration.band_summary(samples) + assert summary.n == 10 + assert summary.mean_algorithm_score == 55.0 + assert summary.mean_ml_score == 57.0 + assert_in_delta summary.alg_ml_rmse, 2.0, 1.0e-9 + assert_in_delta summary.alg_distance_spearman, 1.0, 1.0e-9 + assert_in_delta summary.ml_distance_spearman, 1.0, 1.0e-9 + end + + test "returns a zero-row summary for empty input" do + summary = Calibration.band_summary([]) + assert summary.n == 0 + assert summary.alg_ml_rmse == nil + end + end + + describe "disagreements/2" do + test "returns top N samples by absolute alg-ml score difference, descending" do + samples = [ + %{id: 1, algorithm_score: 50, ml_score: 51}, + %{id: 2, algorithm_score: 50, ml_score: 80}, + %{id: 3, algorithm_score: 50, ml_score: 30}, + %{id: 4, algorithm_score: 50, ml_score: 49}, + %{id: 5, algorithm_score: 50, ml_score: 0} + ] + + top = Calibration.disagreements(samples, 3) + assert Enum.map(top, & &1.id) == [5, 2, 3] + end + end + + describe "detect_drift/3" do + test "flags when current alg_ml_rmse exceeds the drift threshold" do + assert {:drift, msg} = Calibration.detect_drift(:ml_drift, 12.0, threshold: 8.0) + assert msg =~ "ML model has drifted" + end + + test "no drift when within threshold" do + assert :ok = Calibration.detect_drift(:ml_drift, 4.0, threshold: 8.0) + end + + test "flags algorithm drift when current correlation falls below history median by margin" do + history_corrs = [0.42, 0.45, 0.43, 0.44] + + assert {:drift, msg} = + Calibration.detect_drift(:algorithm_drift, 0.25, history: history_corrs, margin: 0.10) + + assert msg =~ "Algorithm correlation dropped" + end + + test "no algorithm drift if history empty" do + assert :ok = Calibration.detect_drift(:algorithm_drift, 0.10, history: [], margin: 0.10) + end + end +end