Retail Demand Forecasting and Statistical Modelling
Dunnhumby The Complete Journey · Independent Forecasting Study
This independent study investigates grocery purchasing histories through data engineering, weekly forecasting and comparative statistical evaluation. Ten forecasting methods are assessed using chronological training, validation and test periods, with future work examining how forecasts inform inventory replenishment decisions.
Work Demonstrated
1. Data Engineering & Analysis
A Python-based data preparation workflow transforms grocery transaction records into daily and weekly SKU-level purchase series, supported by data validation and exploratory analysis of purchasing patterns.
Overview
Transaction records are processed into daily and weekly purchase histories for 12 selected grocery products. Exploratory diagnostics examine trends, autocorrelation and potential seasonal patterns.
Product metadata, data preparation and research-sample provenance.
| Product ID | Commodity | Subcategory | Pack size |
|---|---|---|---|
| Loading product catalogue… | |||
Daily and weekly purchase quantities with chronological training, validation and test partitions.
Investigate underlying trends, autocorrelation and potential recurring patterns using training observations.
2. Forecasting Methodology
Developed a comparative forecasting framework encompassing ten statistical and time-series methods, including naïve forecasting, regression-based approaches, exponential smoothing, and Auto-ARIMA.
Overview
Ten statistical forecasting methods are formulated and implemented, covering naïve approaches, regression, exponential smoothing and ARIMA. Each method produces forecasts across horizons of 1–4 weeks.
Common mathematical notation Shared across the forecasting methods
Equations are written for one SKU at a time; the SKU index is omitted for clarity.
- yt
- Observed purchase quantity in week \(t\).
- ŷt+h|t
- Forecast for week \(t\)+\(h\), issued using observations available through week \(t\).
- t
- Forecast origin: the latest observed week.
- h
- Forecast horizon, \(h\) ∈ {1, 2, 3, 4} weeks.
Direct forecasting through observation repetition, historical averaging, and lag-based selection.
Last Observed
Model description
Uses the most recent observed weekly purchase quantity as the forecast for every horizon. As the rolling origin advances, the forecast is updated with the latest observation.
Mathematical formulation
Model configuration
No parameters require tuning or estimation.
Moving Average
Model description
Calculates the average weekly purchase quantity over the most recent \(k\) weeks and uses it to predict future purchases. As new weekly observations become available, the average is recalculated using the latest \(k\) weeks.
Mathematical formulation
where
- \(k\)
- Number of recent historical weeks included in the average.
- \(i\)
- Historical observation index, running from 0 to \(k\)−1.
- \(y_{t-i}\)
- Observed purchases i weeks before the forecast origin.
Model configuration
Averaging window (\(k\)): Select an integer from 2–72 weeks for each SKU by minimising pooled rolling-validation WAPE across horizons 1–4.
Lag Forecast
Model description
Uses the purchase quantity observed \(k\) weeks earlier to predict future weekly purchases. When the required historical observation is not yet available, recursive forecasting is applied by reusing an earlier prediction.
Mathematical formulation
where
- \(k\)
- Number of weeks separating a target from the observation or prediction copied.
- \(y_{t+h-k}\)
- Observed purchase quantity used when the lag points to known history.
Model configuration
Historical lag (\(k\)): Select an integer from 2–72 weeks for each SKU by minimising pooled rolling-validation WAPE across horizons 1–4.
Trend estimation using linear regression, with optional seasonal adjustment based on historical patterns.
Linear Trend
Model description
Uses linear regression to estimate whether weekly purchase quantities are increasing or decreasing over time. The estimated trend is extended into future weeks to generate forecasts.
Mathematical formulation
where
- \(\hat a_t\)
- Estimated intercept using the current rolling window.
- \(\hat b_t\)
- Estimated slope (change in purchases per week) using the current rolling window.
- \(n_t\)
- Number of observations in the regression window (up to 26); time indices start at zero.
Model configuration
Regression coefficients (\(\hat a_t,\hat b_t\)): Estimate the intercept and slope using ordinary least squares over the latest 26 observed weeks.
Trend + Seasonal Index
Model description
Combines linear trend estimation with seasonal indices to account for both changes over time and recurring purchasing patterns. Future forecasts are generated by extending the estimated trend and adjusting it according to the corresponding seasonal index.
Mathematical formulation
where
- \(\hat a_t,\hat b_t\)
- Intercept and slope estimated from historical weekly purchases at forecast origin \(t\).
- \(\hat S_r\)
- Seasonal index for position \(r\) within the repeating cycle.
- \(r_{t,h}\)
- Seasonal position for the target week, defined by \(r_{t,h}=n_t-1+h-m\left\lfloor\frac{n_t-1+h}{m}\right\rfloor\), taking values from \(0\) to \(m-1\).
- \(m\)
- Cycle length, in weeks.
- \(n_t\)
- Number of observed weeks in the current history; the trend uses zero-based time indices.
Model configuration
Seasonal period (\(m\)): Select an integer from 2–24 weeks for each SKU by minimising pooled rolling-validation WAPE across horizons 1–4.
Trend coefficients (\(\hat a_t,\hat b_t\)) and seasonal indices (\(\hat S_r\)): Estimate from the available historical observations.
Recursive estimation of the underlying level, with extensions for trend, damping, and seasonality.
Simple Exponential Smoothing
Model description
Uses exponential smoothing to estimate the underlying level of weekly purchases, giving greater weight to more recent observations. The estimated level is updated as new observations become available and used to forecast future purchases.
Mathematical formulation
where
- \(\ell_t\)
- Smoothed level after observing week \(t\).
Show level update equation
Model configuration
Level smoothing coefficient (\(\alpha\)): Estimate statistically during model fitting.
Holt Linear Trend
Model description
Extends exponential smoothing by estimating both the underlying purchase level and its trend over time. The estimated trend is projected forward to generate forecasts for future weeks.
Mathematical formulation
where
- \(\ell_t\)
- Estimated level after week \(t\).
- \(b_t\)
- Estimated change per week (trend).
Show level and trend update equations
Model configuration
Level and trend smoothing coefficients (\(\alpha,\beta\)): Estimate statistically during model fitting.
Holt Damped Trend
Model description
Extends Holt’s linear method by introducing a damping parameter that gradually reduces the influence of the estimated trend. This produces forecasts in which the projected trend becomes progressively weaker at longer horizons.
Mathematical formulation
where
- \(\ell_t\)
- Estimated level after week \(t\).
- \(b_t\)
- Estimated damped-trend state.
- \(\phi\)
- Damping factor applied to future trend contributions.
Show level and damped-trend updates
Model configuration
Level and trend smoothing coefficients (\(\alpha,\beta\)) and damping factor (\(\phi\)): Estimate statistically during model fitting.
Holt-Winters
Model description
Extends exponential smoothing by incorporating recurring seasonal patterns alongside the underlying purchase level and an optional trend. The level, trend and seasonal components are updated as new observations become available to generate future forecasts.
Mathematical formulation
where
- \(\ell_t\)
- Estimated level.
- \(b_t\)
- Estimated trend state, when enabled.
- \(s_t\)
- Additive seasonal state.
- \(m\)
- Seasonal period, in weeks.
Show seasonal and trend extensions
This extends the seasonal index across cycles; the trend contribution also covers trend-free and damped specifications.
Model configuration
Seasonal period (\(m\)): Select an integer from 2–24 weeks using pooled rolling-validation WAPE.
Trend specification: Select no trend, linear trend, or damped linear trend.
Smoothing coefficients (\(\alpha,\beta,\gamma\)): Estimate statistically during model fitting, where applicable.
Autoregressive and moving-average structures, optional differencing, and automatic model selection.
Auto-ARIMA
Model description
Uses autoregressive (AR) and moving-average (MA) components, with differencing when required, to model patterns in weekly purchases. The selected model uses relationships between past observations and forecast errors to predict future purchase quantities.
Mathematical formulation
where
- \(p,d,q\)
- Candidate autoregressive, differencing and moving-average orders.
- \(\mathcal C\)
- Set of candidate structures considered by the stepwise procedure.
- \(\operatorname{AICc}\)
- Akaike Information Criterion with small-sample correction.
- \(p^*,d^*,q^*\)
- Orders selected for the SKU.
Model configuration
Model orders (\(p,d,q\)): Select automatically through a stepwise search using the corrected Akaike Information Criterion (AICc).
Model coefficients: Estimate statistically for the selected ARIMA specification.
AR
Model description
Models the current purchase value using its own previous values, with dependence governed by autoregressive coefficients.
Mathematical formulation
where
- \(p\)
- Number of autoregressive lags.
- \(\varphi_i\)
- Coefficient for autoregressive lag i.
- \(c\)
- Constant term.
- \(\varepsilon_t\)
- Model innovation at week \(t\).
- \(i\)
- Autoregressive lag index.
Relationship to Auto-ARIMA
AR is a special case of ARIMA with \(d=0\) and \(q=0\). Auto-ARIMA may select an autoregressive-only specification when it fits a SKU best.
MA
Model description
Models current purchase fluctuations using current and previous innovations (model errors), rather than averaging recent observed purchases.
Mathematical formulation
where
- \(q\)
- Number of innovation lags.
- \(\mu\)
- Mean of a stationary series.
- \(\theta_j\)
- Coefficient of innovation lag j.
- \(\varepsilon_t\)
- Model innovation.
- \(j\)
- Innovation-lag index.
Relationship to Auto-ARIMA
MA is a special case of ARIMA with \(p=0\) and \(d=0\). Here, moving average refers to past model errors, not the historical averaging method in Section 2.1.
ARMA
Model description
Combines autoregressive dependence on previous values with moving-average dependence on past innovations, for stationary time series.
Mathematical formulation
where
- \(p,q\)
- Autoregressive and innovation-lag orders.
- \(\varphi_i,\theta_j\)
- AR and MA coefficients.
- \(c\)
- Model constant.
- \(\varepsilon_t\)
- Innovation.
- \(i,j\)
- Lag indices.
Relationship to Auto-ARIMA
ARMA combines AR and MA components with \(d=0\). It is one of the non-differenced specifications that Auto-ARIMA can consider.
ARIMA
Model description
Extends ARMA by differencing the observed series to address certain non-stationary patterns before modelling autoregressive and moving-average dynamics.
Mathematical formulation
where
- \(p,d,q\)
- AR order, differencing order and MA order.
- \(B\)
- Backshift operator: Byt = yt−1.
- \(\varphi(B)\)
- Autoregressive lag polynomial of degree \(p\).
- \(\theta(B)\)
- Moving-average lag polynomial of degree \(q\).
- \(\varepsilon_t\)
- Model innovation.
- \(c\)
- Constant term, if included.
Show ARIMA component equations
Relationship to Auto-ARIMA
ARIMA includes AR, MA and ARMA as special cases and allows differencing when needed. Auto-ARIMA searches these specifications rather than treating them as separate experimental methods.
3. Hyperparameter Tuning
Implemented SKU-specific hyperparameter tuning using rolling-origin validation to select forecasting model configurations based on predictive accuracy across multiple forecast horizons.
Overview
Forecasting methods are configured individually for each product through validation-based selection, statistical parameter estimation or information-criterion optimisation, as appropriate. The selected configurations are subsequently used for test evaluation.
Method-specific configuration searches, coefficient fitting and model selection.
Each SKU is configured independently. Candidate configurations are evaluated at rolling validation origins over forecast horizons 1–4; coefficients are fitted using historical observations available at each origin. Different methods use different selection procedures, as shown below.
Method-specific tuning
Select a method family to see what is searched, fitted or held fixed.
Last Observed
Fixed ruleUses the latest observed quantity. No configuration search or fitted coefficients.
Moving Average
Validation searchAveraging window (\(k\)): Select an integer from 2–72 weeks by minimising pooled rolling-validation WAPE across horizons 1–4.
Lag Forecast
Validation searchHistorical lag (\(k\)): Select an integer from 2–72 weeks by minimising pooled rolling-validation WAPE across horizons 1–4.
Linear Trend
Statistical fittingRegression coefficients (\(\hat a_t,\hat b_t\)): Estimate the intercept and slope using ordinary least squares over the latest 26 observed weeks. No discrete hyperparameter search.
Trend + Seasonal Index
Validation search + fittingSeasonal period (\(m\)): Select an integer from 2–24 weeks by minimising pooled rolling-validation WAPE across horizons 1–4.
Trend coefficients (\(\hat a_t,\hat b_t\)) and seasonal indices (\(\hat S_r\)): Estimate from historical observations.
Simple Exponential Smoothing
Statistical fittingLevel smoothing coefficient (\(\alpha\)): Estimate statistically during model fitting at each origin. No discrete validation search.
Holt Linear Trend
Statistical fittingLevel and trend smoothing coefficients (\(\alpha,\beta\)): Estimate statistically during model fitting. No discrete validation search.
Holt Damped Trend
Statistical fittingLevel and trend smoothing coefficients (\(\alpha,\beta\)) and damping factor (\(\phi\)): Estimate statistically during model fitting. No discrete validation search.
Holt-Winters
Validation search + fittingSeasonal period (\(m\)): Select an integer from 2–24 weeks using pooled rolling-validation WAPE.
Trend specification: Select no trend, linear trend or damped linear trend.
Smoothing coefficients (\(\alpha,\beta,\gamma\)): Estimate statistically during model fitting, where applicable. The seasonal component is additive.
Nonseasonal Auto-ARIMA
Training AICcModel orders (\(p,d,q\)): Select using stepwise training AICc on Weeks 2–73, subject to \(p\le 3\), \(d\le 2\), \(q\le 3\), and \(p+q\le 5\).
Model coefficients: Re-estimate statistically at rolling origins with the selected orders held fixed. Validation WAPE does not select the orders.
Explore selected settings by product and forecasting method.
Moving Average
4. Forecast Evaluation & Results
Evaluated forecasting performance across twelve grocery SKUs using out-of-sample predictions, multiple accuracy metrics, and horizon-specific comparisons of observed and forecast purchase quantities.
Overview
Forecast accuracy is evaluated across ten methods, 12 products and four prediction horizons. Interactive comparisons and evidence-based interpretations examine predictive performance, generalisation and product-specific differences.
Definitions and interpretation of the forecasting accuracy measures used in the evaluation.
Common Mathematical Notation
- \(y_i\)
- Observed purchase quantity for evaluated prediction \(i\).
- \(\hat y_i\)
- Forecast quantity for the corresponding target week.
- \(N\)
- Number of evaluated forecast–observation pairs.
Weighted Absolute Percentage Error (WAPE)
Metric description
Measures total absolute forecasting error relative to the total observed purchase quantity. A lower value indicates better aggregate accuracy.
Mathematical formulation
Interpretation
Expressed as a percentage; zero is ideal. The denominator must be non-zero.
Mean Absolute Error (MAE)
Metric description
Measures the average absolute difference between predicted and observed purchase quantities, in the same units as the original series.
Mathematical formulation
Interpretation
Lower values indicate smaller average errors; zero is ideal.
Root Mean Squared Error (RMSE)
Metric description
Squares errors before averaging, giving larger errors greater influence than in MAE. The resulting value remains in purchase-quantity units.
Mathematical formulation
Interpretation
Lower values are preferable; RMSE is particularly sensitive to large forecast misses.
Mean Forecast Bias (Mean Bias)
Metric description
Measures the average signed forecast error, retaining its direction rather than cancelling it through absolute values.
Mathematical formulation
Interpretation
Positive values indicate overforecasting; negative values indicate underforecasting. Values nearer zero indicate less average directional error.
Compare WAPE, MAE, RMSE and bias across test and validation horizons.
Method Performance
Compare observed weekly purchases with rolling-origin forecasts across methods and horizons.
Observed vs Forecast
What the forecasting results suggest about generalisation, model complexity and product-level variation.
Generalisation beyond validation
Finding. A method that performs best during configuration selection may not retain its advantage on unseen weeks. Lag Forecast ranks first on validation WAPE at every horizon, but drops to ninth or tenth on the test period.
Interpretation. The reversal is evidence of limited transfer between evaluation periods, not evidence of a particular cause. A lag that represents earlier purchase behaviour well may be less informative later; the independent test period exposes this risk.
Model complexity and forecasting behaviour
Finding. Extra trend or seasonal components do not consistently improve test accuracy. The comparison below places Simple Exponential Smoothing alongside selected trend and seasonal alternatives across all four horizons.
Interpretation. A more elaborate structure is helpful only when its assumptions are informative for future observations. Recent-level smoothing can remain competitive when additional trend or seasonal components do not translate into lower test error. This comparison does not establish that model complexity itself causes the difference.
Product-specific performance
Finding. No pooled ranking can reveal which method is strongest for every product. The heatmap compares methods within each SKU, so differences in purchase scale do not dominate the colours.
1 = lowest WAPE among the ten methods for that SKU; 10 = highest. Hover or focus on a cell to see the underlying WAPE. Ties share a rank.
Interpretation. Changing winners across products indicates that a single pooled recommendation may not serve every SKU. Method selection should consider the specific purchase history and forecast horizon, not just aggregate accuracy.
Understanding the magnitude of forecast errors
What a 30–40% WAPE means. Across the evaluated observations, total absolute forecast error is approximately 30–40% of total observed purchase quantity. It does not mean each weekly forecast misses by that percentage.
What may contribute. The measured series records purchases from a sampled household panel, not total retailer demand. Weekly quantities vary, and abrupt changes can be difficult to anticipate from past observations alone. Because WAPE divides absolute errors by observed purchase volume, similar absolute errors can also produce a larger percentage when purchase totals are lower.
What we can conclude. These are plausible contributing factors, not established explanations for any model's error. Testing them would require additional analysis of weekly variability, prediction residuals, product behaviour and potential shifts between the training and evaluation periods.
5. Future Extensions
Explore possible next steps, including machine-learning forecasting or synthetic inventory replenishment experiments.
Overview
This section is reserved for the next research extension. Potential directions include machine-learning forecasting and a synthetic inventory replenishment environment. Neither has been selected or evaluated as part of the current ten-method comparison.
5.1 Planned Environment
Specify opening stock, review frequency, shelf-life assumptions, lead times and supplier replenishment rules in a fully synthetic but realistic grocery setting.
Compare alternative replenishment logics, such as simple order-up-to policies, forecast-based review rules, and scenario-driven planning responses.
Evaluate service level, stockouts, waste, overstocks and cost implications under consistent demand scenarios.
Use the demand forecasts as one layer inside a broader simulated operational environment, without claiming access to confidential real inventory data.