Retail Demand Forecasting#
Retailers often deal with large datasets that exceed the memory capacity of a single machine. These datasets also exhibit diverse characteristics, particularly when considering the varying behaviors of different stores. Consequently, data scientists face several challenges, including:
Limitations of Pandas/Python in handling big datasets, necessitating the use of Pyspark.
The need to create separate models for different groups in the dataset (e.g., stores or items) to optimize performance.
Despite the dataset’s overall size, each group (items/stores) is small enough to leverage models similar to those provided by
scikit-learn.
This tutorial demonstrates how data scientists can overcome these challenges using ForecastFlowML.
Situation#
Given the considerable size of our sample dataset, it exceeds the memory capacity for direct loading. Consequently, we must leverage PySpark, a distributed processing framework, to effectively manage the data processing operations.
Fortunately, within this extensive dataset, we possess store-specific information.
The individual store data, being relatively smaller in size, can comfortably fit within the memory of a single machine.
Consequently, we can leverage this advantage to employ scikit-learn for model creation, without encountering any memory limitations.
Goal#
Develop independent models for each store in the dataset.
Parallelize the feature engineering, training, and inference steps.
Employ
LightGBMas the machine learning algorithm.Utilize a direct multi-step forecasting approach.
Conduct backtesting to assess the model’s performance.
Import Packages#
from forecastflowml import ForecastFlowML
from forecastflowml import FeatureExtractor
from forecastflowml.data.loader import load_walmart_m5
from lightgbm import LGBMRegressor
from pyspark.sql import SparkSession
import pyspark.sql.functions as F
import plotly.express as px
import plotly.io as pio
import pandas as pd
import sys
import os
os.environ["PYSPARK_PYTHON"] = sys.executable
pd.set_option("display.max_columns", 40)
Initialize Spark#
spark = (
SparkSession.builder.master("local[4]")
.config("spark.driver.memory", "4g")
.config("spark.sql.shuffle.partitions", "4")
.config("spark.sql.execution.pyarrow.enabled", "true")
.getOrCreate()
)
Sample Dataset#
df = load_walmart_m5(spark).localCheckpoint()
df.show(10)
+--------------------+-----------+-------+------+--------+--------+----------+-----+
| id| item_id|dept_id|cat_id|store_id|state_id| date|sales|
+--------------------+-----------+-------+------+--------+--------+----------+-----+
|FOODS_1_002_TX_1_...|FOODS_1_002|FOODS_1| FOODS| TX_1| TX|2015-01-15| 3.0|
|FOODS_1_002_TX_1_...|FOODS_1_002|FOODS_1| FOODS| TX_1| TX|2015-01-16| 0.0|
|FOODS_1_002_TX_1_...|FOODS_1_002|FOODS_1| FOODS| TX_1| TX|2015-01-17| 1.0|
|FOODS_1_002_TX_1_...|FOODS_1_002|FOODS_1| FOODS| TX_1| TX|2015-01-18| 0.0|
|FOODS_1_002_TX_1_...|FOODS_1_002|FOODS_1| FOODS| TX_1| TX|2015-01-19| 0.0|
|FOODS_1_002_TX_1_...|FOODS_1_002|FOODS_1| FOODS| TX_1| TX|2015-01-20| 0.0|
|FOODS_1_002_TX_1_...|FOODS_1_002|FOODS_1| FOODS| TX_1| TX|2015-01-21| 0.0|
|FOODS_1_002_TX_1_...|FOODS_1_002|FOODS_1| FOODS| TX_1| TX|2015-01-22| 0.0|
|FOODS_1_002_TX_1_...|FOODS_1_002|FOODS_1| FOODS| TX_1| TX|2015-01-23| 0.0|
|FOODS_1_002_TX_1_...|FOODS_1_002|FOODS_1| FOODS| TX_1| TX|2015-01-24| 0.0|
+--------------------+-----------+-------+------+--------+--------+----------+-----+
only showing top 10 rows
df.select("store_id").dropDuplicates().show()
+--------+
|store_id|
+--------+
| TX_1|
| WI_1|
| CA_1|
| TX_2|
+--------+
Approach#
There are 4 different stores in our dataset. Since each store might have different behaviours, it would be better to build an independent model for each of them. Moreover, there will be 4 models per store that is going to predict 7 days. In total, we are going to have 4 x 4 = 16 models.
As we are going to build 4 weekly independent direct forecasting models, each of them are going to use different lag features based upon their forecasting horizon.
Feature Engineering#
Let’s create the features that will be used for forecasting.
Our package has a seperate module FeatureExtractor to extract features.
Spesifically, we are going to extract:
Lags
Roling means
Date features
Stockout features (counting consecutive no sales periods)
History length (number of periods passed after the start of the time series)
feature_extractor = FeatureExtractor(
id_col="id",
date_col="date",
target_col="sales",
lag_window_features={
"lag": [7 * (i + 1) for i in range(7)],
"mean": [[window, lag] for lag in [7, 14, 21, 28] for window in [7, 14, 30]],
},
date_features=[
"day_of_month",
"day_of_week",
"week_of_year",
"quarter",
"month",
"year",
],
count_consecutive_values={
"value": 0,
"lags": [7, 14, 21, 28],
},
history_length=True,
)
df_features = feature_extractor.transform(df).localCheckpoint()
df_features.show(10)
+--------------------+-----------+-------+------+--------+--------+----------+-----+-----+------+------+------+------+------+------+-------------------+--------------------+--------------------+--------------------+---------------------+---------------------+--------------------+---------------------+---------------------+--------------------+---------------------+---------------------+-----------------------------+------------------------------+------------------------------+------------------------------+--------------+------------+-----------+------------+-------+-----+----+
| id| item_id|dept_id|cat_id|store_id|state_id| date|sales|lag_7|lag_14|lag_21|lag_28|lag_35|lag_42|lag_49|window_7_lag_7_mean|window_14_lag_7_mean|window_30_lag_7_mean|window_7_lag_14_mean|window_14_lag_14_mean|window_30_lag_14_mean|window_7_lag_21_mean|window_14_lag_21_mean|window_30_lag_21_mean|window_7_lag_28_mean|window_14_lag_28_mean|window_30_lag_28_mean|count_consecutive_value_lag_7|count_consecutive_value_lag_14|count_consecutive_value_lag_21|count_consecutive_value_lag_28|history_length|day_of_month|day_of_week|week_of_year|quarter|month|year|
+--------------------+-----------+-------+------+--------+--------+----------+-----+-----+------+------+------+------+------+------+-------------------+--------------------+--------------------+--------------------+---------------------+---------------------+--------------------+---------------------+---------------------+--------------------+---------------------+---------------------+-----------------------------+------------------------------+------------------------------+------------------------------+--------------+------------+-----------+------------+-------+-----+----+
|FOODS_1_002_TX_1_...|FOODS_1_002|FOODS_1| FOODS| TX_1| TX|2015-01-15| 3.0| null| null| null| null| null| null| null| null| null| null| null| null| null| null| null| null| null| null| null| null| null| null| null| 1| 15| 4| 3| 1| 1|2015|
|FOODS_1_002_TX_1_...|FOODS_1_002|FOODS_1| FOODS| TX_1| TX|2015-01-16| 0.0| null| null| null| null| null| null| null| null| null| null| null| null| null| null| null| null| null| null| null| null| null| null| null| 2| 16| 5| 3| 1| 1|2015|
|FOODS_1_002_TX_1_...|FOODS_1_002|FOODS_1| FOODS| TX_1| TX|2015-01-17| 1.0| null| null| null| null| null| null| null| null| null| null| null| null| null| null| null| null| null| null| null| null| null| null| null| 3| 17| 6| 3| 1| 1|2015|
|FOODS_1_002_TX_1_...|FOODS_1_002|FOODS_1| FOODS| TX_1| TX|2015-01-18| 0.0| null| null| null| null| null| null| null| null| null| null| null| null| null| null| null| null| null| null| null| null| null| null| null| 4| 18| 7| 3| 1| 1|2015|
|FOODS_1_002_TX_1_...|FOODS_1_002|FOODS_1| FOODS| TX_1| TX|2015-01-19| 0.0| null| null| null| null| null| null| null| null| null| null| null| null| null| null| null| null| null| null| null| null| null| null| null| 5| 19| 1| 4| 1| 1|2015|
|FOODS_1_002_TX_1_...|FOODS_1_002|FOODS_1| FOODS| TX_1| TX|2015-01-20| 0.0| null| null| null| null| null| null| null| null| null| null| null| null| null| null| null| null| null| null| null| null| null| null| null| 6| 20| 2| 4| 1| 1|2015|
|FOODS_1_002_TX_1_...|FOODS_1_002|FOODS_1| FOODS| TX_1| TX|2015-01-21| 0.0| null| null| null| null| null| null| null| null| null| null| null| null| null| null| null| null| null| null| null| null| null| null| null| 7| 21| 3| 4| 1| 1|2015|
|FOODS_1_002_TX_1_...|FOODS_1_002|FOODS_1| FOODS| TX_1| TX|2015-01-22| 0.0| 3.0| null| null| null| null| null| null| 3.0| 3.0| 3.0| null| null| null| null| null| null| null| null| null| 0| null| null| null| 8| 22| 4| 4| 1| 1|2015|
|FOODS_1_002_TX_1_...|FOODS_1_002|FOODS_1| FOODS| TX_1| TX|2015-01-23| 0.0| 0.0| null| null| null| null| null| null| 1.5| 1.5| 1.5| null| null| null| null| null| null| null| null| null| 1| null| null| null| 9| 23| 5| 4| 1| 1|2015|
|FOODS_1_002_TX_1_...|FOODS_1_002|FOODS_1| FOODS| TX_1| TX|2015-01-24| 0.0| 1.0| null| null| null| null| null| null| 1.3333333333333333| 1.3333333333333333| 1.3333333333333333| null| null| null| null| null| null| null| null| null| 0| null| null| null| 10| 24| 6| 4| 1| 1|2015|
+--------------------+-----------+-------+------+--------+--------+----------+-----+-----+------+------+------+------+------+------+-------------------+--------------------+--------------------+--------------------+---------------------+---------------------+--------------------+---------------------+---------------------+--------------------+---------------------+---------------------+-----------------------------+------------------------------+------------------------------+------------------------------+--------------+------------+-----------+------------+-------+-----+----+
only showing top 10 rows
Training/Test Dataset#
df_train = df_features.filter(F.col("date") < "2016-04-25")
df_test = df_features.filter(F.col("date") >= "2016-04-25")
Initialize ForecastFlowML#
forecast_flow = ForecastFlowML(
group_col="store_id", # column to slice the dataframe
id_col="id", # time series identifier column
date_col="date", # date column
target_col="sales", # target column
categorical_cols=[
"item_id",
"dept_id",
"cat_id",
], # columns to be treated as categorical
date_frequency="days", # dataset date frequency
model_horizon=7, # single model predicts 7 days
max_forecast_horizon=28, # total forecast horizon will be 28 days
model=LGBMRegressor(), # algorithm
use_lag_range=21, # allow to use 21 days of extra lag features to use based upon the forecast horizon
)
Cross Validation#
First, let’s perform some backtesting to understand the model performance. We are going to predict the last 3x28 days.
cv_forecast = forecast_flow.cross_validate(df_train, n_cv_splits=3).localCheckpoint()
cv_forecast.show(10)
+--------+--------------------+-------------------+---+-----+----------+
|store_id| id| date| cv|sales|prediction|
+--------+--------------------+-------------------+---+-----+----------+
| CA_1|FOODS_1_064_CA_1_...|2016-03-28 00:00:00| 0| 2.0| 0.979212|
| CA_1|FOODS_1_064_CA_1_...|2016-03-29 00:00:00| 0| 0.0| 1.047827|
| CA_1|FOODS_1_064_CA_1_...|2016-03-30 00:00:00| 0| 0.0|0.92554504|
| CA_1|FOODS_1_064_CA_1_...|2016-03-31 00:00:00| 0| 0.0| 0.9589823|
| CA_1|FOODS_1_064_CA_1_...|2016-04-01 00:00:00| 0| 3.0| 1.305723|
| CA_1|FOODS_1_064_CA_1_...|2016-04-02 00:00:00| 0| 0.0| 1.5510178|
| CA_1|FOODS_1_064_CA_1_...|2016-04-03 00:00:00| 0| 3.0| 1.4713657|
| CA_1|FOODS_1_121_CA_1_...|2016-03-28 00:00:00| 0| 0.0| 0.560704|
| CA_1|FOODS_1_121_CA_1_...|2016-03-29 00:00:00| 0| 0.0| 0.5866643|
| CA_1|FOODS_1_121_CA_1_...|2016-03-30 00:00:00| 0| 0.0| 0.5995418|
+--------+--------------------+-------------------+---+-----+----------+
only showing top 10 rows
Calculate RMSE#
(
cv_forecast.withColumn(
"error_squared", F.pow(F.abs(F.col("sales") - F.col("prediction")), 2)
)
.groupBy("cv")
.agg(F.pow(F.avg("error_squared"), 0.5).alias("rmse"))
.orderBy("cv")
.show()
)
+---+------------------+
| cv| rmse|
+---+------------------+
| 0|1.9499041165007749|
| 1|2.0947956593876733|
| 2| 1.925051121492435|
+---+------------------+
Visualize Cross Validation Predictions#
cv_past_future = (
df_train.select("id", "store_id", "date", "sales")
.join(
cv_forecast.select("id", "date", "cv", "prediction"),
on=["id", "date"],
how="left",
)
.groupBy("id", "store_id", "date", "sales")
.pivot("cv")
.sum("prediction")
.groupBy("store_id", "date")
.agg(
F.sum("sales").alias("sales"),
*[F.sum(f"{i}").alias(f"cv_{i}") for i in range(3)],
)
.orderBy("store_id", "date")
).toPandas()
pio.renderers.default = "notebook"
fig = px.line(
cv_past_future,
x="date",
y=["sales", *[f"cv_{i}" for i in range(3)]],
facet_row_spacing=0.04,
facet_col="store_id",
facet_col_wrap=2,
height=700,
width=720,
)
fig.update_layout(
legend=dict(orientation="h", yanchor="top", y=1.09, xanchor="center", x=0.5),
margin=dict(l=0, r=10, t=5, b=5),
legend_title="",
)
fig.update_traces(line=dict(width=1.7))
fig.update_yaxes(matches=None, title="")
fig.update_xaxes(type="date", range=["2015-11-01", "2016-04-24"])
Training#
Here we know that the trained models will be small enough to fit into the memory, therefore; we use the local_result=True which will keep the trained models as an attribute. However, if the trained models are expected to not fit into the memory, then local_result=False (default parameter) will return models in PySpark DataFrame.
forecast_flow.train(df_train, local_result=True)
forecast_flow.model_
| store_id | forecast_horizon | model | start_time | end_time | elapsed_seconds | |
|---|---|---|---|---|---|---|
| 0 | CA_1 | [[1, 2, 3, 4, 5, 6, 7], [8, 9, 10, 11, 12, 13,... | [[128, 4, 149, 236, 1, 0, 0, 0, 0, 0, 0, 140, ... | 19-May-2023 (23:41:17) | 19-May-2023 (23:41:18) | 0.6 |
| 1 | TX_1 | [[1, 2, 3, 4, 5, 6, 7], [8, 9, 10, 11, 12, 13,... | [[128, 4, 149, 236, 1, 0, 0, 0, 0, 0, 0, 140, ... | 19-May-2023 (23:41:18) | 19-May-2023 (23:41:19) | 0.6 |
| 2 | WI_1 | [[1, 2, 3, 4, 5, 6, 7], [8, 9, 10, 11, 12, 13,... | [[128, 4, 149, 236, 1, 0, 0, 0, 0, 0, 0, 140, ... | 19-May-2023 (23:41:19) | 19-May-2023 (23:41:19) | 0.5 |
| 3 | TX_2 | [[1, 2, 3, 4, 5, 6, 7], [8, 9, 10, 11, 12, 13,... | [[128, 4, 149, 236, 1, 0, 0, 0, 0, 0, 0, 140, ... | 19-May-2023 (23:41:16) | 19-May-2023 (23:41:17) | 0.7 |
Prediction#
forecast = forecast_flow.predict(df_test, spark=spark).localCheckpoint()
forecast.show(10)
+--------+--------------------+-------------------+----------+
|store_id| id| date|prediction|
+--------+--------------------+-------------------+----------+
| CA_1|FOODS_1_064_CA_1_...|2016-04-25 00:00:00| 1.1891526|
| CA_1|FOODS_1_064_CA_1_...|2016-04-26 00:00:00| 1.1587629|
| CA_1|FOODS_1_064_CA_1_...|2016-04-27 00:00:00| 1.1314319|
| CA_1|FOODS_1_064_CA_1_...|2016-04-28 00:00:00| 1.1495334|
| CA_1|FOODS_1_064_CA_1_...|2016-04-29 00:00:00| 1.311355|
| CA_1|FOODS_1_064_CA_1_...|2016-04-30 00:00:00| 1.2894475|
| CA_1|FOODS_1_064_CA_1_...|2016-05-01 00:00:00| 1.375367|
| CA_1|FOODS_1_121_CA_1_...|2016-04-25 00:00:00|0.63143337|
| CA_1|FOODS_1_121_CA_1_...|2016-04-26 00:00:00| 0.5990241|
| CA_1|FOODS_1_121_CA_1_...|2016-04-27 00:00:00| 0.6636088|
+--------+--------------------+-------------------+----------+
only showing top 10 rows
Visualize Predictions#
past_future = (
df.select("id", "store_id", "date", "sales")
.join(forecast.select("id", "date", "prediction"), on=["id", "date"], how="left")
.groupBy("store_id", "date")
.agg(
F.sum("sales").alias("sales"),
F.sum("prediction").alias("prediction"),
)
.orderBy("store_id", "date")
.toPandas()
)
pio.renderers.default = "notebook"
fig = px.line(
past_future,
x="date",
y=["sales", "prediction"],
facet_row_spacing=0.04,
facet_col="store_id",
facet_col_wrap=2,
height=700,
width=720,
)
fig.update_layout(
legend=dict(orientation="h", yanchor="top", y=1.09, xanchor="center", x=0.5),
margin=dict(l=0, r=10, t=5, b=5),
legend_title="",
)
fig.update_traces(line=dict(width=1.7))
fig.update_yaxes(matches=None, title="")
fig.update_xaxes(type="date", range=["2015-11-01", "2016-05-22"])