# An Ambit Stochastic Approach to Pricing Electricity Forward Contracts: The Case of the German Energy Market

Authors: Luca Di Persio; Isacco Perin
Journal: Journal of Probability and Statistics
 (2015)

---

## Abstract

We propose an ambit stochastic model to study the electricity forward prices. We provide a detailed analysis of the probabilistic properties of such model, discussing the related martingale conditions and deriving concrete implementation of it for the related underlying spot price. The latter is obtained from the forward model through a limiting argument. Furthermore, we show, also providing a concrete example, that a proper specification of these models is able to effectively forecast prices of forward contracts written on the European Energy Exchange (EEX) AG, or German Energy Exchange, market.


---

## Body

## 1. Introduction

In the last two decades, energy markets have been liberalised in many areas in the world and this has led to the creation of completely new markets as in the case of the Nordic Nord Pool market, the European Energy Exchange (EEX) market, the Italian GME, where GME stands forGestore dei Mercati Energetici, market, and so forth. The latter phenomenon has been underlined by an increasing interest in trading within such commodity frameworks, particularly with respect to electricity, oil, gas, and coal exchange. There is no doubt that such markets will play a vital role in the future given the constant expansion of global demand for energy. From the financial point of view, standard products traded on energy markets are spot prices, forward and futures contracts, and options written on them.In the present work, we focus our attention on economics aspect linked to the electricity commodity. Our choice is mainly due to the fact that electricity cannot be stored; therefore, it has some peculiar characteristics if compared to other commodities like oil or gas. In particular, as soon as electricity has been produced, it has to be delivered to the grid. This leads to very atypical price patterns, such as large price spikes, strong short-term volatility, and negative electricity prices for short, but not zero, periods of time. The latter is the case if the supply considerably exceeds the demand. A previous type of price behaviours suggests both practitioners and academicians to developad hoc techniques for this particular framework, hence trying to get finer approaches to electricity trading than for traditional financial markets.Recent literature has studied, in particular, a class of stochastic models for the electricity price dynamics based onambit processes andambit fields (see, e.g., [1, 2]). Ambit processes are defined as stochastic integrals with respect to a multivariate random measure, where the integrand is given by a product of a deterministickernel function and astochastic volatility field and the integration is carried out over a specific region in space-time, calledambit set. This kind of stochastic process was first developed in the study of turbulence, but, because of their flexible structure, they have been applied in a heterogeneous set of areas (see, e.g., [3, 4]) and, between them, in the financial framework to model dynamic processes (see, e.g., [5]).The wide applicability of ambit processes is due to their flexibility and analytical tractability. As an example, if we consider the energy markets, ambit processes easily incorporateleptokurtic behaviour in returns, stochastic volatility,leverage effects, and the observedSamuelson effect in the volatility, namely, the fact that when the time to maturity approaches zero, the volatility of the forward price increases and converges to the volatility of the underlying spot price, provided the forward price converges to the spot price.In what follows, we develop a class of models for electricity forward prices which will be based on ambit processes and ambit fields. In addition, we show that a proper specification of these models efficiently forecasts the price of the Germanmonthly peak forward contracts within the EEX market scenario.Although many stylized facts concerning energy markets can be incorporated in an ambit framework, one may question whether ambit processes are not in fact too general to be a good building block for financial models. In fact, we always have to find the right tradeoff between flexibility, effectiveness, and concrete application of the proposed model, particularly in the financial framework where decisions have often to be taken rather quickly. Concerning the latter point of view, one has to consider that ambit processes, at least in their most general definition, do not satisfy the martingale property, a fact that can be considered as a particularly relevant drawback within thetraditional approach to financial modeling. Nevertheless, we underline that the martingale framework can be recovered assuming appropriate conditions to hold, paying the price to work with a restricted class of ambit processes. Anyway, it is not obvious that we should adopt the martingale framework if our aim is to model electricity forward contracts. In particular, lines of empirical evidence show that electricity prices are not martingales even at the spot price level, but this does not imply any arbitrage because electricity cannot be stored. The last point is of particular relevance since it ensures that our analysis can be carried out using ambit models, without being worried about arbitrage opportunities. We would also like to mention recent results (see, e.g., [6, 7]), which show that subclasses of nonsemimartingales can be used to model financial assets without necessarily giving rise to arbitrage opportunities in markets exhibiting market frictions, such as, for example, transaction costs.The paper is organized as follows. In Section2, we briefly review the key definitions and concepts on the basis of both ambit fields analysis and ambit processes; in Section 3, we introduce the modeling framework for electricity forward markets, and we study its key properties highlighting the most relevant model specifications; then, we present suitable conditions under which the ambit processes adopted in the presented models are martingales and we also derive a model for the underlying spot price exploiting the forward model through a limiting argument; finally, in Section 4, we concretely apply the approach proposed in Section 3 to forecast the price of a particular forward contract written within the German EEX market.


## 2. Ambit Stochastic Approach

The concepts of ambit fields and ambit processes form the building blocks of the model for the electricity forward price that we will present in Section3. Therefore, in this section, we briefly review the ideas behind ambit stochastic approaches, as they have been introduced by Barndorff-Nielsen and Schmiegel in [3] and then further discussed in [4, 8]; we also refer to [9] and the references therein for further details.


### 2.1. The General Framework

The general background setting for ambit process and ambit field consists of a stochastic fieldY=Yt(x)t∈R, considered in space-time domain D×R, a curve ω(θ)=(x(θ),t(θ)), which also belongs to D×R, and the values of the field along such a curve; namely, Xθ=Yt(θ)(x(θ)); then it follows that X=Xθ is a stochastic process. In most applications, the space D is chosen to be Rn, for n=1,2, or 3. Moreover, the stochastic field Y is supposed to be generated by innovations in space-time and the values Yt(x) are assumed to depend only on those innovations that occur prior to or at time t. More precisely, at each point (x,t), only the innovations in some subset At(x) of D×Rt, where Rt=-∞,t, determine the value of Yt(x). We refer to At(x) as the ambit set associated with (x,t) and to Y and X as an ambit field and an ambit process, respectively.Without further assumptions, nothing interesting can be said about the fieldY and the process X; therefore, we have to specify a suitable mathematical structure in the next subsection; nevertheless, without being rigorous, Yt(x) is defined as a stochastic integral plus a smooth term and the integrand in the stochastic integral consists of a deterministic kernel times a positive random variate which is taken to embody the volatility orintermittency of the field Y. We also underline that the volatility field, denoted by σ, is given as an ambit field and a central issue is what can be learned about σ from observation of Y or X.


### 2.2. Ambit Fields and Ambit Processes

In ambit stochastic approach, an essential role is played by theLévy bases, which we give the definition of in what follows mainly following the extensive and detailed discussion provided by [10, 11]. In addition, the classical concept of subordination of Lévy processes is generalized to subordination of Lévy bases orextended subordination, by introducing the so-calledmeta-times in [10].LetS be aBorel set in Rn and let B(S) be theBorelσ-algebra on S. We denote by Bb(S) the sub-σ-algebra of B(S) consisting of the elements of B(S) bounded with respect to theLebesgue measureμ; that is, BbS=A∈B(S)∣μ(A)<∞.Definition 1 (Lévy basis and homogeneous Lévy basis).
A familyL(A)∣A∈BbS of random vectors in Rn is an Rn-valued Lévy basis onS if the following three properties are satisfied: (i)
The law ofL(A) isinfinitely divisible for all A∈Bb(S).(ii)
For any sequenceA1,…,Ak of disjoint elements of Bb(S), one has that the random variables L(A1),…,L(Ak) areindependent.(iii)
For any sequenceA1,A2,… of disjoint elements of Bb(S) which satisfy ⋃i=1∞Ai∈Bb(S), one has that (1)L⋃i=1∞Ai=∑i=1∞LAi,where the convergence on the right hand side is intended almost surely.
If a Lévy basis has astationary law, then it is called ahomogeneous Lévy basis.Remark 2 (about Lévy bases).
Properties (ii) and (iii) define anindependently scattered random measure; hence, the class of Lévy bases is a subclass of independently scattered random measures.We are now in position to state the definitions of ambit field and ambit process.Definition 3 (ambit field and ambit process).
Anambit fieldY is a random field such that (2)Ytx=η+∫Atxhx,t;ξ,zσzξLdξ,dz+∫Dtxqx,t;ξ,zazξdξdz,provided the integrals exist over ambit sets At(x) and Dt(x), where η is a constant, h and q are deterministic functions, σ≥0 and a are stochastic fields, and L is a Lévy basis. Moreover, we define the correspondingambit processX as the stochastic process, defined on a given probability space Ω,F,P, given by the evaluation of the field Y along a curve ω(θ)=(x(θ),t(θ))⊂D×R; namely, (3)Xθ=Ytθxθ=η+∫Aθhtθ;ξ,zσzξLdξ,dz+∫Dθqtθ;ξ,zazξdξdz,where A(θ)=At(θ)(x(θ)) and D(θ)=Dt(θ)(x(θ)).Remark 4 (particular cases of ambit processes).
Let us specify a particular class of ambit fields that is of particular interest in many applications, which is specified as follows:(4)Ytx=η+∫Atxhx-ξ,t-zσzξLdξ,dz+∫Dtxqx-ξ,t-zazξdξdz,where the ambit sets At(x) and Dt(x) are taken to be homogeneous andnonanticipative; that is, At(x) is of the form At(x)=a+(x,t)∣a∈A, where A only involves negative time coordinates, similar to Dt(x). In fact, this specific framework allows modeling stationary random fields, since Yt(x) in (4) turns out to be stationary when the fields a and σ are stationary and σ is independent of the driving Lévy basis L. Furthermore, other potentially important features of a random field, such as, for example, isotropy or skewness, can be easily modeled with (4) via an appropriate choice of the deterministic kernels h and q and the stochastic fields a and σ.Remark 5 (stochastic integration).
In order to build relevant models based on ambit fields, we need a suitable integration theory. In [12], Rajput and Rosinski have proposed an integration theory for appropriate deterministic integrands with respect to Lévy bases. However, in the context of ambit fields, we need to integrate stochastic integrands. To solve the latter issue, we can count on an alternative integration concept due to Walsh. Such an idea works even if the integrand depends on the particular Lévy basis chosen (see [13]). Unfortunately, the Walsh approach can be used to define ambit fields only when the driving Lévy basis is square-integrable. In the recent work [14], the ideas of Walsh, Rajput, and Rosinski have been combined to propose an integration concept for random integrands and general Lévy bases. Such an approach relies on [15], which is an earlier work by Bichteler and Jacod. We refer to [16] as a complete reference for the comparison of various integration concepts.


## 2.1. The General Framework

The general background setting for ambit process and ambit field consists of a stochastic fieldY=Yt(x)t∈R, considered in space-time domain D×R, a curve ω(θ)=(x(θ),t(θ)), which also belongs to D×R, and the values of the field along such a curve; namely, Xθ=Yt(θ)(x(θ)); then it follows that X=Xθ is a stochastic process. In most applications, the space D is chosen to be Rn, for n=1,2, or 3. Moreover, the stochastic field Y is supposed to be generated by innovations in space-time and the values Yt(x) are assumed to depend only on those innovations that occur prior to or at time t. More precisely, at each point (x,t), only the innovations in some subset At(x) of D×Rt, where Rt=-∞,t, determine the value of Yt(x). We refer to At(x) as the ambit set associated with (x,t) and to Y and X as an ambit field and an ambit process, respectively.Without further assumptions, nothing interesting can be said about the fieldY and the process X; therefore, we have to specify a suitable mathematical structure in the next subsection; nevertheless, without being rigorous, Yt(x) is defined as a stochastic integral plus a smooth term and the integrand in the stochastic integral consists of a deterministic kernel times a positive random variate which is taken to embody the volatility orintermittency of the field Y. We also underline that the volatility field, denoted by σ, is given as an ambit field and a central issue is what can be learned about σ from observation of Y or X.


## 2.2. Ambit Fields and Ambit Processes

In ambit stochastic approach, an essential role is played by theLévy bases, which we give the definition of in what follows mainly following the extensive and detailed discussion provided by [10, 11]. In addition, the classical concept of subordination of Lévy processes is generalized to subordination of Lévy bases orextended subordination, by introducing the so-calledmeta-times in [10].LetS be aBorel set in Rn and let B(S) be theBorelσ-algebra on S. We denote by Bb(S) the sub-σ-algebra of B(S) consisting of the elements of B(S) bounded with respect to theLebesgue measureμ; that is, BbS=A∈B(S)∣μ(A)<∞.Definition 1 (Lévy basis and homogeneous Lévy basis).
A familyL(A)∣A∈BbS of random vectors in Rn is an Rn-valued Lévy basis onS if the following three properties are satisfied: (i)
The law ofL(A) isinfinitely divisible for all A∈Bb(S).(ii)
For any sequenceA1,…,Ak of disjoint elements of Bb(S), one has that the random variables L(A1),…,L(Ak) areindependent.(iii)
For any sequenceA1,A2,… of disjoint elements of Bb(S) which satisfy ⋃i=1∞Ai∈Bb(S), one has that (1)L⋃i=1∞Ai=∑i=1∞LAi,where the convergence on the right hand side is intended almost surely.
If a Lévy basis has astationary law, then it is called ahomogeneous Lévy basis.Remark 2 (about Lévy bases).
Properties (ii) and (iii) define anindependently scattered random measure; hence, the class of Lévy bases is a subclass of independently scattered random measures.We are now in position to state the definitions of ambit field and ambit process.Definition 3 (ambit field and ambit process).
Anambit fieldY is a random field such that (2)Ytx=η+∫Atxhx,t;ξ,zσzξLdξ,dz+∫Dtxqx,t;ξ,zazξdξdz,provided the integrals exist over ambit sets At(x) and Dt(x), where η is a constant, h and q are deterministic functions, σ≥0 and a are stochastic fields, and L is a Lévy basis. Moreover, we define the correspondingambit processX as the stochastic process, defined on a given probability space Ω,F,P, given by the evaluation of the field Y along a curve ω(θ)=(x(θ),t(θ))⊂D×R; namely, (3)Xθ=Ytθxθ=η+∫Aθhtθ;ξ,zσzξLdξ,dz+∫Dθqtθ;ξ,zazξdξdz,where A(θ)=At(θ)(x(θ)) and D(θ)=Dt(θ)(x(θ)).Remark 4 (particular cases of ambit processes).
Let us specify a particular class of ambit fields that is of particular interest in many applications, which is specified as follows:(4)Ytx=η+∫Atxhx-ξ,t-zσzξLdξ,dz+∫Dtxqx-ξ,t-zazξdξdz,where the ambit sets At(x) and Dt(x) are taken to be homogeneous andnonanticipative; that is, At(x) is of the form At(x)=a+(x,t)∣a∈A, where A only involves negative time coordinates, similar to Dt(x). In fact, this specific framework allows modeling stationary random fields, since Yt(x) in (4) turns out to be stationary when the fields a and σ are stationary and σ is independent of the driving Lévy basis L. Furthermore, other potentially important features of a random field, such as, for example, isotropy or skewness, can be easily modeled with (4) via an appropriate choice of the deterministic kernels h and q and the stochastic fields a and σ.Remark 5 (stochastic integration).
In order to build relevant models based on ambit fields, we need a suitable integration theory. In [12], Rajput and Rosinski have proposed an integration theory for appropriate deterministic integrands with respect to Lévy bases. However, in the context of ambit fields, we need to integrate stochastic integrands. To solve the latter issue, we can count on an alternative integration concept due to Walsh. Such an idea works even if the integrand depends on the particular Lévy basis chosen (see [13]). Unfortunately, the Walsh approach can be used to define ambit fields only when the driving Lévy basis is square-integrable. In the recent work [14], the ideas of Walsh, Rajput, and Rosinski have been combined to propose an integration concept for random integrands and general Lévy bases. Such an approach relies on [15], which is an earlier work by Bichteler and Jacod. We refer to [16] as a complete reference for the comparison of various integration concepts.


## 3. Modeling Electricity Forward Price by Ambit Fields

Due to their structure, we can use ambit fields to capture many characteristics of energy markets such as strong seasonal patterns, very pronouncedvolatility clusters, high spikes/jumps, and the Samuelson effect. For these reasons, we review a general model, presented in [5], for electricity forward prices which is based on ambit fields.


### 3.1. A Model for Forward Prices

LetΩ,F,P be a probability space, Lz(A)≔L(A,z)=L(A×0,z) a Lévy basis, σ=σz(A) a stochastic field for A∈B(R+), and z∈R. In addition, we consider the filtration F=Ftt∈R defined by (5)Ft=⋂n=1∞Ft+1/n0,where (6)Ft0=σLA,z,σzA∣A∈BR+,z≤t∨N,and where N denotes the P-null sets of F. Note that Ft is right continuous by construction. Using an ambit field, we model the forward price ft(x) as a stochastic process for each fixed x; namely, (7)ftx≔∫Atxkξ,t-z;xσzξLdξ,dz,where t∈R denotes the current time, σz(ξ)>0 denotes the volatility of the forward market as a whole, At(x) is an ambit set, and k is the kernel function. A complete model description is given specifying an ambit set At(x), a kernel function k, a stochastic volatility field σz(ξ), and a Lévy basis L. Latter choices are often based on qualitative analysis of the markets we want to study as well as concerning the analytical tractability of the resulting model. However, in general, we would like ft(x) to be well defined in the sense of Walsh (see [13]) and stationary in time. In order to have integrability in the sense of Walsh, we assume that the following conditions hold: (i)
The Lévy basisL is square-integrable and has zero-mean.(ii)
The stochastic volatility fieldσ is adapted to the filtration Ftt∈R and independent of the Lévy basis L.(iii)
The kernel functionk is nonnegative and such that k(ξ,u;x)=0 for u<0.(iv)
The convolutionk∗σ is integrable with respect to L.Moreover, in order to ensure thatft(x) is stationary in time, we assume that the following conditions hold:(v)
The stochastic volatility fieldσz(ξ) is stationary in z.(vi)
The ambit setAt(x) is such that (8)Atx=a+0,t∣a∈A0x,whereA0(x)=(ξ,z)∣ξ≥0,z≤t.Throughout this section, we suppose thatft(x) is defined as in (7) and satisfies conditions (i)–(vi).Remark 6 (whenx is a function of t).
The aforementioned conditions onft(x) ensure having stationarity in time; however, as soon as we replace x by a function of t, that is, x=x(t), ft(x(t)) is not in general stationary any more. For instance, if we want to reflect the fact that the forward price depends also on the time to maturity, let us indicate it by T>0; placing x(t)=T-t, we lose the stationarity of ft(x(t)), with respect to time.Remark 7 (aboutf0(x)).
If we consider (7), we have that the forward price at time 0 is given by (9)f0x=∫A0xkξ,-z;xσzξLdξ,dz,where in the right hand side in (9) we perform integration with respect to a random measure. It follows that the forward price at time 0 is thevalue of the random variable f0(x) given in (9), contrary to most other models where f0(x) is considered as a deterministic quantity equal to the observed price.Under previous conditions (i)–(vi), the ambit field specification provided by (7) is highly analytically tractable and itsconditional cumulant function is given by the following proposition.Proposition 8 (conditional cumulant function offt(x)).
LetL be a homogeneous Lévy basis; hence, its control measure is proportional to the Lebesgue measure and one can assume, without loss of generality, that the proportionality constant equals 1. Then, the cumulant function of ft(x), conditioned to σ, is given by (10)Cσζ‡ftx=∫-∞t∫0∞Cζkξ,t-z;xσzξ‡L′dξdz,where L′ is the Lévy seed associated with L; see Subsection 3.2 of [5] for a rigorous definition of Lévy seed associated with a Lévy basis. Moreover, if L is Gaussian, then (11)Cζkξ,t-z;xσzξ‡L′=-12ζ2k2ξ,t-z;xσz2ξ.Proof.
The proposition is an immediate consequence of Proposition2.6 in [12].


### 3.2. Examples of Model Specifications

The forward modelft(x) based on an ambit field has a very general structure; consequently, it is possible to concretely specify exploiting several choices a peculiarity which is highly useful in many applications. There are three components of the model to specify, namely, the Lévy basis L, the kernel function k, and the volatility σ. In the following, we present some examples of specification for theparametersL and k, while we discuss volatility modulation in the next subsection.


#### 3.2.1. Specification of the Lévy Basis

Our model is based on the zero-mean and square-integrable assumptions; hence, we can choose any infinitely divisible distribution satisfying these two hypotheses. A natural choice is to consider a Gaussian Lévy basis which implies a smooth random field. Alternatively, we can take anormal inverse Gaussian Lévy basis. We could need to relax the zero-mean assumption for the Lévy basis, in order to choose, for example, aGamma or aninverse Gaussian Lévy basis.


#### 3.2.2. Specification of the Kernel Function

The kernel functionk plays a key role in our setting since(i)
the kernel function completely determines the space-timeautocorrelation structure of a zero-mean ambit field (see Section 3.4);(ii)
the kernel function characterises the Samuelson effect as we will see in Theorem28;(iii)
the kernel function determines whether the forward price is a martingale or not (see Theorem19 and Corollary 21).Recall that the kernelk is a function of the three variables ξ,t-z, and x, where t-z refers to the temporal component and ξ,x refer to the spatial dimension. A rather natural approach for specifying a kernel function is to assume a factorisation property that we will present by two different realizations each of which has its own relevance in a specific framework. First, we study factorisation into a temporal and a spatial kernel. In particular, we assume that the kernel function factorises as follows: (12)kξ,t-z;x=ϕξ,xψt-z,for suitable functions ψ and ϕ representing the temporal part and the spatial part, respectively. Equation (12) allows us to study specifications of ϕ and ψ separately. In empirical work, it will be particularly interesting to focus in more detail on the question of how to model the spatial kernel function ϕ, which determines the correlation between various forward contracts. Even if we are allowed to choose similar or identical types of functions for both the temporal and the spatial dimension, nevertheless, we will see in Section 3.6 that particular choices for ϕ lead to a rather natural relation between forward and implied spot prices. Let us briefly study an example which is included in the presented modeling framework.Example 9 (exponential kernel function).
LetL be a homogeneous symmetric normal inverse Gaussian Lévy basis NIG(α,0,0,δ); namely, it has Lévy seed L′ with density (13)π-1δαx-1Kαx,where K denotes the modified Bessel function of the second kind, while α, δ>0. Then, (14)Cζ‡L′=δα-δα2+ζ2.If the kernel function k factorises as in (12), and σz(ξ)≡1, then (15)log⁡Eivftx=∫AtxCvkξ,t-z;x‡L′dξdz=∫Atxδα-δα2+ϕξ,xψt-z2dξdz, and the latter integral can be computed explicitly for suitable kernel functions, for example, for ϕ(ξ,x)=e-λ(ξ+x) and ψ(t-z)=e-λ(t-z) for λ>0. Then, (16)kξ,t-z;x=e-λξ+xe-λt-z,and, taking x(t)=T-t as in Remark 6, we have (17)kξ,t-z;x=e-λξe-λT-z.Therefore, the cumulant function is given by (18)log⁡Eivftx=δα∫-∞t∫0∞1-1+c2e-2λξ-sdξdz,for c=ve-2λT/α, and the integral in (18) can be easily computed.Alternative factorisation of the kernel function is given as follows:(19)kξ,t-z;x=ΦξΨt-z,x,for suitable functions Φ and Ψ. Although the latter does not look very natural at first sight, it is very important since it naturally includes cases where t is removed, for example, choosing x(t)=T-t and Ψ(t-z,x)=Ψ(t-z+x); in fact, in this case, we have Ψ(t-z,x)=Ψ(t-z+x)=Ψ(T-z). The latter case is crucial when we want to formulate martingale conditions for the forward price (see Section 3.5). Let us consider two particular realizations for (19).Example 10 (exponentialΨ).
Miming the Ornstein-Uhlenbeck case, forα>0, we can set (20)kt-z,x=Ψt-z,x=e-αt-z+x,hence taking Φ(ξ)≡1 from factorisation (19).Example 11 (hyperbolic kernel function).
In [17], the authors propose a model for electricity forward prices with a Lévy basis obtained from a standard Brownian motion and kernel function given by (21)kz,T=aT-z+b+c,where a, b, and c are positive constants. The authors argued that the Samuelson effect in electricity markets is much steeper than in other commodity markets, defending the choice of a hyperbolic function rather than an exponential one. Note that we can obtain this kernel function from factorisation (19) with Φ(ξ)≡1 and Ψ(t-z,x)=a/(x+t-z+b)+c when x=T-t.


### 3.3. Space-Time Stochastic Volatility

Although a variety of purely temporal stochastic volatility models can be found in the literature, suitable space-time stochastic volatility models still need to be developed. Until now, volatility modulation within the framework of an ambit field can be achieved by four complementary methods: by introducing a stochastic integrand, by extended subordination, byprobability mixing, or byLévy mixing.


#### 3.3.1. Volatility Modulation via a Stochastic Integrand

This method has already been suggested in the initial definition of ambit fields (see (2)). In that case, there are essentially two approaches that can be used to construct a stochastic volatility field. In fact, we can specify the stochastic field directly as a random field, for example, as another ambit field, or we can start with a purely temporal or purely spatial stochastic volatility process, respectively, and then suitably generalise it to a random field. In the following, we will present examples of both types of construction. First, we focus on the modeling approach where we directly specify a random field for the volatility field. A natural starting point for modeling the volatility is to combine kernel smoothing of a homogeneous Lévy basis with (nonlinear) transformation to ensure positivity. For instance, let (22)σt2x=V∫R×RnKx,t-z;ξLσdξ,dz,where Lσ is a homogeneous Lévy basis independent of L, K:Rn×R×Rn→R+ is an integrable kernel function satisfying K(x,u;ξ)=0 for u<0, and V:R→R+ is a continuous, nonnegative function. Clearly, the kernel function K determines the space-time autocorrelation structure of the volatility field σ2.Remark 12 (about stochastic volatility).
Note thatσ2 defined by (22) is stationary in the temporal dimension and if we have (23)Kx,t-z;ξ=K∗x-ξ,t-z,for some suitable function K∗, then the stochastic volatility is both stationary in time and homogeneous in space.Example 13 (X2 distribution).
Consider stochastic volatility as in (22), let Lσ be a standard normal Lévy basis, and let V(x)=x2. Then, σz2(ξ) is clearly positive and has a pointwise X2 distribution with one degree of freedom.In what follows, we show how to construct a stochastic volatility field by extending a purely temporal stochastic process to a random field. Note that our objective is to construct a stochastic volatility field which is stationary, at least in the temporal dimension. There are many possibilities to reach the latter goal, but in the following we decide to focus our attention on a particularly relevant one, namely, theOrnstein-Uhlenbeck type volatility field (OUTVF). The choice of an Ornstein-Uhlenbeck process as the stationary base component is motivated by its analytical tractability and also because it tends to perform well in practice, at least in the purely temporal case. Without loss of generality, we restrict our attention to the case n=1; that is, we consider one spatial dimension. In order to present this type of volatility fields, we have to define theLévy supraprocesses.Definition 14 (Lévy supraprocess generated byX^).
Suppose thatX^=X^tt∈R is a stationary, positive, and infinitely divisible process on R. Let Z·=Zx·x∈R+ be a family of stationary processes such that Z· has independent increments and for each x the cumulant function of Zx· is given by (24)Cm‡Zx·=xCm‡X^,where (25)Cm‡X^=log⁡EeimX^,with (26)mX^=∫X^zmdz,for an arbitrary signed measure m on R. Then, for any fixed t^∈R, Zt^=Zxt^x∈R+ is a Lévy process calledLévy supraprocess generated byX^.Definition 15 (OUTVF).
Suppose thatY^ is a positive Ornstein-Uhlenbeck process with rate parameter λ>0 and such that it is generated by a Lévy subordinator Y; that is, (27)Y^t=∫-∞te-λt-zdYz.A stochastic volatility field σt2(x) on R×R is an OUTVF if it is defined as (28)τtx=σt2x=e-ηxY^t+∫0xe-ηx-ξdZξt,where η>0 is the spatial rate parameter and Z=Z·tt∈R+ is a family of Lévy supraprocesses.Note that, in the above construction, we start from an Ornstein-Uhlenbeck process in time,τt(0)=Y^t, and the spatial structure is then introduced through two steps. First, we multiply the process with an exponential weight e-ηx dependent on a spatial component, which reaches its maximum for x=0, decaying as it moves away from the purely temporal case. Second, we add an integral which is similar to an Ornstein-Uhlenbeck process in the spatial variable x. However, although the stochastic volatility field τ is stationary in time, note here that the integration starts from 0 rather than from -∞; hence, the resulting component is not stationary in the spatial variable x, even if this could be changed without affecting the structure if such a property is required in a particular application. Finally, note that the process τt(x) is in general not predictable, which is disadvantageous given that we want to construct stochastic integrals in the sense of Walsh. However, if we choose X^ as an Ornstein-Uhlenbeck process, we obtain a predictable stochastic volatility process.Example 16 (covariance function ofτ when X^ is an Ornstein-Uhlenbeck process).
Suppose thatX^ is an Ornstein-Uhlenbeck process with rate parameter κ>0 and generated by a Lévy process X. Then, using the notations in Definition 15, we have (29)Covτtx,τt~x~=12VarY1λ-1e-λt-t~-ηx+x~+VarX^0η-1e-κt-t~-ηx-x~-VarX^0η-1e-κt-t~-ηx+x~.Furthermore, if Var(Y1)=Var(X^0) and κ=λ=η, then for fixed x and x~ the covariance function of τ is (30)Covτtx,τt~x~=e-κt-t~+x-x~.


#### 3.3.2. Volatility Modulation via Extended Subordination

An alternative way of volatility modulation is by means of extended subordination. In this modulation, the volatilityσt2(x)=τt(x) is incorporated in the modeling through a meta-time T associated with an absolutely continuous measure T on S; namely, (31)TA=∫Aτtxdxdt,where A∈Bb(S). A natural choice of meta-time is T(A)=T(x,t)∣(x,t)∈A, where T is the mapping given by T(x,t)=x,τt+(x), while (32)τt+x=∫0tτzxdz.The above construction of a meta-time in a space-time model is very general, since it allows defining a variety of models for the random field τt(x), leading to new model specifications. For instance, one could model τt(x) by an OUTVF or any other model discussed previously. It follows that this way to modulate the volatility is more general than the previous one.


#### 3.3.3. Volatility Modulation via Probability Mixing and Lévy Mixing

Volatility modulation can also be obtained through a probability mixing approach which is realized defining new distributions by randomising a parameter from a given parametric distribution, as in the following classical example.Example 17 (normal variance-mean mixture distribution).
Consider as base model(33)Ytx=∫Atxhx,t;ξ,zLdξ,dz,assuming that the corresponding Lévy basis L is homogeneous and Gaussian; namely, the corresponding Lévy seed is given by L′~Nη+βσ2,σ2 with η,β,σ∈R and σ2>0. Then, we use probability mixing supposing that σ2 is random. Hence, the conditional law of the Lévy seed is given by L′∣σ~Nη+βσ2,σ2. Due to the scaling property of the Gaussian distribution, such a model can be represented as in (2); therefore, in this particular case, probability mixing and volatility modulation via a stationary stochastic integrand produce the same result. Such a construction falls into the class ofnormal variance-mean mixtures.In this context, it is important to note that probability mixing does not generally lead to infinitely divisible distributions. In order to overcome this problem, we present the concept of Lévy mixing. LetL be afactorisable Lévy basis on Rn withcharacteristic quadruplet(a~,m~,ν(dx,·),c) (see Subsection 2.2 of [10] for a rigorous definition of characteristic quadruplet of a Lévy basis and factorisable Lévy basis). Suppose that ν(dx,·) depends on the parameter θ∈Π, where Π denotes the space of the parameter θ. Then, theLévy measure of L is given by ν(dx,θ)c(ds). Let γ be a measure on Π and define (34)n^dx,ds=∫Πνdx,θγdθcds,where we assume that (35)∫R1∧x2n^dx,ds<∞;then there exists a Lévy basis L^ which has n^ as its Lévy measure. We define L^ as theLévy basis obtained by Lévy mixing ofLwith the measureγ. Let us study a concrete example of Lévy mixing.Example 18 (superposition of Ornstein-Uhlenbeck type process).
Let us consider the example of superposition of Ornstein-Uhlenbeck type process. LetL be a subordinator without drift and with Lévy measure νL, and consider an Ornstein-Uhlenbeck process (36)Yt=∫-∞te-θt-zdLz,for θ>0. By a simple calculation, we obtain the following expression for the cumulant function of Yt: (37)Cζ‡Yt=∫0∞eiζx-1νdx,θ,for ζ∈R, where (38)νdx,θ=∫0∞νLeθudxduis a mixture of νL with the Lebesgue measure. Lévy mixing can be carried out with respect to the parameter θ; namely, (39)ν~dx=∫0∞νdx,θγdθ,where γ is a measure on 0,∞ satisfying(40)∫0∞xν~dx<∞.Taking L~ to be the Lévy basis with Lévy measure νL(dx)duγ(dθ) and defining the superposition of Ornstein-Uhlenbeck type process Y~t with respect to L~ by (41)Y~t=∫0∞∫-∞te-λt-zL~dz,dλ,where λ>0, we have that the cumulant function of Y~t is given by (42)Cζ‡Y~t=∫0∞∫-∞tCζe-θt-z‡L1dzγdθ=∫0∞∫0∞eiζe-θu-1νLdxduγdθ=∫0∞eiζx-1ν~dx,hence showing that superposition of Ornstein-Uhlenbeck type process can be obtained from an Ornstein-Uhlenbeck process through Lévy mixing.


### 3.4. Autocorrelation and Cross-Correlation

We would like to underline that the proposed approach allows modeling the entire forward curve. Hence, it is interesting to study the correlation structure between various forward contracts implied by the presented modeling framework. Starting from the definition of covariance function of an ambit field (see SectionA.4 of [5]), we deduce theautocorrelation function of a forward contract ft(x). More precisely, assuming that the Lévy basis is homogeneous, and taking t, h∈R, h≥0, and x, x~≥0, we have that the ambit set defined in (8) can be expressed as At(x)∩At+h(x~)=A-t(x)≔A-t; then we define the autocorrelation function as (43)Corftx,ft+hx~≔∫0,∞×0,∞kξ,u;xkξ,u+h;x~Eσ02ξdξdu∫0,∞×0,∞k2ξ,u;xEσ02ξdξdu∫0,∞×0,∞k2ξ,u′;x~Eσ02ξdξdu′.Note that the correlation structure is principally determined by three factors, which are the intersection of the corresponding ambit sets, the kernel function, and the autocorrelation structure of the stochastic volatility field. Furthermore, we could model various commodity forward contracts, such as electricity and natural gas futures, simultaneously. It follows that we can specify different ambit sets, kernel functions, stochastic volatility fields, and Lévy bases in order to obtain a rather flexible correlation structure. In such a situation, it becomes even more clear how flexible the ambit setting is. The details of these multivariate extensions can be found in Section A.5 of [5].


### 3.5. Martingale Condition

In what follows, we give sufficient conditions for the presented model to be a martingale, even if such property may be too restrictive if we are interested in modeling electricity forward contracts. In fact, there are at least two arguments that support the choice of more general classes of stochastic processes than martingales. First, in the energy context, it might not be as crucial thatft(T-t) is a martingale as it is in the context of modeling interest rates. In fact, in many emerging electricity markets, a seller may not be able to find any buyer to get rid of a forward; oppositely, it may also happen that there are no sellers when one wants to buy a contract. Hence, the illiquidity prevents possible arbitrage opportunities from being exercised and we need a martingale condition only in the absence of arbitrage setting (see [18]). Second, independently of the particular structure of energy markets, the recent literature in mathematical finance (see, e.g., [6, 18]) has highlighted that some classes of nonsemimartingales, in particular stochastic processes with conditional full support, do not necessarily give rise to arbitrage opportunities reflecting more realistic characteristics of the market. In the null-spatial setting, Pakkanen proves thatBrownian semistationary processes have in fact conditional full support (see [7]). In future research, it will be interesting to study extensions of the latter result to more general ambit fields. However, the question of establishing martingale conditions for ambit fields is interesting in order to point out a subset of the class we are considering, which can be exploited to model general forward prices, hence without restricting ourselves to the electricity case. In what follows, we formulate the martingale condition for ambit fields more general than those defined by (7), where the ambit set At(x) is chosen as in (8), also showing that such condition simplifies in the modeling framework described in Section 3.1.Theorem 19 (martingale condition for general ambit fields).
Letx=T-t for some T>0 and, for fixed t∈R, let (44)Ytx≔YtT-t=∫Atxhξ,z;T-t,tσzξLdξ,dz,where At(x) is defined as in (8), h is a deterministic kernel function, σ is an adapted, nonnegative random field, and L is a Lévy basis satisfying both (i) and (iv). Further, σ and L are assumed to be independent. Then, Yt(T-t)t∈R is a martingale with respect to Ftt∈R defined as in (5) if and only if, for all ξ>0,z≤t≤T, it holds (45)hξ,z;T-t,t=h~ξ,z;T,for some deterministic kernel function h~.Proof.
Measurability and integrability properties are trivially verified according to what we have seen so far; therefore, we just have to show that(46)EYtT-t∣Ft~=Yt~T-t~,for all t~≤t. Note that, for t~≤t, we have that At~(x)⊂At(x). Using the independence property of σ and L and the fact that L is a zero-mean process, we find (47)EYtT-t∣Ft~=E∫At~xhξ,z;T-t,tσzξLdξ,dz+∫Atx∖At~xhξ,z;T-t,tσzξLdξ,dz∣Ft~=∫At~xhξ,z;t,T-tσzξLdξ,dz=Yt~T-t~+It~T-t~, where (48)It~T-t~=∫At~xhξ,z;T-t,t-hξ,z;T-t~,t~σzξLdξ,dz.Without loss of generality, we can assume that Var(L)=1. Since L is a Lévy basis with zero-mean, we know that EIt~(T-t~)=0 and from theItô isometry (see Theorem 2.5 in [13]) we therefore get that (49)VarIt~T-t~=∫At~xhξ,z;T-t,t-hξ,z;T-t~,t~2Eσz2ξc~dξ,dz,where c~ is thecovariance measure induced by L; see, for example, [13], for a rigorous definition of covariance measure induced by a Lévy basis. Thus, in order to obtain It~(T-t~)=0, we need that, for all ξ>0,z≤t~≤t≤T, (50)hξ,z;T-t,t=hξ,z;T-t~,t~.Note that there is only one class of functions which satisfy the above condition, namely, functions of the form (51)hξ,z;T-t,t=h~ξ,z;T,for all ξ>0, z≤t~≤t≤T and for some deterministic kernel function h~.Remark 20 (Lévy basis with mean different from zero).
If we would like to work with Lévy basesL which do not have zero-mean, then the martingale condition has to be extended by an additional drift condition.Corollary 21 (martingale condition forft(x)).
One gets thatft(T-t)t∈R is a martingale with respect to Ftt∈R if and only if, for all ξ>0, z≤t≤T, the following holds: (52)kξ,t-z;T-t=k~ξ,T-z,for some deterministic kernel function k~; hence, the latter is a special case of the factorisation given in (19).Proof.
This result follows directly from Theorem19 considering Yt(T-t)=ft(T-t).Clearly, the martingale condition is rather strong and it is necessary to check whether there are actually any relevant cases, not excluded by condition (52). Therefore, let us study some examples.Example 22 (exponential kernel function).
If we choosek to be of the form (53)kξ,t-z;x=kξ,t-z;T-t=e-αξ+x+t-z=e-αξ+T-z,for some α>0, then the above martingale condition is clearly satisfied. Note that this choice of the kernel function belongs to both factorisation classes (12) and (19).Example 23 (factorisation (19)).
Let us consider kernel functionsk which factorise as in (19); that is, (54)kξ,t-z;x=ΦξΨt-z,x;then the choice of the function Φ does not have any impact concerning whether the ambit field is a martingale; indeed such a property is solely determined by the function Ψ. It follows that every choice of the form (55)Ψt-z,T-t=Ψ~T-zsatisfies the martingale condition.


### 3.6. Derivation of a Spot Model from the Forward One

Under the absence of arbitrage assumption, the forward price of a contract tends to be equal to the spot price of the underlying as the time to maturity tends to zero. In this subsection, we investigate in detail the nature of the spot price model implied by the above modeling framework for the forward price. Note that this study should be understood as a theoretical exercise, since we do not observe convergence of the electricity forward price to the electricity spot price in real markets.Letft(x) be the forward price and let St be the spot price of the underlying. When the time to maturity x=T-t tends to 0, we would have that ft(x) tends to St; hence, we define the spot price as (56)St≔∫-∞t∫0∞kξ,t-z;0σzξLdξ,dz.Lemma 24 (convergence offt(x) to St).
Suppose that(57)limx↓0⁡∫-∞t∫0∞kξ,t-z;x-kξ,t-z;02Eσz2ξdξdz=0;then ft(x)→St in L2Ω,F,P∗ as time to maturity x tends to zero, where P∗ denotes the risk neutral probability measure.Proof.
The result follows immediately by applying the Itô isometry, as proved by Theorem2.5 in [13].The previous lemma gives us that the forward price will tend continuously in variance to the spot price as time to maturity decreases to zero. Note also that wheneverσz(ξ) is a stationary field, the condition in the above lemma is translated to convergence of k(·,·;x) to k(·,·;0) in L2R+2.Using Proposition8 for the forward price, we can easily derive the conditional cumulant function for the implied spot price.Proposition 25 (conditional cumulant function ofSt).
LetL be a homogeneous Lévy basis. Then, for St defined as in (56), the cumulant function conditioned to σ is given by (58)Cσζ‡St=∫-∞t∫0∞Cζkξ,t-z;0σzξ‡L′dξdz,where L′ is the Lévy seed associated with L.Proof.
This result is an immediate consequence of Proposition8.A case of some special interest is given by considering the driving Lévy basisL of the ambit field to be a homogeneous Gaussian Lévy basis; in fact, in this case, we get the following result.Corollary 26 (spot price whenL is a Gaussian Lévy basis).
LetSt be as in (56) and let L be a homogeneous Gaussian Lévy basis with (59)Cζ‡L′=-12ζ2,where L′ is the Lévy seed associated with L. Then, the cumulant function of St conditioned to σ is given by (60)Cσζ‡St=-12ζ2∫-∞t∫0∞k2ξ,t-z;0σz2ξdξdz.If the kernel function k can be factorized as in (12), then (61)Cσζ‡St=-12ζ2∫-∞tψ2t-zωz2ξdz,where (62)ωz2=∫0∞ϕ2ξ,0σz2ξdξ,implying that St is equal in law to (63)St=∫-∞tψ2t-zωz-ξdWz,where W is the standard Brownian motion.Proof.
The result follows by straightforward calculations, exploiting Proposition25 and the definition of Gaussian Lévy basis.Remark 27 (modeling spot price).
Note thatSt in (63) is a Brownian semistationary process, which has been used as a model for energy spot prices in [1].


#### 3.6.1. The Samuelson Effect

Recall that the Samuelson effect describes the empirical fact that the volatility of the forward price increasingly converges to the volatility of the underlying spot price when the time to maturity tends to zero. Here, we show that this effect is naturally included in our modeling framework.Theorem 28 (nondecreasing of the conditional variance offt(x)).
Assume that the functionx↦k(ξ,u;x) is monotonically nondecreasing in x≥0 for every (ξ,u)∈R+2. The conditional variance of the forward price ft(x) given by (64)vtx≔c∫-∞t∫0∞k2ξ,t-z;xσz2ξdξdzis monotonically nondecreasing in x, for t≥0, where the term c is a suitable constant; see Section A.2 of [5] for a rigorous definition of c.Proof.
Forx≥0, we have (65)vtx=c∫-∞t∫0∞k2ξ,t-z;xσz2ξdξdz,and, taking 0≤x≤x~, then (66)vtx~-vtx=∫-∞t∫0∞k2ξ,t-z;x~-k2ξ,t-z;xσz2ξdξdz≥0,since k(ξ,t-z;x) is nondecreasing in x.Remark 29 (Samuelson effect).
The conditional variance of the spot price is given byvt(0) and from the above theorem it follows that vt(x)≥vt(0) for x≥0. Since we can consider vt(x) as a monotonically increasing sequence bounded from below by vt(0), there exists the limit limx↓0vt(x)=l. Under the condition in Lemma 24, this limit is given by l=vt(0); hence, it equals the conditional variance of spot price, which implies that we have a Samuelson effect.


### 3.7. Simulation Algorithm

In what follows, we present a possible algorithm, taken from [5], to simulate an ambit field Yt(x) of the following type: (67)Ytx=∫-∞t∫0∞kξ,t-z;xσzξLdξ,dz.A problem that we immediately note is that we perform integration over an ambit set At(x)=-∞,t×0,∞, which is unbounded; hence, we have to truncate the ambit set before performing any computations. Therefore, we introduce the ambit field (68)Y-tx=∫M1t∫0M2kξ,t-z;xσzξLdξ,dz,where M1<t and M2>0. Obviously, if M1→-∞ and M2→∞, then Y-t(x) converges to Yt(x).Takingn,m∈N+, we construct a grid for M1,t×0,M2 by dividing the interval M1,t into n equidistant intervals of length (t-M1)/n, where t=t1>t2>⋯>tn=M1, and by dividing 0,M2 into m equidistant intervals of length M2/m, where 0=x1<x2<⋯<xm=M2. Next, we simulate an ambit field through the following steps:(1)
Simulate the stochastic field on the grid points(ti,xj) in order to obtain the values of stochastic volatility σti(xj) for i=1,…,n and j=1,…,m. In case we have no stochastic volatility, we set σti(xj)=1, for all i=1,…,n and for all j=1,…,m.(2)
Simulaten·m random variables Zi+j~L(Δ), where L is the Lévy basis computed on (69)Δ=Δn,m,M1,M2=t-M1nM2m,fori=1,…,n and j=1,…,m.(3)
Compute the approximated ambit field in this way:(70)Y^tx=∑i=1n-1∑j=1m-1kxj,t-ti;xσtixjZi+j.We note that the previous algorithm can be used also to simulate stochastic volatility. In the last step of the simulation algorithm, we use the definition of stochastic integral in the sense of Walsh for simple processes, which ensures the convergence of (70) to the actual stochastic integral, when Δ↓0.


## 3.1. A Model for Forward Prices

LetΩ,F,P be a probability space, Lz(A)≔L(A,z)=L(A×0,z) a Lévy basis, σ=σz(A) a stochastic field for A∈B(R+), and z∈R. In addition, we consider the filtration F=Ftt∈R defined by (5)Ft=⋂n=1∞Ft+1/n0,where (6)Ft0=σLA,z,σzA∣A∈BR+,z≤t∨N,and where N denotes the P-null sets of F. Note that Ft is right continuous by construction. Using an ambit field, we model the forward price ft(x) as a stochastic process for each fixed x; namely, (7)ftx≔∫Atxkξ,t-z;xσzξLdξ,dz,where t∈R denotes the current time, σz(ξ)>0 denotes the volatility of the forward market as a whole, At(x) is an ambit set, and k is the kernel function. A complete model description is given specifying an ambit set At(x), a kernel function k, a stochastic volatility field σz(ξ), and a Lévy basis L. Latter choices are often based on qualitative analysis of the markets we want to study as well as concerning the analytical tractability of the resulting model. However, in general, we would like ft(x) to be well defined in the sense of Walsh (see [13]) and stationary in time. In order to have integrability in the sense of Walsh, we assume that the following conditions hold: (i)
The Lévy basisL is square-integrable and has zero-mean.(ii)
The stochastic volatility fieldσ is adapted to the filtration Ftt∈R and independent of the Lévy basis L.(iii)
The kernel functionk is nonnegative and such that k(ξ,u;x)=0 for u<0.(iv)
The convolutionk∗σ is integrable with respect to L.Moreover, in order to ensure thatft(x) is stationary in time, we assume that the following conditions hold:(v)
The stochastic volatility fieldσz(ξ) is stationary in z.(vi)
The ambit setAt(x) is such that (8)Atx=a+0,t∣a∈A0x,whereA0(x)=(ξ,z)∣ξ≥0,z≤t.Throughout this section, we suppose thatft(x) is defined as in (7) and satisfies conditions (i)–(vi).Remark 6 (whenx is a function of t).
The aforementioned conditions onft(x) ensure having stationarity in time; however, as soon as we replace x by a function of t, that is, x=x(t), ft(x(t)) is not in general stationary any more. For instance, if we want to reflect the fact that the forward price depends also on the time to maturity, let us indicate it by T>0; placing x(t)=T-t, we lose the stationarity of ft(x(t)), with respect to time.Remark 7 (aboutf0(x)).
If we consider (7), we have that the forward price at time 0 is given by (9)f0x=∫A0xkξ,-z;xσzξLdξ,dz,where in the right hand side in (9) we perform integration with respect to a random measure. It follows that the forward price at time 0 is thevalue of the random variable f0(x) given in (9), contrary to most other models where f0(x) is considered as a deterministic quantity equal to the observed price.Under previous conditions (i)–(vi), the ambit field specification provided by (7) is highly analytically tractable and itsconditional cumulant function is given by the following proposition.Proposition 8 (conditional cumulant function offt(x)).
LetL be a homogeneous Lévy basis; hence, its control measure is proportional to the Lebesgue measure and one can assume, without loss of generality, that the proportionality constant equals 1. Then, the cumulant function of ft(x), conditioned to σ, is given by (10)Cσζ‡ftx=∫-∞t∫0∞Cζkξ,t-z;xσzξ‡L′dξdz,where L′ is the Lévy seed associated with L; see Subsection 3.2 of [5] for a rigorous definition of Lévy seed associated with a Lévy basis. Moreover, if L is Gaussian, then (11)Cζkξ,t-z;xσzξ‡L′=-12ζ2k2ξ,t-z;xσz2ξ.Proof.
The proposition is an immediate consequence of Proposition2.6 in [12].


## 3.2. Examples of Model Specifications

The forward modelft(x) based on an ambit field has a very general structure; consequently, it is possible to concretely specify exploiting several choices a peculiarity which is highly useful in many applications. There are three components of the model to specify, namely, the Lévy basis L, the kernel function k, and the volatility σ. In the following, we present some examples of specification for theparametersL and k, while we discuss volatility modulation in the next subsection.


### 3.2.1. Specification of the Lévy Basis

Our model is based on the zero-mean and square-integrable assumptions; hence, we can choose any infinitely divisible distribution satisfying these two hypotheses. A natural choice is to consider a Gaussian Lévy basis which implies a smooth random field. Alternatively, we can take anormal inverse Gaussian Lévy basis. We could need to relax the zero-mean assumption for the Lévy basis, in order to choose, for example, aGamma or aninverse Gaussian Lévy basis.


### 3.2.2. Specification of the Kernel Function

The kernel functionk plays a key role in our setting since(i)
the kernel function completely determines the space-timeautocorrelation structure of a zero-mean ambit field (see Section 3.4);(ii)
the kernel function characterises the Samuelson effect as we will see in Theorem28;(iii)
the kernel function determines whether the forward price is a martingale or not (see Theorem19 and Corollary 21).Recall that the kernelk is a function of the three variables ξ,t-z, and x, where t-z refers to the temporal component and ξ,x refer to the spatial dimension. A rather natural approach for specifying a kernel function is to assume a factorisation property that we will present by two different realizations each of which has its own relevance in a specific framework. First, we study factorisation into a temporal and a spatial kernel. In particular, we assume that the kernel function factorises as follows: (12)kξ,t-z;x=ϕξ,xψt-z,for suitable functions ψ and ϕ representing the temporal part and the spatial part, respectively. Equation (12) allows us to study specifications of ϕ and ψ separately. In empirical work, it will be particularly interesting to focus in more detail on the question of how to model the spatial kernel function ϕ, which determines the correlation between various forward contracts. Even if we are allowed to choose similar or identical types of functions for both the temporal and the spatial dimension, nevertheless, we will see in Section 3.6 that particular choices for ϕ lead to a rather natural relation between forward and implied spot prices. Let us briefly study an example which is included in the presented modeling framework.Example 9 (exponential kernel function).
LetL be a homogeneous symmetric normal inverse Gaussian Lévy basis NIG(α,0,0,δ); namely, it has Lévy seed L′ with density (13)π-1δαx-1Kαx,where K denotes the modified Bessel function of the second kind, while α, δ>0. Then, (14)Cζ‡L′=δα-δα2+ζ2.If the kernel function k factorises as in (12), and σz(ξ)≡1, then (15)log⁡Eivftx=∫AtxCvkξ,t-z;x‡L′dξdz=∫Atxδα-δα2+ϕξ,xψt-z2dξdz, and the latter integral can be computed explicitly for suitable kernel functions, for example, for ϕ(ξ,x)=e-λ(ξ+x) and ψ(t-z)=e-λ(t-z) for λ>0. Then, (16)kξ,t-z;x=e-λξ+xe-λt-z,and, taking x(t)=T-t as in Remark 6, we have (17)kξ,t-z;x=e-λξe-λT-z.Therefore, the cumulant function is given by (18)log⁡Eivftx=δα∫-∞t∫0∞1-1+c2e-2λξ-sdξdz,for c=ve-2λT/α, and the integral in (18) can be easily computed.Alternative factorisation of the kernel function is given as follows:(19)kξ,t-z;x=ΦξΨt-z,x,for suitable functions Φ and Ψ. Although the latter does not look very natural at first sight, it is very important since it naturally includes cases where t is removed, for example, choosing x(t)=T-t and Ψ(t-z,x)=Ψ(t-z+x); in fact, in this case, we have Ψ(t-z,x)=Ψ(t-z+x)=Ψ(T-z). The latter case is crucial when we want to formulate martingale conditions for the forward price (see Section 3.5). Let us consider two particular realizations for (19).Example 10 (exponentialΨ).
Miming the Ornstein-Uhlenbeck case, forα>0, we can set (20)kt-z,x=Ψt-z,x=e-αt-z+x,hence taking Φ(ξ)≡1 from factorisation (19).Example 11 (hyperbolic kernel function).
In [17], the authors propose a model for electricity forward prices with a Lévy basis obtained from a standard Brownian motion and kernel function given by (21)kz,T=aT-z+b+c,where a, b, and c are positive constants. The authors argued that the Samuelson effect in electricity markets is much steeper than in other commodity markets, defending the choice of a hyperbolic function rather than an exponential one. Note that we can obtain this kernel function from factorisation (19) with Φ(ξ)≡1 and Ψ(t-z,x)=a/(x+t-z+b)+c when x=T-t.


## 3.2.1. Specification of the Lévy Basis

Our model is based on the zero-mean and square-integrable assumptions; hence, we can choose any infinitely divisible distribution satisfying these two hypotheses. A natural choice is to consider a Gaussian Lévy basis which implies a smooth random field. Alternatively, we can take anormal inverse Gaussian Lévy basis. We could need to relax the zero-mean assumption for the Lévy basis, in order to choose, for example, aGamma or aninverse Gaussian Lévy basis.


## 3.2.2. Specification of the Kernel Function

The kernel functionk plays a key role in our setting since(i)
the kernel function completely determines the space-timeautocorrelation structure of a zero-mean ambit field (see Section 3.4);(ii)
the kernel function characterises the Samuelson effect as we will see in Theorem28;(iii)
the kernel function determines whether the forward price is a martingale or not (see Theorem19 and Corollary 21).Recall that the kernelk is a function of the three variables ξ,t-z, and x, where t-z refers to the temporal component and ξ,x refer to the spatial dimension. A rather natural approach for specifying a kernel function is to assume a factorisation property that we will present by two different realizations each of which has its own relevance in a specific framework. First, we study factorisation into a temporal and a spatial kernel. In particular, we assume that the kernel function factorises as follows: (12)kξ,t-z;x=ϕξ,xψt-z,for suitable functions ψ and ϕ representing the temporal part and the spatial part, respectively. Equation (12) allows us to study specifications of ϕ and ψ separately. In empirical work, it will be particularly interesting to focus in more detail on the question of how to model the spatial kernel function ϕ, which determines the correlation between various forward contracts. Even if we are allowed to choose similar or identical types of functions for both the temporal and the spatial dimension, nevertheless, we will see in Section 3.6 that particular choices for ϕ lead to a rather natural relation between forward and implied spot prices. Let us briefly study an example which is included in the presented modeling framework.Example 9 (exponential kernel function).
LetL be a homogeneous symmetric normal inverse Gaussian Lévy basis NIG(α,0,0,δ); namely, it has Lévy seed L′ with density (13)π-1δαx-1Kαx,where K denotes the modified Bessel function of the second kind, while α, δ>0. Then, (14)Cζ‡L′=δα-δα2+ζ2.If the kernel function k factorises as in (12), and σz(ξ)≡1, then (15)log⁡Eivftx=∫AtxCvkξ,t-z;x‡L′dξdz=∫Atxδα-δα2+ϕξ,xψt-z2dξdz, and the latter integral can be computed explicitly for suitable kernel functions, for example, for ϕ(ξ,x)=e-λ(ξ+x) and ψ(t-z)=e-λ(t-z) for λ>0. Then, (16)kξ,t-z;x=e-λξ+xe-λt-z,and, taking x(t)=T-t as in Remark 6, we have (17)kξ,t-z;x=e-λξe-λT-z.Therefore, the cumulant function is given by (18)log⁡Eivftx=δα∫-∞t∫0∞1-1+c2e-2λξ-sdξdz,for c=ve-2λT/α, and the integral in (18) can be easily computed.Alternative factorisation of the kernel function is given as follows:(19)kξ,t-z;x=ΦξΨt-z,x,for suitable functions Φ and Ψ. Although the latter does not look very natural at first sight, it is very important since it naturally includes cases where t is removed, for example, choosing x(t)=T-t and Ψ(t-z,x)=Ψ(t-z+x); in fact, in this case, we have Ψ(t-z,x)=Ψ(t-z+x)=Ψ(T-z). The latter case is crucial when we want to formulate martingale conditions for the forward price (see Section 3.5). Let us consider two particular realizations for (19).Example 10 (exponentialΨ).
Miming the Ornstein-Uhlenbeck case, forα>0, we can set (20)kt-z,x=Ψt-z,x=e-αt-z+x,hence taking Φ(ξ)≡1 from factorisation (19).Example 11 (hyperbolic kernel function).
In [17], the authors propose a model for electricity forward prices with a Lévy basis obtained from a standard Brownian motion and kernel function given by (21)kz,T=aT-z+b+c,where a, b, and c are positive constants. The authors argued that the Samuelson effect in electricity markets is much steeper than in other commodity markets, defending the choice of a hyperbolic function rather than an exponential one. Note that we can obtain this kernel function from factorisation (19) with Φ(ξ)≡1 and Ψ(t-z,x)=a/(x+t-z+b)+c when x=T-t.


## 3.3. Space-Time Stochastic Volatility

Although a variety of purely temporal stochastic volatility models can be found in the literature, suitable space-time stochastic volatility models still need to be developed. Until now, volatility modulation within the framework of an ambit field can be achieved by four complementary methods: by introducing a stochastic integrand, by extended subordination, byprobability mixing, or byLévy mixing.


### 3.3.1. Volatility Modulation via a Stochastic Integrand

This method has already been suggested in the initial definition of ambit fields (see (2)). In that case, there are essentially two approaches that can be used to construct a stochastic volatility field. In fact, we can specify the stochastic field directly as a random field, for example, as another ambit field, or we can start with a purely temporal or purely spatial stochastic volatility process, respectively, and then suitably generalise it to a random field. In the following, we will present examples of both types of construction. First, we focus on the modeling approach where we directly specify a random field for the volatility field. A natural starting point for modeling the volatility is to combine kernel smoothing of a homogeneous Lévy basis with (nonlinear) transformation to ensure positivity. For instance, let (22)σt2x=V∫R×RnKx,t-z;ξLσdξ,dz,where Lσ is a homogeneous Lévy basis independent of L, K:Rn×R×Rn→R+ is an integrable kernel function satisfying K(x,u;ξ)=0 for u<0, and V:R→R+ is a continuous, nonnegative function. Clearly, the kernel function K determines the space-time autocorrelation structure of the volatility field σ2.Remark 12 (about stochastic volatility).
Note thatσ2 defined by (22) is stationary in the temporal dimension and if we have (23)Kx,t-z;ξ=K∗x-ξ,t-z,for some suitable function K∗, then the stochastic volatility is both stationary in time and homogeneous in space.Example 13 (X2 distribution).
Consider stochastic volatility as in (22), let Lσ be a standard normal Lévy basis, and let V(x)=x2. Then, σz2(ξ) is clearly positive and has a pointwise X2 distribution with one degree of freedom.In what follows, we show how to construct a stochastic volatility field by extending a purely temporal stochastic process to a random field. Note that our objective is to construct a stochastic volatility field which is stationary, at least in the temporal dimension. There are many possibilities to reach the latter goal, but in the following we decide to focus our attention on a particularly relevant one, namely, theOrnstein-Uhlenbeck type volatility field (OUTVF). The choice of an Ornstein-Uhlenbeck process as the stationary base component is motivated by its analytical tractability and also because it tends to perform well in practice, at least in the purely temporal case. Without loss of generality, we restrict our attention to the case n=1; that is, we consider one spatial dimension. In order to present this type of volatility fields, we have to define theLévy supraprocesses.Definition 14 (Lévy supraprocess generated byX^).
Suppose thatX^=X^tt∈R is a stationary, positive, and infinitely divisible process on R. Let Z·=Zx·x∈R+ be a family of stationary processes such that Z· has independent increments and for each x the cumulant function of Zx· is given by (24)Cm‡Zx·=xCm‡X^,where (25)Cm‡X^=log⁡EeimX^,with (26)mX^=∫X^zmdz,for an arbitrary signed measure m on R. Then, for any fixed t^∈R, Zt^=Zxt^x∈R+ is a Lévy process calledLévy supraprocess generated byX^.Definition 15 (OUTVF).
Suppose thatY^ is a positive Ornstein-Uhlenbeck process with rate parameter λ>0 and such that it is generated by a Lévy subordinator Y; that is, (27)Y^t=∫-∞te-λt-zdYz.A stochastic volatility field σt2(x) on R×R is an OUTVF if it is defined as (28)τtx=σt2x=e-ηxY^t+∫0xe-ηx-ξdZξt,where η>0 is the spatial rate parameter and Z=Z·tt∈R+ is a family of Lévy supraprocesses.Note that, in the above construction, we start from an Ornstein-Uhlenbeck process in time,τt(0)=Y^t, and the spatial structure is then introduced through two steps. First, we multiply the process with an exponential weight e-ηx dependent on a spatial component, which reaches its maximum for x=0, decaying as it moves away from the purely temporal case. Second, we add an integral which is similar to an Ornstein-Uhlenbeck process in the spatial variable x. However, although the stochastic volatility field τ is stationary in time, note here that the integration starts from 0 rather than from -∞; hence, the resulting component is not stationary in the spatial variable x, even if this could be changed without affecting the structure if such a property is required in a particular application. Finally, note that the process τt(x) is in general not predictable, which is disadvantageous given that we want to construct stochastic integrals in the sense of Walsh. However, if we choose X^ as an Ornstein-Uhlenbeck process, we obtain a predictable stochastic volatility process.Example 16 (covariance function ofτ when X^ is an Ornstein-Uhlenbeck process).
Suppose thatX^ is an Ornstein-Uhlenbeck process with rate parameter κ>0 and generated by a Lévy process X. Then, using the notations in Definition 15, we have (29)Covτtx,τt~x~=12VarY1λ-1e-λt-t~-ηx+x~+VarX^0η-1e-κt-t~-ηx-x~-VarX^0η-1e-κt-t~-ηx+x~.Furthermore, if Var(Y1)=Var(X^0) and κ=λ=η, then for fixed x and x~ the covariance function of τ is (30)Covτtx,τt~x~=e-κt-t~+x-x~.


### 3.3.2. Volatility Modulation via Extended Subordination

An alternative way of volatility modulation is by means of extended subordination. In this modulation, the volatilityσt2(x)=τt(x) is incorporated in the modeling through a meta-time T associated with an absolutely continuous measure T on S; namely, (31)TA=∫Aτtxdxdt,where A∈Bb(S). A natural choice of meta-time is T(A)=T(x,t)∣(x,t)∈A, where T is the mapping given by T(x,t)=x,τt+(x), while (32)τt+x=∫0tτzxdz.The above construction of a meta-time in a space-time model is very general, since it allows defining a variety of models for the random field τt(x), leading to new model specifications. For instance, one could model τt(x) by an OUTVF or any other model discussed previously. It follows that this way to modulate the volatility is more general than the previous one.


### 3.3.3. Volatility Modulation via Probability Mixing and Lévy Mixing

Volatility modulation can also be obtained through a probability mixing approach which is realized defining new distributions by randomising a parameter from a given parametric distribution, as in the following classical example.Example 17 (normal variance-mean mixture distribution).
Consider as base model(33)Ytx=∫Atxhx,t;ξ,zLdξ,dz,assuming that the corresponding Lévy basis L is homogeneous and Gaussian; namely, the corresponding Lévy seed is given by L′~Nη+βσ2,σ2 with η,β,σ∈R and σ2>0. Then, we use probability mixing supposing that σ2 is random. Hence, the conditional law of the Lévy seed is given by L′∣σ~Nη+βσ2,σ2. Due to the scaling property of the Gaussian distribution, such a model can be represented as in (2); therefore, in this particular case, probability mixing and volatility modulation via a stationary stochastic integrand produce the same result. Such a construction falls into the class ofnormal variance-mean mixtures.In this context, it is important to note that probability mixing does not generally lead to infinitely divisible distributions. In order to overcome this problem, we present the concept of Lévy mixing. LetL be afactorisable Lévy basis on Rn withcharacteristic quadruplet(a~,m~,ν(dx,·),c) (see Subsection 2.2 of [10] for a rigorous definition of characteristic quadruplet of a Lévy basis and factorisable Lévy basis). Suppose that ν(dx,·) depends on the parameter θ∈Π, where Π denotes the space of the parameter θ. Then, theLévy measure of L is given by ν(dx,θ)c(ds). Let γ be a measure on Π and define (34)n^dx,ds=∫Πνdx,θγdθcds,where we assume that (35)∫R1∧x2n^dx,ds<∞;then there exists a Lévy basis L^ which has n^ as its Lévy measure. We define L^ as theLévy basis obtained by Lévy mixing ofLwith the measureγ. Let us study a concrete example of Lévy mixing.Example 18 (superposition of Ornstein-Uhlenbeck type process).
Let us consider the example of superposition of Ornstein-Uhlenbeck type process. LetL be a subordinator without drift and with Lévy measure νL, and consider an Ornstein-Uhlenbeck process (36)Yt=∫-∞te-θt-zdLz,for θ>0. By a simple calculation, we obtain the following expression for the cumulant function of Yt: (37)Cζ‡Yt=∫0∞eiζx-1νdx,θ,for ζ∈R, where (38)νdx,θ=∫0∞νLeθudxduis a mixture of νL with the Lebesgue measure. Lévy mixing can be carried out with respect to the parameter θ; namely, (39)ν~dx=∫0∞νdx,θγdθ,where γ is a measure on 0,∞ satisfying(40)∫0∞xν~dx<∞.Taking L~ to be the Lévy basis with Lévy measure νL(dx)duγ(dθ) and defining the superposition of Ornstein-Uhlenbeck type process Y~t with respect to L~ by (41)Y~t=∫0∞∫-∞te-λt-zL~dz,dλ,where λ>0, we have that the cumulant function of Y~t is given by (42)Cζ‡Y~t=∫0∞∫-∞tCζe-θt-z‡L1dzγdθ=∫0∞∫0∞eiζe-θu-1νLdxduγdθ=∫0∞eiζx-1ν~dx,hence showing that superposition of Ornstein-Uhlenbeck type process can be obtained from an Ornstein-Uhlenbeck process through Lévy mixing.


## 3.3.1. Volatility Modulation via a Stochastic Integrand

This method has already been suggested in the initial definition of ambit fields (see (2)). In that case, there are essentially two approaches that can be used to construct a stochastic volatility field. In fact, we can specify the stochastic field directly as a random field, for example, as another ambit field, or we can start with a purely temporal or purely spatial stochastic volatility process, respectively, and then suitably generalise it to a random field. In the following, we will present examples of both types of construction. First, we focus on the modeling approach where we directly specify a random field for the volatility field. A natural starting point for modeling the volatility is to combine kernel smoothing of a homogeneous Lévy basis with (nonlinear) transformation to ensure positivity. For instance, let (22)σt2x=V∫R×RnKx,t-z;ξLσdξ,dz,where Lσ is a homogeneous Lévy basis independent of L, K:Rn×R×Rn→R+ is an integrable kernel function satisfying K(x,u;ξ)=0 for u<0, and V:R→R+ is a continuous, nonnegative function. Clearly, the kernel function K determines the space-time autocorrelation structure of the volatility field σ2.Remark 12 (about stochastic volatility).
Note thatσ2 defined by (22) is stationary in the temporal dimension and if we have (23)Kx,t-z;ξ=K∗x-ξ,t-z,for some suitable function K∗, then the stochastic volatility is both stationary in time and homogeneous in space.Example 13 (X2 distribution).
Consider stochastic volatility as in (22), let Lσ be a standard normal Lévy basis, and let V(x)=x2. Then, σz2(ξ) is clearly positive and has a pointwise X2 distribution with one degree of freedom.In what follows, we show how to construct a stochastic volatility field by extending a purely temporal stochastic process to a random field. Note that our objective is to construct a stochastic volatility field which is stationary, at least in the temporal dimension. There are many possibilities to reach the latter goal, but in the following we decide to focus our attention on a particularly relevant one, namely, theOrnstein-Uhlenbeck type volatility field (OUTVF). The choice of an Ornstein-Uhlenbeck process as the stationary base component is motivated by its analytical tractability and also because it tends to perform well in practice, at least in the purely temporal case. Without loss of generality, we restrict our attention to the case n=1; that is, we consider one spatial dimension. In order to present this type of volatility fields, we have to define theLévy supraprocesses.Definition 14 (Lévy supraprocess generated byX^).
Suppose thatX^=X^tt∈R is a stationary, positive, and infinitely divisible process on R. Let Z·=Zx·x∈R+ be a family of stationary processes such that Z· has independent increments and for each x the cumulant function of Zx· is given by (24)Cm‡Zx·=xCm‡X^,where (25)Cm‡X^=log⁡EeimX^,with (26)mX^=∫X^zmdz,for an arbitrary signed measure m on R. Then, for any fixed t^∈R, Zt^=Zxt^x∈R+ is a Lévy process calledLévy supraprocess generated byX^.Definition 15 (OUTVF).
Suppose thatY^ is a positive Ornstein-Uhlenbeck process with rate parameter λ>0 and such that it is generated by a Lévy subordinator Y; that is, (27)Y^t=∫-∞te-λt-zdYz.A stochastic volatility field σt2(x) on R×R is an OUTVF if it is defined as (28)τtx=σt2x=e-ηxY^t+∫0xe-ηx-ξdZξt,where η>0 is the spatial rate parameter and Z=Z·tt∈R+ is a family of Lévy supraprocesses.Note that, in the above construction, we start from an Ornstein-Uhlenbeck process in time,τt(0)=Y^t, and the spatial structure is then introduced through two steps. First, we multiply the process with an exponential weight e-ηx dependent on a spatial component, which reaches its maximum for x=0, decaying as it moves away from the purely temporal case. Second, we add an integral which is similar to an Ornstein-Uhlenbeck process in the spatial variable x. However, although the stochastic volatility field τ is stationary in time, note here that the integration starts from 0 rather than from -∞; hence, the resulting component is not stationary in the spatial variable x, even if this could be changed without affecting the structure if such a property is required in a particular application. Finally, note that the process τt(x) is in general not predictable, which is disadvantageous given that we want to construct stochastic integrals in the sense of Walsh. However, if we choose X^ as an Ornstein-Uhlenbeck process, we obtain a predictable stochastic volatility process.Example 16 (covariance function ofτ when X^ is an Ornstein-Uhlenbeck process).
Suppose thatX^ is an Ornstein-Uhlenbeck process with rate parameter κ>0 and generated by a Lévy process X. Then, using the notations in Definition 15, we have (29)Covτtx,τt~x~=12VarY1λ-1e-λt-t~-ηx+x~+VarX^0η-1e-κt-t~-ηx-x~-VarX^0η-1e-κt-t~-ηx+x~.Furthermore, if Var(Y1)=Var(X^0) and κ=λ=η, then for fixed x and x~ the covariance function of τ is (30)Covτtx,τt~x~=e-κt-t~+x-x~.


## 3.3.2. Volatility Modulation via Extended Subordination

An alternative way of volatility modulation is by means of extended subordination. In this modulation, the volatilityσt2(x)=τt(x) is incorporated in the modeling through a meta-time T associated with an absolutely continuous measure T on S; namely, (31)TA=∫Aτtxdxdt,where A∈Bb(S). A natural choice of meta-time is T(A)=T(x,t)∣(x,t)∈A, where T is the mapping given by T(x,t)=x,τt+(x), while (32)τt+x=∫0tτzxdz.The above construction of a meta-time in a space-time model is very general, since it allows defining a variety of models for the random field τt(x), leading to new model specifications. For instance, one could model τt(x) by an OUTVF or any other model discussed previously. It follows that this way to modulate the volatility is more general than the previous one.


## 3.3.3. Volatility Modulation via Probability Mixing and Lévy Mixing

Volatility modulation can also be obtained through a probability mixing approach which is realized defining new distributions by randomising a parameter from a given parametric distribution, as in the following classical example.Example 17 (normal variance-mean mixture distribution).
Consider as base model(33)Ytx=∫Atxhx,t;ξ,zLdξ,dz,assuming that the corresponding Lévy basis L is homogeneous and Gaussian; namely, the corresponding Lévy seed is given by L′~Nη+βσ2,σ2 with η,β,σ∈R and σ2>0. Then, we use probability mixing supposing that σ2 is random. Hence, the conditional law of the Lévy seed is given by L′∣σ~Nη+βσ2,σ2. Due to the scaling property of the Gaussian distribution, such a model can be represented as in (2); therefore, in this particular case, probability mixing and volatility modulation via a stationary stochastic integrand produce the same result. Such a construction falls into the class ofnormal variance-mean mixtures.In this context, it is important to note that probability mixing does not generally lead to infinitely divisible distributions. In order to overcome this problem, we present the concept of Lévy mixing. LetL be afactorisable Lévy basis on Rn withcharacteristic quadruplet(a~,m~,ν(dx,·),c) (see Subsection 2.2 of [10] for a rigorous definition of characteristic quadruplet of a Lévy basis and factorisable Lévy basis). Suppose that ν(dx,·) depends on the parameter θ∈Π, where Π denotes the space of the parameter θ. Then, theLévy measure of L is given by ν(dx,θ)c(ds). Let γ be a measure on Π and define (34)n^dx,ds=∫Πνdx,θγdθcds,where we assume that (35)∫R1∧x2n^dx,ds<∞;then there exists a Lévy basis L^ which has n^ as its Lévy measure. We define L^ as theLévy basis obtained by Lévy mixing ofLwith the measureγ. Let us study a concrete example of Lévy mixing.Example 18 (superposition of Ornstein-Uhlenbeck type process).
Let us consider the example of superposition of Ornstein-Uhlenbeck type process. LetL be a subordinator without drift and with Lévy measure νL, and consider an Ornstein-Uhlenbeck process (36)Yt=∫-∞te-θt-zdLz,for θ>0. By a simple calculation, we obtain the following expression for the cumulant function of Yt: (37)Cζ‡Yt=∫0∞eiζx-1νdx,θ,for ζ∈R, where (38)νdx,θ=∫0∞νLeθudxduis a mixture of νL with the Lebesgue measure. Lévy mixing can be carried out with respect to the parameter θ; namely, (39)ν~dx=∫0∞νdx,θγdθ,where γ is a measure on 0,∞ satisfying(40)∫0∞xν~dx<∞.Taking L~ to be the Lévy basis with Lévy measure νL(dx)duγ(dθ) and defining the superposition of Ornstein-Uhlenbeck type process Y~t with respect to L~ by (41)Y~t=∫0∞∫-∞te-λt-zL~dz,dλ,where λ>0, we have that the cumulant function of Y~t is given by (42)Cζ‡Y~t=∫0∞∫-∞tCζe-θt-z‡L1dzγdθ=∫0∞∫0∞eiζe-θu-1νLdxduγdθ=∫0∞eiζx-1ν~dx,hence showing that superposition of Ornstein-Uhlenbeck type process can be obtained from an Ornstein-Uhlenbeck process through Lévy mixing.


## 3.4. Autocorrelation and Cross-Correlation

We would like to underline that the proposed approach allows modeling the entire forward curve. Hence, it is interesting to study the correlation structure between various forward contracts implied by the presented modeling framework. Starting from the definition of covariance function of an ambit field (see SectionA.4 of [5]), we deduce theautocorrelation function of a forward contract ft(x). More precisely, assuming that the Lévy basis is homogeneous, and taking t, h∈R, h≥0, and x, x~≥0, we have that the ambit set defined in (8) can be expressed as At(x)∩At+h(x~)=A-t(x)≔A-t; then we define the autocorrelation function as (43)Corftx,ft+hx~≔∫0,∞×0,∞kξ,u;xkξ,u+h;x~Eσ02ξdξdu∫0,∞×0,∞k2ξ,u;xEσ02ξdξdu∫0,∞×0,∞k2ξ,u′;x~Eσ02ξdξdu′.Note that the correlation structure is principally determined by three factors, which are the intersection of the corresponding ambit sets, the kernel function, and the autocorrelation structure of the stochastic volatility field. Furthermore, we could model various commodity forward contracts, such as electricity and natural gas futures, simultaneously. It follows that we can specify different ambit sets, kernel functions, stochastic volatility fields, and Lévy bases in order to obtain a rather flexible correlation structure. In such a situation, it becomes even more clear how flexible the ambit setting is. The details of these multivariate extensions can be found in Section A.5 of [5].


## 3.5. Martingale Condition

In what follows, we give sufficient conditions for the presented model to be a martingale, even if such property may be too restrictive if we are interested in modeling electricity forward contracts. In fact, there are at least two arguments that support the choice of more general classes of stochastic processes than martingales. First, in the energy context, it might not be as crucial thatft(T-t) is a martingale as it is in the context of modeling interest rates. In fact, in many emerging electricity markets, a seller may not be able to find any buyer to get rid of a forward; oppositely, it may also happen that there are no sellers when one wants to buy a contract. Hence, the illiquidity prevents possible arbitrage opportunities from being exercised and we need a martingale condition only in the absence of arbitrage setting (see [18]). Second, independently of the particular structure of energy markets, the recent literature in mathematical finance (see, e.g., [6, 18]) has highlighted that some classes of nonsemimartingales, in particular stochastic processes with conditional full support, do not necessarily give rise to arbitrage opportunities reflecting more realistic characteristics of the market. In the null-spatial setting, Pakkanen proves thatBrownian semistationary processes have in fact conditional full support (see [7]). In future research, it will be interesting to study extensions of the latter result to more general ambit fields. However, the question of establishing martingale conditions for ambit fields is interesting in order to point out a subset of the class we are considering, which can be exploited to model general forward prices, hence without restricting ourselves to the electricity case. In what follows, we formulate the martingale condition for ambit fields more general than those defined by (7), where the ambit set At(x) is chosen as in (8), also showing that such condition simplifies in the modeling framework described in Section 3.1.Theorem 19 (martingale condition for general ambit fields).
Letx=T-t for some T>0 and, for fixed t∈R, let (44)Ytx≔YtT-t=∫Atxhξ,z;T-t,tσzξLdξ,dz,where At(x) is defined as in (8), h is a deterministic kernel function, σ is an adapted, nonnegative random field, and L is a Lévy basis satisfying both (i) and (iv). Further, σ and L are assumed to be independent. Then, Yt(T-t)t∈R is a martingale with respect to Ftt∈R defined as in (5) if and only if, for all ξ>0,z≤t≤T, it holds (45)hξ,z;T-t,t=h~ξ,z;T,for some deterministic kernel function h~.Proof.
Measurability and integrability properties are trivially verified according to what we have seen so far; therefore, we just have to show that(46)EYtT-t∣Ft~=Yt~T-t~,for all t~≤t. Note that, for t~≤t, we have that At~(x)⊂At(x). Using the independence property of σ and L and the fact that L is a zero-mean process, we find (47)EYtT-t∣Ft~=E∫At~xhξ,z;T-t,tσzξLdξ,dz+∫Atx∖At~xhξ,z;T-t,tσzξLdξ,dz∣Ft~=∫At~xhξ,z;t,T-tσzξLdξ,dz=Yt~T-t~+It~T-t~, where (48)It~T-t~=∫At~xhξ,z;T-t,t-hξ,z;T-t~,t~σzξLdξ,dz.Without loss of generality, we can assume that Var(L)=1. Since L is a Lévy basis with zero-mean, we know that EIt~(T-t~)=0 and from theItô isometry (see Theorem 2.5 in [13]) we therefore get that (49)VarIt~T-t~=∫At~xhξ,z;T-t,t-hξ,z;T-t~,t~2Eσz2ξc~dξ,dz,where c~ is thecovariance measure induced by L; see, for example, [13], for a rigorous definition of covariance measure induced by a Lévy basis. Thus, in order to obtain It~(T-t~)=0, we need that, for all ξ>0,z≤t~≤t≤T, (50)hξ,z;T-t,t=hξ,z;T-t~,t~.Note that there is only one class of functions which satisfy the above condition, namely, functions of the form (51)hξ,z;T-t,t=h~ξ,z;T,for all ξ>0, z≤t~≤t≤T and for some deterministic kernel function h~.Remark 20 (Lévy basis with mean different from zero).
If we would like to work with Lévy basesL which do not have zero-mean, then the martingale condition has to be extended by an additional drift condition.Corollary 21 (martingale condition forft(x)).
One gets thatft(T-t)t∈R is a martingale with respect to Ftt∈R if and only if, for all ξ>0, z≤t≤T, the following holds: (52)kξ,t-z;T-t=k~ξ,T-z,for some deterministic kernel function k~; hence, the latter is a special case of the factorisation given in (19).Proof.
This result follows directly from Theorem19 considering Yt(T-t)=ft(T-t).Clearly, the martingale condition is rather strong and it is necessary to check whether there are actually any relevant cases, not excluded by condition (52). Therefore, let us study some examples.Example 22 (exponential kernel function).
If we choosek to be of the form (53)kξ,t-z;x=kξ,t-z;T-t=e-αξ+x+t-z=e-αξ+T-z,for some α>0, then the above martingale condition is clearly satisfied. Note that this choice of the kernel function belongs to both factorisation classes (12) and (19).Example 23 (factorisation (19)).
Let us consider kernel functionsk which factorise as in (19); that is, (54)kξ,t-z;x=ΦξΨt-z,x;then the choice of the function Φ does not have any impact concerning whether the ambit field is a martingale; indeed such a property is solely determined by the function Ψ. It follows that every choice of the form (55)Ψt-z,T-t=Ψ~T-zsatisfies the martingale condition.


## 3.6. Derivation of a Spot Model from the Forward One

Under the absence of arbitrage assumption, the forward price of a contract tends to be equal to the spot price of the underlying as the time to maturity tends to zero. In this subsection, we investigate in detail the nature of the spot price model implied by the above modeling framework for the forward price. Note that this study should be understood as a theoretical exercise, since we do not observe convergence of the electricity forward price to the electricity spot price in real markets.Letft(x) be the forward price and let St be the spot price of the underlying. When the time to maturity x=T-t tends to 0, we would have that ft(x) tends to St; hence, we define the spot price as (56)St≔∫-∞t∫0∞kξ,t-z;0σzξLdξ,dz.Lemma 24 (convergence offt(x) to St).
Suppose that(57)limx↓0⁡∫-∞t∫0∞kξ,t-z;x-kξ,t-z;02Eσz2ξdξdz=0;then ft(x)→St in L2Ω,F,P∗ as time to maturity x tends to zero, where P∗ denotes the risk neutral probability measure.Proof.
The result follows immediately by applying the Itô isometry, as proved by Theorem2.5 in [13].The previous lemma gives us that the forward price will tend continuously in variance to the spot price as time to maturity decreases to zero. Note also that wheneverσz(ξ) is a stationary field, the condition in the above lemma is translated to convergence of k(·,·;x) to k(·,·;0) in L2R+2.Using Proposition8 for the forward price, we can easily derive the conditional cumulant function for the implied spot price.Proposition 25 (conditional cumulant function ofSt).
LetL be a homogeneous Lévy basis. Then, for St defined as in (56), the cumulant function conditioned to σ is given by (58)Cσζ‡St=∫-∞t∫0∞Cζkξ,t-z;0σzξ‡L′dξdz,where L′ is the Lévy seed associated with L.Proof.
This result is an immediate consequence of Proposition8.A case of some special interest is given by considering the driving Lévy basisL of the ambit field to be a homogeneous Gaussian Lévy basis; in fact, in this case, we get the following result.Corollary 26 (spot price whenL is a Gaussian Lévy basis).
LetSt be as in (56) and let L be a homogeneous Gaussian Lévy basis with (59)Cζ‡L′=-12ζ2,where L′ is the Lévy seed associated with L. Then, the cumulant function of St conditioned to σ is given by (60)Cσζ‡St=-12ζ2∫-∞t∫0∞k2ξ,t-z;0σz2ξdξdz.If the kernel function k can be factorized as in (12), then (61)Cσζ‡St=-12ζ2∫-∞tψ2t-zωz2ξdz,where (62)ωz2=∫0∞ϕ2ξ,0σz2ξdξ,implying that St is equal in law to (63)St=∫-∞tψ2t-zωz-ξdWz,where W is the standard Brownian motion.Proof.
The result follows by straightforward calculations, exploiting Proposition25 and the definition of Gaussian Lévy basis.Remark 27 (modeling spot price).
Note thatSt in (63) is a Brownian semistationary process, which has been used as a model for energy spot prices in [1].


### 3.6.1. The Samuelson Effect

Recall that the Samuelson effect describes the empirical fact that the volatility of the forward price increasingly converges to the volatility of the underlying spot price when the time to maturity tends to zero. Here, we show that this effect is naturally included in our modeling framework.Theorem 28 (nondecreasing of the conditional variance offt(x)).
Assume that the functionx↦k(ξ,u;x) is monotonically nondecreasing in x≥0 for every (ξ,u)∈R+2. The conditional variance of the forward price ft(x) given by (64)vtx≔c∫-∞t∫0∞k2ξ,t-z;xσz2ξdξdzis monotonically nondecreasing in x, for t≥0, where the term c is a suitable constant; see Section A.2 of [5] for a rigorous definition of c.Proof.
Forx≥0, we have (65)vtx=c∫-∞t∫0∞k2ξ,t-z;xσz2ξdξdz,and, taking 0≤x≤x~, then (66)vtx~-vtx=∫-∞t∫0∞k2ξ,t-z;x~-k2ξ,t-z;xσz2ξdξdz≥0,since k(ξ,t-z;x) is nondecreasing in x.Remark 29 (Samuelson effect).
The conditional variance of the spot price is given byvt(0) and from the above theorem it follows that vt(x)≥vt(0) for x≥0. Since we can consider vt(x) as a monotonically increasing sequence bounded from below by vt(0), there exists the limit limx↓0vt(x)=l. Under the condition in Lemma 24, this limit is given by l=vt(0); hence, it equals the conditional variance of spot price, which implies that we have a Samuelson effect.


## 3.6.1. The Samuelson Effect

Recall that the Samuelson effect describes the empirical fact that the volatility of the forward price increasingly converges to the volatility of the underlying spot price when the time to maturity tends to zero. Here, we show that this effect is naturally included in our modeling framework.Theorem 28 (nondecreasing of the conditional variance offt(x)).
Assume that the functionx↦k(ξ,u;x) is monotonically nondecreasing in x≥0 for every (ξ,u)∈R+2. The conditional variance of the forward price ft(x) given by (64)vtx≔c∫-∞t∫0∞k2ξ,t-z;xσz2ξdξdzis monotonically nondecreasing in x, for t≥0, where the term c is a suitable constant; see Section A.2 of [5] for a rigorous definition of c.Proof.
Forx≥0, we have (65)vtx=c∫-∞t∫0∞k2ξ,t-z;xσz2ξdξdz,and, taking 0≤x≤x~, then (66)vtx~-vtx=∫-∞t∫0∞k2ξ,t-z;x~-k2ξ,t-z;xσz2ξdξdz≥0,since k(ξ,t-z;x) is nondecreasing in x.Remark 29 (Samuelson effect).
The conditional variance of the spot price is given byvt(0) and from the above theorem it follows that vt(x)≥vt(0) for x≥0. Since we can consider vt(x) as a monotonically increasing sequence bounded from below by vt(0), there exists the limit limx↓0vt(x)=l. Under the condition in Lemma 24, this limit is given by l=vt(0); hence, it equals the conditional variance of spot price, which implies that we have a Samuelson effect.


## 3.7. Simulation Algorithm

In what follows, we present a possible algorithm, taken from [5], to simulate an ambit field Yt(x) of the following type: (67)Ytx=∫-∞t∫0∞kξ,t-z;xσzξLdξ,dz.A problem that we immediately note is that we perform integration over an ambit set At(x)=-∞,t×0,∞, which is unbounded; hence, we have to truncate the ambit set before performing any computations. Therefore, we introduce the ambit field (68)Y-tx=∫M1t∫0M2kξ,t-z;xσzξLdξ,dz,where M1<t and M2>0. Obviously, if M1→-∞ and M2→∞, then Y-t(x) converges to Yt(x).Takingn,m∈N+, we construct a grid for M1,t×0,M2 by dividing the interval M1,t into n equidistant intervals of length (t-M1)/n, where t=t1>t2>⋯>tn=M1, and by dividing 0,M2 into m equidistant intervals of length M2/m, where 0=x1<x2<⋯<xm=M2. Next, we simulate an ambit field through the following steps:(1)
Simulate the stochastic field on the grid points(ti,xj) in order to obtain the values of stochastic volatility σti(xj) for i=1,…,n and j=1,…,m. In case we have no stochastic volatility, we set σti(xj)=1, for all i=1,…,n and for all j=1,…,m.(2)
Simulaten·m random variables Zi+j~L(Δ), where L is the Lévy basis computed on (69)Δ=Δn,m,M1,M2=t-M1nM2m,fori=1,…,n and j=1,…,m.(3)
Compute the approximated ambit field in this way:(70)Y^tx=∑i=1n-1∑j=1m-1kxj,t-ti;xσtixjZi+j.We note that the previous algorithm can be used also to simulate stochastic volatility. In the last step of the simulation algorithm, we use the definition of stochastic integral in the sense of Walsh for simple processes, which ensures the convergence of (70) to the actual stochastic integral, when Δ↓0.


## 4. Empirical Study

In what follows, we exploit results and techniques developed throughout previous sections to give a detailed analysis as well as a forecast analysis concerning forward contracts within the German energy markets, namely, within the framework of the EEX. The latter ones are characterized by various types of forward contracts whose structure mainly depends on the characterizing time slots each of them takes into consideration, the duration of the contract, whether there are some lower or upper bounds in the quantity of electricity exchanged, and so forth. We focus our attention on a widely used forward contract that assures the supply of electricity at a fixed price for a month from 9 a.m. to 8 p.m., but only during working days. The latter type of contract is traded in theover-the-counter market and it is called a monthly peak forward contract. Trades for such contracts take place in the three months before the effective supply of requested electricity. In particular, we analyze twelve peak contracts, one for each month of the year, from April 2012 to March 2013.


### 4.1. Model Specification and Parameters Calibration

Following the approach presented in Section3.1, we model the forward price ft(x) as an ambit process for each fixed x; namely, (71)ftx=∫-∞t∫Iξkξ,t-z;xσzξLdξ,dz,where t∈R denotes the current time, I(ξ)⊂0,∞ is an interval that depends on the peak contract considered, σz(ξ)>0 denotes the volatility of the forward market as a whole, L is a Lévy basis, and k is the kernel function. We also suppose that ft(x) satisfies conditions (i)–(vi) stated in Section 3.1, and T>0 being the maturity time, we model the fact that the forward price depends also on time to maturity x≔T-t; hence, we use thespatial variable x to define a second time variable, namely, thetime to maturity one.There are three components of the model which we have to specify, namely, the Lévy basisL, the kernel function k, and the volatility σ. The choices we made are mainly based on qualitative analysis of the market of interest, but also taking into consideration analytical as well as numerical tractability of the resulting model and, last but not least, the empirical evidence.


#### 4.1.1. Kernel Function

Every demand in general, and electricity consumption in particular, is affected by seasonality, which is also the case for forward prices. In fact, they depend on the period of the year; indeed it is easy to recognize a cyclic behaviour for such values. In order to depict this kind of periodicity, in (71), we consider a kernel function of the following type: (72)kξ,t-z;T-t=f0+ΦξΨt-z,T-t,where f0 is the price of the forward contract at time 0 taken from the data, Φ takes into account the seasonal influences that affect electricity markets, and Ψ models the average price trend over the peak contracts considered, which are relative to one year. Clearly, the periodicity acts solely on the average trend of forward contracts and this is the reason why Φ multiplies only Ψ and not also f0 in (72).Exploiting related time series data and the fact that in winter the demand for electricity is greater than in summer, we choose the following bell-shaped function to model the effect of seasonality on forward prices:(73)Φξ=1+Kμ+αe-ξ-β2/γ2,where μ=0.8, α=2.4, β=8.5, and γ=1.5, while K is a constant that depends only on the length of the forecast. Short-term forecasts are less influenced by seasonality than long-term ones; therefore, the parameter K takes care of such phenomenon. In particular, we choose K=l/T, where l is the length of the forecast, so that the periodicity function Φ is significant only in medium- or long-term previsions. Moreover, the variable ξ spans the different months of the year; therefore, the integration interval I(ξ) is different for every peak contract considered (see Table 1).Table 1
In the first column, we have the months of electricity supply of the peak contracts considered, while in the second one we have the intervalI(ξ) relative to each contract. Note that a contract traded for 3 months is integrated in an interval of length 3.
Month
I
(
ξ
)
April
0,3
May
1,4
June
2,5
July
3,6
August
4,7
September
5,8
October
6,9
November
7,10
December
8,11
January
9,12
February
10,13
March
11,14Before entering into details concerning the specification ofΨ, let us note that, as a rather general rule, forward contracts tend to depreciate getting close to maturity time. Even if the behaviour of peak contracts considered is not the same, based on the latter ones, it seems reasonable to choose a concave function Ψ. In particular, we choose the following exponential function: (74)Ψt-z,T-t=ηe-λT-z,ifz≥0,0,ifz<0,where η is a constant computed as the average daily price variation of all peak contracts, which is equal to -0.063, times the length of the forecast l and λ is a constant that models the price fall. We have estimated this parameter studying the average fall in price over the peak contracts and the trends of the latter ones, obtaining λ=0.04. This choice of kernel function is also motivated by the empirical studies in [1], where the authors show that such a specification in a model based onLévy semistationary processes effectively forecasts spot price dynamics on the German EEX market.


#### 4.1.2. Lévy Basis

Several publications (see, e.g., [2, 5, 19, 20]) have shown the applicability of the generalised hyperbolic distribution, in particular the subfamily consisting of the normal inverse Gaussian distributions, to describe appropriately financial datasets. Exploiting these results, we choose a normal inverse Gaussian Lévy basis plugging it in the model described by (71). This kind of Lévy basis can be obtained starting from a normal inverse Gaussian Lévy process of parameters θ, k, and ω, and the latter is the result of subordination of an inverse Gaussian Lévy process, which has variance k, with a Brownian motion with volatility ω and drift θ. In order to use the right parameters, note that the prices of peak contracts are quite stable, being characterized by small daily variation; hence, in general, there are no fat tail phenomena to be taken into account, at least at such a time scale. Latter observation suggests taking small values of ω and k. In particular, also using lines of empirical evidence coming from the analysis of related time series, we take ω=0.01 and k=0.2, also setting θ=1, since we model the entire forward curve.


#### 4.1.3. Stochastic Volatility

The time series of monthly peak forward contracts we have taken into consideration show a rather regular trend. Therefore, there is no empirical evidence that causes us to consider this term. Furthermore, the complexity of the resulting model increases by introducing stochastic volatility; hence, we have followed what is suggested in [5], and we have not introduced it, instead fixing σz(ξ)≡1 in (71).


### 4.2. Short-Term Analysis

In this short-term analysis, we take the price of the peak contract at a day to forecast the value of the contract from1 to 5 days ahead. We perform this study simulating numerically the trajectories of (71) using the algorithm presented in Section 3.7. The results are listed in Table 2 and are discussed in Section 4.5.Table 2
In the first column, we have the months of electricity supply of the peak contracts considered, while in the next five columns we have the mean absolute percentage errors (MAPE) obtained in the forecasts of the contracts value at1,2,3,4, and 5 days ahead, respectively.
Month
1
2
3
4
5
April 2012
1.10%
1.81%
2.12%
2.42%
2.87%
May 2012
0.84%
1.28%
1.62%
1.90%
2.01%
June 2012
0.72%
0.97%
1.13%
1.24%
1.33%
July 2012
0.69%
1.05%
1.25%
1.48%
1.54%
August 2012
0.62%
0.84%
1.13%
1.35%
1.51%
September 2012
0.73%
1.13%
1.42%
1.60%
1.70%
October 2012
0.44%
0.62%
0.79%
0.88%
0.95%
November 2012
0.59%
0.82%
0.92%
1.00%
1.07%
December 2012
0.60%
0.80%
0.99%
1.11%
1.17%
January 2013
0.62%
0.88%
1.09%
1.25%
1.39%
February 2013
0.68%
1.09%
1.54%
1.89%
2.19%
March 2013
0.61%
0.82%
1.09%
1.30%
1.61%


### 4.3. Long-Term Analysis

In the short-term analysis, we have taken the price at a day to forecast the value of the contract some days ahead; instead, in order to study the long-term behaviour, we use the first observed price of each peak contract to predict the whole trend of it. In this case, a very important role is played by seasonality. We perform this study simulating numerically the trajectories of (71), using the algorithm presented in Section 3.7. The results are listed in Tables 3 and 4 and are discussed in Section 4.5.Table 3
In the first column, we have the months of electricity supply of the peak contracts considered, while in the second column we have the mean absolute percentage errors obtained in the long-term forecasts of the contracts value.
Month
MAPE
April 2012
2.96%
May 2012
1.94%
June 2012
1.73%
July 2012
1.88%
August 2012
2.86%
September 2012
2.36%
October 2012
1.97%
November 2012
1.99%
December 2012
1.79%
January 2013
2.79%
February 2013
4.15%
March 2013
4.95%Table 4
In the first column, we have the months of electricity supply of the peak contracts considered, while in the next three columns we have the mean absolute percentage errors obtained in the long-term forecasts of the contracts value during the first, second, and third month of trading, respectively.
Month
1st month
2nd month
3rd month
April 2012
3.98%
2.25%
2.65%
May 2012
2.10%
1.30%
2.39%
June 2012
0.97%
2.42%
1.79%
July 2012
0.68%
0.75%
3.99%
August 2012
2.24%
3.90%
2.44%
September 2012
1.40%
0.75%
4.69%
October 2012
0.62%
3.25%
2.02%
November 2012
1.44%
0.97%
3.49%
December 2012
0.83%
2.49%
2.03%
January 2013
0.70%
1.87%
5.79%
February 2013
2.20%
4.99%
5.15%
March 2013
0.74%
4.39%
9.48%


### 4.4. Ambit Approach versus a Standard Approach

In order to assess whether the ambit approach proposed is convenient compared to a standard approach, we have performed the same analyses proposed in Sections4.2 and 4.3 modeling the peak contracts with ageometric Brownian motion, that is, the stochastic process St that satisfies the following stochastic differential equation: (75)dSt=μStdt+σStdWt,where Wt is a standard Brownian motion, μ is a constant calledpercentage drift, and σ is a constant calledpercentage volatility. For an arbitrary initial value S0, (75) has the following analytic solution: (76)St=S0eμ-σ2/2t+σWt.In these analyses, we have estimated the parameters μ and σ using themethod of moments estimation. The results are listed in Tables 5, 6, and 7 and are discussed in Section 4.5.Table 5
In the first column, we have the months of electricity supply of the peak contracts considered, while in the next five columns we have the mean absolute percentage errors (MAPE) obtained in the forecasts of the contracts value at1,2,3,4, and 5 days ahead, respectively, using a geometric Brownian motion.
Month
1
2
3
4
5
April 2012
4,00%
3,96%
3,98%
3,98%
4,04%
May 2012
3,04%
3,03%
3,02%
3,05%
3,09%
June 2012
2,61%
2,62%
2,63%
2,66%
2,68%
July 2012
2,04%
2,05%
2,04%
2,10%
2,08%
August 2012
1,97%
1,98%
1,98%
2,01%
2,02%
September 2012
3,62%
3,66%
3,67%
3,69%
3,71%
October 2012
2,03%
2,03%
2,04%
2,05%
2,10%
November 2012
3,13%
3,15%
3,18%
3,23%
3,26%
December 2012
2,15%
2,17%
2,17%
2,19%
2,22%
January 2013
2,24%
2,26%
2,29%
2,32%
2,33%
February 2013
5,05%
5,04%
5,10%
5,14%
5,21%
March 2013
5,21%
5,28%
5,32%
5,41%
5,52%Table 6
In the first column, we have the months of electricity supply of the peak contracts considered, while in the second column we have the mean absolute percentage errors obtained in the long-term forecasts of the contracts value using a geometric Brownian motion.
Month
MAPE
April 2012
4,60%
May 2012
4,22%
June 2012
3,76%
July 2012
2,12%
August 2012
2,44%
September 2012
4,27%
October 2012
2,21%
November 2012
4,25%
December 2012
3,10%
January 2013
2,07%
February 2013
5,23%
March 2013
7,83%Table 7
In the first column, we have the months of electricity supply of the peak contracts considered, while in the next three columns we have the mean absolute percentage errors obtained in the long-term forecasts of the contracts value during the first, second, and third month of trading, respectively, using a geometric Brownian motion.
Month
1st month
2nd month
3rd month
April 2012
4,13%
3,98%
5,69%
May 2012
4,33%
3,58%
4,72%
June 2012
4,16%
2,87%
4,21%
July 2012
1,94%
1,91%
2,48%
August 2012
2,07%
2,94%
2,32%
September 2012
4,12%
3,51%
5,11%
October 2012
1,81%
2,42%
2,40%
November 2012
3,95%
3,41%
5,35%
December 2012
3,49%
2,54%
3,27%
January 2013
1,61%
1,78%
2,84%
February 2013
3,69%
4,22%
7,53%
March 2013
7,50%
6,09%
9,79%


### 4.5. Conclusions and Remarks

In the short-term framework, the forecasts are very little affected by seasonality and the model presented in this section reflects this fact. Looking closely at Table2, we note that in the short-term prediction the error tends to grow when the length of the forecast increases, as we can expect, but it remains very small. The forecasts obtained for the first day ahead case have always an error less than 1.1% and forecasts on 5 days ahead only once overtake a 2.5% error. Such good performances gained by our approach are slightly facilitated by the small variance characterizing the chosen Lévy basis. We also note that in such short-term analysis the trend of the simulated trajectories is very close to the one of real data (see Figures 1 and 2).Figure 1
In this figure, we present the price curves obtained by simulating the peak contract of June2012 using (71) with the specification above stated. The graphics in the figure present forecasts of the contract value at 1,2,3, and 4 days ahead.Figure 2
In this figure, we present the price curves obtained by simulating the peak contract of November2012 using (71) with the specification above stated. The graphics in the figure present forecasts of the contract value at 1,2,3, and 4 days ahead.Usually, the error increases enlarging the period to forecast, but in the above long-term analysis the mean absolute percentage errors are all less than 5% (see Table3), hence showing that our model provides robust and effective results in the presence of seasonality. In fact, the model often forecasts the behavior of the price some week ahead, even if not always the trend presented by the model is strictly close to the one shown by real data (see Figure 3 and Table 4). The latter result is not surprising since it is related to a “three months from now” forecast, a framework in which energy markets tend to be characterized by small sensitivity from initial data. However, the contracts with a bigger error are the ones which are more complicated to forecast such as, for example, that of February 2013.Figure 3
In this figure, we present the price curves obtained by simulating all the peak contracts using (71) with the specification above stated. The graphics in the figure present long-term forecasts of the contracts value and the month written above each graph refers to electricity supply of the contract considered.Comparing the errors obtained using the proposed model with those obtained exploiting the geometric Brownian motion, we notice that the first are significantly lower over all the performed analysis, as shown in Tables2, 3, 4, 5, 6, and 7. Accordingly, the ambit approach turns out to give better results than thestandard one. We would like to underline that the results obtained by our model provide both good approximations and effective forecasts concerning real electricity forward prices. Moreover, compared with standard approaches, our model, due to the high flexibility of the ambit stochastic setting, allows considering more easily the characteristics of the energy markets.


## 4.1. Model Specification and Parameters Calibration

Following the approach presented in Section3.1, we model the forward price ft(x) as an ambit process for each fixed x; namely, (71)ftx=∫-∞t∫Iξkξ,t-z;xσzξLdξ,dz,where t∈R denotes the current time, I(ξ)⊂0,∞ is an interval that depends on the peak contract considered, σz(ξ)>0 denotes the volatility of the forward market as a whole, L is a Lévy basis, and k is the kernel function. We also suppose that ft(x) satisfies conditions (i)–(vi) stated in Section 3.1, and T>0 being the maturity time, we model the fact that the forward price depends also on time to maturity x≔T-t; hence, we use thespatial variable x to define a second time variable, namely, thetime to maturity one.There are three components of the model which we have to specify, namely, the Lévy basisL, the kernel function k, and the volatility σ. The choices we made are mainly based on qualitative analysis of the market of interest, but also taking into consideration analytical as well as numerical tractability of the resulting model and, last but not least, the empirical evidence.


### 4.1.1. Kernel Function

Every demand in general, and electricity consumption in particular, is affected by seasonality, which is also the case for forward prices. In fact, they depend on the period of the year; indeed it is easy to recognize a cyclic behaviour for such values. In order to depict this kind of periodicity, in (71), we consider a kernel function of the following type: (72)kξ,t-z;T-t=f0+ΦξΨt-z,T-t,where f0 is the price of the forward contract at time 0 taken from the data, Φ takes into account the seasonal influences that affect electricity markets, and Ψ models the average price trend over the peak contracts considered, which are relative to one year. Clearly, the periodicity acts solely on the average trend of forward contracts and this is the reason why Φ multiplies only Ψ and not also f0 in (72).Exploiting related time series data and the fact that in winter the demand for electricity is greater than in summer, we choose the following bell-shaped function to model the effect of seasonality on forward prices:(73)Φξ=1+Kμ+αe-ξ-β2/γ2,where μ=0.8, α=2.4, β=8.5, and γ=1.5, while K is a constant that depends only on the length of the forecast. Short-term forecasts are less influenced by seasonality than long-term ones; therefore, the parameter K takes care of such phenomenon. In particular, we choose K=l/T, where l is the length of the forecast, so that the periodicity function Φ is significant only in medium- or long-term previsions. Moreover, the variable ξ spans the different months of the year; therefore, the integration interval I(ξ) is different for every peak contract considered (see Table 1).Table 1
In the first column, we have the months of electricity supply of the peak contracts considered, while in the second one we have the intervalI(ξ) relative to each contract. Note that a contract traded for 3 months is integrated in an interval of length 3.
Month
I
(
ξ
)
April
0,3
May
1,4
June
2,5
July
3,6
August
4,7
September
5,8
October
6,9
November
7,10
December
8,11
January
9,12
February
10,13
March
11,14Before entering into details concerning the specification ofΨ, let us note that, as a rather general rule, forward contracts tend to depreciate getting close to maturity time. Even if the behaviour of peak contracts considered is not the same, based on the latter ones, it seems reasonable to choose a concave function Ψ. In particular, we choose the following exponential function: (74)Ψt-z,T-t=ηe-λT-z,ifz≥0,0,ifz<0,where η is a constant computed as the average daily price variation of all peak contracts, which is equal to -0.063, times the length of the forecast l and λ is a constant that models the price fall. We have estimated this parameter studying the average fall in price over the peak contracts and the trends of the latter ones, obtaining λ=0.04. This choice of kernel function is also motivated by the empirical studies in [1], where the authors show that such a specification in a model based onLévy semistationary processes effectively forecasts spot price dynamics on the German EEX market.


### 4.1.2. Lévy Basis

Several publications (see, e.g., [2, 5, 19, 20]) have shown the applicability of the generalised hyperbolic distribution, in particular the subfamily consisting of the normal inverse Gaussian distributions, to describe appropriately financial datasets. Exploiting these results, we choose a normal inverse Gaussian Lévy basis plugging it in the model described by (71). This kind of Lévy basis can be obtained starting from a normal inverse Gaussian Lévy process of parameters θ, k, and ω, and the latter is the result of subordination of an inverse Gaussian Lévy process, which has variance k, with a Brownian motion with volatility ω and drift θ. In order to use the right parameters, note that the prices of peak contracts are quite stable, being characterized by small daily variation; hence, in general, there are no fat tail phenomena to be taken into account, at least at such a time scale. Latter observation suggests taking small values of ω and k. In particular, also using lines of empirical evidence coming from the analysis of related time series, we take ω=0.01 and k=0.2, also setting θ=1, since we model the entire forward curve.


### 4.1.3. Stochastic Volatility

The time series of monthly peak forward contracts we have taken into consideration show a rather regular trend. Therefore, there is no empirical evidence that causes us to consider this term. Furthermore, the complexity of the resulting model increases by introducing stochastic volatility; hence, we have followed what is suggested in [5], and we have not introduced it, instead fixing σz(ξ)≡1 in (71).


## 4.1.1. Kernel Function

Every demand in general, and electricity consumption in particular, is affected by seasonality, which is also the case for forward prices. In fact, they depend on the period of the year; indeed it is easy to recognize a cyclic behaviour for such values. In order to depict this kind of periodicity, in (71), we consider a kernel function of the following type: (72)kξ,t-z;T-t=f0+ΦξΨt-z,T-t,where f0 is the price of the forward contract at time 0 taken from the data, Φ takes into account the seasonal influences that affect electricity markets, and Ψ models the average price trend over the peak contracts considered, which are relative to one year. Clearly, the periodicity acts solely on the average trend of forward contracts and this is the reason why Φ multiplies only Ψ and not also f0 in (72).Exploiting related time series data and the fact that in winter the demand for electricity is greater than in summer, we choose the following bell-shaped function to model the effect of seasonality on forward prices:(73)Φξ=1+Kμ+αe-ξ-β2/γ2,where μ=0.8, α=2.4, β=8.5, and γ=1.5, while K is a constant that depends only on the length of the forecast. Short-term forecasts are less influenced by seasonality than long-term ones; therefore, the parameter K takes care of such phenomenon. In particular, we choose K=l/T, where l is the length of the forecast, so that the periodicity function Φ is significant only in medium- or long-term previsions. Moreover, the variable ξ spans the different months of the year; therefore, the integration interval I(ξ) is different for every peak contract considered (see Table 1).Table 1
In the first column, we have the months of electricity supply of the peak contracts considered, while in the second one we have the intervalI(ξ) relative to each contract. Note that a contract traded for 3 months is integrated in an interval of length 3.
Month
I
(
ξ
)
April
0,3
May
1,4
June
2,5
July
3,6
August
4,7
September
5,8
October
6,9
November
7,10
December
8,11
January
9,12
February
10,13
March
11,14Before entering into details concerning the specification ofΨ, let us note that, as a rather general rule, forward contracts tend to depreciate getting close to maturity time. Even if the behaviour of peak contracts considered is not the same, based on the latter ones, it seems reasonable to choose a concave function Ψ. In particular, we choose the following exponential function: (74)Ψt-z,T-t=ηe-λT-z,ifz≥0,0,ifz<0,where η is a constant computed as the average daily price variation of all peak contracts, which is equal to -0.063, times the length of the forecast l and λ is a constant that models the price fall. We have estimated this parameter studying the average fall in price over the peak contracts and the trends of the latter ones, obtaining λ=0.04. This choice of kernel function is also motivated by the empirical studies in [1], where the authors show that such a specification in a model based onLévy semistationary processes effectively forecasts spot price dynamics on the German EEX market.


## 4.1.2. Lévy Basis

Several publications (see, e.g., [2, 5, 19, 20]) have shown the applicability of the generalised hyperbolic distribution, in particular the subfamily consisting of the normal inverse Gaussian distributions, to describe appropriately financial datasets. Exploiting these results, we choose a normal inverse Gaussian Lévy basis plugging it in the model described by (71). This kind of Lévy basis can be obtained starting from a normal inverse Gaussian Lévy process of parameters θ, k, and ω, and the latter is the result of subordination of an inverse Gaussian Lévy process, which has variance k, with a Brownian motion with volatility ω and drift θ. In order to use the right parameters, note that the prices of peak contracts are quite stable, being characterized by small daily variation; hence, in general, there are no fat tail phenomena to be taken into account, at least at such a time scale. Latter observation suggests taking small values of ω and k. In particular, also using lines of empirical evidence coming from the analysis of related time series, we take ω=0.01 and k=0.2, also setting θ=1, since we model the entire forward curve.


## 4.1.3. Stochastic Volatility

The time series of monthly peak forward contracts we have taken into consideration show a rather regular trend. Therefore, there is no empirical evidence that causes us to consider this term. Furthermore, the complexity of the resulting model increases by introducing stochastic volatility; hence, we have followed what is suggested in [5], and we have not introduced it, instead fixing σz(ξ)≡1 in (71).


## 4.2. Short-Term Analysis

In this short-term analysis, we take the price of the peak contract at a day to forecast the value of the contract from1 to 5 days ahead. We perform this study simulating numerically the trajectories of (71) using the algorithm presented in Section 3.7. The results are listed in Table 2 and are discussed in Section 4.5.Table 2
In the first column, we have the months of electricity supply of the peak contracts considered, while in the next five columns we have the mean absolute percentage errors (MAPE) obtained in the forecasts of the contracts value at1,2,3,4, and 5 days ahead, respectively.
Month
1
2
3
4
5
April 2012
1.10%
1.81%
2.12%
2.42%
2.87%
May 2012
0.84%
1.28%
1.62%
1.90%
2.01%
June 2012
0.72%
0.97%
1.13%
1.24%
1.33%
July 2012
0.69%
1.05%
1.25%
1.48%
1.54%
August 2012
0.62%
0.84%
1.13%
1.35%
1.51%
September 2012
0.73%
1.13%
1.42%
1.60%
1.70%
October 2012
0.44%
0.62%
0.79%
0.88%
0.95%
November 2012
0.59%
0.82%
0.92%
1.00%
1.07%
December 2012
0.60%
0.80%
0.99%
1.11%
1.17%
January 2013
0.62%
0.88%
1.09%
1.25%
1.39%
February 2013
0.68%
1.09%
1.54%
1.89%
2.19%
March 2013
0.61%
0.82%
1.09%
1.30%
1.61%


## 4.3. Long-Term Analysis

In the short-term analysis, we have taken the price at a day to forecast the value of the contract some days ahead; instead, in order to study the long-term behaviour, we use the first observed price of each peak contract to predict the whole trend of it. In this case, a very important role is played by seasonality. We perform this study simulating numerically the trajectories of (71), using the algorithm presented in Section 3.7. The results are listed in Tables 3 and 4 and are discussed in Section 4.5.Table 3
In the first column, we have the months of electricity supply of the peak contracts considered, while in the second column we have the mean absolute percentage errors obtained in the long-term forecasts of the contracts value.
Month
MAPE
April 2012
2.96%
May 2012
1.94%
June 2012
1.73%
July 2012
1.88%
August 2012
2.86%
September 2012
2.36%
October 2012
1.97%
November 2012
1.99%
December 2012
1.79%
January 2013
2.79%
February 2013
4.15%
March 2013
4.95%Table 4
In the first column, we have the months of electricity supply of the peak contracts considered, while in the next three columns we have the mean absolute percentage errors obtained in the long-term forecasts of the contracts value during the first, second, and third month of trading, respectively.
Month
1st month
2nd month
3rd month
April 2012
3.98%
2.25%
2.65%
May 2012
2.10%
1.30%
2.39%
June 2012
0.97%
2.42%
1.79%
July 2012
0.68%
0.75%
3.99%
August 2012
2.24%
3.90%
2.44%
September 2012
1.40%
0.75%
4.69%
October 2012
0.62%
3.25%
2.02%
November 2012
1.44%
0.97%
3.49%
December 2012
0.83%
2.49%
2.03%
January 2013
0.70%
1.87%
5.79%
February 2013
2.20%
4.99%
5.15%
March 2013
0.74%
4.39%
9.48%


## 4.4. Ambit Approach versus a Standard Approach

In order to assess whether the ambit approach proposed is convenient compared to a standard approach, we have performed the same analyses proposed in Sections4.2 and 4.3 modeling the peak contracts with ageometric Brownian motion, that is, the stochastic process St that satisfies the following stochastic differential equation: (75)dSt=μStdt+σStdWt,where Wt is a standard Brownian motion, μ is a constant calledpercentage drift, and σ is a constant calledpercentage volatility. For an arbitrary initial value S0, (75) has the following analytic solution: (76)St=S0eμ-σ2/2t+σWt.In these analyses, we have estimated the parameters μ and σ using themethod of moments estimation. The results are listed in Tables 5, 6, and 7 and are discussed in Section 4.5.Table 5
In the first column, we have the months of electricity supply of the peak contracts considered, while in the next five columns we have the mean absolute percentage errors (MAPE) obtained in the forecasts of the contracts value at1,2,3,4, and 5 days ahead, respectively, using a geometric Brownian motion.
Month
1
2
3
4
5
April 2012
4,00%
3,96%
3,98%
3,98%
4,04%
May 2012
3,04%
3,03%
3,02%
3,05%
3,09%
June 2012
2,61%
2,62%
2,63%
2,66%
2,68%
July 2012
2,04%
2,05%
2,04%
2,10%
2,08%
August 2012
1,97%
1,98%
1,98%
2,01%
2,02%
September 2012
3,62%
3,66%
3,67%
3,69%
3,71%
October 2012
2,03%
2,03%
2,04%
2,05%
2,10%
November 2012
3,13%
3,15%
3,18%
3,23%
3,26%
December 2012
2,15%
2,17%
2,17%
2,19%
2,22%
January 2013
2,24%
2,26%
2,29%
2,32%
2,33%
February 2013
5,05%
5,04%
5,10%
5,14%
5,21%
March 2013
5,21%
5,28%
5,32%
5,41%
5,52%Table 6
In the first column, we have the months of electricity supply of the peak contracts considered, while in the second column we have the mean absolute percentage errors obtained in the long-term forecasts of the contracts value using a geometric Brownian motion.
Month
MAPE
April 2012
4,60%
May 2012
4,22%
June 2012
3,76%
July 2012
2,12%
August 2012
2,44%
September 2012
4,27%
October 2012
2,21%
November 2012
4,25%
December 2012
3,10%
January 2013
2,07%
February 2013
5,23%
March 2013
7,83%Table 7
In the first column, we have the months of electricity supply of the peak contracts considered, while in the next three columns we have the mean absolute percentage errors obtained in the long-term forecasts of the contracts value during the first, second, and third month of trading, respectively, using a geometric Brownian motion.
Month
1st month
2nd month
3rd month
April 2012
4,13%
3,98%
5,69%
May 2012
4,33%
3,58%
4,72%
June 2012
4,16%
2,87%
4,21%
July 2012
1,94%
1,91%
2,48%
August 2012
2,07%
2,94%
2,32%
September 2012
4,12%
3,51%
5,11%
October 2012
1,81%
2,42%
2,40%
November 2012
3,95%
3,41%
5,35%
December 2012
3,49%
2,54%
3,27%
January 2013
1,61%
1,78%
2,84%
February 2013
3,69%
4,22%
7,53%
March 2013
7,50%
6,09%
9,79%


## 4.5. Conclusions and Remarks

In the short-term framework, the forecasts are very little affected by seasonality and the model presented in this section reflects this fact. Looking closely at Table2, we note that in the short-term prediction the error tends to grow when the length of the forecast increases, as we can expect, but it remains very small. The forecasts obtained for the first day ahead case have always an error less than 1.1% and forecasts on 5 days ahead only once overtake a 2.5% error. Such good performances gained by our approach are slightly facilitated by the small variance characterizing the chosen Lévy basis. We also note that in such short-term analysis the trend of the simulated trajectories is very close to the one of real data (see Figures 1 and 2).Figure 1
In this figure, we present the price curves obtained by simulating the peak contract of June2012 using (71) with the specification above stated. The graphics in the figure present forecasts of the contract value at 1,2,3, and 4 days ahead.Figure 2
In this figure, we present the price curves obtained by simulating the peak contract of November2012 using (71) with the specification above stated. The graphics in the figure present forecasts of the contract value at 1,2,3, and 4 days ahead.Usually, the error increases enlarging the period to forecast, but in the above long-term analysis the mean absolute percentage errors are all less than 5% (see Table3), hence showing that our model provides robust and effective results in the presence of seasonality. In fact, the model often forecasts the behavior of the price some week ahead, even if not always the trend presented by the model is strictly close to the one shown by real data (see Figure 3 and Table 4). The latter result is not surprising since it is related to a “three months from now” forecast, a framework in which energy markets tend to be characterized by small sensitivity from initial data. However, the contracts with a bigger error are the ones which are more complicated to forecast such as, for example, that of February 2013.Figure 3
In this figure, we present the price curves obtained by simulating all the peak contracts using (71) with the specification above stated. The graphics in the figure present long-term forecasts of the contracts value and the month written above each graph refers to electricity supply of the contract considered.Comparing the errors obtained using the proposed model with those obtained exploiting the geometric Brownian motion, we notice that the first are significantly lower over all the performed analysis, as shown in Tables2, 3, 4, 5, 6, and 7. Accordingly, the ambit approach turns out to give better results than thestandard one. We would like to underline that the results obtained by our model provide both good approximations and effective forecasts concerning real electricity forward prices. Moreover, compared with standard approaches, our model, due to the high flexibility of the ambit stochastic setting, allows considering more easily the characteristics of the energy markets.



---

*Source: 626020-2015-10-27.xml*