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.
74,290 characters · 13 sections · 65 citation commands
Multivariate Simulation-based Forecasting for Intraday Power Markets: Modelling Cross-Product Price Effects
Keywords: Intraday Electricity Markets, Electricity Price Forecasting, Volatility Forecasting, Copula, Probabilistic Forecasting, Monte-Carlo Methods, SIDC \newline Acknowledgements: Simon Hirsch is employed at Statkraft Trading GmbH and gratefully acknowledges support through Statkraft (\url{https://www.statkraft.com/}). This work contains the author’s opinion and does not necessarily reflect Statkraft’s position. The authors declare no conflict of interest. \newline
Intraday electricity markets are used to balance the short-term intermittency of renewable energy assets. The increasing penetration of wind and solar yields the need for accurate probabilistic price forecasts. This need is underscored by the heightened volatility in light of the European energy crisis in 2022/23. To date, the literature on (probabilistic) electricity price forecasting on intraday markets remains scarce narajewski2020ensemble, hirsch2022simulation, uniejewski2019understanding, narajewski2020econometric and is exclusively focused on univariate approaches. However, in most European countries, the intraday market is a continuous forward market, where all delivery periods for a delivery day $d$ are traded in parallel, thus univariate approaches neglect the complex dependency structure in these markets. As products close with the physical delivery of electricity, it is not possible to “glue” subsequent trading sessions together, as it is commonly done with equity markets. Additionally, while the spot market is driven by the absolute level of fundamentals such as wind and solar production, the intraday prices are influenced more by the changes in forecasts ziel2017modeling. Our work focuses on these challenges in modelling the dependency structure and fills the according gap in the literature. Based on the work of narajewski2020ensemble, hirsch2022simulation on the marginal distribution of the price process, we employ copulas to model the time-dependent correlation. We model the dependency parameter as latent variable using beta regression. We validate our results in a forecasting study for the German intraday electricity market. Our results indicate that modelling the dependence structure between the different trading session improves the forecasting performance. Additionally, we provide evidence on market efficiency in the German short-term market. Our fundamental and parametric approach allows us to shed further light on the impact of the cross-border shared order books of the single intraday coupling (SIDC) and the driving factors of the distribution parameters.
Let us give an illustrative example for the trading schedule of the German intraday electricity market. Trading for physical delivery on $d$ starts on the previous day $d-1$, 15:00 hours and lasts till few minutes before the start of physical delivery of electricity. For example, trading for the delivery on $d$, 18:00-19:00 started on $d-1$ at 15:00 hours and closes at $d$, 17:55 hours. During the beginning of this trading window, traders will be able to trade power in many neighbouring delivery hours, from delivery at $d$ 16:00-17:00 (which closes at $d-1$ 15:55), to all delivery periods for the next delivery day $d+1$ (which start trading at $d$, 15:00 hours). Figure (ref) shows all trades in the two trading sessions for delivery on 2022-12-01 and 2022-12-02 for all hourly products on the German intraday market.
The literature on intraday electricity price forecasting can be divided into three main groups: $(1)$ paper which treat the intraday market in a similar fashion as the day-ahead market and predict (index) prices along the delivery time line uniejewski2019understanding, cramer2022multivariate, kath2021conformal and $(2)$ paper which predict prices along the trading time for single delivery periods serafin2022trading, hirsch2022simulation, narajewski2020econometric, narajewski2020ensemble, janke2019forecasting, marcjasz2020beating. The correlation between different trading windows has, to the best knowledge of the authors, not been investigated so far, while the correlation between day-ahead and intraday (index) prices has been the subject of studies such as andrade2017probabilistic, cramer2022multivariate. Lastly, $(3)$, empirical, in-sample studies on price formation in intraday electricity markets have been conducted by kiesel2017econometric, kremer2021econometric, kath2019modeling.
In light of the aforementioned challenges and the gap identified in the literature, our contributions are:
Our main strategy is the canonical inference-for-margins patton2012review approach commonly used for copula modelling and can be summarized as follows: We use the probabilistic models developed by narajewski2020ensemble, hirsch2022simulation as a starting point to gaussianize the intraday market observations. We use the pseudo-Gaussian observations to fit the dependency structure. For forecasting, we simulate (multivariate) Gaussian random variates and use the inverse probability integral transformation to receive samples in the desired marginal distribution. Within the energy markets literature, similar approaches have been used for simulating wind and load forecast deviations tastu2015space, carmona2021glasso, carmona2022joint, day-ahead electricity prices manner2019forecasting, pircalabu2017regime, berrisch2023modeling and the design of hedging strategies pircalabu2017mixed. Our results show that modelling the dependence structure leads to significantly improved forecasting performance compared to univariate approaches. However, we find that time-dependent modelling of the dependence structure is of little added value compared to constant dependence. Additionally, we provide new insight on the effects of the opening and closing of the cross-border order books during SIDC. We interpret our results in light of the market efficiency hypothesis and discuss reflections on modelling already highly volatile prices during a period of increased uncertainty. Figure (ref) gives an illustrative example of our 24-dimensional forecast
Our results offer multiple avenues for further research. First, we restrict ourselves to Gaussian dependence structures. Future work might improve our methods by using copulas that reflect possible tail dependence effects or use Vine-copulas to better approximate the dependence structure. Secondly, we note that the zero-inflated Johnson's $S_U$ distribution offers a good, but not yet perfect fit for the marginal distribution in the intraday market and hence further research in the marginal distribution of intraday electricity prices is required. Third, in light of the continued development of SIDC and the planned introduction of intraday auctions within the SIDC system and the possible integration of interconnector cables in the SIDC system entsoe2023single, our results on the impact of SIDC on the trading activity provide a fruitful starting point for further research kath2019modeling.
What is more, our results are also relevant for researchers and practitioners working on stochastic optimization of bidding strategies for the intraday markets. For storage assets such as batteries and pumped-hydro, modelling the dependency structure is important as charging (pumping) positions some periods depends on the ability to discharge (generate) later during the day and therefore depends on the dependency structure. Recent works as boukas2021deep, nolzen2022market, finhold2023optimizing use sampled paths for the intraday market, but model the dependency structure only implicitly, if at all, and thus might produce too optimistic results.
The remainder of this paper is structured as follows: the following Section (ref) gives a detailed introduction to the German short-term electricity market. Section (ref) introduces our data set, preparation and summary statistics. We present our modelling approach in Section (ref). Section (ref) describes the forecasting study design and scoring rules. Finally, Section (ref) and (ref) present and discuss our results and conclude this paper.
Electricity markets are structured as forward markets. The German short-term electricity market consists of three major parts: (1) the daily spot auction, (2) the continuous intraday market and (3) the balancing market. The following discussion focuses on the spot and intraday markets. The spot market is the main electricity market in Germany. It is organized as a daily auction at noon on which electricity for all 24 delivery hours on the following day is traded. The following intraday market is used to balance deviations in forecasts after the day-ahead market. It is organized as continuous trading similar to equity markets. The continuous trading starts at $d-1$, 15:00 hours and closes shortly before delivery. After the delivery period ends, remaining imbalances between the traded position and the actual production are settled in the balancing market with the TSO. However, strict market regulation in Germany prohibit explicit active position taking in the balancing market. Figure (ref) depicts the daily procedure for a single delivery hour.
Let us introduce some nomenclature to ease the following discussion of the market structure of the spot and intraday market. We refer to the delivery time as the time of actual production power, while the trading time refers to the time at which a trade for a certain delivery period is conducted. As a general rule, we try to denote delivery time in superscript, while we denote trading time in subscript. We refer to a trading session on the intraday market as the time window between market opening on $d-1$, 15:00 hours to gate closure. Note that for different delivery periods, trading sessions are of different length.
The spot market in Germany is organized by EPEX Spot and Nordpool AS with shared order books. On the spot market, electricity for all 24 hours for the following day is traded. The market is organized as a pay-as-cleared auction and its order book closes at $d-1$, 12:00 hours. Results are published at $d-1$, 12:42 hours. The minimum price is currently set to -500 EUR/MWh and the maximum price is set to 3000 EUR/MWh.
The intraday market is structured as continuous forward market. For all hourly delivery periods with delivery on $d$, trading starts at $d-1$, 15:00 hours and ends 5 minutes before the actual start of delivery. Within the Single Intraday Coupling (SIDC), the order books of all major continuous intraday markets across Europe are coupled, as long as there is sufficient cross-border transmission capacity. The coupling proceeds in two waves: first, at $d-1$, 18:00 hours, the order books of Germany, Denmark, Sweden, Poland and Norway and Netherlands\footnote{Norway and Netherlands are coupled through the NorNed high-voltage submarine cable. Within the central European Core region, the Netherlands is coupled to the geographic neighbours at 22:00 hours.} are coupled. At $d-1$, 22:00, France, Netherlands, Belgium, Austria, the Czech Republic, Hungary, Romania follow nemo2021single. All shared order books close 60 minutes before the start of physical delivery and trading resumes with Germany wide delivery. 30 minutes before the start of physical delivery, the Germany wide delivery and trading resumes on a TSO/grid zone-level. Finally, 5 minutes before delivery, the grid-zone trading closes as well. Note that all hourly (and also half-hourly and quarter-hourly) delivery periods for a delivery day $d$ are traded in parallel, as it shown in Figure (ref). A more detailed overview on intraday electricity markets can be found in shinde2019literature and viehmann2017state.
The following chapter gives a brief introduction of the data used in this paper and the required pre-processing. We use intraday transaction data from EPEX, the anonymous day-ahead spot auction bid curves from EPEX, and wind, solar and demand forecasts from SMARD respectively ENTSO-E. We also provide summary statistics. As a rule of thumb, superscript indices denote delivery periods while subscript indices denote trading time. We hope this makes the forward market structure of the short-term electricity markets more clear.
Intraday transactions are individual trades conducted on the continuous intraday market. We use only trades conducted with either the buy- or sell leg in one of the 4 German grid zones. As trading happens continuously in the intraday market, transactions are irregular spaced in time and need to be aggregated. In line with hirsch2022simulation, narajewski2020ensemble, narajewski2020econometric, serafin2022trading, we aggregate all trades for delivery period $d, h$ on a 15-minute equidistant grid along the trading time (denoted with $t$), where $t=0$ denotes the first 15 minutes after trading start. Let $P_t^{d,h}$ denote the volume-weighted average price of all trades with delivery period $d, h$ belonging to bucket $t$ and $\alpha_t^{d,h}$ denote a Boolean indicator whether there was at least one trade. The spot price for each delivery period is denote as $P_\text{Spot}^{d,h}$. We drop all trades in the local trading phase in the last 30 minutes to the start of physical delivery. The full aggregation process can be seen in Figure (ref). Summary statistics for the price differences are given in Table (ref). We especially note an increasing volatility in year 2022 driven by the Russian invasion in Ukraine and the energy crisis in Europe. Additionally, we note that trading activity in general increases as the share of no-trade event decreases. Figure (ref) gives the pairwise correlation between price changes for all $24 \times 24$ intraday delivery periods during the training set. We already note a cluster of high correlation for the night hours and for the afternoon peak hours.
We use wind on- and offshore, solar and demand forecasts from ENTSO-E. The data is aggregated to hourly resolution using a simple arithmetic average. Forecasts are generated by the transmission system operator for each delivery period $d, h$ and available at the day-ahead stage, i.e. latest at $d-1$, 12:00 o'clock. We denote the forecasts as $\ensuremath{\text{WindOn}^{d,h}}, \ensuremath{\text{WindOff}^{d,h}}, \ensuremath{\text{Solar}^{d,h}}$ and $\ensuremath{\text{Load}^{d,h}}$. For all data, we adjust the daylight saving times by (back-) filling the missing hour in spring and averaging the double hour in autumn as it is standard in the electricity price forecasting literature.
Additionally, we employ a metric for the merit-order regime. In the classical model, the merit-order is defined as the supply side of the electricity market, sorted by the marginal production costs. The intersection of the supply and demand curve gives the market price. Depending on the slope of the merit-order, changes in the supply or demand have different impacts on the market price. kremer2021econometric and hirsch2022simulation have shown that the different merit-order regimes explain the size and volatility of price changes in the intraday markets. However, modelling the merit-order in short-term power markets is not straight forward and several approaches have proposed. We follow the approach of hirsch2022simulation by using anonymous bid and offer curves from the day-ahead market to model the intraday merit order. The curves are published by EPEX Spot around 14:30 and are therefore available at the time of forecasting.
The overview of the strategy for the calculation of the merit-order regime coefficient $\text{MO}^{d,h}$ is given in Figure (ref). We take the anonymous supply and demand curves from the auction and transform these into an elastic supply curve and an inelastic demand curve using the same transformation as in hirsch2022simulation, kulakov2021impact. This transformation is based on the idea that buying 50 MW up to a price of 100 EUR/MWh is the same as buying 50 MW at any price and placing a sell order at 100.01 EUR/MWh. Intuitively, we can therefore create a price independent buy curve and move the price-dependent bids to the sell side. The resulting demand and supply curves are depicted in the middle panel of Figure (ref). Note that the equilibrium price at the intersection between supply and demand curve is unchanged. Lastly, we calculate the slope around the equilibrium price as finite difference quotient. For the exact calculation, we refer the reader to hirsch2022simulation.
Our approach follows the widely used inference for margins approach for copula models. We first use a univariate model to estimate the time-varying conditional marginal distribution, apply the probability integral transform and estimate the copula distribution respectively the time-varying dependence parameter in the second step. The following two subsections (ref) and (ref) describe our modelling in more detail. Section (ref) describes the general simulation set-up. Subsection (ref) describes the different (nested) models for our forecasting study and our benchmark models.
Let us remark that strictly, Sklar's theorem sklar1973random is only valid for continuous marginal distributions, as $\mathcal{F}_X(X) \sim \mathcal{U}_{[0, 1]}$ is not true for discrete distributions $\mathcal{F}$. This leads to the issue that the copula distribution is not necessarily identifiable with discrete or mixed discrete-continuous marginals. A practical approach to alleviate this issue is to fill the gaps induced through the discreteness of the marginal distribution by some uniform 'jitters', which results in the so-called checkerboard copula, a strategy we employ in this paper as well geenens2020copula, faugeras2017inference.
We model the distribution of the price changes $\Delta P_t^{d,h}$ as a mixture distribution $\mathcal{D}^{d,h}_t$ to account for the zero-inflation in price changes.
where $\mathcal{D}^{d,h}_t$ is a mixture distribution from a continuous distribution and the Dirac distribution with an atom at 0, denoted as $\delta_0$.
where $\alpha_t^{d,h} = 1$ indicates a trade event and is modelled as an binomial variable $\alpha^{d,h}_t \sim \mathcal{B}^{d,h}_t(\pi_t^{d,h})$.
A similar approach is taken by narajewski2020ensemble, hirsch2022simulation, who estimate the mixture distribution in a two-stage procedure. Owing to stylized facts on intraday electricity prices, narajewski2020ensemble, hirsch2022simulation assume $\mathcal{F}$ to follow the (skewed) Student-$t$ or Johnson's $S_U$ distribution. With regards to hirsch2022simulation remarks on the issues with estimation stability using the skew $t$-distribution, we generally use Johnson's $S_U$ distribution in this work. We use a two-step estimation procedure for the estimation. First, we estimate the conditional mixing probability $\pi_t^{d,h}$ using regularized logistic regression. Subsequently, we estimate the conditional distribution parameters on all non-zero price changes $\Delta P_t^{d,h} \mid \alpha_t^{d,h} = 1$ using maximum likelihood. The remainder of this section describes our feature engineering, the exact specifications for our models for the conditional mixing probabilities and the conditional distribution parameters as well as our hyperparameter tuning scheme.\footnote{We have also experimented with estimating the zero-inflation in a single estimation step. However, there are two kinds of zero-inflation present in the intraday trade data: We have periods without trades and periods where trading happens, but the trading does not lead to a change in the price. A joint estimation does not differentiate between both effects and provided inferior results in initial testing. Additionally, the estimation of discrete-continuous mixture distributions is not straight-forward using maximum-likelihood as the probability mass and likelihood functions live on different scales.}
Our predictive variables can be grouped into five groups:
ReLU splines are piecewise linear approximations on the domain of an explanatory variable. For an arbitrary explanatory variable $x$ we define a set of thresholds $\tau \in \mathcal{T}_x$ and we define the ReLU splines as
thereby the clipped values on the domain can contribute with individual slope coefficients. This allows an efficient approximation of non-linear functional relationships while preserving linearity in the coefficients nair2010rectified. ReLU splines are especially useful for variables with known and bounded domain such as the trading time $t$ in our application, but can in theory be used for any continuous variable. Even though the thresholds $\mathcal{T}_x$ can be chosen arbitrarily, we generally use an equidistant grid.
The logistic regression model used to estimate the conditional mixing probabilities $\ensuremath{\pi_t^{d,h}}$ is defined as
where $\beta_0$ denotes the intercept, $\beta_1$ to $\beta_4$ are the coefficients for the fundamental wind, solar and load forecasts, $\beta_5$ to $\beta_9$ are the coefficients for the single intraday coupling dummies $\text{SIDC}_{\{\text{O}, 1, 2, \text{C}, \text{L}\}}(d,h,t)$. The dummies denote whether the cross-border order books are open, the start of the 1st and 2nd wave, the closing of the cross-border order books 1 hour before the gate closure and the last periods of (local) trading. Table (ref) gives the specification of the SIDC related dummies. We include day-of-the-week dummies for the delivery day $d$ with the coefficients $\beta_{10}$ to $\beta_{15}$. The last term denotes a ReLU spline for the trading time $t$ on an equidistant grid $\mathcal{T}_t$, with $\otimes$ denoting the Kronecker product. We conduct a grid search on the validation data set to select the step size of the grid $\mathcal{T}_t$ and the regularization parameter based on the validation accuracy. The step size is constant across all $h$. An exemplary fit for the ReLU spline along the trading time $t$ can be seen in Figure (ref). The estimation is regularized using the $L_2$ norm, which we choose from an exponential grid on 0.001 to 1000. We also evaluate the Bayesian Information Criterion (BIC) for each combination of step size and $L_2$ regularization and find that the results align.
The following paragraph describes the probabilistic model for the marginal distribution $\mathcal{F}$ with time-varying location, scale and shape parameters. Our general framework follows the generalized additive models for location, scale and shape introduced by rigby2005generalized. Let $Y = (Y_1, Y_2, ..., Y_n)$ be a vector of $n$ independent observations $Y_i$ and $Y_i$ have the probability (density) function
where each distribution parameter can be a smooth function of explanatory variables. Denote $\theta_i^k = (\mu_i, \sigma_i, \nu_i, \tau_i) = (\theta^1, \theta^2, \theta^3, \theta^k)$ as the $n \times k$ parameter vector with $k$ location, scale and shape parameters. We have
Let $g_k(\cdot)$ be a known monotonic link function for each distribution parameter relating $\theta^k$ to explanatory variables through an additive model
where $\boldsymbol{X}_k$ is a $n \times J_k$ known design matrix of $J_k$ exogenous regressors. Additionally, the model can consists of non-linear effects as in classical generalized additive models. Note that each distribution parameter can have an individual design matrix $\boldsymbol{X}_k$.
As noted, we assume Johnson's $S_U$ distribution for the marginal distribution of $\Delta P_t^{d,h}$. The probability density function of Johnson's $S_U$ distribution is defined as follows:
where $-\infty \leq \mu \leq \infty, 0 < \sigma \leq \infty, 0 < \nu \leq \infty, - \infty \leq \tau \leq \infty$ represent the location, scale, tail and skewness parameters and their respective domains. To ensure all distribution parameters are within their domain, we employ the link functions
where the link functions for $\mu$ and $\tau$ are the identity and the link functions for $\sigma$ and $\nu$ are known as Softplus link function with constants $\epsilon = 1^{-3}$ and $\gamma = 0.1$ to improve numerical stability sonnenschein2022probabilistic. Formally, we define the model as follows:
We model the conditional location, scale and shape parameters based on fundamental variables. We guide our selection from the literature hirsch2022simulation, narajewski2020ensemble, janke2019forecasting, narajewski2020econometric and regularize the estimation to avoid overfitting.
As the modelling of the distributions higher moments is dependent on the quality of the lower moments' models ziel2022m5, we refrain from making the models for $\nu$ and $\tau$ as complex as the model for the scale parameter $\sigma$ and decrease complexity accordingly. For the expected value of intraday price changes, many studies have shown indications of weak-form market efficiency hirsch2022simulation, narajewski2020econometric, narajewski2020ensemble, uniejewski2019understanding, janke2019forecasting, lohndorf2023value, hence we keep the model simple as possible without comprising the quality of the estimation for $\sigma , \nu$ and $\tau$.
Our estimations are implemented using the scikit-learn scikit-learn and the tensorflow-probability packages dillon2017tensorflow, abadi2016tensorflow in python. We automatically tune the hyperparameters of our model using the well-known optuna library, a Bayesian framework specifically tailored towards the optimization of hyperparameters of machine learning models optuna_2019. The sampling space for the hyperparameter tuning framework can be found in Table (ref). We initialize the coefficients for the estimation of the conditional location parameter $\mu$ as 0, while the remaining coefficient vectors are initialized uniformly sampled. We tune the $L_1$-regularization for all distribution parameters, the learning rate and introduce a dropout layer during the training process to reduce the risk of overfitting bengio2012practical, srivastava2014dropout, ziel2022m5. During model fitting, we reserve 25% of our training data as validation set to employ early stopping if the training-validation loss does not improve after 25 epochs ziel2022m5, marcjasz2022distributional. We run 250 iterations of the optuna algorithm and observe the best trial at iteration 48. Diagnostic plots can be found in the Appendix (see Figure (ref)).
Remember that we simulate the $T \times H$ price difference vector $\Delta P^d_t = (\Delta P^{d,1}_t, \Delta P^{d,2}_t, ..., \Delta P^{d,23}_t)$, where the different delivery periods can be correlated. Our strategy follows the general setup of carmona2021glasso, carmona2022joint. We estimate the dependence as the covariance after we gaussianize our data in the spirit of Sklar's theorem.
In our copula-based modelling approach, we employ three different structures for the dependence with increasing complexity. The simplest model assumes that all delivery periods are independent and is denoted as Mix.Ind. The next complex model assumes a constant cross-product dependence across the trading window and is denoted as Mix.CD. Lastly, we estimate a time-varying dependence parameter across the trading time $t$. The resulting most complex model is denoted as Mix.TD.
We denote the dependence parameter between two delivery hours $A, B$ as $\rho^{A,B}$. For the constant dependence model, we estimate the pairwise dependence parameter as:
The pair-wise estimation is necessary as not all delivery periods have the same trading length. For our time-varying dependence model Mix.TD, we fit $\rho$ as a latent variable using the beta-regression.
To estimate the beta-regression, we use the correlation coefficient on the pseudo-Gaussian observations, which we map to the (0, 1) space required by the beta-regression model. We use the BIC to decide on the width of the grid of the spline on $t$ and transform the estimated values back to the covariance using the inverse mapping and the empirical standard deviation of the pseudo-Gaussian observations.
We simulate $m = 1, ..., M$ paths for the $T \times H$-dimensional intraday price path vector. Denote with $\ensuremath{\Delta P_{t}^{d,h, [m]}}$ the simulated price change at trading time $t$ for delivery period $d, h$. The price path can be seen as the cumulative sum of all price changes and the initial price $P_0^{d,h, [m]}$. In line with previous works nolzen2022market, lohndorf2023value we initialise the intraday price as the spot price $P_{0}^{d,h,[m]} = P_\text{Spot}^{d,h}$. Formally, we have:
where it is important to note that the end of trading, $T$, depends on $h$. Therefore, the vector is not square but has the asymmetric trapezoid shape already visible in Figure (ref) and Figure (ref).
We employ three benchmark models from the literature. First, we employ the well performing naive model introduced by narajewski2020ensemble, hirsch2022simulation. The model re-uses past price trajectories by sampling. We employ the model in two versions, first assuming independence between different delivery hours (i.e sampling each delivery hour $h = 0, ..., H$ individually) and second by sampling the full $T \times H$ path vector. The third benchmark model is a arithmetic random walk model introduced by lohndorf2023value drawing from the empirical price difference distribution. The following paragraphs introduce the benchmark models formally.
The Naive.Ind and Naive.Dep model is defined as
where $d'$ is a random day sampled from the training data set $d' \sim \mathcal{U}(\{0, ... , d-1\})$. Note that for the Naive.Ind we sample $d'$ independent for each delivery hour $h = 0, ..., H$ and for the Naive.Dep, we use the same $d'$ for all $h$. hirsch2022simulation, narajewski2020ensemble have shown that this type of benchmark models provides very good point and probabilistic forecasting performance in ensemble forecasting settings.
The RW.Emp has been used by lohndorf2023value to generate paths for the intraday market in an application study for grid-scale storage optimization. The price process is defined as random walk, where the innovations are drawn from a discrete distribution of the centered empirical price changes. It is defined as:
where again $d'$ is a random day sampled from the training data set $d' \sim \mathcal{U}(\{0, ... , d-1\})$ and the second term centers all residuals to mean zero, ensuring the price process to be a martingale.
For all three benchmark models, the size of the training data $d = 0, ..., d-1$ is a tuning parameter, which we optimize through a grid search.
We employ the well-known rolling window forecasting study design. Our data comprises of the three years between 2020-01-01 and 2023-01-01. We split our data set in a training, validation and test data set. Our initial training set is 2020-01-01 to 2021-07-01, the validation set contains 2021-07-01 to 2022-01-01 and our test set contains the final year 2022-01-01 to 2023-01-01. The split is also depicted in Figure (ref). We use the initial training set to develop and estimate our models and the validation data set to calibrate hyperparameters. We use the test set to run a rolling window forecasting study with monthly re-training of all models using the most recent two years of data. This is due to the high computational burden of training and tuning the probabilistic models in tensorflow-probability.
We evaluate our simulations from a point and probabilistic forecasting perspective using strictly proper scoring rules gneiting2007strictly. We evaluate the mean and median simulation paths using the well known root mean squared error (RMSE) and median absolute error (MAE). The marginal fit is evaluated using the continuous ranked probability score (CRPS), which we approximate on a dense grid of quantiles using the pinball score (PB). We evaluate the scenario paths using the energy score (ES). The energy score is the multivariate generalization of the continuous ranked probability score, taking into account the correlation structure. We establish significance using the Diebold-Mariano test diebold2002comparing, diebold2015comparing for comparing predictive accuracy. The following few paragraphs introduce our metrics formally.
The root mean squared error is defined as:
The RMSE is a strictly proper scoring rule for the expected value. We also report the RMSE averaged additionally across the delivery hours $h$ and delivery days $d$.
Denote the median trajectory over all $M$ simulated trajectories as $\operatorname{med}(P_t^{d,h, [m]})$. The mean absolute error is defined as:
We evaluate the probabilistic forecasting performance using the continuous ranked probability score (CRPS) and the energy score (ES). Both are strictly proper scoring rules for the marginal distribution respectively the multidimensional predictive distribution gneiting2007strictly, ziel2019multivariate and routinely used in probabilistic energy forecasting hirsch2022simulation, narajewski2020ensemble, berrisch2023modeling, nowotarski2018recent.
The CRPS is approximated using the pinball score on a dense grid of quantiles of our scenarios. Let $Q_{\tau,t}^{d,h}(P_t^{d,h,[m]})$ denote the $\tau$-quantile of our simulation paths $P_t^{d,h,[m]}$. The $\tau$-PB can be defined as:
The CRPS is the average across all quantile levels $\mathcal{T} = \{ 0.01, .., 0.99 \}$ of length 99:
For the energy score, we implement the $K$-band estimator as given in ziel2019multivariate:
for an integer $1 \leq K \leq M$ and where we set $P_t^{d,h,[M+k]} = P_t^{d,h,[k]}$. $\left\lVert \cdot \right\rVert_{2}$ denotes the $L_2$ or Euclidean norm. We use $K = 10$ as trade-off between computational complexity and estimation accuracy. We evaluate the energy score for the full scenario trajectory and for the last three hours of before the start of physical delivery for each model. The later evaluation acknowledges the importance of the last hours of trading and appeals to practitioners in the field.
We use the Diebold-Mariano (DM) test to compare the predictive accuracy of the forecasts nowotarski2018recent, diebold2002comparing, diebold2015comparing. Intuitively, the DM-test evaluates the null hypothesis ($H_0$) that the difference in means between the loss series of two models is statistically significantly different from zero. Formally, for two models $A$ and $B$, let $L_A$ and $L_B$ the loss series. The loss differential $\Delta_{A,B}$ is defined as
for the $i$-norm. For each model pair, we test two one-sided tests for the null hypothesis (1) $\operatorname{E}\left[ \Delta_{A,B}^i \right] > 0$ and (2) $\operatorname{\mathbb{E}}\left[ \Delta_{A,B}^i \right] < 0$, i.e. (1) the forecasts of model $B$ outperform the forecasts of model $A$ and (2) the forecasts of model $A$ outperform the forecasts of model $B$. These tests are complimentary. Note that the Diebold-Mariano test assumes the loss differential series to be stationary. We test this assumption using the augmented Dickey-Fuller test dickey1979distribution.
Our parametric approach to modelling the distribution parameters allows us to derive some fundamental insight in the driving factors of the intraday price process. The following section presents first presents some in-sample results from our modelling and subsequently presents the results from our forecasting study.
Our parametric modelling set-up allows us to analyse the influence of the driving factors for the location, shape and scale parameters of the price distribution. Our focus here is on the impact of the SIDC on the trading activity and the impact of fundamental variables on the volatility.
The evolution of the trading activity respectively the share of no-trade events can be seen in Figure (ref). We note that in the first hours of trading and that trading activity rises non-linearly towards gate closure. The opening of SIDC induces to spikes in trading activity when the cross-border order books are coupled. At this point, orders in markets with previously different price levels are instantly matched, leading high trading activity for short periods of time. After the SIDC coupling, trading activity increases towards the gate closure. The last trading periods have (almost) no no-trade events.
We show the contribution of individual groups of regressors on the scale parameter in Figure (ref). Recall that we model the volatility by four main groups of regressors: fundamental forecasts, time-derived variables, SIDC-related variables, and variables related to the trading activity. Generally, we note a pattern of higher volatility in the beginning of the trading session, followed by decrease and an increase closer to delivery. We note that SIDC has a distinct impact on the scale parameter: The opening of the cross-border order books at 18:00 and 22:00 for the leads to clearly visible spikes in the volatility. Subsequently, the phase during which the order-books are coupled is characterized by lower volatility. The closing of the cross-border order books shortly before delivery leads to a spike in volatility. The dampening effect of the open SIDC order books on volatility is likely due to the increased liquidity available to market participants, while the volatility spikes during the opening at 18:00 and 22:00 hours can be explained by matching the order books at different price levels. Overall, the effects contradict the results of kath2019modeling, who finds no effects of SIDC, but align with narajewski2020ensemble, hirsch2022simulation on the effect of the SIDC closing period. Comparing the different delivery hours $h = 0, 6, 12$ and $18$, we note the different impact of the time to delivery and trading time splines. The time-to-delivery spline kicks in early, but keeps rather constant shortly before the delivery. On the other side, the trading time $t$ spline rises with trading time.
The influence of fundamental variables like wind, solar and demand forecasts is small. One reason for this might be, that these forecasts are generated at the day-ahead stage and are not updated throughout the trading window, as previous works kremer2021econometric, ziel2017modeling, hirsch2022simulation have shown that the intraday market is more impacted by forecast changes compared, while the level of forecasts is less important. Additionally, the conclusion that fundamental variables seem not to improve forecasts supports the notion of market efficiency as already indicated by narajewski2020econometric, narajewski2020ensemble, hirsch2022simulation.
This section presents the results of the out-of-sample forecasting study. Aggregate error metrics are given in Table (ref). Figure (ref) presents the error metrics by the delivery hour $h$. Figure (ref) gives the results of our pairwise Diebold-Mariano tests.
Remember that we have both, the naive model and the mixture model in in at least two versions: one version assuming independence between different delivery hours and at least one version that considers the correlation structure. Additionally, the RW.Emp considers the correlation structure implicitly. The results for the energy score show that considering the correlation structure leads to (significantly, see Figure (ref)) better forecasts than assuming independence if the remaining model structure is unchanged. This holds for moving from the Naive.Ind to the Naive.Dep and for moving from the Mix.Ind to Mix.CD. Interestingly, modelling the dependence structure in a potentially time-varying fashion does not improve the forecasting performance. For the energy score, the difference in ordering for the models between the full path and the last three hours of trading. Note however, that the scale of both is not directly comparable. This is an important result for modelling the intraday market in applications such as battery/storage optimisation, where the neglecting the correlation structure can therefore lead to too optimistic results nolzen2022market, lohndorf2023value.
On an aggregate level, we see that the RW.Emp yields the best point forecasting performance (MAE and RMSE) and the Naive.Dep yields the best probabilistic forecasting performance. The Diebold-Mariano test confirms the statistical significance of superior probabilistic forecasting forecasts for the Naive.Dep. We note that the forecasting performance of the mixture models is not as good as the rather simple benchmark models. Additionally, we note that the Mix.Ind exhibits a lower CRPS than the Mix.CD and Mix.TD, which suggests some some cross-propagation of errors. On the other hand, the \textbf{Mix.Ind} yields worse scores for the energy score than its sister models including a dependence structure.
The hourly shape of forecasts errors throughout the day is depicted in Figure (ref). It follows the typical shape of prices in electricity prices, we see low forecast errors in the morning and higher errors through the day. Let us note that the benchmark models exhibit a slightly higher MAE during the afternoon peak, but show lower RMSE and ES during the same periods. Figure (ref) presents the PB across the quantile range and the delivery hour. We can see that the highest errors occur during the afternoon peak hours and in the higher distribution quantiles. Note that the pinball loss only measures the marginal fit to the distribution. In an analysis of the in-sample, transformed observations we also note that throughout the rolling window study, the calibration decreases and we experience an underdispersed forecast (see Figure (ref) in the Appendix).
Our results also emphasize the importance of robust approaches to modelling and forecasting in periods of high volatility and black swan events such as the Russian invasion of Ukraine and the subsequent energy crisis in Europe 2022/23. The Naive.Dep and Naive.Ind model already have shown very good probabilistic forecasting performance in the respective studies of hirsch2022simulation, narajewski2020ensemble and are robust to extreme events, as out-of-support situations with extreme prices cannot happen for these models. Additionally, our results can be viewed in the light of the market efficiency hypothesis. Our result indicate that including more data, especially data that does not change throughout the trading period (such as day-ahead forecasts) does not improve the modelling or the price. The superior performance of the benchmark models, especially of the RW.Emp, which ensures the Martingale assumption, underscores this notion. Similar results with respect to market efficiency have been found by narajewski2020econometric.
This paper presents a a forecasting study for multivariate, simulation-based forecasting for intraday electricity markets. We provide insight in the dependence structure of short-term electricity markets and extend previous works of narajewski2020ensemble, hirsch2022simulation to include cross-product price effects.
We develop a probabilistic model for the marginal distribution of the intraday price path, accounting for the impact of fundamental driving variables on the location, scale and shape parameters. As novelty, we employ copulas to model the (time-dependent) dependency between different delivery periods and the according parallel trading sessions and allow the dependency parameter to be time-dependent. We validate our results in a forecasting study for the German intraday electricity market. Our results indicate that modelling the dependence structure between the different trading session improves the forecasting performance. Additionally, we provide evidence on market efficiency in the German short-term market. Our fundamental and parametric approach allows us to shed further light on the impact of the cross-border shared order books of the single intraday coupling (SIDC) and the driving factors of the distribution parameters. Our case study employs data from the German intraday electricity market, but our method is directly transferable to other European electricity markets.
While we are able to show that modelling the dependence structure improves the forecasting performance, our methods can surely be improved and our results offer a multitude of further research areas: First, a further investigation of the dependency structure seems worthwhile. However, the improved modelling of the correlation structure is dependent on having a suitable probabilistic marginal model at hand. We note that the proposed mixture of Johnson's $S_U$ distribution still does not provide a perfect fit and struggles to cope with periods of extreme volatility. Hence, the modelling of the price distribution is an important field for further research. Third, we provide new evidence on the impact of the single intraday coupling (SIDC) on the price distribution and trading activity. As the SIDC system is dynamically changing and little researched, this might be an interesting direction for further research into the micro-structure of intraday electricity markets. The literature on intraday markets is still scarce compared to the fast growth of renewable energy sources and intraday electricity markets.