# Journal of Energy Markets

**ISSN:**

1756-3607 (print)

1756-3615 (online)

**Editor-in-chief:** Derek W. Bunn

# Optimal intraday power trading with a Gaussian additive process

####
Need to know

- It is in principle possible to exploit the characteristics of intra-day power prices to define a fruitful trading strategy.
- By modeling intra-day prices as mean-reverting processes and maximizing a CRRA utility, we have an explicit optimal strategy, whose risky position increases as maturity approaches.

####
Abstract

Trading activity in intraday (ID) electricity markets has increased significantly over the last few years. We study the problem of a financial agent wishing to maximize a constant relative risk-aversion expected utility of their terminal wealth while operating in an ID market. Assuming that the price of traded hours follows an additive Ornstein– Uhlenbeck process, we derive the optimal strategy via the Hamilton–Jacobi–Bellman equation. The optimal portfolio in the log case is totally myopic with respect to time to maturity, while in the power case it becomes more and more risky as final maturity approaches. In order to implement our strategy, it is necessary to estimate the model parameters. One cannot resort to known results, as it is typical for time series to be unevenly time spaced, with more and more transactions as maturity approaches. Thus, we present an estimation procedure for unevenly spaced observations, based on maximum likelihood estimation and a bootstrap bias correction, in order to compensate for having few observations at the beginning of the observation frame. Finally, we backtest our method and conclude.

####
Introduction

1 Introduction

Trading activity in intraday (ID) electricity markets has increased significantly over the last few years. This is partially due to the growing penetration of nonprogrammable renewable energy sources (such as wind and solar) in the power mix, necessitating ID adjustments with respect to the day-ahead (DA) schedule. However, since ID markets are more liquid than in the past, “purely financial” traders are also beginning to enter into these markets in search of mere financial gains. The aim of this paper is to study exactly how a purely financial trader could invest optimally in a continuous-time ID market.

In more detail, we analyze the case of a financial trader operating in markets organized by the European Power Exchange (EPEX) SPOT market (http://www.epexspot.com). EPEX SPOT is the power market of Central Western Europe (France, Germany, Austria, Switzerland, the Netherlands, Belgium and Luxembourg). Like most power markets, it is mainly divided into two sections. The DA market is the main one: it is auction based, occurs the day before the physical delivery of the energy and is based on twenty-four system marginal price auctions, one for each hour of the following day. The ID market allows agents to change their schedules as defined by the end of the DA market, following on from the DA for the same day. Unlike the DA, the ID is a continuous-time market. For each delivery hour, it starts at 15:00 on the day $D-1$ (where $D$ is the date when delivery takes place) and closes thirty minutes before the beginning of each delivery hour on the day $D$.

By exploiting the continuous-time nature of the ID market, this paper attempts to establish whether it is possible for a “small” financial agent to employ a fruitful pure trading strategy. This means that our agent does not have significant price impact (in contrast with Aïd et al (2016) or Tan and Tankov (2016), for instance) and would not wish to physically buy or sell power at the end of the market session; they would only be aiming to trade in order to make financial gains before the session. A significant example of this kind of behavior is “proprietary trading”, where a producer wants to sell the power that they have produced but is entitled to trade into power markets where they could, in principle, make mere financial gains. Obviously, the stochastic nature of energy price dynamics induces some level of market risk for every trading strategy. For this reason, it is necessary to define a criterion that allows us to keep market risk under control. To do this, we choose the classical criterion of maximizing a utility function of the final profit for our agent. A similar problem of maximizing wealth due to ID trading has been treated in Aïd et al (2016) and Tan and Tankov (2016). These works explicitly model possible price impact; however, they are also risk-neutral with respect to gains in the ID market, reserving risk aversion for positions in the balancing market. By contrast, we assume that our investor trades in a risk-averse way only in the ID market, ending up with a null position in the balancing market.

In order to properly define an optimal strategy, we must first choose a model. We choose a parsimonious yet quite realistic model by assuming that the price of each traded hour evolves as an additive Ornstein–Uhlenbeck (OU) model driven by a Brownian motion. This choice is in line with many models used in the literature for spot (Benth et al 2007, 2012; Meyer-Brandis and Tankov 2008) and ID prices (Aïd et al 2016; Kiesel and Paraschiv 2017; Tan and Tankov 2016). It also acknowledges the possibility of negative prices, observed more and more in power markets in recent years (EPEX 2015).

Inspired by Kiesel and Paraschiv (2017), we choose an OU process as a natural version in continuous time of the model therein. Based upon this model, we present the solution of the optimal portfolio problem: this is done via the standard tools of dynamic programming and the Hamilton–Jacobi–Bellman (HJB) equation. The optimal portfolio in the log case turns out to be totally myopic with respect to time to maturity. In the power case, however, the optimal strategy is that of the log case multiplied by a decreasing function of time to maturity; hence, the investor takes more and more risky positions as the final maturity approaches, as in other models in the literature (Battauz et al 2014; Kim and Omberg 1996). The derivation of this result had to be explicitly carried out in our paper, as the literature only provides solutions for the case of asset prices following geometric models, while our model is arithmetic.

In order to implement the optimal strategy that we find, it is necessary to estimate the model parameters based on the time series of each traded hour. Again, this is not a standard task. In fact, parameter estimations for OU processes can be easily found in cases where observations are equispaced in time. Unfortunately, this is not the case here; in fact, for each traded hour, ID time series prices are most active in the last two to three hours before maturity. During the first hours, few transactions occur. Thus, we are forced to develop a technique for estimating the parameters from nonequispaced observations. We do this via a maximum likelihood (ML) estimation, with a bootstrap bias correction inspired by Tang and Chen (2009) to compensate for the fact that little data may be observed during the first hours. Finally, we backtest our model and estimators, in order to assess whether it is possible to make a financial profit on the ID market.

We present a brief outline of this paper. In Section 2, we define a Gaussian additive model for ID energy prices. In Section 3, we formalize the utility maximization problem and solve it using the dynamic programming technique. In Section 4, we present a procedure for estimating the parameters of our model. Finally, in Section 5 we present the results of a backtest and conclude.

## 2 The ID price dynamics

We now specify a parsimonious, yet quite realistic, model for the price dynamics of traded hours in the ID market.

We start from a data set of high-frequency daily data from January 1 to July 31, 2014 from the EPEX German market. The EPEX German market is one of the most liquid intraday markets working in continuous time. In our data set, a total of 767 651 transactions on single hours are present. Transactions per hourly block or quarter, though present in the EPEX market, have not been included in our database. In this sample, for each hourly contract we observe in mean

- •
150 transactions,

- •
14.58 MWh of exchanged average volume per transaction.

As we can see from Figure 1 in greater detail, peak hours are the most liquidly traded. For example, the maximum average occurs during hour 14, with 242 daily transactions and 16.05 MWh exchanged.

On finer analysis, we note the following important ID energy price features for each traded hour:

- •
liquidity grows during the last two to three hours of transactions (see Figure 2 for typical price paths);

- •
there is high volatility;

- •
there are local trends, ie, trends which appear only in part of the price path (typically in the last two to three hours of transactions);

- •
sometimes it is possible to observe negative prices (EPEX 2015).

A negative price is a signal that occurs when high inflexible power generation meets with either low demand or unforeseen high production from renewables. (According to the Erneuerbare Energien Gesetz (EEG), German electricity network operators are obligated to accept the entire electrical output of renewable plants (see Mitchell et al 2006).) Power generation units are inflexible and cannot be shut down and restarted in a quick, cost-efficient manner (as is the case with renewable energy plants). In this scenario, producers usually deem it more convenient to keep their power plants online and pay an additional cost, rather than wind them down and restart them later.

The possibility of negative prices has also already been acknowledged by more recent spot and ID price models. Benth et al (2007, 2012) and Meyer-Brandis and Tankov (2008) propose detailed models for the spot price of electricity that are additive and driven by general Lévy processes. If the driving process is a subordinator (ie, a process with only positive jumps), then the spot price remains positive. However, due to negative spot prices being observed more and more in power markets over the past few years (EPEX 2015), a choice of driving process that also comprises negative increments could be viable. In fact, the ID price models by Aïd et al (2016) and Tan and Tankov (2016) are both drifted Brownian motions, while Kiesel and Paraschiv (2017) propose an alternative kind that is autoregressive and driven by Gaussian noise, with a regime switch driven by a latent variable linked with demand. However, the technique used by Kiesel and Paraschiv (2017) to estimate this latent variable is quite time consuming and needs the entire price path of the day; of course, this makes the model quite unsuitable for any investor implementing a dynamic trading strategy, as this investor could only observe the price history up to the present time.

For the same kinds of reasons, we did not include any dependence on weather forecast. First, forecasts are not given in continuous time but only at fixed dates in time; thus, their impact upon different traded hours may be different.^{1}^{1}Tan and Tankov (2016) present a partial solution to this, but at the cost of adding one more state variable that would complicate our model. However, the main reason why we did not explicitly include forecasts is because we assume that they are broadcast publicly and known to every market participant; thus their information can be absorbed by the prices’ dynamics as soon as they are revealed. This will be reflected in the iterated estimation procedure in Section 5. As soon as a new price is present in the market, the parameters will be estimated again, thus incorporating the new information. Considering these behaviors for the price of each traded hour, we choose an arithmetic dynamic for ID energy price ${({S}_{t})}_{t}$. In particular, we assume an OU process to describe the ID price dynamics of any given traded hour:

$$\begin{array}{cccc}\hfill \mathrm{d}{S}_{t}& =(\mu -\lambda {S}_{t})\mathrm{d}t+\sigma \mathrm{d}{W}_{t},\hfill & \hfill t& \in (0,T],\hfill \\ \hfill {S}_{0}& ={s}_{0},\hfill & \hfill t& =0.\hfill \end{array}\}$$ | (2.1) |

This specification of an OU process driven by a Brownian motion may seem an excessive simplification. However, it is in line with many models used in the literature for ID prices (see Kiesel and Paraschiv (2017) and the references therein), which are additive and driven by Gaussian noise. Unlike Kiesel and Paraschiv (2017), we choose not to introduce regime switches because the techniques needed to assess a regime change are too computationally intensive and typically require knowledge of the time series in the future, making them hard to use for dynamic trading strategies.

## 3 The optimization problem

To define an optimal strategy for the ID market, it is necessary to deal with the risk arising from the stochastic price movements included in every strategy. In principle, the agent could decide to trade in the ID market in all twenty-four hours (or in another product, such as quarters or blocks of hours). Thus, in full generality, this would require a multidimensional model. However, such a model would be much more complicated to specify and calibrate for the reasons already seen in Section 2.

More specifically, the first reason is that, for each of the twenty-four hours, trading stops thirty minutes before the physical delivery, causing all delivery hours to differ from one another. For example, the ID market for hour 1 starts at 15:00 on the day before the delivery date $D-1$ and ends at 00:30 on the delivery date $D$, while the ID market for the hour 24 starts at 15:00 on day $D-1$ and ends at 23:30 on day $D$. Thus, there is no natural final time horizon for our trader.

The second reason is that, for each traded hour, transactions are concentrated during the last two to three hours prior to delivery. As a consequence, when the price of hour 1 moves, the price of hour 24 stays mostly still; conversely, when the price of hour 24 is moving, the price of hour 1 no longer exists. This makes the estimation of correlations between hourly prices very difficult, as one never observes contemporary changes in price. This could lead to a bias toward zero known as the Epps effect taking place in the covariance estimators (Epps 1979).

For these reasons, we assume that the prices of individual traded hours are independent of one another. This is equivalent, in our formulation, to assuming that the various Brownian motions in (2.1) are independent. Moreover, we make the naive choice of representing our optimal investment problem with a single-asset optimal portfolio problem, in the sense that the agent divides their initial capital between the hours and performs an optimal portfolio problem for each given hour that is independent of all other hours.

Thus, for each hour we frame our problem as a classical utility maximization problem

$$\underset{\pi}{\mathrm{max}}?[U({X}_{T})],$$ | (3.1) |

where the asset’s dynamics are given by (2.1); the portfolio value’s dynamics are

$$\begin{array}{cccc}\hfill \mathrm{d}{X}_{t}& ={\pi}_{t}\mathrm{d}{S}_{t},\hfill & \hfill t& \in (0,T],\hfill \\ \hfill {X}_{0}& ={x}_{0},\hfill & \hfill t& =0;\hfill \end{array}\}$$ | (3.2) |

the trading strategy $\pi ={({\pi}_{t})}_{t\in [0,T]}$ is such that ${\pi}_{t}\in \mathbb{R}$ for all $t\in [0,T)$; and

$${\pi}_{T}=0\mathit{\hspace{1em}}\text{for}t=T$$ | (3.3) |

is imposed. This signifies that we want to neither buy nor sell power at the end of the ID market; instead we are in a pure trading context. Here ${s}_{0}$ is the DA energy price, ${x}_{0}$ is our initial risk capital and $T$ is the closure time of the ID market.

###### Remark 3.1.

Imposing (3.3) requires that the market is liquidly accepting all positions at all times $t$, including $t=T$. Should this condition not be satisfied, the trader would pay a cost in the form of a penalty in the so-called balancing market. This would also be the case for a physical producer wanting to impose ${\mathrm{\Pi}}_{T}$ equal to a certain nonnegative quantity sold by them. Should this quantity come from a nonprogrammable renewable source (such as solar or wind), there would be a positive probability that, even if the financial position ${\pi}_{T}$ is met, the final quantity produced would be different from that. In this case, the producer would also end up with a position in the balancing market. We leave this more structured problem for future research. From now on we assume that (3.3) can be satisfied, ie, the market is perfectly liquid at all times $t\le T$.

We may take into account the risk associated with every ID market strategy with the use of a utility function $U$. Further, we consider the classical choice of a constant relative risk-aversion (CRRA) utility function; this could be either the power utility function $U(x)={x}^{\gamma}/\gamma $, with $\gamma \in (0,1)$, or the logarithmic utility function $U(x)=\mathrm{log}(x)$, which we will indicate conventionally as “the case $\gamma =0$”.

###### Definition 3.2 (Admissible controls, objective function and value function).

We denote by $?([t,T])$ the set of admissible controls, ie, predictable processes $\pi ={({\pi}_{u})}_{u}$ on $[t,T]$, such that (3.2) has a unique strong solution ${X}^{t,x,s;\pi}$ for each initial condition $(x,s)$ at time $t$, and ${\mathrm{\Pi}}_{T}=0$.

We call objective function $J$ the function defined by

$$J(t,x,s,\pi )=?[U({X}_{T}^{t,x,s;\pi})],$$ | (3.4) |

where $({X}^{t,x,s;\pi},{S}^{t,x,s;\pi})$ is the two-dimensional controlled Markov process starting from $(x,s)$ at time $t$ with the dynamics defined by (2.1) and (3.2) with the control $\pi \in ?([t,T])$.

Finally, we call value function $V$ the function defined by

$$V(t,x,s)=\underset{\pi \in ?([t,T])}{sup}J(t,x,s,\pi ).$$ | (3.5) |

Now we can write the HJB equation associated with our optimization problem as

${V}_{t}+\underset{\pi}{sup}{A}^{\pi}V$ | $=0,$ | $\mathrm{\hspace{1em}}(t,s,x)$ | $\in [0,T)\times \mathbb{R}\times \mathbb{R},$ | (3.6) | ||

$V(T,s,x)$ | $=U(x),$ | $\mathrm{\hspace{1em}}(s,x)$ | $\in \mathbb{R}\times \mathbb{R},$ | (3.7) |

with the infinitesimal generator ${A}^{\pi}$ given by

${A}^{\pi}V(t,x,s)$ | $=({V}_{x}(t,x,s)\pi +{V}_{s}(t,x,s))(\mu -\lambda s)+\frac{1}{2}{\pi}^{2}{\sigma}^{2}{V}_{xx}(t,x,s)$ | |||

$\mathrm{\hspace{1em}\hspace{1em}}+\pi {\sigma}^{2}{V}_{sx}(t,x,s)+\frac{1}{2}{\sigma}^{2}{V}_{ss}(t,x,s).$ | (3.8) |

First-order conditions on (3.6) to (3.8) imply that a candidate optimal control $\pi $ has the particular form

$$\widehat{\pi}:=\mathrm{arg}\underset{\pi}{\mathrm{max}}{A}^{\pi}V=-\frac{{V}_{x}(\mu -\lambda s)+{\sigma}^{2}{V}_{sx}}{{\sigma}^{2}{V}_{xx}},$$ | (3.9) |

and the HJB equation (3.6) can be written as

$${V}_{t}+{V}_{s}(\mu -\lambda s)+\frac{1}{2}{\sigma}^{2}{V}_{ss}-\frac{1}{2}{\widehat{\pi}}^{2}{\sigma}^{2}{V}_{xx}=0,(t,s,x)\in [0,T)\times \mathbb{R}\times \mathbb{R}.$$ | (3.10) |

The usual procedure for solving our optimization problem is to guess at a particular solution to the HJB equation and for a given utility function $U$. We then test whether this particular candidate satisfies the hypothesis of a so-called verification theorem. To carry out this plan, as previously stated, we consider the case of CRRA utility functions, which can be either the power utility function $U(x)={x}^{\gamma}/\gamma $, with $\gamma \in (0,1)$, or the logarithmic utility function $U(x)=\mathrm{log}(x)$.

###### Lemma 3.3 (CRRA utility function).

Let $U\mathit{}\mathrm{(}x\mathrm{)}\mathrm{=}{x}^{\gamma}\mathrm{/}\gamma $, with $$, $\gamma \mathrm{\ne}\mathrm{0}$, or $U\mathit{}\mathrm{(}x\mathrm{)}\mathrm{=}\mathrm{log}\mathit{}\mathrm{(}x\mathrm{)}$. Then the following holds.

- (1)
The system of ordinary differential equations

${a}_{2}^{\prime}$ $=-{\displaystyle \frac{2\gamma {\sigma}^{2}}{1-\gamma}}{a}_{2}^{2}+{\displaystyle \frac{2\lambda}{1-\gamma}}{a}_{2}-{\displaystyle \frac{{\lambda}^{2}}{2(1-\gamma ){\sigma}^{2}}},$ (3.11) ${a}_{1}^{\prime}$ $={a}_{1}\left({\displaystyle \frac{\lambda}{1-\gamma}}-{\displaystyle \frac{2\gamma {\sigma}^{2}}{1-\gamma}}{a}_{2}\right)-{\displaystyle \frac{2\mu {a}_{2}}{1-\gamma}}+{\displaystyle \frac{\mu \lambda}{(1-\gamma ){\sigma}^{2}}},$ (3.12) ${a}_{0}^{\prime}$ $=-{\displaystyle \frac{\mu}{1-\gamma}}{a}_{1}-{\sigma}^{2}{a}_{2}-{\displaystyle \frac{\gamma {\sigma}^{2}}{2(1-\gamma )}}{a}_{1}^{2}-{\displaystyle \frac{{\mu}^{2}}{2(1-\gamma ){\sigma}^{2}}},$ (3.13) with terminal conditions ${a}_{2}(T)={a}_{1}(T)={a}_{0}(T)=0$ , admits a unique solution for $$ , such that

${a}_{2}(t)$ $={\displaystyle \frac{\lambda}{2{\sigma}^{2}}}{\displaystyle \frac{\mathrm{sinh}m(T-t)}{\mathrm{sinh}m(T-t)+\sqrt{1-\gamma}\mathrm{cosh}m(T-t)}}$ (3.14) $={\displaystyle \frac{\lambda}{2{\sigma}^{2}}}{\displaystyle \frac{{\mathrm{e}}^{2m(T-t)}-1}{(1+\sqrt{1-\gamma}){\mathrm{e}}^{2m(T-t)}+\sqrt{1-\gamma}-1}},$ ${a}_{1}(t)$ $={\displaystyle \frac{\mu}{{\sigma}^{2}}}{\displaystyle \frac{1-{\mathrm{e}}^{2m(T-t)}}{(1+\sqrt{1-\gamma}){\mathrm{e}}^{2m(T-t)}+\sqrt{1-\gamma}-1}}=-{\displaystyle \frac{2\mu}{\lambda}}{a}_{2}(t),$ (3.15) with $m:=\lambda /\sqrt{1-\gamma}$ . Moreover, ${a}_{2}\ge 0$ for all $t\le T$ .

- (2)
- (3)
The portfolio strategy ${({\pi}_{u}^{*})}_{u\ge t}$ such that ${\pi}_{u}^{*}:=\widehat{\pi}(u,{S}_{u},{X}_{u})$ for $u\in [t,T)$ and ${\pi}_{T}^{*}=0$ , with

$$\widehat{\pi}(t,s,x)=\frac{(\mu -\lambda s)x}{(1-\gamma ){\sigma}^{2}}+\frac{\gamma}{1-\gamma}{g}_{s}(t,s)x=\frac{(\mu -\lambda s)}{(1-\gamma ){\sigma}^{2}}xh(T-t)$$ (3.16) with

$$h(u)=\frac{1+\sqrt{1-\gamma}\mathrm{tanh}mu}{1-\gamma +\sqrt{1-\gamma}\mathrm{tanh}mu}$$ (3.17) defines an admissible control in the sense of Definition 3.2 .

###### Proof.

See online Appendix A. ∎

We now present a verification theorem summarizing the different choices for $U$ (log or power case).

###### Theorem 3.4 (Verification theorem).

Let $U\mathit{}\mathrm{(}x\mathrm{)}$ be equal either to $\mathrm{log}\mathit{}\mathrm{(}x\mathrm{)}$, letting $\gamma \mathrm{=}\mathrm{0}$, or to $U\mathit{}\mathrm{(}x\mathrm{)}\mathrm{=}{x}^{\gamma}\mathrm{/}\gamma $, with $\gamma \mathrm{\in}\mathrm{(}\mathrm{0}\mathrm{,}\mathrm{1}\mathrm{)}$. Further, let $H$ be as in Lemma 3.3(ii). Then, $H\mathrm{=}V$, and the portfolio strategy ${\mathrm{(}{\pi}_{u}^{\mathrm{*}}\mathrm{)}}_{u\mathrm{\ge}t}$, such that ${\pi}_{u}^{\mathrm{*}}\mathrm{:=}\widehat{\pi}\mathit{}\mathrm{(}u\mathrm{,}{S}_{u}\mathrm{,}{X}_{u}\mathrm{)}$ for $u\mathrm{\in}\mathrm{[}t\mathrm{,}T\mathrm{)}$ and ${\pi}_{T}^{\mathrm{*}}\mathrm{=}\mathrm{0}$, defines an admissible optimal control, with $\widehat{\pi}$ defined by (3.16) and (3.17).

###### Proof.

What remains to be proved is that the value function satisfies good integrability conditions such that, when applying the Ito formula to the process ${(H(u,{X}_{,}{S}_{u}))}_{u\in [t,T]}$, the stochastic integral appearing in its stochastic differential is a proper martingale. For the logarithmic case, this is quite straightforward, and the final result follows from standard versions of the verification theorem (see, for example, Bjork 2003, Theorem 14.6). The power case is more complicated, but it can be treated with techniques similar to those of Benth and Karlsen (2005), Propositions 4.1 and 4.2. The proof of the latter case is long and quite technical, and we do not provide it here. However, interested readers may obtain it from the authors upon request. ∎

###### Remark 3.5.

The optimal control turns out to be proportional to the optimal one for the log case, which is totally myopic with respect to $t$, multiplied by the proportionality factor $h(T-t)$ in (3.17). This factor tends uniformly to 1 as $\gamma \to 0$, increases when $\gamma $ increases (ie, when the risk-aversion coefficient $1-\gamma $ decreases) and explodes to $+\mathrm{\infty}$ as $\gamma \to 1$. This means that the exposure associated with the risky asset increases as the risk-aversion coefficient $1-\gamma $ decreases. Moreover, this factor depends on $t$ only in the time leading up to maturity $u=T-t$, as $t\to T$ increases up to the finite limit

$$\underset{t\to T}{lim}h(T-t)=\underset{u\to 0}{lim}h(u)\frac{1}{1-\gamma}$$ |

for all $\gamma \in (0,1)$. A sample of this behavior can be seen in Figure 3. This means that the investor is more and more nonmyopic with respect to time as $\gamma $ increases. For this reason, investors tend to assume more extreme positions as time $t$ approaches maturity $T$. However, these positions do not explode as in Battauz et al (2014); instead they tend to a Merton-like position (Sharpe ratio divided by $1-\gamma $). This phenomenon has already been observed in other models (see, for example, Battauz et al 2014; Kim and Omberg 1996).

## 4 Parameter estimation

We now present a procedure for estimating the parameters of our model, in order to implement the optimal portfolio strategy found in Section 3. As seen in Section 2, a typical time series for the price of each traded hour will have unevenly spaced observations in time, since transactions become more and more frequent as the delivery time approaches. Indeed, most transactions happen in the last two to three hours. For this reason, we cannot rely upon known parameter estimation techniques designed to work on evenly spaced observations. Instead, we must develop an estimation procedure for observations that are unevenly spaced in time.

Let us suppose that we have $n+1$ unevenly time-spaced observations ${\{{s}_{{t}_{k}}\}}_{k=0,\mathrm{\dots},n}$ of the process $S$ at times ${t}_{0}=0,{t}_{1},\mathrm{\dots},{t}_{n}$. We use the notation ${\delta}_{k}={t}_{k}-{t}_{k-1}$ and ${s}_{k}:={s}_{{t}_{k}}$ for all $k=1,\mathrm{\dots},n$. To estimate the parameters of the OU process we need to produce maximum likelihood estimators (MLEs). For this purpose, first of all, we have to write the loglikelihood function of the process (2.1). It is easy to see that the loglikelihood function of the unevenly spaced observations of the OU process (2.1) at times ${({t}_{k})}_{k=0,\mathrm{\dots},n}$ is

$\mathcal{L}(\lambda ,\mu ,{\sigma}^{2};{s}_{0},\mathrm{\dots},{s}_{n})$ | $=-{\displaystyle \frac{n}{2}}\mathrm{log}(2\pi {\sigma}^{2})-{\displaystyle \frac{1}{2}}{\displaystyle \sum _{k=1}^{n}}\mathrm{log}\left({\displaystyle \frac{(1-{\mathrm{e}}^{-2\lambda {\delta}_{k}})}{\lambda}}\right)$ | |||

$\mathrm{\hspace{1em}\hspace{1em}}-{\displaystyle \frac{\lambda}{{\sigma}^{2}}}{\displaystyle \sum _{k=1}^{n}}{\displaystyle \frac{{(({s}_{k}-\mu /\lambda )-({s}_{k-1}-\mu /\lambda ){\mathrm{e}}^{-\lambda {\delta}_{k}})}^{2}}{1-{\mathrm{e}}^{-2\lambda {\delta}_{k}}}}$ | (4.1) |

(see also Franco 2003). To obtain MLEs for the parameters of the OU process, we solve the system

$$\nabla \mathcal{L}=0,$$ | (4.2) |

where $\mathcal{L}$ is given by (4.1). In this case, it turns out to be impossible to solve (4.2) in closed form (as we would with evenly time-spaced observations). Thus, it is not possible to obtain closed formulas for the MLEs of $\lambda $, $\mu $ and ${\sigma}^{2}$. Instead, we apply the method outlined in Franco (2003), which is similar to the one used by Edoli et al (2013). First, we express $\mu $ and $\sigma $ as functions of $\lambda $, for $\lambda >0$, as

$\widehat{\mu}=\mu (\lambda )$ | $=\lambda \left({\displaystyle \sum _{i}}{\displaystyle \frac{{s}_{i}-{s}_{i-1}{\mathrm{e}}^{-\lambda {\delta}_{i}}}{1+{\mathrm{e}}^{-\lambda {\delta}_{i}}}}\right){\left({\displaystyle \sum _{i}}{\displaystyle \frac{1-{\mathrm{e}}^{-\lambda {\delta}_{i}}}{1+{\mathrm{e}}^{-\lambda {\delta}_{i}}}}\right)}^{-1},$ | ||

${\widehat{\sigma}}^{2}={\sigma}^{2}(\lambda ,\mu (\lambda ))$ | $={\displaystyle \frac{2\lambda}{n}}{\displaystyle \sum _{i}}{\displaystyle \frac{{({s}_{i}-\mu (\lambda )/\lambda -({s}_{i-1}-\mu (\lambda )/\lambda ){\mathrm{e}}^{-\lambda {\delta}_{i}})}^{2}}{1-{\mathrm{e}}^{-2\lambda {\delta}_{i}}}}$ |

and then we perform a numerical maximization of $\mathcal{L}$ in the unique real variable $\lambda $. This is a much faster and more efficient technique for numerically maximizing the three real variables $\lambda $, $\mu $ and ${\sigma}^{2}$. Their MLEs are

$\widehat{\lambda}$ | $=\underset{\lambda >0}{\mathrm{max}}\mathcal{L}(\lambda ,\mu (\lambda ),{\sigma}^{2}(\lambda ,\mu (\lambda ));{s}_{0},\mathrm{\dots},{s}_{n}),$ | ||

$\widehat{\mu}$ | $=\mu (\widehat{\lambda}),$ | ||

${\widehat{\sigma}}^{2}$ | $={\sigma}^{2}(\widehat{\lambda},\mu (\widehat{\lambda})).$ |

If, as we saw, $\lambda =0$, the process (2.1) reduces to a simple arithmetic Brownian motion. In this case we can write two closed formulas for the MLEs of $\mu $ and ${\sigma}^{2}$, namely

$\widehat{\mu}$ | $={\displaystyle \frac{{S}_{T}-{S}_{0}}{T}},$ | (4.3) | ||

${\widehat{\sigma}}^{2}$ | $={\displaystyle \frac{1}{n}}{\displaystyle \sum _{i=1}^{n}}{\displaystyle \frac{{({S}_{i}-{S}_{i-1}-\widehat{\mu}{\delta}_{i})}^{2}}{{\delta}_{i}}}.$ | (4.4) |

At this point, we can define the following estimation procedure for the process (2.1).

- (1)
Solve the maximization problem described above and obtain the values of $\widehat{\lambda}$, $\widehat{\mu}$ and ${\widehat{\sigma}}^{2}$.

- (2)
Perform a hypothesis test on $\widehat{\lambda}$ to assess whether it is statistically different from zero. If this is the case, then the calibration procedure is complete; if it is not, we progress to (3).

- (3)
- (4)
If $\widehat{\mu}$ is statistically different from zero, then the procedure is complete; otherwise, define $\widehat{\sigma}$ using (4.4) with $\mu =0$.

An important problem (that we need to solve) concerns the bias of the MLEs for parameters ${\theta}_{\mathrm{d}}$ of the drift function of a general diffusion process such as

$$\mathrm{d}{X}_{t}=\mu ({X}_{t},{\theta}_{\mathrm{d}})\mathrm{d}t+\sigma ({X}_{t},{\theta}_{\mathrm{v}})\mathrm{d}{W}_{t}.$$ | (4.5) |

Indeed, the bias of the ML estimator ${\widehat{\theta}}_{\mathrm{d}}$ of $\theta $ is caused by the discretization of the continuous-time process and becomes particularly pronounced when we have few observed data in a given time interval (see, for example, Bonest and Yu 2006; Franco 2003). To solve this problem, we follow Tang and Chen (2009) and use a bootstrap bias correction. First of all, we recall the following.

###### Definition 4.1.

We define the bias function of an estimator $\widehat{\theta}$ of a parameter $\theta $ as

$$b(\widehat{\theta},\theta )={?}_{\theta}[\widehat{\theta}]-\theta .$$ |

The idea of the bootstrap bias correction procedure is to obtain an evaluation $\stackrel{~}{b}(\widehat{\theta},\theta )$ of the bias function $b(\widehat{\theta},\theta )$ using Monte Carlo techniques. We may then define the unbiased new estimator of $\theta $ as

$$\stackrel{~}{\theta}=\widehat{\theta}-\stackrel{~}{b}(\widehat{\theta},\theta ).$$ |

Suppose that we observe the process (4.5) at times $$, and let us call ${x}_{t}$ the value of ${X}_{t}$ observed at time $t$. The procedure is as follows.

- (1)
Find the MLEs, ${\widehat{\theta}}_{\mathrm{d}}$ of ${\theta}_{\mathrm{d}}$ and ${\widehat{\theta}}_{\mathrm{v}}$ of ${\theta}_{\mathrm{v}}$, using the observed path ${\{{x}_{t}\}}_{t}$.

- (2)
Simulate $M$ paths ${\{{x}^{i}\}}_{i=1,\mathrm{\dots},M}$ using the dynamics

$$\mathrm{d}{X}_{t}=\mu ({X}_{t},{\widehat{\theta}}_{\mathrm{d}})\mathrm{d}t+\sigma ({X}_{t},{\widehat{\theta}}_{\mathrm{v}})\mathrm{d}{W}_{t}$$ as well as the same instants $$ for the time grid of these simulations.

- (3)
For each path in (2), find the MLE ${\widehat{\mathrm{\Theta}}}_{\mathrm{d}}^{i}$ of ${\theta}_{\mathrm{d}}$.

- (4)
Define the new correct estimator of ${\theta}_{\mathrm{d}}$ as

$${\stackrel{~}{\theta}}_{\mathrm{d}}=2{\widehat{\theta}}_{\mathrm{d}}-\frac{{\sum}_{j=1}^{M}{\widehat{\mathrm{\Theta}}}_{\mathrm{d}}^{j}}{M}.$$

For the OU process we apply this procedure to the parameters $\lambda $ and $\mu $ only, since the MLE for ${\sigma}^{2}$ is already unbiased (see Tang and Chen 2009). It is also worth noting that the bootstrap procedure described above can be iterated several times.

To test this procedure, we simulate $N=1000$ paths with three sets of known parameters. We then estimate the parameters from the simulated data. Table 1 compares the results of different methods of numerical estimation: without the bootstrap procedure (NL); with the bootstrap procedure applied only once (NB); and with the bootstrap procedure iterated twice (N2B). In this test, we define $\theta :=\mu /\lambda $, that is, the long-term mean of the OU process.

(a) NL | |||||||||

$?$ | $?$ | $?$ | $?$ | $?$ | $?$ | $?$ | $?$ | $?$ | |

Real | 0.851 | 8.510 | 10.000 | 0.251 | 2.510 | 10.000 | 0.080 | 0.800 | 10.000 |

Mean | 0.233 | 2.345 | 10.067 | 0.064 | 0.660 | 10.324 | 0.030 | 0.324 | 10.750 |

MAE | 0.618 | 6.165 | 0.067 | 0.187 | 1.850 | 0.324 | 0.050 | 0.476 | 0.750 |

MAPE (%) | 73 | 72 | 0.7 | 74 | 74 | 3.2 | 62 | 60 | 7.5 |

RMSE | 0.382 | 38.003 | 0.005 | 0.035 | 3.421 | 0.105 | 0.002 | 0.227 | 0.562 |

(b) NB | |||||||||

$?$ | $?$ | $?$ | $?$ | $?$ | $?$ | $?$ | $?$ | $?$ | |

Real | 0.851 | 8.510 | 10.000 | 0.251 | 2.510 | 10.000 | 0.080 | 0.800 | 10.000 |

Mean | 0.420 | 4.200 | 10.006 | 0.117 | 1.180 | 10.059 | 0.053 | 0.540 | 10.144 |

MAE | 0.431 | 4.310 | 0.006 | 0.134 | 1.330 | 0.059 | 0.027 | 0.260 | 0.144 |

MAPE (%) | 51 | 51 | 0 | 53 | 53 | 0.6 | 33 | 32 | 1.4 |

RMSE | 0.186 | 18.578 | 0.000 | 0.018 | 1.768 | 0.003 | 0.001 | 0.068 | 0.021 |

(c) N2B | |||||||||

$?$ | $?$ | $?$ | $?$ | $?$ | $?$ | $?$ | $?$ | $?$ | |

Real | 0.851 | 8.510 | 10.000 | 0.251 | 2.510 | 10.000 | 0.080 | 0.800 | 10.000 |

Mean | 0.733 | 7.302 | 9.966 | 0.214 | 2.131 | 9.942 | 0.094 | 0.931 | 9.854 |

MAE | 0.118 | 1.208 | 0.034 | 0.037 | 0.379 | 0.058 | 0.014 | 0.131 | 0.146 |

MAPE (%) | 14 | 14 | 0.3 | 15 | 15 | 0.6 | 18 | 16 | 1.5 |

RMSE | 0.014 | 1.459 | 0.001 | 0.001 | 0.144 | 0.003 | 0.000 | 0.017 | 0.021 |

It is interesting to note that, within Table 1, the estimation of the parameter $\theta $ is very precise from the beginning, even without the bootstrap bias correction procedure. However, the numerical loglikelihood estimation of both $\lambda $ and $\mu $ improves significantly with the application of the iterated bootstrap bias correction procedure.

## 5 Backtest and results

We perform a backtest to evaluate the effectiveness of our solutions to the optimization problem as well as our calibration process for the various classes of CRRA utility function used therein. We compare the results obtained with the utility functions for different values of the risk parameter. Each day, we start our exercise with a risk capital of ${x}_{0}=$ €1000, shared equally among the hours in the day, so that each hour has a starting capital of €41.66. This can easily be rescaled to suit any initial amount by adjusting the properties of CRRA utility functions. This simply entails multiplying the optimal strategies by the correct factor, taking into account that the liquidity of the ID market is limited, as shown in Figure 1. Moreover, we need to make the following simplifying hypotheses.

- •
Market liquidity (ie, every offer of the optimal quantity $\widehat{\pi}\in \mathbb{R}$) will always be accepted; in reality, the minimum tradable quantity is $0.1$ MWh. This gives quite a fine granularity, and short sellings are allowed, provided we end with a nonnegative position at time $T$.

- •
Our financial guarantees always allow us to implement our optimal strategy.

- •
Transaction costs do not apply.

- •
When we reach time ${t}_{k}$ and observe the price ${s}_{k}$, we may buy or sell the quantity ${\widehat{\pi}}_{k}$ at the same price ${s}_{k}$ exactly at that time. This, in some sense, corresponds to the fact that our agent is “small” and has a negligible market impact.

Finally, before the backtest, we replace simultaneous transactions with a volume-weighted average price. We consider a three-month time period (from March to May in 2014) and proceed as follows, for each trading day and traded hour.

- (1)
We observe forty transactions:

- –
should we observe fewer than forty transactions, we take no position;

- –
once we have observed at least forty transactions, we go ahead to point (2).

- –
- (2)
We apply the calibration process described in Section 2. To estimate the parameters, we always use the last forty observed prices combined with a double bootstrap procedure.

- (3)
We calculate the optimal strategy $\widehat{\pi}$ (see (3.16)):

- –
if $\widehat{\pi}$ is statistically different from zero, we implement the optimal strategy;

- –
if $\widehat{\pi}$ is statistically equal to zero, we take no position but go ahead to (5).

- –
- (4)
If $\gamma =0$ (as with log utility), we go directly to the next point; otherwise, in the case where no transactions occur within one minute, we restart the bidding procedure as in (2).

- (5)
When we observe a new transaction, we restart the bidding procedure from (2).

In the above list, (4) is peculiar to power utility functions. As noted in Remark 3.5, it is the case with power utility that the optimal position $\widehat{\pi}$ is influenced by time to maturity $T-t$, while in the case of logarithmic utility the optimal strategy is myopic with respect to time. We add an upper bound to our strategies to ensure that the hypothesis of not being market maker is fulfilled. In particular, all our strategies are less than or equal to $50$ MWh. We define the performance of hour $h$ as ${X}_{T}^{h}-{x}_{0}^{h}$, where ${X}_{T}^{h}$ is the portfolio’s value for hour $h$, the final hour in the trading day, and ${x}_{0}^{h}$ is the hourly starting risk capital. Finally, we compare the results obtained by the different utility functions. Figure 4 shows the results of our backtest for several CRRA utility functions, with different values of the parameter $\gamma $.

We note that both the average performance and the average standard deviation grow along with $\gamma $, ie, as risk aversion $1-\gamma $ decreases. We pass from the case of the logarithmic utility function (a), with an average performance of €3.51 across all twenty-four-hour and single-hour performances ranging from $-$€2.01 to €17.17, to the case of $\gamma =0.8$ (e), which shows an average performance of €20.53 and single-hour performances ranging from $-$€8.76 to more than €65.28. Standard deviation scales vary accordingly; in particular, they increase along with the parameter $\gamma $. In all utility functions, the negative performances recorded for a few hours are counterbalanced by much higher positive performances in various other hours.

Further, since we do not take transaction costs into account in this work, we define a simple indicator of the number of transactions made by our trading algorithm. To this purpose, we examine the ratio of the number of bids made to the number of ticks observed; we call this the “activity level”. Figure 5 shows the activity level of our bidding procedure during the backtesting period.

The average activity level during the backtest period is $42\%$. This is obviously very high, meaning that our activities may not be profitable in the presence of transaction costs. However, this also means that our bidding process may lead to some interesting trading chances in the ID power market.

## 6 Conclusions

This work shows that, in principle, it is possible to exploit the characteristics of ID power prices to define a fruitful trading strategy. We modeled ID prices as additive Gaussian mean-reverting processes (capable of assuming the negative values increasingly observed in today’s markets) and formulated the portfolio selection problem by maximizing a CRRA utility function. Our problem turned out to have an explicit solution, which has the structure of the optimal strategy for a logarithmic utility, multiplied by a correction factor that increases with time to maturity, up to a finite final value. The parameter estimation needed to implement the model is nontrivial and relies on a numerical technique based on MLE. This technique gives the mean-reversion level as well as the volatility as functions of the mean-reversion speed, to be estimated numerically. We further improve the numerical results using a bootstrap bias correction technique. Finally, we apply this technique to a data set from January 1 to July 31, 2014, taken from Germany’s EPEX market. We find that, on average, our optimal strategies give a significant gain, resulting in both the mean performance and the standard deviation increasing with respect to the risk-aversion coefficient. Negative performances are present in certain trading hours, but they are counterbalanced by much more positive performances in other hours.

Future work can be developed in three significant directions. First, there is the question of the correlations between various traded hours. These are very difficult to estimate with standard techniques, since very few hours are traded simultaneously. Second, there is the prospect of realistic modeling of bid–ask spreads. In this way, we might assess whether fruitful ID trading strategies could be built for more realistic situations. The third and final direction is the development of ID market models in which prices could possibly exhibit jumps, even in a time-inhomogeneous way, with more and more activity observed as maturity approaches. This is the topic of a work in progress by one of the authors (Piccirilli and Vargiolu 2017).

## Declaration of interest

The authors report no conflicts of interest. The authors alone are responsible for the content and writing of the paper.

## Acknowledgements

This work was partly supported by the grant CPDA158845-2015 of the University of Padova “Multidimensional polynomial processes and applications to new challenges in mathematical finance and in energy markets”. Moreover, the authors wish to thank Roberto Baviera, Stefano Fiorenzani, Paolo Foschi, Juri Hinz, Rüdiger Kiesel, Alfred Müller, Florentina Paraschiv and the participants of the conferences Energy Finance Italia in Camerino (2015), Workshop in Quantitative Finance in Pisa (2016), Energy and Commodity Finance in Paris (2016), and the XL AMASES conference in Catania (2016), together with anonymous referees, for valuable comments.

## References

Aïd, R., Gruet, P., and Pham, H. (2016). An optimal trading problem in intraday electricity markets. Mathematics and Financial Economics 10(1), 49–85.

Battauz, A., De Donno, M., and Sbuelz, A. (2014). Reaching Nirvana with a defaultable asset? Working Paper, Bocconi University.

Benth, F. E., and Karlsen, K. H. (2005). A note on Merton’s portfolio selection problem for the Schwartz mean-reversion model. Stochastic Analysis and Applications 23(4), 687–704.

Benth, F. E., Kallsen, J., and Meyer-Brandis, T. (2007). A non-Gaussian Ornstein–Uhlenbeck process for electricity spot price modelling and derivatives pricing. Applied Mathematical Finance 14(2), 153–169.

Benth, F. E., Kiesel, R., and Nazarova, A. (2012). A critical empirical study of three electricity spot price models. Energy Economics 34(5), 1589–1616.

Bjork, T. (2003). Arbitrage Theory in Continuous Time, 2nd edn. Oxford Univesity Press.

Bonest, P. P. C., and Yu, J. (2006). Maximum likelihood and Gaussian estimation of continuous time models in finance. In Handbook of Financial Time Series, Mikosch, T., Kreiß, J.-P., Davis, R. A., and Andersen, T. G. (eds), pp. 497–530. Springer.

Edoli, E., Tasinato, D., and Vargiolu, T. (2013). Calibration of a multifactor model for the forward markets of several commodities. Optimization 62(11), 1553–1574.

European Power Exchange (2015). Negative prices Q&A: how they occur, what they mean. Report, EPEX SPOT. URL: http://bit.ly/2A0Wqgw.

Epps, T. (1979). Comovements in stock prices in the very short run. Journal of the American Statistical Association 74(366), 291–298.

Franco, J. C. G. (2003). Maximum likelihood estimation of mean reverting processes. Real Option Practice. Working Paper, Onward, Inc. URL: http://bit.ly/2A1ItSF.

Kiesel, R., and Paraschiv, F. (2017). Econometric analysis of 15-minute intraday electricity prices. Energy Economics 64, 77–90 (http://doi.org/gbns77).

Kim, T. S., and Omberg, E.(1996). Dynamic nonmyopic portfolio behavior. Review of Financial Studies 9(1), 141–161.

Meyer-Brandis, T., and Tankov, P. (2008). Multi-factor jump diffusion models of electricity prices. International Journal of Theoretical and Applied Finance 11(5), 503–528.

Mitchell, C., Bauknecht, D., and Connor, P. M. (2006). Effectiveness through risk reduction: a comparison of the renewable obligation in England and Wales and the feed-in system in Germany. Energy Policy 34(3), 297–305.

Piccirilli, M., and Vargiolu, T. (2017). Optimal portfolio in intraday electricity markets modelled by Lèvy-Ornstein–Uhlenbeck processes. In preparation.

Tan, Z., and Tankov, P. (2016). Optimal trading policies for wind energy producer. Working Paper, Université Paris-Diderot.

Tang, C. Y., and Chen, S. X. (2009). Parameter estimation and bias correction for diffusion processes. Journal of Econometrics 149(1), 65–81.

Only users who have a paid subscription or are part of a corporate subscription are able to print or copy content.

To access these options, along with all other subscription benefits, please contact info@risk.net or view our subscription options here: http://subscriptions.risk.net/subscribe

You are currently unable to print this content. Please contact info@risk.net to find out more.

You are currently unable to copy this content. Please contact info@risk.net to find out more.

Copyright Infopro Digital Limited. All rights reserved.

As outlined in our terms and conditions, https://www.infopro-digital.com/terms-and-conditions/subscriptions/ (point 2.4), printing is limited to a single copy.

If you would like to purchase additional rights please email info@risk.net

Copyright Infopro Digital Limited. All rights reserved.

You may share this content using our article tools. As outlined in our terms and conditions, https://www.infopro-digital.com/terms-and-conditions/subscriptions/ (clause 2.4), an Authorised User may only make one copy of the materials for their own personal use. You must also comply with the restrictions in clause 2.5.

If you would like to purchase additional rights please email info@risk.net