Extracted main text — title through conclusion, appendix excluded. This is what our citation measures are computed over, published so the extraction can be checked by eye.
44,712 characters · 17 sections · 48 citation commands
Stealing accuracy: Predicting day-ahead electricity prices with temporal hierarchy forecasting (THieF)
Operational decisions often require tailored short-term forecasts that focus on different levels of detail and granularity pet:etal:22. For instance, models for hourly products in wholesale electricity markets can use different information sets than those for average daily prices wer:14,mac:wer:16. These forecasts may not align, which can lead to suboptimal decisions. To cope with this, the forecasts from each temporal level of the hierarchy should be reconciled to be coherent.
The last decade has seen an unprecedented growth in interest in forecast reconciliation ath:hyn:kou:pan:24:rev and the introduction of temporal hierarchy forecasting ath:hyn:kou:pet:17. The latter can be applied to any time series by means of non-overlapping temporal aggregation -- the predictions computed at all levels of the hierarchy are combined to yield temporally reconciled, accurate and robust forecasts. The concept is new in energy forecasting and there are only a handful of publications on predicting electric load nys:lin:pin:mad:20,mol:nys:mad:24 or wind jeo:pan:pet:19,eng:abo:24,sha:bha:jai:24 and solar yan:qua:dis:r-g:17,dif:gir:23 generation. More importantly, temporal hierarchies have not yet been applied to electricity price forecasting (EPF). The only related study predicts aggregated demand and supply curves by exploiting their intrinsic hierarchical structure ghe:zie:24, but is not concerned with temporal hierarchy forecasting nor EPF.
To fill this gap and provide market participants with a universal tool to improve the predictions of hourly, block and baseload prices (i.e., daily average prices), we conduct an extensive study involving four classes of models and two 4-year test periods (2021-2024) from two major European power markets: Germany and Spain. We use a variety of base forecasting methods with increasing complexity and diverse modeling assumptions, ranging from a parsimonious expert model built on classical linear regression zie:wer:18,ser:wer:25 to Amazon's state-of-the-art tabular foundation model based on a 12-layer 72 million-parameter transformer architecture zha:rob:25.
To our best knowledge, we are the first to use pretrained foundation models not only for EPF, but also for temporal hierarchies. As such, we complement the literature that has focused on evaluating the benefits of THieF with conventional, relatively small, base forecasting methods. We find that temporal reconciliation leads to significant improvements in terms of the mean absolute error (MAE) and the root mean squared error (RMSE) across all classes of models and hierarchy levels considered: up to 5.5% for hourly and 13.4% for baseload prices.
The remainder of the paper is structured as follows. In Section (ref) we briefly review the THieF approach and discuss its application to electricity price forecasting. Next, in Section (ref) we introduce the datasets and in Section (ref) we explain how the base forecasts for the 24 hourly prices are determined using four classes of models: linear regression ($\rightarrow$ Section (ref)), committee machines of shallow feedforward neural networks ($\rightarrow$ Section (ref)), gradient-boosted decision trees ($\rightarrow$ Section (ref)), and a state-of-the-art, pretrained transformer ($\rightarrow$ Section (ref)). In Section (ref) we discuss the results obtained. Finally, in Section (ref), we conclude and speculate on how our approach can be extended in future studies.
THieF can be understood in three steps: (i) construct non-overlapping temporally aggregate time series (levels), (ii) generate base forecasts at these levels independently, (iii) reconcile the base forecasts to combine the diverse information ath:hyn:kou:pet:17. Temporal aggregation is a moving average that filters and strengthens different aspects of the original signal, therefore, the combination of the base forecasts helps recover information that would otherwise be difficult to estimate kou:pet:tra:14.
The benefits of THieF primarily stem from three elements. First, it independently models different views of the time series, which is particularly advantageous in the presence of increased model misspecification ath:hyn:kou:pet:17. Second, it relies on an implicit forecast combination, leveraging the relationships among temporally aggregated series kou:ath:19. Third, it uses restricted estimation of the implied combination weights, which reduces estimation uncertainty pri:sve:kou:21.
In Figure (ref) we plot the out-of-sample base and reconciled forecasts at four aggregation levels for a sample day, market, and model. THieF blends the information of the base forecasts, and the relatively low bias of the average daily price forecast helps to reduce the bias at all other levels. The combination achieved by THieF is affine, but not convex (weights sum to unity but can be negative). This can help to further improve forecasts.
Let $m$ be the frequency of a time series, for instance, $m=24$ for hourly data. Let $k = \{k_i\}_{i=1}^p$ be the divisors of $m$, ordered such that $m = k_p > k_{p-1} > \ldots > k_1 = 1$, and $p$ is the number of factors. For example, for hourly data we have $k_i \in \{24, 12, 8, 6, 4, 3, 2, 1\}$ and $p=8$. We define the summing matrix of size $(\sum_{i=1}^p{k_i}) \times m$ as:
where $\mathbf{S}_{k}^T$ is the transpose of $\mathbf{S}_{k}$. Each $\mathbf{S}_{k}$ is given by:
i.e., the weight $s_k = \frac{1}{k}$ multiplied by the Kronecker product between the identity matrix $\mathbf{I}_{m/k}$ of size $\frac{m}{k} \times \frac{m}{k}$ and a row vector $\mathbf{J}_{1,k}$ of size $1 \times k$ with all its elements equal to 1. For instance, for hourly data we get:
This is illustrated in Figure (ref). The top row $\mathbf{S}_{24}$ represents the baseload block, i.e., the average of the 24 hourly prices. The two rows of $\mathbf{S}_{12}$ correspond to the 12-hour blocks, respectively for hours 1-12 and hours 13-24. Analogously, the three rows of $\mathbf{S}_{8}$ map the hourly prices to the average price for the first, middle and last 8-hourly blocks. The bottom-level $\mathbf{S}_{1}$ is a $24\times24$ identity matrix and corresponds to forecasts for the individual hours.
The temporal aggregation levels prescribed by $k$ are non-overlapping for two reasons. First, overlapping temporal aggregation introduces serial autocorrelation that would need to be captured by the base forecasting methods. This not only increases model complexity but can also obscure identification of the underlying signal. Second, under the reasonable assumption that the base models may be misspecified, such serial autocorrelation can further complicate reconciliation in THieF.
Let $\mathbf{p}_d$ be a column vector that contains the hourly observations for day $d$, then $\mathbf{P}_d=\mathbf{S}\mathbf{p}_d$ implements step (i). Further, let $\mathbf{\hat{P}_{d}}$ be a vector of unreconciled base forecasts $\rightarrow$ step (ii). Following ath:hyn:kou:pet:17, THieF reconciliation is performed each day by computing:
where $\mathbf{W}$ is the covariance matrix of the base forecast errors, estimated on the errors from the training sample $\rightarrow$ step (iii). Since the covariance matrix for our hierarchical forecasts is $60\times60$, the total number of elements to estimate (1800) is larger than the sample size (1085; see Section (ref)). Therefore, we resort to regularized covariance matrix estimators.
Note that our formulation differs from ath:hyn:kou:pet:17 by using the mean for the aggregation, as in kou:pet:tra:14. This is due to the effect of weights $s_k$ in the calculation of each $\mathbf{S}_{k}$ in Eq. (ref). Effectively, this normalizes $\mathbf{W}$.
To see this, let $\mathbf{S}_\Sigma$ be the summing matrix as in ath:hyn:kou:pet:17, which is constructed using Eq. (ref), but with all $s_k \equiv 1$. Further, we set $\mathbf{\Lambda} = diag(1/\mathbf{S}_{\Sigma}\mathbf{1})$ and use it to rescale sum-aggregates into mean-aggregates. Since $\mathbf{S} = \mathbf{\Lambda} \mathbf{S}_{\Sigma}$ and $\mathbf{\Lambda}^T = \mathbf{\Lambda}$, we can rewrite Eq. (ref) as:
Given that $\mathbf{W}_{\Sigma}^{-1} = (\mathbf{\Lambda}^{-1} \mathbf{W} \mathbf{\Lambda}^{-1})^{-1} = \mathbf{\Lambda} \mathbf{W}^{-1} \mathbf{\Lambda}$ and multiplying both sides by $\mathbf{\Lambda}^{-1}$ we obtain:
The terms within the square brackets are identical to the reconciliation using $\mathbf{S}_{\Sigma}$ in ath:hyn:kou:pet:17. Therefore, the use of $\mathbf{S}$ instead of $\mathbf{S}_{\Sigma}$ normalizes $\mathbf{W}$, but otherwise results in the equivalent sum-aggregate forecasts.
We considered three shrinkage approaches for estimating the covariance matrix $\mathbf{W}$, the constant correlation model led:wol:04, the diagonal shrinkage target sch:str:2005, which has been shown to perform well in forecast reconciliation wic:etal:2019, and the so-called `variance scaling' ath:hyn:kou:pet:17 that sets all off-diagonal elements to zero and corresponds to the Weighted Least Squares (WLS) estimator. We found that the two first approaches performed similarly well. The variance scaling also provided improvements, but more mediocre than the alternatives. In Section (ref), we report the results for the less conservative, constant correlation shrinkage of led:wol:04. The results for the remaining approaches are provided in (ref).
To ensure a sound assessment of the THieF approach, we consider two major European power markets: EPEX-DE (Germany) and OMIE (Spain). Both datasets span seven years (05.01.2018-31.12.2024, see Figure (ref)) and include publicly available data:
Note that the EPEX-DE dataset is the same as used by ser:wer:25.
The first 1092 days (until 31.12.2020) are the initial training window. Note that after accounting for the lags in prices, see Eq. (ref) below, this corresponds to $1092-7=1085$ data points for training. Each day the training window is rolled forward by 24 hours.
The remaining 4-year period (starting 01.01.2021) is a challenging test set that includes the COVID-19 pandemic (lockdowns in Germany and Spain began in March 2020, with 15.12.2020 being the first day of hard lockdown in Germany), the Russian invasion of Ukraine (24.02.2022), and the soaring natural gas prices (Q4 2021-Q4 2022). This period is also marked by the appearance of negative price spikes in Germany. For instance, on Sunday 02.07.2023 at hour 15 the price dropped to $-$500 EUR/MWh due to low demand and high renewable generation. On the other hand, negative prices were not seen in the Spanish market before 2024, and even then, prices dropped only a few EUR below zero. However, near-zero prices were reached about 5% of the time in 2024.
To demonstrate the versatility of our approach, we generate base forecasts using four models with four distinct architectures: linear regression, a shallow neural network, gradient boosting, and a transformer. The first two are often used as benchmarks in EPF uni:wer:zie:18,zie:wer:18,hub:mar:wer:19,mac:nit:wer:21,jan:woj:22,bil:gia:del:rav:23,uni:mac:23,ghe:zie:25,ser:wer:25. The third has been reported to perform well in numerous forecasting competitions, including M5 jan:etal:22. The fourth uses a transformer architecture, which is still rare in EPF bot:etal:23. It is AutoGluon's new tabular foundation model called Mitra, which excels on datasets with less than 5,000 samples and 100 features zha:rob:25.
All four models compute the electricity price forecast $\hat{p}_{d, h}$ for day $d$ and block $h$ as a function of 20 features, as in ser:wer:25:
where
All models are trained independently for each block $h$ using a 3-year window of past values; each day the window is rolled forward by one day. Overall, we consider 60 blocks: 24 $\times$ 1H (i.e., 24 one-hour blocks), 12 $\times$ 2H, 8 $\times$ 3H, 6 $\times$ 4H, 4 $\times$ 6H, 3 $\times$ 8H, 2 $\times$ 12H, and 1 $\times$ 24H (i.e., baseload).
For more robust parameter estimation, following uni:wer:zie:18 and lag:mar:sch:wer:21, before training the models, we preprocess the inputs using the area hyperbolic sine transformation: $\text{asinh}\left((y - \hat{\mu}_y)/\hat{\sigma}_y\right)$, where $\hat{\mu}_y$ and $\hat{\sigma}_y$ are the sample mean and sample standard deviation of $y$ estimated on the training set. More precisely, $p_{d-i, h}$ for all $i=1,...,7$ and each selected block $h$ are transformed with a common $\hat{\mu}_y$ and $\hat{\sigma}_y$ calculated for the vector $[p_{d-1092, h}, ..., p_{d-1, h}]$. The extremes $p^{min}_{d-1}$ and $p^{max}_{d-1}$ as well as the exogenous variables $\hat{L}_{d, h}$, $\hat{W}_{d, h}$, $\text{API}_{d-2}$ and $\text{TTF}_{d-2}$ are transformed independently, also using 3-year vectors of past values; weekday dummies are not transformed.
This expert -- in the sense of zie:wer:18 -- model is estimated via ordinary least squares (OLS) and uses the same inputs as ser:wer:25 to allow for direct comparisons. More precisely, the ARX model is given by:
The nonlinear counterpart of ARX approximates $f(\cdot)$ in Eq. (ref) with a shallow feedforward neural network using the series-parallel architecture xie:tan:lia:09. As in hub:mar:wer:19, the hidden layer consists of 5 neurons and uses hyperbolic tangent activation, with a linear function in the output layer. The weights are calculated in Matlab R2025a using the Levenberg-Marquardt algorithm with early stopping based on a 10% validation set. To mitigate the effects of parameter estimation uncertainty, we clip the network output to $[-3, 3]$ before applying the inverse transformation (i.e., the hyperbolic sine, see Section (ref)) and use a committee machine of 10 networks mar:uni:wer:19. More precisely, the final $\hat{p}_{d, h}$ is obtained by training the network 10 times for each $d$ and $h$, and averaging the 10 price forecasts; no hyperparameter optimization is performed.
The third model is an ensemble of gradient-boosted decision trees has:tib:fri:09. We use eXtreme GBDT implemented in the Python package XGBoost v3.0 che:gue:16, with the mean squared error (MSE) as the loss function for training and hyperparameter tuning. The latter is carried out independently for each of the 60 blocks (see Section (ref)) using Bayesian sequential optimization optuna_2019 on the initial 3-year training window (05.01.2018-31.12.2024; see Figure (ref)) with a random 10% left out for validation. We select from the following set of hyperparameters: maximum depth of a single tree (between 2 and 10), learning rate (between $10^{-4}$ and $1$ on log scale), subsample ratio of the training instances (between $\frac12$ and 1), minimum sum of instance weight required for partitioning a leaf node (between 0 and 10), minimum loss reduction required for partioning a leaf node (between 0 and $\frac12$), strengths of the L1 and L2 regularization on model weights (each between $10^{-3}$ and $10$ on log scale). Additionally, we apply early stopping to limit the number of boosting rounds (capped at 1000). Since the hyperparameter selection process is stochastic, we perform it 10 times and the final $\hat{p}_{d, h}$ is obtained by averaging the forecasts generated from 10 XGB instances with different hyperparameter sets.
The final forecaster is Amazon's state-of-the-art tabular foundation model called Mitra, released with AutoGluon v1.4 in July 2025, based on a 12-layer 72 million-parameter transformer architecture zha:rob:25. Since it is pre-trained on purely synthetic data, evaluating its performance on historic time series does not pose data contamination issues. Although Mitra is not specifically tailored to time series, our forecasting task requires only one-step-ahead predictions, with Eq. (ref) being equivalent to a tabular regression problem.
We generate forecasts in zero-shot mode, i.e., do not optimize the hyperparameters and rely on a pre-trained model specialized for regression tasks. We fit Mitra using AutoGluon's TabularPredictor class with finetuning turned off and other settings kept default. Mitra is based on the in-context learning paradigm, hence the training set acts as support examples that condition the prediction. A single forecast $\hat{p}_{d, h}$ is generated by passing input data, see Eq. (ref), as a query. As with the other models considered in our study, we use Mitra in a rolling window approach and train one model for each day $d$ and block $h$ in the test period.
To measure the differences in predictive performance, we compute the mean absolute errors (MAE) and the root mean squared errors (RMSE) of the four models presented in Section (ref) for both markets and over the entire 4-year test period. The errors and the respective percentage gains from forecast reconciliation are reported in Table (ref). Additionally, the gains are plotted in Figure (ref). Clearly, THieF improves forecast accuracy across all aggregation levels, for every considered base model and for both markets.
The relative gains in the German (EPEX-DE) market are in general greater than in the Spanish (OMIE) market; The highest accuracy improvements in MAE are obtained for the XGB model, with a 6.7% average gain in Germany and a 4.8% in Spain. XGB also enjoys the highest average reconciliation gains in terms of RMSE, with improvements of 7.6% and 5.0%, respectively. Additionally, gains are typically higher for larger blocks, which is particularly relevant given the quadrupling of the number of products in the European day-ahead market and the resulting importance of block trading EC:25:96.
Importantly, all reconciled forecasts are significantly more accurate than the corresponding unreconciled ones, as measured by the multivariate variant of the Diebold-Mariano test zie:wer:18,lag:mar:sch:wer:21 based on the L2-norm of forecast errors. For the EPEX-DE market, the test indicated significance at the 1% level for each aggregation level (1H, 2H, 3H, 4H, 6H, 8H, 12H, 24H) and model. The same holds for the OMIE market with the exception of three isolated cases: 1H NARX, 1H Mitra and 24H NARX, for which the differences in unconditional predictive ability were significant at the 5% level.
The fact that the results are consistent across architectures ranging from parsimonious linear regression with 20 parameters to a state-of-the-art transformer with up to 72 million parameters is a strong argument in favor of using the THieF approach in day-ahead EPF. Moreover, the computational overhead of implementing THieF in a forecasting pipeline scales linearly with the baseline. The cost of the reconciliation itself is negligible, hence the only relevant cost is estimating models at additional aggregation levels. Overall, THieF increases the computational cost to about 2.5 times that of the 24-hour baseline. This increase can be considered a minor overhead in the context of day-ahead electricity price forecasting, since the per-day runtime for a single ARX, NARX, or XGB model is sub-second on a consumer-grade CPU, whereas Mitra requires roughly one second on an NVIDIA RTX A600 GPU.
Additionally, to analyze the reconciliation gains not only across different models but also across market regimes, Figure (ref) provides a yearly breakdown of the RMSE-based gains reported in Table (ref). Notably, THieF delivers accuracy gains in the German EPEX market across all years and blocks for ARX, XGB, and Mitra. For NARX, there are a few instances in which the reconciled forecasts are less accurate than the unreconciled ones, in particular for hourly forecasts in 2023-2024 and for 3-, 4-, and 8-hour blocks in 2024.
The results for the Spanish OMIE are less robust. While XGB is consistently improved by THieF in each year considered, ARX does not benefit in 2022 and loses accuracy in 2023. Similarly, the advantage of THieF is less evident for NARX in 2022 and for Mitra in 2023-2024. These findings may reflect how relevant the historical price dynamics from the 2022 energy crisis are for forecasting future prices in these two markets. Note that Spanish electricity prices exhibit lower volatility from 2022Q3 onward: after the significant spikes in the first half of 2022, prices remained below 300 EUR/MWh, with only a few spikes above 200 EUR/MWh in 2023-2024. After the extreme volatility during the crisis, the Spanish market entered a relatively calmer regime, whereas large price spikes were still observed in the EPEX market even in the final days of 2024.
In this study, we demonstrated the benefits of temporal hierarchy forecasting (THieF) for predicting electricity prices. We showed this for the German and Spanish power markets, using a variety of base forecasting methods with increasing complexity and diverse modeling assumptions. We found a consistent and a statistically significant improvement from reconciliation in all setups -- up to 5.5% for hourly and 13.4% for baseload prices -- with increasing gains for higher aggregation levels ($\rightarrow$ larger blocks).
The latter is particularly relevant given the quadrupling of the number of products EC:25:96. On 30 September 2025, the European day-ahead electricity markets moved from hourly ($\rightarrow$ 24 products per day) to 15-minute trading intervals ($\rightarrow$ 96 products per day). Although this finer granularity offers a more accurate representation of expected generation and demand, it also adds complexity to electricity trading, thereby increasing the importance of block trading and temporal reconciliation.
In principle, if the base forecasts would have no misspecification, there would be limited space for THieF to improve accuracy. However, in practice, any forecasting method will have some degree of misspecification from the unknown underlying electricity price generating process. THieF augments these base forecasts by enabling them to better capture the different short and long-term dynamics present in the time series. An interesting finding is that the improvements are consistent for all base forecasts, including Mitra. The foundation model's misspecification is arguably of a different nature because it is used in a zero-shot, pre-trained mode that does not involve estimating model parameters. Moreover, our evaluation provides additional evidence of the usefulness of THieF when base forecasts use covariates, with the majority of applications in the literature focusing on the univariate case.
There are various extensions that stem from this work. We applied the same base forecast to all aggregation levels, however, more diverse model selection could be beneficial, better leveraging the properties of each level, and thus minimizing computational overheads. Likewise, considering a richer set of base models would be interesting, including the well-performing in EPF applications LASSO-Estimated AutoRegressions lag:mar:sch:wer:21 and Distributional Deep Neural Networks mar:nar:wer:zie:23. Furthermore, investigating the impact of temporal reconciliation on probabilistic EPF is a useful extension. Beyond any accuracy gains, THieF provides reconciled forecasts across all levels, which can support the alignment of different decisions and potentially new strategies in energy markets. Finally, the temporal hierarchy was restricted to using aggregation levels based solely on factors of $m$. Unequal aggregation blocks are relevant for decisions supported by EPF, and embedding these in THieF can be beneficial both for forecasting performance and decisions. These are fruitful avenues for future research.
Conceptualization -- NK, AL, RW; Data curation -- AL; Formal analysis -- NK, AL; Funding acquisition -- NK, RW; Investigation -- KB, AL; Methodology -- NK, AL, RW; Software -- KB, AL; Supervision -- RW; Validation -- NK, AL, RW; Visualization -- AL; Writing (original draft) -- KB, NK, AL, RW; Writing (review & editing) -- NK, RW.
The datasets and sources are described in Section (ref). The replication package is publicly available at \url{https://github.com/lipiecki/thief}.
The study was partially supported by the National Science Center (NCN, Poland) through grants no.\ 2018/30/A/HS4/00444 (to AL and KB) and 2021/43/I/HS4/02578 (to RW), and by Riksbankens Jubileumsfond (Sweden) through grant no.\ SAB22-0073 (to NK).