V. Forecasting Models¶

import pandas as pd
from darts import TimeSeries
from darts.models import ARIMA
from darts.metrics import mape
import warnings
import matplotlib as mpl

AR 1¶

AR1 8-week horizon

# Optional: suppress convergence warnings
warnings.filterwarnings("ignore")

# Ensure consistent lengths
min_len = min(len(time_index), len(energy_consumption_data))
time_index = time_index[:min_len]
energy_consumption_data = energy_consumption_data[:min_len]

# --- Convert to Darts TimeSeries (WEEKLY) ---
series = TimeSeries.from_times_and_values(
    pd.DatetimeIndex(time_index),
    energy_consumption_data,
    fill_missing_dates=True,
)

# --- Fit ARIMA and Multi-Step Forecast (Horizon = 8 weeks) ---
model = ARIMA(1, 0, 0)
horizon = 8
start = len(series) - 52 - horizon + 1

forecast = model.historical_forecasts(
    series=series,
    start=start,
    forecast_horizon=horizon,
    stride=1,
    retrain=True,
    train_length=104,
    last_points_only=True,
    show_warnings=False
)

# --- Align actual values ---
actual = series.slice(forecast.start_time(), forecast.end_time())

# Convert to Pandas
actual_df = actual.to_series().to_frame(name="Actual")
forecast_df = forecast.to_series().to_frame(name="Forecast")

# --- Set matplotlib dark mode ---
mpl.rcParams.update({
    "figure.facecolor": "#121212",
    "axes.facecolor": "#121212",
    "axes.edgecolor": "white",
    "axes.labelcolor": "white",
    "xtick.color": "white",
    "ytick.color": "white",
    "text.color": "white",
    "legend.facecolor": "#1e1e1e",
    "legend.edgecolor": "white",
    "grid.color": "gray"
})

# --- Plot using matplotlib ---
plt.figure(figsize=(12, 6))

# Plot actual values (blue)
plt.plot(actual_df.index, actual_df["Actual"], label="Actual Weekly", color='blue', marker='o', linestyle='-')

# Plot forecast values (red)
plt.plot(forecast_df.index, forecast_df["Forecast"], label="Rolling Forecasted Weekly (AR(1,0,0), 8-week horizon)", color='red', marker='o', linestyle='-')

# Title and labels
plt.title("AR(1,0,0) Rolling Forecast (8-week Horizon, Weekly)")
plt.xlabel("Date")
plt.ylabel("Energy Consumption (kWh)")
plt.legend()
plt.grid(True, linestyle='--', alpha=0.5)

# Show plot
plt.tight_layout()
plt.show()

# --- Compute MAPE ---
if len(actual) == len(forecast):
    mape_score = mape(actual, forecast)
    print(f"MAPE (Weekly, 8-week Horizon): {mape_score:.2f}%")
else:
    print("Forecast and actual series do not match in length.")

No description has been provided for this image

MAPE (Weekly, 8-week Horizon): 8.39%

ARMA¶

ARMA 8-week horizon

# --- Optional: suppress convergence warnings ---
warnings.filterwarnings("ignore")

# --- Ensure consistent lengths ---
min_len = min(len(time_index), len(energy_consumption_data))
time_index = time_index[:min_len]
energy_consumption_data = energy_consumption_data[:min_len]

# --- Convert to Darts TimeSeries (WEEKLY) ---
series = TimeSeries.from_times_and_values(
    pd.DatetimeIndex(time_index),
    energy_consumption_data,
    fill_missing_dates=True,
)

# --- Model settings ---
horizon = 8                    # Forecast horizon
train_len = 104                # Training length (2 years)
start_idx = len(series) - 52  # Start index for historical forecasts

# --- Initialize and run ARIMA(1,0,1) model ---
model = ARIMA(1, 0, 1)
start = len(series) - 52 - horizon + 1

forecast = model.historical_forecasts(
    series=series,
    start=start,
    forecast_horizon=horizon,
    stride=1,
    retrain=True,
    train_length=train_len,
    last_points_only=True,
    show_warnings=False
)

# --- Align actual values ---
actual = series.slice(forecast.start_time(), forecast.end_time())

# --- Convert to Pandas ---
actual_df = actual.to_series().to_frame(name="Actual")
forecast_df = forecast.to_series().to_frame(name="Forecast")

# --- Plot ---
plt.figure(figsize=(12, 6))

# Actual values in blue
plt.plot(actual_df.index, actual_df["Actual"], label="Actual Weekly", color="blue", marker="o", linestyle='-')

# Forecast values in red
plt.plot(forecast_df.index, forecast_df["Forecast"], label="Autoregressive Moving Average Model (1,0,1)", color="red", marker="o", linestyle='-')

# Title and labels
plt.title("Autoregressive Moving Average Model (1,0,1), 8-Week Horizon")
plt.xlabel("Date")
plt.ylabel("Energy Consumption (kWh)")
plt.legend()
plt.grid(True, linestyle='--', alpha=0.5)

# Show plot
plt.tight_layout()
plt.show()

# --- Compute MAPE ---
if len(actual) == len(forecast):
    mape_score = mape(actual, forecast)
    print(f"MAPE (Weekly, 8-week Horizon): {mape_score:.2f}%")
else:
    print("Forecast and actual series do not match in length.")

MAPE (Weekly, 8-week Horizon): 8.48%

ARIMA¶

# --- Ensure consistent lengths ---
min_len = min(len(time_index), len(energy_consumption_data))
time_index = time_index[:min_len]
energy_consumption_data = energy_consumption_data[:min_len]

# --- Convert to Darts TimeSeries (WEEKLY) ---
series = TimeSeries.from_times_and_values(
    pd.DatetimeIndex(time_index),
    energy_consumption_data,
    fill_missing_dates=True,
)

# --- ARIMA(1,1,0) model settings ---
model = ARIMA(1, 1, 0)
horizon = 8
start = len(series) - 52 - horizon + 1

forecast = model.historical_forecasts(
    series=series,
    start=start,
    forecast_horizon=horizon,
    stride=1,
    retrain=True,
    train_length=104,
    last_points_only=True,
    show_warnings=False
)

# --- Align actual values ---
actual = series.slice(forecast.start_time(), forecast.end_time())

# --- Convert to Pandas ---
actual_df = actual.to_series().to_frame(name="Actual")
forecast_df = forecast.to_series().to_frame(name="Forecast")


# --- Plot using matplotlib ---
plt.figure(figsize=(12, 6))

# Actual in blue
plt.plot(actual_df.index, actual_df["Actual"], label="Actual Weekly", color="blue", marker="o", linestyle='-')

# Forecast in red
plt.plot(forecast_df.index, forecast_df["Forecast"], label="Autoregressive Integrated Moving Average Model (1,1,0), 8-week horizon", color="red", marker="o", linestyle='-')

# Titles and labels
plt.title("Autoregressive Integrated Moving Average Model (1,1,0), 8-week Horizon")
plt.xlabel("Date")
plt.ylabel("Energy Consumption (kWh)")
plt.legend()
plt.grid(True, linestyle='--', alpha=0.5)

# --- Show plot ---
plt.tight_layout()
plt.show()

# --- Compute MAPE ---
if len(actual) == len(forecast):
    mape_score = mape(actual, forecast)
    print(f"MAPE (Weekly, 8-week Horizon): {mape_score:.2f}%")
    # Save MAPE to file
    with open("/Users/omarabdesslem/Documents/Projects/Thesis-Website/public/ARIMA/arima_mape.txt", "w") as f:
        f.write(f"{mape_score:.2f}")
else:
    print("Forecast and actual series do not match in length.")

MAPE (Weekly, 8-week Horizon): 8.42%

SARIMA¶

8 week horizon

# --- Ensure consistent lengths ---
min_len = min(len(time_index), len(energy_consumption_data))
time_index = time_index[:min_len]
energy_consumption_data = energy_consumption_data[:min_len]

# --- Convert to Darts TimeSeries (WEEKLY) ---
series = TimeSeries.from_times_and_values(
    pd.DatetimeIndex(time_index),
    energy_consumption_data,
    fill_missing_dates=True,
)

# --- Seasonal ARIMA model: ARIMA(1,0,0)(1,0,0)[52] ---
model = ARIMA(1, 0, 0, seasonal_order=(1, 0, 0, 52))
horizon = 8
start = len(series) - 52

forecast = model.historical_forecasts(
    series=series,
    start=start,
    forecast_horizon=horizon,
    stride=1,
    retrain=True,
    train_length=104,
    last_points_only=True
)

# --- Align actual values ---
actual = series.slice(forecast.start_time(), forecast.end_time())

# --- Convert to Pandas ---
actual_df = actual.to_series().to_frame(name="Actual")
forecast_df = forecast.to_series().to_frame(name="Forecast")


# --- Plot using matplotlib ---
plt.figure(figsize=(12, 6))

# Actual values (blue)
plt.plot(actual_df.index, actual_df["Actual"], label="Actual Weekly", color="blue", marker="o", linestyle="-")

# Forecast values (red)
plt.plot(forecast_df.index, forecast_df["Forecast"], label="SARIMA(1,0,0)(1,0,0)[52], 8-week horizon", color="red", marker="o", linestyle="-")

# Titles and labels
plt.title("Seasonal ARIMA Model SARIMA(1,0,0)(1,0,0)[52], 8-Week Horizon")
plt.xlabel("Date")
plt.ylabel("Energy Consumption (kWh)")
plt.legend()
plt.grid(True, linestyle="--", alpha=0.5)

# --- Show plot ---
plt.tight_layout()
plt.show()

# --- Compute and display MAPE ---
if len(actual) == len(forecast):
    mape_score = mape(actual, forecast)
    print(f"MAPE (Weekly, 8-week Horizon): {mape_score:.2f}%")
    # Save MAPE to file
    with open("/Users/omarabdesslem/Documents/Projects/Thesis-Website/public/SARIMA/sarima_mape.txt", "w") as f:
        f.write(f"{mape_score:.2f}")
else:
    print("Forecast and actual series do not match in length.")

MAPE (Weekly, 8-week Horizon): 6.60%

SARIMAX¶

Seasonal ARIMA with eXogenous variables forecasts time series data that shows trends, seasonality, and is affected by external variables.

In our setup:

Non-seasonal order: $(1, 0, 0)$
Seasonal order: $(1, 0, 0, 52)$
Exogenous input: weekly temperature values

The model looks like:

$y_t = c + \phi_1 y_{t-1} + \Phi_1 y_{t-52} + \beta X_t + \varepsilon_t$

Where:

$y_t$ is the value at time $t$ (energy consumption)
$\phi_1$ is the non-seasonal autoregressive term
$\Phi_1$ is the seasonal autoregressive term (lag of 52 weeks)
$\theta_1$ is the moving average term
$\varepsilon_t$ is the error (white noise)
$X_t$ is the exogenous variable (weekly temperature)

# --- Step 1: Preprocessing and Alignment ---
weekly_temp = weekly_temp.dropna()

# Ensure same length
min_len = min(len(time_index), len(energy_consumption_data), len(weekly_temp))
time_index = time_index[:min_len]
energy_consumption_data = energy_consumption_data[:min_len]
weekly_temp = weekly_temp[:min_len]

# Create TimeSeries objects
series = TimeSeries.from_times_and_values(
    pd.DatetimeIndex(time_index),
    energy_consumption_data,
    fill_missing_dates=True,
)

weekly_temp.index = pd.DatetimeIndex(weekly_temp.index)
temp_cov = TimeSeries.from_series(weekly_temp)

# Slice to common intersection
series = series.slice_intersect(temp_cov)
temp_cov = temp_cov.slice_intersect(series)

# --- Step 2: SARIMAX Model with Exogenous Temp (Weekly Temp) ---
model = ARIMA(1, 0, 0, seasonal_order=(1, 0, 0, 52))

# Rolling 8-week forecasts across the last 52 weeks
forecast_list = model.historical_forecasts(
    series=series,
    future_covariates=temp_cov,
    start=len(series) - 52,
    forecast_horizon=8,
    stride=1,
    retrain=True,
    train_length=104,
    last_points_only=True,
    show_warnings=False
)

# Combine forecasts
forecast = forecast_list[0]
for f in forecast_list[1:]:
    forecast = forecast.append(f)

# Match actuals
actual = series.slice(forecast.start_time(), forecast.end_time())

# Convert to Pandas
actual_df = actual.to_series().to_frame(name="Actual")
forecast_df = forecast.to_series().to_frame(name="Forecast")


# --- Step 4: Plot with Matplotlib ---
plt.figure(figsize=(12, 6))

# Plot actual
plt.plot(actual_df.index, actual_df["Actual"], label="Actual Weekly", color="blue", marker="o", linestyle='-')

# Plot forecast
plt.plot(forecast_df.index, forecast_df["Forecast"], label="SARIMAX(1,0,0)(1,0,0)[52] + Weekly Temp, 8-week horizon", color="red", marker="o", linestyle='-')

# Titles and labels
plt.title("SARIMAX(1,0,0)(1,0,0)[52] + Weekly Temp | 8-Week Horizon")
plt.xlabel("Date")
plt.ylabel("Energy Consumption (kWh)")
plt.legend()
plt.grid(True, linestyle='--', alpha=0.5)

# --- Show plot ---
plt.tight_layout()
plt.show()

# --- Step 6: MAPE ---
if len(actual) == len(forecast):
    mape_score = mape(actual, forecast)
    print(f"MAPE (SARIMAX, Rolling 8-week Forecasts): {mape_score:.2f}%")
    with open("/Users/omarabdesslem/Documents/Projects/Thesis-Website/public/SARIMAX/sarimax_mape.txt", "w") as f:
        f.write(f"{mape_score:.2f}")
else:
    print("Forecast and actual series do not match in length.")

MAPE (SARIMAX, Rolling 8-week Forecasts): 4.83%