Backtest one-day-ahead demand forecasts
Construct past-only lag features, refit expanding windows and compare ridge with two transparent persistence baselines.
Goal and how to run
This lesson assumes basic Python knowledge. Read each method and caption before changing parameters. AI exercises are optional; no AI account is needed for the core workflow. Data lives in notebook metadata, so there is no companion data download. Keep the notebook saved in the active working folder.
Your goal is to understand the research question, execute transparent calculations, check the results, and state the limits of the conclusion.
Setup and parameters
from pathlib import Path
import base64, gzip, hashlib, json, platform
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
from IPython.display import display, HTML
from datetime import datetime
from IPython.display import display as _display
def display(value):
if isinstance(value,pd.DataFrame):
table=value.to_html(max_rows=12,escape=True)
style="<style>.vibeit-readable-table{max-width:100%;overflow-x:auto}.vibeit-readable-table table{width:max-content;max-width:none;border-collapse:collapse}.vibeit-readable-table td,.vibeit-readable-table th{white-space:nowrap!important;word-break:normal!important;overflow-wrap:normal!important}</style>"
return _display(HTML(style+'<div class="vibeit-readable-table">'+table+'</div>'))
return _display(value)
RECIPE_ID='ds06-bike-backtesting'
OUTPUT_DIR=Path.cwd()/"results"/RECIPE_ID
plt.rcParams.update({"figure.dpi":120,"font.size":10,"axes.spines.top":False,"axes.spines.right":False,
"axes.prop_cycle":plt.cycler(color=["#18785f","#416aa6","#bc6541","#9974a5"])})
print("Python:",platform.python_version(),"; NumPy:",np.__version__,"; Pandas:",pd.__version__)
print("Working folder:",Path.cwd())
print("Core data mode: embedded snapshot; no network requests")
Python: 3.13.13 ; NumPy: 2.5.3 ; Pandas: 3.0.6 Working folder: /Users/zhiluo/Developer/vibeit-cookbook/cookbook/downloads Core data mode: embedded snapshot; no network requests
Load the embedded snapshot
def load_snapshot(recipe_id, expected_sha256, snapshot_key):
"""Read embedded data from this notebook, including a renamed copy."""
for path in sorted(Path.cwd().glob("*.ipynb")):
try:
notebook = json.loads(path.read_text(encoding="utf-8"))
metadata = notebook.get("metadata", {}).get("vibeit_cookbook", {})
blob = metadata.get("snapshots", {}).get(snapshot_key)
if metadata.get("recipe_id") != recipe_id or not blob:
continue
if blob.get("sha256") != expected_sha256:
continue
raw = gzip.decompress(base64.b64decode(blob["gzip_base64"], validate=True))
if hashlib.sha256(raw).hexdigest() != expected_sha256:
raise ValueError("Embedded snapshot checksum mismatch: " + snapshot_key)
return json.loads(raw)
except (OSError, json.JSONDecodeError):
continue
raise FileNotFoundError(
"Save/open the complete .ipynb in the active working folder. "
"The embedded snapshot could not be located for " + recipe_id)
SNAPSHOT_HASHES={'bike': '99d075d1b2d3c2c9ce73f256dcc370081a51d1df60347f59674a8ed3edf92f92'}
def snapshot(key):
return load_snapshot(RECIPE_ID,SNAPSHOT_HASHES[key],key)
print("Embedded snapshots:",list(SNAPSHOT_HASHES))
Embedded snapshots: ['bike']
Read the method functions
def ridge_fit_predict(train_values,train_target,test_values,penalty=1.0):
"""Scale using training statistics and solve ridge by an augmented least-squares system."""
train_values=np.asarray(train_values,dtype=float);test_values=np.asarray(test_values,dtype=float)
target=np.asarray(train_target,dtype=float)
if train_values.ndim!=2 or len(train_values)==0 or target.shape!=(len(train_values),) or penalty<0:raise ValueError('Invalid training data or penalty')
if test_values.ndim!=2 or train_values.shape[1]!=test_values.shape[1]:raise ValueError('Training/test feature dimensions differ')
if not np.isfinite(train_values).all() or not np.isfinite(test_values).all() or not np.isfinite(target).all():raise ValueError('Non-finite training data')
mean=train_values.mean(axis=0);scale=train_values.std(axis=0);scale=np.where(scale==0,1,scale)
design=np.column_stack([np.ones(len(train_values)),(train_values-mean)/scale])
test_design=np.column_stack([np.ones(len(test_values)),(test_values-mean)/scale])
regularizer=np.sqrt(penalty)*np.eye(design.shape[1]);regularizer[0,0]=0
coefficients=np.linalg.lstsq(np.vstack([design,regularizer]),np.r_[target,np.zeros(design.shape[1])],rcond=None)[0]
return dict(predictions=test_design@coefficients,coefficients=coefficients,mean=mean,scale=scale,train_predictions=design@coefficients)
print("Method ready: ridge_fit_predict")
Method ready: ridge_fit_predict
def regression_metrics(truth,predictions):
truth=np.asarray(truth,dtype=float);predictions=np.asarray(predictions,dtype=float)
if truth.shape!=predictions.shape or truth.size==0 or not np.isfinite(truth-predictions).all():raise ValueError('Invalid regression comparison')
error=predictions-truth;denominator=np.sum((truth-truth.mean())**2)
return dict(mae=float(np.mean(np.abs(error))),rmse=float(np.sqrt(np.mean(error**2))),r_squared=float(1-np.sum(error**2)/denominator) if denominator>0 else float('nan'))
print("Method ready: regression_metrics")
Method ready: regression_metrics
def iso_datetimes(values):
"""Parse source ISO dates with the standard library; keep tables as ISO strings."""
return [datetime.fromisoformat(str(value)) for value in values]
print("Method ready: iso_datetimes")
Method ready: iso_datetimes
Data and model boundary
The source contains one daily rental count per date in 2011–2012. These are rental events, not unique people. Casual and registered counts sum to the response and are excluded from predictors. Source-normalized daily weather fields are not silently relabelled as degrees Celsius.
Steps
1. Make strictly past-only features
Shift counts before creating the seven-day average. Each row predicts its own day using previous observations and known calendar information. This is sequential one-day-ahead forecasting; observed test outcomes become legitimate inputs for later days, not for their own prediction.
bike=snapshot("bike")
frame=pd.DataFrame(bike["rows"],columns=bike["columns"])
frame["date"]=frame["dteday"].astype(str)
frame=frame.sort_values("date").reset_index(drop=True)
assert len(frame)==731 and frame["date"].is_unique
dates=iso_datetimes(frame.date)
assert all((b-a).total_seconds()==86400 for a,b in zip(dates,dates[1:]))
assert (frame["casual"]+frame["registered"]).eq(frame["cnt"]).all()
assert frame["cnt"].ge(0).all()
print("Complete daily observations:",len(frame),"; target: rentals/day")
display(frame[["date","cnt","temp","hum","windspeed"]].head(7))
frame["lag1"]=frame.cnt.shift(1);frame["lag7"]=frame.cnt.shift(7)
frame["past7_mean"]=frame.cnt.shift(1).rolling(7).mean()
frame["weekday_sin"]=np.sin(2*np.pi*np.array([value.weekday() for value in dates])/7)
frame["weekday_cos"]=np.cos(2*np.pi*np.array([value.weekday() for value in dates])/7)
model_frame=frame.dropna(subset=["lag1","lag7","past7_mean"]).reset_index(drop=True)
features=["lag1","lag7","past7_mean","weekday_sin","weekday_cos"]
assert np.isclose(model_frame.iloc[0].past7_mean,frame.cnt.iloc[:7].mean())
display(model_frame[["date","cnt"]+features].head(6))
Complete daily observations: 731 ; target: rentals/day
| date | cnt | temp | hum | windspeed | |
|---|---|---|---|---|---|
| 0 | 2011-01-01 | 985 | 0.344167 | 0.805833 | 0.160446 |
| 1 | 2011-01-02 | 801 | 0.363478 | 0.696087 | 0.248539 |
| 2 | 2011-01-03 | 1349 | 0.196364 | 0.437273 | 0.248309 |
| 3 | 2011-01-04 | 1562 | 0.200000 | 0.590435 | 0.160296 |
| 4 | 2011-01-05 | 1600 | 0.226957 | 0.436957 | 0.186900 |
| 5 | 2011-01-06 | 1606 | 0.204348 | 0.518261 | 0.089565 |
| 6 | 2011-01-07 | 1510 | 0.196522 | 0.498696 | 0.168726 |
| date | cnt | lag1 | lag7 | past7_mean | weekday_sin | weekday_cos | |
|---|---|---|---|---|---|---|---|
| 0 | 2011-01-08 | 959 | 1510.0 | 985.0 | 1344.714286 | -0.974928 | -0.222521 |
| 1 | 2011-01-09 | 822 | 959.0 | 801.0 | 1341.000000 | -0.781831 | 0.623490 |
| 2 | 2011-01-10 | 1321 | 822.0 | 1349.0 | 1344.000000 | 0.000000 | 1.000000 |
| 3 | 2011-01-11 | 1263 | 1321.0 | 1562.0 | 1340.000000 | 0.781831 | 0.623490 |
| 4 | 2011-01-12 | 1162 | 1263.0 | 1600.0 | 1297.285714 | 0.974928 | -0.222521 |
| 5 | 2011-01-13 | 1406 | 1162.0 | 1606.0 | 1234.714286 | 0.433884 | -0.900969 |
2. Refit expanding training windows
For each 2012 day, fit only rows with earlier dates. Penalty=10 is fixed for this exercise. Persistence uses yesterday’s demand; seasonal persistence uses demand seven days earlier. Neither baseline nor model can read the same-day target.
start=int(np.flatnonzero(np.array([value.year==2012 for value in iso_datetimes(model_frame.date)]))[0])
forecasts=[]
for position in range(start,len(model_frame)):
history=model_frame.iloc[:position];current=model_frame.iloc[position:position+1]
fitted=ridge_fit_predict(history[features].to_numpy(),history.cnt.to_numpy(),current[features].to_numpy(),penalty=10.0)
forecasts.append({"date":current.date.iloc[0],"training_end":history.date.iloc[-1],"observed":float(current.cnt.iloc[0]),"ridge":float(fitted["predictions"][0]),"yesterday":float(current.lag1.iloc[0]),"last_week":float(current.lag7.iloc[0])})
result=pd.DataFrame(forecasts)
assert result.training_end.lt(result.date).all()
scores=pd.DataFrame({name:regression_metrics(result.observed,result[name]) for name in ["ridge","yesterday","last_week"]}).T
display(scores)
print("Expanding-window predictions:",len(result),"; all training dates precede forecast dates")
| mae | rmse | r_squared | |
|---|---|---|---|
| ridge | 797.233632 | 1093.985277 | 0.624896 |
| yesterday | 870.174863 | 1246.368239 | 0.513120 |
| last_week | 1110.396175 | 1561.849883 | 0.235447 |
Expanding-window predictions: 366 ; all training dates precede forecast dates
3. Read a held-out forecast interval
The first 90 test days show observed rentals and three forecasts in rentals/day. Adjacent points belong to a time series, so errors are not automatically independent. A historical winner in this interval is not guaranteed to remain best.
fig,ax=plt.subplots(figsize=(9,4));part=result.iloc[:90]
ax.plot(iso_datetimes(part.date),part.observed,color="black",label="Observed",lw=1.5)
for name in ["ridge","yesterday","last_week"]:ax.plot(iso_datetimes(part.date),part[name],label=name,alpha=.75,lw=1)
ax.set(xlabel="Forecast date",ylabel="Rentals/day",title="One-day-ahead forecasts; first 90 test days");ax.legend(ncol=2);plt.tight_layout();plt.show()
4. Compare errors by month
Monthly MAE is the mean absolute daily error within each month. Every method is evaluated on the same dates. Separate months expose seasonal changes that an overall average can hide; this retrospective breakdown is not a new training rule.
monthly=pd.DataFrame({name:(result[name]-result.observed).abs().groupby([value.month for value in iso_datetimes(result.date)]).mean() for name in ["ridge","yesterday","last_week"]})
fig,ax=plt.subplots(figsize=(9,4));monthly.plot.bar(ax=ax)
ax.set(xlabel="Test month in 2012",ylabel="MAE (rentals/day)",title="Forecast errors on common dates");plt.tight_layout();plt.show()
5. Inspect residual dependence
Lag correlations compare paired residuals separated by 1–30 days. They are descriptive coefficients, not significance tests. Serial dependence motivates time-respecting validation and block-based uncertainty methods rather than a random row split.
residual=(result.ridge-result.observed).to_numpy()
lags=np.arange(1,31);autocorrelation=np.array([np.corrcoef(residual[:-lag],residual[lag:])[0,1] for lag in lags])
fig,ax=plt.subplots(figsize=(8,4));ax.stem(lags,autocorrelation)
ax.set(xlabel="Lag (days)",ylabel="Residual correlation",ylim=(-1,1),title="Descriptive dependence of held-out errors");plt.tight_layout();plt.show()
6. Audit information time and export the backtest
Reconcile baseline predictions with shifted original observations, check every training cutoff and export the full backtest. Resetting saved outputs before rerunning must reproduce the same predictions within numerical precision.
assert len(result)==366 and result.date.is_unique
assert np.allclose(result.yesterday,model_frame.iloc[start:].lag1)
assert np.allclose(result.last_week,model_frame.iloc[start:].lag7)
assert np.isclose(scores.loc["ridge","mae"],np.mean(np.abs(residual)))
OUTPUT_DIR.mkdir(parents=True,exist_ok=True)
result.to_csv(OUTPUT_DIR/"one-day-backtest.csv",index=False);scores.to_csv(OUTPUT_DIR/"forecast-metrics.csv")
monthly.to_csv(OUTPUT_DIR/"monthly-mae.csv")
(OUTPUT_DIR/"protocol.json").write_text(json.dumps({"horizon_days":1,"features":features,"penalty":10.0,"training":"all earlier dates","test_year":2012,"observed_test_outcomes":"available only for subsequent forecast days"},indent=2))
print("Information-time checks passed; exports:",OUTPUT_DIR)
Information-time checks passed; exports: /Users/zhiluo/Developer/vibeit-cookbook/cookbook/downloads/results/ds06-bike-backtesting
Exercises and reference explanation
- Explain why rolling(7).mean() must follow shift(1).
- Explain how a sequential one-day forecast differs from a fixed multi-day forecast.
Reference explanation. Without the shift, the current target enters its own feature. Sequential evaluation can use newly observed outcomes for subsequent days; a fixed multi-day forecast cannot use future realized counts.
Keep a copy before changing parameters. Compare old and new figures and record which inputs, calculations and conclusions changed.
Optional AI coding exercises
Prompt 1
Using NumPy/Pandas/Matplotlib, add a past-14-day mean and retain strict earlier-date training cutoffs.
Prompt 2
Using the same packages, compare expanding windows with a 180-day rolling training window, keeping the same test dates and declared penalties.
Human acceptance. Each training_end precedes its forecast date, lag means equal explicit past slices, 366 common forecasts are exported and no same-day target enters a predictor.
Review the edited source, then manually run affected cells and subsequent checks. Ask the assistant to explain the purpose, packages and data provenance of its changes, and reconcile that explanation with actual output.
Optional online extension
RUN_ONLINE_UPDATE=False
if RUN_ONLINE_UPDATE:
import requests
try:
response=requests.get('https://archive.ics.uci.edu/dataset/275',timeout=20)
response.raise_for_status()
print("Source:",response.url)
print(response.text[:700])
except requests.RequestException as error:
print("Online extension unavailable:",type(error).__name__,str(error)[:160])
else:
print("Online extension skipped; embedded core data unchanged")
Online extension skipped; embedded core data unchanged
Sources, snapshots and next steps
- UCI bike — CC BY 4.0; credit UCI and original dataset authors; retrieved 2026-10-09. SHA-256
b70182d0d0508e9abbb79306ce5c0cec34869000f8220175ac83d11dbe845401. Extract day.csv from original download; preserve all rows and declared units; no imputation. Fanaee-T, H. (2013). Bike Sharing. UCI Machine Learning Repository. DOI: 10.24432/C5W894.
Next steps. Transfer these methods to your own public or authorized data while retaining input checks, recorded parameters and result validation. Document the new experimental design and method limits; this example does not validate a new experiment. Full data licenses and generation instructions accompany the cookbook.