回测一步向前的需求预测

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

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']

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

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

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")

Expanding-window predictions: 366 ; all training dates precede forecast dates

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()

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()

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()

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

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