Forecasting energy demand using Skforecast

Tour of the project

Let’s start with a diagram summarizing the interactions of the various modules of the codebase. I will go through most of them, describe how I implemented them and how they fit into the bigger picture.

CLI Interface

All the functionalities of the program are accessible through a command-line interface, built using the argparse module, in main.py. This makes it very easy to then automate behaviour in a CI/CD workflow. To use the various functions of the interface, you can simply call python main.py command,³ where command can be one of {download, merge, train, predict}.

You can provide an optional starting and end date arguments to the download command if you’d like to download a specific period, otherwise it will default to the period since the last download. The merge command merges all the fragmented timeseries into a single .csv file. It’s automatically applied during the download operation, but if you modify any of the raw files manually (or even delete some), it can be useful to call it manually.

When you call the train command, the program will create a new model and tune its hyperparameters using a portion of the dataset—this can take a couple of minutes so feel free to go grab a coffee. The model is only retrained if more than a week has passed since the end of the last training.

Finally, once the model is fully tuned and trained, the predict command will infer the timeseries value on the future (or test) part of the dataset. You can also provide the optional argument --plot to this command if you want to create a plot of the result and have it displayed on the index.html page.

This entire workflow is scheduled through Github actions, every Monday-Wednesday-Friday, using the following .yaml workflow configuration:

name: Download and Train Model
on:
    schedule:
        - cron: "0 0 * * MON,WED,FRI" # Scheduled three times a week
    workflow_dispatch: # Allows for manual triggering
jobs:
    download-tune:
        runs-on: ubuntu-latest
        steps:
            - name: Checkout repo content
              uses: actions/checkout@v4
              with:
                token: ${{ secrets.PERSONAL_ACCESS_TOKEN }}  # Use the PAT instead of the default GITHUB_TOKEN
            - name: Install uv
              uses: astral-sh/setup-uv@v5
              with:
                # Install a specific version of uv.
                version: "0.5.1"
                enable-cache: true
                cache-dependency-glob: "uv.lock"
            - name: Setup python
              uses: actions/setup-python@v5
              with:
                python-version-file: "pyproject.toml"
                cache: 'pip'
            - name: Install the project
              run: uv sync --all-extras --dev
            - name: Make envfile
              uses: SpicyPizza/create-envfile@v2.0
              with:
                envkey_API_KEY: ${{ secrets.API_KEY }} # The API key is stored in the .env file
                file_name: .env
                fail_on_empty: false
                sort_keys: false
            - name: Run download # Downlaod new data
              run: uv run main.py download
            - name: Run training # Train model
              run: uv run main.py train
            - name: Make demo page # Build demo page with prediction
              run: uv run main.py predict --plot
            - name: Check for changes # create env variable indicating if any changes were made
              id: git-check
              run: |
                git config user.name 'github-actions'
                git config user.email 'github-actions@github.com'
                git add .
                git diff --staged --quiet || echo "changes=true" >> $GITHUB_ENV 
            - name: Commit and push if changes
              if: env.changes == 'true' # if changes made push new data to repo
              run: |
                git pull
                git commit -m "downloaded new data and trained model"
                git push

Data

The module data.py is responsible for fetching and processing the data. I initially used a custom function to query the API and parse its response. But in the few months I let the pipeline run by itself, I quickly ran into issues where the API would return unexpected responses, breaking my parser and requiring a lot of error handling. After spending a bit too much time debugging, I decided to switch to the entsoe-py Python library. It was a game-changer! This library is specifically designed to handle Entso-e’s API, so it’s better equipped to handle changes in the response or fetching errors, and in the long run, it will save me a lot of time.

The downloaded data sits inside the data/raw folder as .csv files. Ideally, a more complex project would involve some form of data versioning while the data would be stored on some external storage (like AWS or GCP). In this particular case, the timeseries consist of two float64 columns and a datetime index, with about 9k rows per year. This amounts to less than 400kB of data and doesn’t require anything too complicated.

Models

To build my models and handle predictions, I have chosen to use the Skforecast library, a scikit-learn compatible library which abstracts away a lot of the specificities of timeseries while allowing for ample configuration and customization of the models. In this project, I have made use of a lot of these features so I hope this little tour of the project can also serve as a short overview of the library. For more details about both Skforecast and generally how to work with timeseries, I whole-heartedly recommend their documentation which features an excellent introduction to forecasting, in-depth user guides, and examples of how to use the library for common use-cases. I learnt a lot from going through this documentation and I can only recommend it to you too!

Now let’s dive into the code. In the module models.py, I have created a ForecasterRecursiveModel class which wraps around the ForecasterRecursive class from Skforecast. The forecast classes from Skforecast take care of creating lag and window features during predictions, and on top of this, the ForecasterRecursive class specializes in making multi-step predictions, using each prediction as input for the next one.⁴ If you choose a regressor (say LGBM), you can then define your forecast as:

from lightgbm import LGBMRegressor
from skforecast.recursive import ForecasterRecursive

forecaster = ForecasterRecursive(
    regressor=LGBMRegressor(random_state=1234, n_jobs=-1, verbose=-1),
    lags=12,
)

The goal of the ForecasterRecursiveModel class is to automatically handle hyperparameter tuning, predictions, and backtesting of the ForecasterRecursive model.

from scipy.signal import periodogram

frequencies, spectrum = periodogram(
    y, # Timeseries
    window="boxcar",
    scaling="spectrum",
)

from skforecast.model_selection import TimeSeriesFold

cv = TimeSeriesFold(
    steps=24, # Number of steps to predict
    refit=False, # Whether to refit between each fold. If an int is given refit intermittently
    initial_train_size=24*365*3, # Train over 3 years of data
    allow_incomplete_fold=True, # If the last fold is not complete, still use it
)

from skforecast.model_selection import bayesian_search_forecaster
from optuna.trial import Trial

def search_space_lgbm(trial: Trial) -> dict:
  """Function describing the search space for hyperparameters of an LGBM model."""
  search_space = {
      "num_leaves": trial.suggest_int("num_leaves", 8, 256),
      "max_depth": trial.suggest_int("max_depth", 3, 16),
      "learning_rate": trial.suggest_float("learning_rate", 0.001, 0.2, log=True),
      "n_estimators": trial.suggest_int("n_estimators", 50, 1000, log=True),
      "bagging_fraction": trial.suggest_float("bagging_fraction", 0.5, 1),
      "feature_fraction": trial.suggest_float("feature_fraction", 0.5, 1),
      "reg_alpha": trial.suggest_float("reg_alpha", 0.01, 100),
      "reg_lambda": trial.suggest_float("reg_lambda", 0.01, 100),
      "lags": trial.suggest_categorical("lags", c.lags_consider),
  }
  return search_space

# The Bayesian search forecaster uses Optuna
results, _ = bayesian_search_forecaster(
    forecaster=forecaster, # The forecaster object
    y=y, # Timeseries data of shape (N,1)
    cv=cv, # The cross-validation folds
    search_space=search_space_lgbm, # Range of parameters to consider
    metric="mean_absolute_error",
    exog=X, # Exogeneous features of shape (N, n_features)
    return_best=False, # If True, the model will be retrained on the best parameters found
    random_state=c.random_state,
    verbose=False,
    n_trials=10, # Number of trials to find the right parameters
)

print(f"Best parameters:\n {results.iloc[0]}")

Feature importance: Can we say anything more about this model?

An important part of any ML project is explainability: it’s all good and well to make predictive models, but it’s important to be able to explain how your model works and why it’s making some predictions. Thankfully, the shap library comes to the rescue! Based on game theory,⁶ this library brings a simple way to evaluate which feature matters the most to make predictions. In the case of an LGBM model trained over the last 3 years, these are the top features in terms of averaged shap values:

Let’s take a closer look at this plot and figure out why each of these features contribute to the model in this way:

• the first lag: this is usually an important feature of timeseries as the next time step tends to be strongly correlated to the current time step.
• The daily, weekly and yearly cyclical features. As I’ve mentioned before, energy demand tends to trend strongly with seasons and with daily usage so it’s not surprising to see these features matter a lot here.
• Rolling window features: these features tend to average seasonal effects and instead encode for evolving trends.
• Other lag features: most timeseries tend to depend on higher-order lag features. Putting my physicist cap on, a timeseries which would only depend on the first lag would most likely be some exponentially damped system,⁷ while a second-order lag dependency would be reminiscent of periodic systems, like harmonic oscillators. More complicated systems will either involve more complex couplings or higher-order dependencies.

Since I am using LGBM as a regressor, I can also peek into the output of its native get_feature_importance() method and see the top 10 most important features:

While the exact number itself isn’t that important, the order and relative magnitude gives an idea as to which features matter most. It’s a good sign that this list agrees with the SHAP summary plot and tells the same story.

Revisiting the project: key lessons

One of the most interesting takeaways from this project was how tough it can be to improve upon an existing project. There’s a fine line between over-engineering a codebase by trying to foresee future refactoring that won’t happen, and accumulating technical debt by not writing maintainable and extendable code. I personally tend to fall on the side of the former, creating more and more abstraction than is generally needed. Surprisingly when I came back to this codebase, I realized I was instead fighting with most of the decisions I had made early on.

This experience also reaffirmed something I had learnt before: it’s easy to think that improving a project means adding new features, but it’s just as important to take a step back and simplify things. After refactoring the code, I realized that simplifying the structure and focusing on the essentials helped me achieve better results with less complexity.