buzeraandClaude Opus 5.5
26c415fa6c
The Forecast Error sensors replayed the model with the relay, electrical heat and window history that actually followed, so they measured the model, not a forecast anyone could have made. They are renamed Model Error (unique IDs unchanged) and gain RMSE, effective sample count and a quality label; their state stays unknown until there are at least 5 roughly independent spans. A new Forecast Error sensor scores forecasts as published, with the plan and assumptions of the moment they were made. A forecast is scored only when the relay followed its plan; otherwise it is counted as counterfactual, not as error. Statistics are saved with the model and dropped with a zone's recalibration. Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
209 lines
11 KiB
Python
209 lines
11 KiB
Python
"""Walk-forward forecast evaluation, without Home Assistant."""
|
|
import sys
|
|
import types
|
|
import unittest
|
|
from datetime import datetime, timedelta, timezone
|
|
from pathlib import Path
|
|
|
|
COMPONENTS = Path(__file__).resolve().parents[1] / 'custom_components' / 'heating_mpc'
|
|
# Fake only the package container, to avoid importing Home Assistant from __init__.
|
|
sys.modules.setdefault('heating_mpc', types.ModuleType('heating_mpc'))
|
|
sys.modules['heating_mpc'].__path__ = [str(COMPONENTS)]
|
|
from heating_mpc import evaluation as ev # noqa: E402
|
|
|
|
THETA = (0.004, 0.12, 0.0005)
|
|
START = datetime(2026, 10, 26, tzinfo=timezone.utc)
|
|
|
|
|
|
def plant(polls, theta=THETA, model_theta=None, outside=5.0, lag_at=lambda k: 1.0 if k % 30 < 12 else 0.0,
|
|
evaluator=None, uncertainty_at=lambda k: 0.0, start=START, step=300):
|
|
"""Feed polls of a room that follows the 5-minute model exactly with `theta`."""
|
|
a, b, c = theta
|
|
evaluator = evaluator or ev.ForecastEvaluator(0.1)
|
|
temp = 20.0
|
|
for k in range(polls):
|
|
lag = lag_at(k)
|
|
if k:
|
|
temp += a * (outside - temp) + b * lag + c
|
|
evaluator.observe(start + timedelta(seconds=k * step), {'sensor.room': temp}, outside, lag,
|
|
uncertainty_at(k), {'sensor.room': model_theta or theta})
|
|
return evaluator
|
|
|
|
|
|
class EvaluationTests(unittest.TestCase):
|
|
def test_exact_model_scores_zero_error_and_beats_persistence(self):
|
|
summary = plant(60).summary(['sensor.room'])
|
|
for horizon, expected_samples in ((15, 57), (30, 54), (60, 48), (120, 36)):
|
|
with self.subTest(horizon=horizon):
|
|
stats = summary[horizon]
|
|
self.assertEqual(stats['samples'], expected_samples)
|
|
self.assertEqual(stats['rooms']['sensor.room']['samples'], expected_samples)
|
|
self.assertAlmostEqual(stats['mae'], 0, places=6)
|
|
self.assertGreater(stats['persistence_mae'], 0.05)
|
|
self.assertEqual(stats['skill'], 1)
|
|
|
|
def test_parameters_from_the_forecast_start_are_used(self):
|
|
# A model without heat gain under-predicts warming: negative bias.
|
|
stats = plant(10, model_theta=(0.004, 0.0, 0.0005), lag_at=lambda k: 1.0).summary(['sensor.room'])
|
|
self.assertLess(stats[15]['bias'], -0.3)
|
|
self.assertEqual(stats[15]['bias'], -stats[15]['mae'])
|
|
|
|
def test_replay_matches_hand_calculation(self):
|
|
temp = 20.0
|
|
for _ in range(3):
|
|
temp += 0.004 * (5 - temp) + 0.12 * 0.5 + 0.0005
|
|
self.assertAlmostEqual(ev.replay(20.0, 5.0, THETA, [(0.5, 300)] * 3), temp)
|
|
# Jitter is normalized like learning: 330 s counts as 1.1 steps.
|
|
self.assertAlmostEqual(ev.replay(20.0, 5.0, (0, 0.1, 0), [(1.0, 330)]), 20.11)
|
|
|
|
def test_breaks_uncertain_lag_and_missing_rooms_are_not_scored(self):
|
|
evaluator = plant(4)
|
|
self.assertEqual(evaluator.summary(['sensor.room'])[15]['samples'], 1)
|
|
evaluator.reset_window()
|
|
plant(3, evaluator=evaluator, start=START + timedelta(minutes=20))
|
|
self.assertEqual(evaluator.summary(['sensor.room'])[15]['samples'], 1)
|
|
# Irregular spacing restarts the window as well.
|
|
jittery = plant(8, step=400)
|
|
self.assertEqual(jittery.summary(['sensor.room'])[15]['samples'], 0)
|
|
uncertain = plant(8, uncertainty_at=lambda k: 0.5 if k == 4 else 0.0)
|
|
self.assertEqual(uncertain.summary(['sensor.room'])[15]['samples'], 2) # k = 3 and 7 only
|
|
self.assertEqual(plant(8).summary(['sensor.other'])[15]['samples'], 0)
|
|
|
|
def test_rmse_effective_samples_and_quality(self):
|
|
stats = plant(10, model_theta=(0.004, 0.0, 0.0005), lag_at=lambda k: 1.0).summary(['sensor.room'])
|
|
s15 = stats[15]
|
|
self.assertGreaterEqual(s15['rmse'], s15['mae'])
|
|
# 7 overlapping 15-minute forecasts started 5 minutes apart: ~2.3 independent ones.
|
|
self.assertEqual((s15['samples'], s15['effective_samples'], s15['quality']), (7, 2.3, 'insufficient'))
|
|
self.assertEqual(ev.quality(5), 'indicative')
|
|
self.assertEqual(ev.quality(30), 'ok')
|
|
self.assertEqual(ev.effective_samples(120, 60, zones=2), 5.0)
|
|
|
|
def test_window_keeps_seven_days(self):
|
|
evaluator = plant(10)
|
|
before = evaluator.summary(['sensor.room'])[15]['samples']
|
|
plant(10, evaluator=evaluator, start=START + timedelta(days=7))
|
|
self.assertEqual(evaluator.summary(['sensor.room'])[15]['samples'], 7)
|
|
self.assertEqual(before, 7)
|
|
|
|
def test_full_model_replays_actual_inputs_and_is_compared_with_base(self):
|
|
base = (0.004, 0.12, 0.0005)
|
|
full = {'base': base, 'uses_extended': True,
|
|
'extended': (['loss', 'heat', 'bias', 'internal'], base + (0.15,)),
|
|
'ventilation': (['vent_natural'], (0.01,))}
|
|
evaluator, temp = ev.ForecastEvaluator(0.1), 20.0
|
|
for k in range(40):
|
|
internal = 0.4 if (k // 8) % 2 else 0.0 # gaming on and off
|
|
opening = 1.0 if k % 10 in (3, 4) else 0.0 # short ventilation
|
|
if k:
|
|
temp += 0.004 * (5 - temp) + 0.12 * 0.5 + 0.0005 + 0.15 * internal + 0.01 * opening * (5 - temp)
|
|
evaluator.observe(START + timedelta(seconds=300 * k), {'r': temp}, 5.0, 0.5, 0.0, {'r': base},
|
|
models={'r': full},
|
|
inputs={'r': {'internal': internal, 'ventilation': {'opening': opening}}})
|
|
stats = evaluator.summary(['r'])[60]
|
|
self.assertAlmostEqual(stats['mae'], 0, places=6)
|
|
self.assertGreater(stats['base_model_mae'], 0.05)
|
|
self.assertEqual(stats['improvement_vs_base'], 1.0)
|
|
|
|
def test_saved_statistics_round_trip_and_reject_garbage(self):
|
|
evaluator = plant(10)
|
|
restored = ev.ForecastEvaluator(0.1)
|
|
restored.load(evaluator.as_dict())
|
|
self.assertEqual(restored.summary(['sensor.room']), evaluator.summary(['sensor.room']))
|
|
day = START.date().isoformat()
|
|
restored.load({'days': {
|
|
'not-a-day': {'15': {'sensor.room': [1, 0.0, 0.1, 0.1]}},
|
|
day: {'45': {'sensor.room': [1, 0.0, 0.1, 0.1]},
|
|
'15': {'sensor.a': [1, 0.0, -0.1, 0.1], 'sensor.b': [True, 0, 0, 0],
|
|
'sensor.c': [1, float('nan'), 0, 0], 'sensor.d': [2, 0.1, 0.2, 0.4],
|
|
'sensor.e': [2, 0.1, 0.2, 0.4, 0.3], 'sensor.f': [2, 0.1, 0.2, 0.4, -1]}}}})
|
|
# A 0.3 bucket has no base-model column; its model was the base model then.
|
|
# Neither 0.3 nor 0.4 buckets have squared errors: they do not count for RMSE.
|
|
self.assertEqual(restored.as_dict(), {'days': {day: {'15': {
|
|
'sensor.d': [2, 0.1, 0.2, 0.4, 0.2, 0.0, 0], 'sensor.e': [2, 0.1, 0.2, 0.4, 0.3, 0.0, 0]}}}})
|
|
self.assertIsNone(restored.summary(['sensor.d'])[15]['rmse'])
|
|
for garbage in (None, [], {'days': []}):
|
|
restored.load(garbage)
|
|
self.assertEqual(restored.as_dict(), {'days': {}})
|
|
|
|
|
|
T0 = START.timestamp()
|
|
|
|
|
|
def run_ledger(plan, relay_on, temps, trajectory, start=20.0, polls=25, unknown_relay_at=None, skip=()):
|
|
"""One forecast at T0, then a poll every 5 minutes with the given relay state per poll."""
|
|
ledger = ev.OperationalForecastLedger()
|
|
ledger.observe(T0, {'r': start}, None, None)
|
|
ledger.record(T0, {'r': start}, {'r': trajectory}, plan, 'rev1', {'strategy': 'mpc_plan'})
|
|
for k in range(1, polls):
|
|
if k in skip:
|
|
continue
|
|
on = None if k == unknown_relay_at else (300.0 if relay_on(k) else 0.0)
|
|
ledger.observe(T0 + 300 * k, {'r': temps(k)}, on, 300.0)
|
|
return ledger.summary(['r'])
|
|
|
|
|
|
class OperationalForecastTests(unittest.TestCase):
|
|
PLAN_OFF = [0] * 8
|
|
TRAJECTORY = [19.9, 19.8, 19.7, 19.6, 19.5, 19.4, 19.3, 19.2] # every 15 minutes
|
|
|
|
def test_followed_plan_is_scored_at_every_horizon_against_persistence(self):
|
|
stats = run_ledger(self.PLAN_OFF, lambda k: False, lambda k: 20.0 - k / 30, self.TRAJECTORY)
|
|
for horizon, actual in ((15, 19.9), (30, 19.8), (60, 19.6), (120, 19.2)):
|
|
with self.subTest(horizon=horizon):
|
|
s = stats[horizon]
|
|
self.assertEqual((s['samples'], s['counterfactual'], s['unverifiable']), (1, 0, 0))
|
|
self.assertAlmostEqual(s['mae'], 0, places=3)
|
|
self.assertAlmostEqual(s['persistence_mae'], 20.0 - actual, places=3)
|
|
self.assertEqual(s['improvement_vs_persistence'], 1.0)
|
|
self.assertEqual(s['quality'], 'insufficient') # one forecast proves little
|
|
|
|
def test_plan_not_followed_is_counterfactual_not_scored(self):
|
|
stats = run_ledger([1] * 8, lambda k: False, lambda k: 20.0, self.TRAJECTORY)
|
|
for horizon in (15, 30, 60, 120):
|
|
self.assertEqual((stats[horizon]['samples'], stats[horizon]['counterfactual']), (0, 1))
|
|
self.assertIsNone(stats[horizon]['mae'])
|
|
|
|
def test_plan_followed_for_a_while(self):
|
|
# The thermostat switches the boiler on in the 4th block (after 45 minutes).
|
|
stats = run_ledger(self.PLAN_OFF, lambda k: k > 9, lambda k: 20.0, self.TRAJECTORY)
|
|
self.assertEqual([stats[h]['samples'] for h in (15, 30, 60, 120)], [1, 1, 0, 0])
|
|
self.assertEqual([stats[h]['counterfactual'] for h in (15, 30, 60, 120)], [0, 0, 1, 1])
|
|
|
|
def test_unknown_relay_history_or_missing_observations(self):
|
|
stats = run_ledger(self.PLAN_OFF, lambda k: False, lambda k: 20.0, self.TRAJECTORY, unknown_relay_at=2)
|
|
self.assertEqual((stats[15]['samples'], stats[15]['unverifiable'], stats[60]['unverifiable']), (0, 1, 1))
|
|
missing = run_ledger(self.PLAN_OFF, lambda k: False, lambda k: None if k == 12 else 20.0, self.TRAJECTORY)
|
|
self.assertEqual((missing[60]['samples'], missing[60]['missing'], missing[30]['samples']), (0, 1, 1))
|
|
# No poll near the due time (polls 11 to 13 missing): the 60-minute point cannot be scored.
|
|
gap = run_ledger(self.PLAN_OFF, lambda k: False, lambda k: 20.0, self.TRAJECTORY, skip=(11, 12, 13))
|
|
self.assertEqual((gap[60]['samples'], gap[60]['missing']), (0, 1))
|
|
|
|
def test_due_time_alignment(self):
|
|
ledger = ev.OperationalForecastLedger()
|
|
ledger.record(T0, {'r': 20.0}, {'r': self.TRAJECTORY}, self.PLAN_OFF, 'rev1', {})
|
|
for offset in (300, 600):
|
|
ledger.observe(T0 + offset, {'r': 20.0}, 0.0, 300.0)
|
|
ledger.observe(T0 + 900 - 200, {'r': 20.0}, 0.0, 100.0) # too early for +15 min
|
|
self.assertEqual(ledger.summary(['r'])[15]['samples'], 0)
|
|
ledger.observe(T0 + 900 + 100, {'r': 19.9}, 0.0, 300.0) # within half a poll of the due time
|
|
self.assertEqual(ledger.summary(['r'])[15]['samples'], 1)
|
|
|
|
def test_round_trip_forget_and_garbage(self):
|
|
stats_ledger = ev.OperationalForecastLedger()
|
|
stats_ledger.record(T0, {'r': 20.0}, {'r': self.TRAJECTORY}, [1] * 8, 'rev1', {})
|
|
stats_ledger.record(T0, {'r': 20.0}, {'r': self.TRAJECTORY}, self.PLAN_OFF, 'rev1', {})
|
|
for k in range(1, 4):
|
|
stats_ledger.observe(T0 + 300 * k, {'r': 19.9}, 0.0, 300.0)
|
|
restored = ev.OperationalForecastLedger(stats_ledger.as_dict())
|
|
self.assertEqual(restored.summary(['r']), stats_ledger.summary(['r']))
|
|
self.assertEqual(restored.summary(['r'])[15]['counterfactual'], 1)
|
|
restored.forget('r')
|
|
self.assertEqual(restored.summary(['r'])[15]['samples'], 0)
|
|
for garbage in (None, [], {'days': {'x': {}}, 'verdicts': {'2026-10-27': {'15': {'bogus': 3}}}}):
|
|
self.assertEqual(ev.OperationalForecastLedger(garbage).as_dict(), {'days': {}, 'verdicts': {}})
|
|
|
|
|
|
if __name__ == '__main__':
|
|
unittest.main()
|