Distributionally robust optimization with multiple time scales: valuation of a thermal power plant

The valuation of a real option is preferably done with the inclusion of uncertainties in the model, since the value depends on future costs and revenues, which are not perfectly known today. The usual value of the option is defined as the maximal expected (discounted) profit one may achieve under optimal management of the operation. However, also this approach has its limitations, since quite often the models for costs and revenues are subject to model error. Under a prudent valuation, the possible model error should be incorporated into the calculation. In this paper, we consider the valuation of a power plant under ambiguity of probability models for costs and revenues. The valuation is done by stochastic dynamic programming and on top of it, we use a dynamic ambiguity model for obtaining the prudent minimax valuation. For the valuation of the power plant under model ambiguity we introduce a distance based on the Wasserstein distance. Another highlight of this paper is the multiscale approach, since decision stages are defined on a weekly basis, while the random costs and revenues appear on a much finer scale. The idea of bridging stochastic processes is used to link the weekly decision scale with the finer simulation scale. The applicability of the introduced concepts is broad and not limited to the motivating valuation problem.


Introduction
Since the deregulation of the energy market, the question of how to determine the value of a power plant can be asked. The traditional approach of valuing it within a given portfolio of other assets in a coordinated way against one's customer load is one possibility. A second approach is to adopt the ideas of real option pricing in finance. In the first case one ends up with models resembling unit commitment (e.g., van Ackooij et al. 2018) but at a long time scale. Although the actual operation of the power plant can be presented in great detail, it will be harder to incorporate other features in the model. This will typically be the case for uncertainty, where one ends up with multi-stage mixed-integer programs which are not easily solved. One can also argue that it is unreasonable to model the system as fully coordinated. In contrast, when modelling the power plant as a real option, thus operating it in the face of a set of market signals, the setting becomes that of perfect competition. Uncertainty is also naturally modelled, but it comes at the expense of modelling the plant as an independent production unit and thus with less realism in that sense.
However, the price of the real option may well serve as a financial reference base between two parties. For example between the power plant owner and a trading entity actually operating on the market. Taking the option pricing perspective, it must be emphasized that energy markets are by far not as "granular" as the equity markets. For instance, on the electricity market one cannot buy a contract of delivery for a given hour 6 months from now. The classical pricing-hedging duality argument is thus not feasible. Moreover, when operating a power plant, generation will be bound locally by a given power output level. This can be either the result of ramping conditions or minimum up/down times. It is therefore reasonable to try to model the power plant with sufficient realism for the above discrepancies to be minimal. This is the stance that we have taken in the current work. It will lead us to consider a multiscale stochastic program in the sense of , i.e., a multistage stochastic optimization problem where each stage itself is subdivided into a given set of time instants.
To account for uncertainty, we start out with a set of typical stochastic models for underlying prices, which are based on multi-factor models (e.g., Clelow and Strickland 2000) driven by Brownian motions. Clearly, such (commonly used) idealized modelling assumptions are rather unrealistic. It is thus the aim and the core part of the present paper to relax such strong assumptions by computing distributionally robust solutions to the studied operational problem and to investigate how the resulting valuation deviates when considering model ambiguity. Distributionally robust optimization is a field which has recently gained a lot of popularity in the literature (see (Pflug and Pichler 2014, pp. 232-233) for a review of different approaches). In particular, ambiguity sets based on distance concepts between probability measures (such as the Wasserstein distance) are well-supported by theory and frequently applied (e.g., Pflug and Wozabal 2007;Esfahani and Kuhn 2018;Duan et al. 2018;Gao and Kleywegt 2016). However, to the best of our knowledge, the effects of distributional robustness in (especially multistage) real-world applications, have not been investigated yet.
In order to solve the formulated problem numerically, the given uncertainty model will be discretized on a scenario lattice. The multiscale structure could then simply mean that uncertainty is lost within a given stage (cf., e.g., Moriggia et al. 2018). More advanced approaches do consider some uncertainty [e.g., the so-called multihorizon approach originally suggested in Kaut et al. (2014) and subsequently studied and applied in Seljom and Tomasgard (2017), Skar et al. (2016), Werner et al. (2013), Zhonghua et al. (2015), Maggioni et al. (2019)], but the resulting paths do not necessarily connect with subsequent elements in the scenario tree/lattice. Hence, the multi-horizon approach is not appropriate for the present problem, as the key requirement of two time scales which may be assumed to run completely independent from each other, is not given. Indeed, we deal with two different granularities associated with one and the same stochastic process reflecting the evolution of the underlying market prices. A framework for such situations, where sub-stage paths in the lattice are carefully connected, has recently been proposed in . We test the multiscale stochastic programming approach suggested in  in the context of the present real-world application.
Although the resulting ideas will be illustrated through the power plant real option framework, their potential usage is readily seen to be beyond this specific application. In terms of contributions we can can state: • For the application of real option pricing, we investigate more reasonable exercise patterns. In order to keep computational burden low, this naturally leads to multiscale stochastic programs. We also consider model ambiguity to mitigate the fairly ideal models for market prices. From a high-level perspective, we thus extend the literature on real-world applications of dealing with two fundamental problems in stochastic programming, namely the problem of time scales with multiple granularities as well as the problem of model ambiguity. • With respect to multistage model ambiguity, we propose a new concept based on the Wasserstein distance. It is tailored with a computational intention, namely in such a way that (on a discrete scenario tree/lattice) the applicability of a classical backward dynamic programming recursion can be maintained. In particular, the suggested framework leads to solutions that are robust w.r.t. model misspecification in a ball around each conditional transition probability distribution. The size of these balls may be controlled uniformly by a single input parameter. We also link the concept to the nested distance in such a way that it inherits a favourable stability property of the latter. • In the context of Wasserstein ambiguity sets, we propose a state-dependent metric as a basis for the Wasserstein distance. Thereby we account for more realistic worst-case scenarios. We discuss that the well-appreciated statistical motivation for using Wasserstein balls is not invalidated by doing so.
The paper is organized as follows. Section 2 describes the valuation model and the uncertainty model. As typical for real-world energy applications, a sound mathematical framework reflecting all peculiarities of the problem requires carefulness in all details. The underlying uncertainty factors are modelled by a continuous time stochastic process. However, in the light of the nature of the decision problem, we will eventually apply a stochastic dynamic programming algorithm which operates backwards in (discrete) time. To prepare for the computational solution, we therefore discuss all discretization steps required by the multiscale stochastic programming framework that we adopt. Section 3 is dedicated to model ambiguity. We introduce and discuss a new concept which is tailor-made for incorporating model ambiguity into dynamic stochastic optimization models on discrete structures. All numerical experiments and aspects of the computational solution algorithm are given in Sect. 4. Section 5 concludes. Some technical details and examples are deferred to the "Appendix".

The model
Our valuation problem belongs to the class of discrete time sequential decision problems with finite horizon T , decisions u t , state variables z t , and a Markovian driving process ξ t : Here T is the number of decision stages and g t (z t , u t , ξ t ) is the state transition function. The driving stochastic process ξ t is assumed to belong to L 1 ( t , F t ; R m ) and the feasible decision variables at stage t are defined by the set U t (z t ) ⊆ R m . The set of all reachable state variables is denoted by Z t ⊆ R d 1 . The stage-wise profit function h t : R d 1 × R m × R d 2 → R is continuous and satisfies the following growth condition: for all (z, u, x) ∈ R d 1 +m+d 2 and some constant K . We choose the discount factor β t = β t for some constant β ∈ (0, 1] throughout the paper. Any decision u t to be made at time t may only depend on the current state z t and the most recent observation of exogenous information ξ t−1 . This is the non-anticipativity condition. The initial conditions for the random process ξ and the state vector z are that ξ 0 and z 0 are assumed to be constant. In our application, the decisions u t represent the weekly electricity production plan for a thermal power plant. The latter is characterized by many technical constraints, such as minimum up/down times or ramping constraints. Fine grain constraints can be incorporated into the model by increasing the dimension of the state vector and accounting for the number of hours the plant has been offline/online. Such state-representations of constraints on generation assets have received attention in the literature (see, e.g., Martinéz et al. 2008;Frangioni and Gentile 2006;Frangioni et al. 2008 and the references therein). Finer granularity of the time dimension and/or the state variable would result in a significant increase of time steps T (and reduction of the time step size t) as well as an increase in the dimension of z t . For this reason, we introduce here the idea of a multiscale model: While the production decisions u t are made on a weekly scale, the production costs and revenues are calculated on a finer time scale. To make the dynamic optimization algorithm tractable, we make the assumption that the decisions, i.e. the production profiles, must be chosen from a pre-specified set with finite cardinality. The profiles are set up such that they reflect realistic operating conditions and key choices, such as generating at minimal stable generation (MSG) at off peak hours.
Just prior to presenting the specific instantiation of (1), let us emphasize once more that the idea of subdividing a "stage" to mitigate (the curse of) dimensionality goes largely beyond the presented application. Typical other energy problems with similar mechanisms are cascaded reservoir management problems [e.g., see the extensive discussion in van Ackooij et al. (2014) as well as Escudero et al. (1996Escudero et al. ( , 1999, Zéphyr et al. (2015), Cervellera et al. (2006), Aasgård et al. (2014), Séguin et al. (2017), Fleten et al. (2011)].

Instantiation of the problem: the valuation model
In our instantiation of problem (1) the time horizon is spanned by T weeks. Each week t = 1, ..., T is subdivided into S equally sized blocks of hours. With respect to our earlier introduced notation, we now present the following specific versions: • the price process ξ t,s = (ξ e t,s , ξ f t,s , ξ c t,s ) ∈ R 3 represents the electricity price in GBP per megawatt hour (£/MWh), the fuel price in USD per tonne ($/tonne(fuel)) and the CO 2 allowances price in EUR per tonne (e/tonne(carbon)), for each block s within week t, for s = 0, . . . , S. The information up to stage t is the information up to ξ t,0 . The values within the weeks are ξ t,s for s = 1, . . . , S with the convention that ξ t,S = ξ t+1,0 coincides with the initial prices of the next stage. With this convention we ensure continuity of prices in between weeks. In this way, the information up to stage t is to be understood as the information up to the value ξ t,0 .
• the control u t = {u t,s } S−1 s=0 ∈ U ⊆ R S + represents the production profile vector for week t, where u t,s is given in megawatt (MW) and denotes the production at block s. Before the beginning of intermediate values of week t, we determine u t . Then, u t is ξ t,0 − measurable.
• the state vector z t is two-dimensional, i.e., z t = (x t , y t ) with x t ∈ R + representing the amount of CO 2 allowances (measured in tonnes of carbon), that are left for week t. y t ∈ Z + representing the number of hours the power plant was offline before the beginning of week t.
The objective function h t : R × Z + × R S + × R + → R is given by The profit at each block s within week t is defined as follows: whereū t := S−1 s=0 u t,s . Costs incurred are based on the following component functions: • f CO 2 : R × R + × R + → R + gives the cost of buying more CO 2 allowances at the beginning of week t (before the values within week t are known); • f start : Z 2 + × R S + × R 2 + → R + gives the start-up cost if the power plant has been offline prior to (one of its arguments) block s; • f tr : R S + → R + represents fuel transportation costs linked to a selected production profile at the beginning of week t. Table 1 summarizes constants used above or in the sequel. The way in which each state variable is updated will be described now. First we will focus on the variables regarding the CO 2 allowances. Although in the past, a given set of allowances was allocated for free, in principle, they are now obtained from a non modelled auction process. Within our model, the variable I t will represent the number of additional CO 2 allowances received from the regulator at the beginning of week t (measured in tonnes of carbon). Note that this variable will typically be equal to zero but sometimes it will take a relatively high value. The latter happens exactly at the rare events when new allowances are obtained. Now, H 4ūt s is the amount of generated CO 2 during stage t. Hence, together with x t being the remaining stock level and I t the "inflows", the amount of allowances one needs to buy at stage t is α t = [x t + I t − H 4ūt s] − . 1 In the case where α t is positive, we follow a procurement strategy based on a low/middle/high price range partition resulting from some pre-market analysis. Prices in the interval [b, b], with 0 ≤ b < b, are considered middle range. This is formalized as follows: where C α is a constant that determines the size of the extra amount to be bought. Recall, our implicit assumption is that new allowances are always bought before the prices within week t are known. The cost of buying more certificates for week t is then given by The amount x t+1 of remaining allowances after the previous purchase, is updated as follows: The second state variable accounts for start-ups and related costs. The latter depend on the amount of time the power plant was offline. In our model this time frame will be partitioned into C different intervals of hours denoted by (c j , c j+1 ], (c C , ∞), for j = 1, . . . , C − 1, c 1 = 0, over which the start-up costs are assumed to be constant. The associated costs are in terms of power, fuel burnt and extra costs. Depending on y t and the chosen profile u t , one can readily figure out in which interval each start-up of u t falls.
The induced start-up costs at block s within week t are given by: where • W s (y t , u t ) is the amount of works power (MWh) for a start-up at s; • B s (y t , u t ) is the amount of solid fuel burnt (GJ) during a start-up at s; • E s (y t , u t ) denotes engineering and imbalance costs (£) during a start-up at s.
The updated state y t+1 is given by: As a further cost factor, we account for the fuel transportation costs associated to each profile: where the unit transportation cost (in £ per tonne of fuel) is given by the constant factor C tr .

Underlying price processes
To model the underlying uncertainties, i.e., the stochastic price-evolution of electricity, fuel and CO 2 allowances, we postulate a version of a classical two-factor model. The latter are commonly used for the modeling of commodity markets (cf. Clelow and Strickland 2000;Ewald et al. 2018;Ribeiro and Hodges 2004;Farkas et al. 2017). More specifically, in our model the electricity price behaviour is governed by a long term and a short term factor, whereas fuel and CO 2 allowances prices evolve according to a one factor model. In summary, we get a three-dimensional geometric Brownian motion model driven by four correlated one-dimensional Brownian components B e,sh , B e,lo , B f , B c . In particular, the dynamics of the underlying stochastic process F are described by the SDE where the superscripts sh and lo refer to short-term and long-term, respectively. Volatility is allowed to be time-dependent but deterministic. The double-index notation F t,t expresses the fact that we model the spot price as a special case of the forward price.
In particular, the forward price F 0,t (as observed in the market at time 0) will enter the solution of (6) at time t. In this way, we account for the well-known seasonality (peak-hours and off-peak-hours) inherent in electricity prices. Figure 1 visualizes this typical effect. To avoid notational clutter, we henceforth write ξ t with one index for the spot price as a short hand for F t,t . Note that we are dealing with a continuous time stochastic model here. The index notation should not be confused with the discretetime multiscale indexes used in the valuation model; the context will always make clear what is meant.
Regarding the dependence structure between the underlying assets, we allow for a time dependent correlation matrix Using the (lower triangular) matrix L t resulting from a Cholesky decomposition of ρ t , we may replace the Brownian factors t ] , such that the underlying prices are driven by independent Wiener processes W s 1 , W 2 , W 3 , W l 1 . Multiplying the volatility matrix in (6) with L t , we can write the model in the form The non-zero components of the above coefficient matrix involve nasty terms with combinations of the various correlations. The precise parameters can be found in the "Appendix". The solution of SDEs of such a form as in (7) is well known to be of the geometric Brownian motion type (e.g., see (Oksendal 2000, p. 62)). In particular, the random vector ξ t = [ξ e t , ξ f t , ξ c t ] follows a three-dimensional log-normal distribution. The corresponding parameters can again be found in the "Appendix".

Discretization and the associated bridge process
For our numerical solution framework, which is discussed in detail in Sect. 4, the process ξ will first be discretized in all decision stages. Then, an approximate solution of the problem will be computed by stochastic dynamic programming with a backward recursion. In each decision stage, the algorithm relies on the expected profit/loss associated with any decision to be made for the upcoming observation blocks of the following week. To compute such values, we will exploit the structure of the valuation model, the uncertainty model and the backwards recursion. In particular, we are able to compute the expected profits by an analytical formula.
Let us start with the discretization step. To account for the two different time scales explained in Sect. 2.1 above, namely the weekly decision scale and the much finer intra-week observation scale, we use the notation ξ t,s , where t = 0, . . . , T runs in weeks and s = 0, . . . , S in hour-blocks of equal size s (such that t + S · s = t + 1).
Considering the fact that intra-week data of fuel and CO 2 allowances prices typically show -if even available -a rather stable evolution with low fluctuations, we assume those prices to be constant from Monday to Sunday. On the contrary, the electricity price dimension is truly stochastic even on a fine time-scale. Looking at the expected profit function (2) at some block s during a week t, it turns out that (on the basis of our assumptions) the problem boils down to the expected value of the electricity price ξ e t,s given both the initial value ξ e t,0 as well as the final value ξ e t,S of week t. This is due to the fact that the function h t (·) is linear in ξ e t,s . Mathematically speaking, we are left with the computation of the conditional expected value at time t + s · s of the stochastic bridge process linking the values ξ e t,0 and ξ e t,S , for all s = 1, . . . , S. All other parts can be computed in a straightforward way.
The one-dimensional process ξ e t,s follows a univariate lognormal distribution. Thus, its transition density δ is available in analytical terms and the transition density of the associated bridge process can be computed explicitly. Let an initial value η 1 of the process at the beginning of some week t and a final value η 2 at the end of that week be given (i.e., ξ e t,0 = η 1 and ξ e t,S = η 2 ). Then, the bridge process transition density, at time s ∈ [0, S], is given by In particular, we get for the conditional expectation Let us emphasize that the above analytical tractability is not due to our restriction of the intra-week stochasticity to one dimension (see Glanzer and Pflug 2019 for a treatment of the more general multi-dimensional case). This restriction is purely motivated by data.

Ambiguity for dynamic stochastic optimization models
It is an application of classical stochastic dynamic programming theory to solve (1) backwards in time on the basis of the following recursion scheme: where V T (z T , ξ T ) ≡ 0, z 0 and ξ 0 are given.
Let ξ be a Markovian process defined on a finite state space 0 ×· · ·× T , where on each t there is a distance d t . Let the cardinality of t be N t with N 0 = 1 (typically nondecreasing in t). Then the transition matrices P t , t = 0, . . . , T − 1 are of the form N t × N t+1 , where the i−th row of the matrix P t is denoted by p t (i) , for all i = 1, . . . , N t . Notice that each row p t (i) describes a probability measure on the metric space ( t+1 , d t+1 ).
Let ξ i t ∈ t be given. Then, the conditional probability to transition to ξ j t+1 ∈ t+1 is given by the jth element of the row vector p t (i), denoted as p t (i, j) for j = 1, . . . , N t+1 and i = 1, . . . , N t . In this discrete case, the objective of the recursion in (9) can be written as

A new concept: uniform Wasserstein distances
In order to consider model ambiguity, we look for alternative transition matrices Q t , which are close to a given matrix P t . Let us first recall the general definition of the Wasserstein distance for discrete models.
Definition 3.1 Let P = n i=1 P i δ ξ i and Q = ñ j=1 Q j δξ j be two discrete measures sitting on the points {ξ 1 , . . . , ξ n } ⊂ and {ξ 1 , . . . ,ξñ} ⊂˜ , respectively. Then, the Wasserstein distance between P and Q is defined as where π = i, j π i, j · δ ξ i ,ξ j is a probability measure on ×˜ with marginals P and Q, and where D i j is a distance between the resp. atoms.
We will now measure the closeness between (discrete) multistage models P and Q by a uniform Wasserstein distance concept. The rows of alternative matrices Q t are denoted by q t (i), i = 1, . . . , N t . The measure q t (i) is sitting on at most N t+1 points in t+1 and is such that supp(q t (i)) ⊆ supp( p t (i)). Then, we define the distance which can be interpreted as a uniform version of scenariowise Wasserstein distances. An ε ball around P is characterized by the fact that all members Q satisfy for all t = 0, . . . , T − 1 and all i = 1, . . . , N t . When introducing ambiguity into the model, we would like to solve problem (1) wherein the objective function is replaced with: where E Q denotes the expectation with respect to the measure Q.
Our choice of the multistage distance makes it possible to keep the decomposed structure of the backward recursion, which reads: Hence, the ambiguity approach just extends the max model to a maximin model. The ambiguous model can also be seen as a risk adverse model in contrast to the basic risk neutral model. If the distance is not of the decomposable form, then the backward recursion does not decompose scenariowise and one has to find all optimal decisions in one very big stagewise but not scenariowise decomposed algorithm. However, decomposability is the key feature of successful methods for dynamic decision problems. Hence, our concept is strongly motivated by its favourable computational properties. However, as we will discuss now, under a mild regularity condition (in the sense that it will always hold for discrete models, which are the basis of the whole computational framework) it can still be shown that optimal solutions are close if the underlying models are close w.r.t. the uniform Wasserstein distance. The general distance concept for stochastic processes (including their discrete representation in the form of scenario trees) is the nested distance introduced in Pflug (2010), Pflug and Pichler (2012) as a multistage generalization of the classical Wasserstein distance. 2 In our case we have a Markov process which can be seen as a lattice process. Notice that a lattice can be interpreted as a compressed form of a tree. It can always be "unfolded" to a tree representing the same filtration structure, by splitting each node according to the number of incoming arcs. Thus, all results applying for trees do hold for lattices as well. The uniform Wasserstein distance introduced above is given by the maximum Wasserstein distance over all conditional transitions. The subsequent stability result holds.
Proposition 3.1 Let P andP be two discrete Markovian probability models defined on the filtered space ( , σ (ξ )). Assume the following Lipschitz condition regarding P to hold for all t = 0, . . . , T − 1 and all values ξ i t , ξ j t , where i, j = 1, . . . , N t : for K t ∈ R. Consider the generic multistage stochastic optimization problem

where the (nonanticipative) decisions x lie in some convex set and where the function c(·, ·) is convex in x and 1-Lipschitz w.r.t. ξ . Then the relation v(P) − v(P) ≤ dI(P,P) ≤ K · W ∞ (P,P)
holds, where dI(·, ·) denotes the nested distance, and where the constant K is given by Proof The first inequality is a well-known result from (Pflug and Pichler 2012, Th. 11). The statement then follows readily from (Pflug and Pichler 2014, Lem. 4.27), by using W ∞ (P,P) as a uniform bound for W(P t (·|ξ t ),P t (·|ξ t )), over all t.

Remark
Notice that for discrete Markov chain models the assumption in Proposition 3.1 always holds, as one can simply choose the ergodic coefficient

Remark
In the above construction, all models contained in the ambiguity set share exactly the same tree structure and node values. Thus, one might conjecture at a first glance that it would be possible to bound the nested distance by a simple sum of the stagewise maximum of conditional Wasserstein distances, weighted by the number of subtrees at the respective stage. A simple example in the "Appendix" shows that such a construction does not work in general.

State-dependent distances
In practice, the worst-case model for an upcoming period may often depend on the current state. In the model considered in the present paper (cf. Sect. 2.1), we decide only at the beginning of each stage about the procurement of additional CO 2 allowances. In particular, we restrict ourselves not to buy any if the current stock is sufficient for whatever we may do during the subsequent week; regardless of their market price. If we neglect this consideration when searching for the optimal distributionally robust production profile, the worst case may reflect a variation in the CO 2 allowances price dimension which in fact will not have an impact on our optimal decision. Thus, we modify the distance on the underlying three-dimensional space by projecting to the electricity price and the fuel price dimension only, given that our stock of CO 2 allowances is sufficient. Otherwise, we keep the usual L 1 norm. More formally, we define with α t defined in Sect. 2.1 and positive weights w e , w f , w c . Notice that D is not a distance, as it does not separate points. However, this fact does not entail any restrictions for our considerations.
When basing the uncertainty model on historical observations, there is a strong statistical argument for using balls w.r.t. the Wasserstein distance as ambiguity sets (cf. Esfahani and Kuhn 2018). In particular, large deviations results are available (see Bolley et al. 2007;Fournier and Guillin 2015 for the case of the Wasserstein distance and  for the case of the nested distance) which provide probabilistic confidence bounds for the true model being contained in the ambiguity set around the (smoothed) empirical measure. Observe that such results are not invalidated by the state-dependency that we introduce: it is evident that a given confidence bound is directly inherited if one neglects some dimension. Notice however that a general state-dependent weighting of the dimensions would require a more careful treatment.

A case study
In the following, we will test the framework elaborated in Sects. 2.1 and 3 for a specific power plant. For the present application, each week t is subdivided in S = 42 blocks, where each block has 4 h. We solve the problem for a quarter ahead, thus we take the horizon to be T = 13 weeks.
In this case, the control variables are given vectors of dimension S = 42. The set of production profiles we use for this case study consists of 10 different production schedules. This set is denoted by U = {u (i) } 10 i=1 (see Fig. 3) and it will remain constant for every stage.
The discrete evolution of prices is given as a lattice process (ξ ) defined on 0 ×· · ·× T , where each space t has N t elements correspondent to the number of nodes, i.e., . . , N t and all t = 0, . . . , T . As explained in Sect. 3, P t will denote the probability transition matrices from stage t to stage t + 1 of dimensions N t × N t+1 . The description of the lattice construction will be explained in Sect. 4.1.1.
The profit function h t at stage t is defined by the expected profit during the upcoming week (see (2)). Profits at every block s in the week are quantified by the functions f s , for s = 0, . . . , S − 1. Costs of buying additional allowances and transportation costs are quantified only at the beginning of each week, while start-up costs need to be assigned at each block s. We proceed to the description of the state variables. The costs of buying new allowances depend on the strategy A defined in (3). For its computation we consider [b, b] = [4.4, 9.6], where the latter values were obtained by applying a simple quantile rule to the available data set. Moreover, we set C α = 2. For the partition of the amount of available CO 2 allowances x t , we consider different possible values from 0 tonne(carbon) to 10 5 tonne(carbon), and we also take into account the allowances needed for each profile. All in all, the partition of state x t has 16 different elements.
For every state in the partition, an example of the procurement strategy is shown in Fig. 4 for different prices. Note that we illustrate this example when we choose full production, i.e., u (10) . The strategy when choosing a different profile is similar, the only change is that we do not need to buy as many allowances as with u (10) . The second state variable y t describes the hours the power plant was off since the last time it was on. The costs associated with restarting the production depend on y t and the chosen profile u t . The cost function is a step function given in Table 2.
Once a profile is chosen for any week, the first time the profile is different than zeros is where we consider initial start up costs. Then, for the initial start up costs we consider the hours the power plant was offline before week t starts (i.e., y t ) in addition to the hours the chosen profile is off before it starts producing. For the rest of the blocks the costs will only depend on the profile. As for the notation of the elements in f start , see (4), E s is the sum of the last two columns of Table 2. We illustrate the initial start up costs for a profile that is on in the beginning of the week for all possible values of y t (see Fig. 5) .
Given that we have a finite set of profiles, we can calculate the value of y t for each profile. The partition of y t will have these values and the limit numbers of the classes in Table 2.

The solution algorithm
We numerically solve the power plant valuation problem by a stochastic dynamic programming algorithm. A lattice structure is used as a discrete representation of the uncertainty model in all decision stages.

On the lattice construction
The state-of-the art approach for the construction of scenario lattices is based on optimal quantization techniques (cf. Bally and Pagès 2003;Löhndorf and Wozabal 2018). For a given number of discretization points, such methods select the optimal locations as well as the associated probabilities in such a way that the Wasserstein distance (or some other distance concept for probability measures) with respect to a (continuous) target distribution is minimized. For the present study, we have implemented a stochastic approximation algorithm for the quantization task, following (Pflug and Pichler 2014, Algorithm 4.5). Referring to the latter algorithm, we first applied the iteration step (ii) in order to find the atoms of all marginal distributions, separately for each stage. We then formed a lattice out of these sets of points by fixing the structure of allowed transitions and then applied step (iv) to determine all conditional transition probabilities. Eventually, this also determines the absolute probabilities of each node in the lattice. The Wasserstein distance of order two has been used as a target measure for the minimization. We use a ternary lattice, i.e., each node has (at most) three successors with a positive transition probability.
As for the notation, recall that we distinguish between stages (weeks) and intra-week blocks. Decisions are taken only at each stage but the profit will be calculated taking into account the random evolution of the prices during the entire week. The discretized process in each stage t and node i is denoted as takes the value of node j in the next stage with probability p t (i, j), for j = 1, . . . , N t+1 . Note that if the lattice structure does not make a link between two nodes in consecutive stages, then the probability of such a transition will be zero.

Computing the expected profit between two decisions
The valuation of the power plant is obtained by solving (10) backwards in time from t = T to t = 0. In this section, we specify how to solve (10) and its robust version (12) at any stage t. We specifically concentrate on the calculation of the expected profit within each week given the current node and a successor node.
As discussed in Sect. 2.2.1, electricity prices within weeks will be modeled by the bridge process that we described. Fuel and CO 2 allowances prices are assumed to remain constant between stages.
We start now with the computation of the expected profit, as defined in (2). In classical stochastic dynamic programming problems, h t exclusively depends on the values observed at time t. In contrast, in Sect. 2.1 we instantiated (1) in such a way that the function h t is defined as an expected value of the the random profits within week t. Hence, given the values of node i at stage t (at block s = 0), as well as initial states (x t , y t ); the weekly profit h t will be calculated as follows We define At stage T we set the terminal condition V T = 0. Given an initial state (x 0 , y 0 ), going backwards in time from t = T − 1 to t = 0, we obtain the power plant value at t = 0. The latter is calculated with respect to the baseline multistage model P and will be denoted as ν 0 (P). The policy associated with ν 0 (P) is denoted with u * P . It can be represented as a probabilistic tree of profiles.
If we incorporate ambiguity in the lattice process using the uniform Wasserstein distance, for all 0 ≤ t ≤ T − 1 we solve: The optimal value is reached at t = 0 and it will be denoted as ν 0 (Q ε ), where Q ε denotes the ambiguity set defined as A worst-case model Q ε * is any multistage probability model contained in Q ε such that ν 0 (Q ε ) = ν 0 (Q ε * ). More concrete, the optimal value is reached at a saddle point (u * Q ε * , Q ε * ) where u * Q ε * is the policy associated with the worst-case model. At each node i, the objective function of the minimization problem is linear in q t (i) under linear constraints. Define . . . , N t+1 . Then, the minimization problem can be written as l t+1 ] is the distance between nodes k and l at stage t + 1, as defined in (13) with specific weights w e = 1, w f = H 1 and w c = H 3 · H 4 .

The value of the power plant
We describe the optimal valuation of the power plant when we compute the iterative system of backward equations in (10) and (12). We assume that the initial state x 0 provides enough allowances to execute any of the profiles and the power plant was not offline before we start, i.e., y 0 = 0. Moreover, the terminal condition for both problems is set to be V T = 0. With the baseline model we obtain an expected profit of approximately ν 0 (P) = 2.3 · 10 6 (£). The optimal decision at t = 0 is to turn off the power plant by choosing u (1) . The valuation of the power plant including ambiguity is obtained for different radii ε ∈ [0, 2]. The different optimal values ν 0 (Q ε ) are Fig. 6 The impact of the ambiguity radius ε on the optimal profit over 13 weeks  Fig. 6 with respect to the ambiguity level. We observe the valuation of the power plant decreases when the ambiguity radius is higher. For ε = 2, the valuation of the power plant decreases to ν 0 (Q ε ) = 6.8 · 10 5 . In order to get an insight into the change of prices in the ambiguity model, we report the changes of electricity prices for the worst case models Q ε * in Table 3. Let B 0 be a vector containing the stagewise expectations of electricity prices with respect to the baseline model and let B ε be the vector containing the expectations with respect to each worst-case model Q ε * . The percentage of change at each stage t in the prices are denoted with the parameter θ t , such that B ε (t) = (1 − θ t )B 0 (t). We observe that the largest change in prices is given in stage 12, where electricity prices decrease up to 30% for ε = 2.

Forward in time
With the iterative solution of the backwards equations we eventually obtain an initial optimal profile at t = 0, namely u * 0 = u (1) . With this initial decision we go forward in time and create a probabilistic tree u * P of the optimal decisions together with their profits. Starting with the given states and the optimal profile at t = 0, the updated states at stage t = 1 are completely determined by the knowledge of x 0 , y 0 and u * 0 . The choice of the optimal profile in t = 1, for each node i = 1, . . . , N 1 , will be made by looking at the nearest location of the updated states in the grid and taking the correspondent profile chosen in the backwards algorithm. We proceed in this way until we obtain all the optimal profiles at stage T − 1. Eventually, we obtain a probabilistic tree with 3 T possible paths. Following the same procedure, we calculate Fig. 7 Stagewise distribution of the optimal profiles u * P Fig. 8 Stage distribution of the optimal profile with respect the ambiguity radius. For each ε we plot the distribution of u * Q ε * the probabilistic tree of optimal profiles u * Q ε * for each worst-case model Q ε * , and the corresponding profits under the worst-case models for different radii ε.
Starting with u (1) is optimal for all models at t = 0. For the subsequent stages the choices of optimal profiles change. Figure 7 shows the stagewise distribution of the optimal profiles chosen with the baseline model P. Figure 8 shows the changes of profile choices when we incorporate ambiguity in the model.
With no ambiguity there is a probability greater than 0 to choose full production in stages 6, 7, 8, 10 and 12. When we start increasing the radius of ambiguity these chances drop to 0. The larger the ambiguity radius is the more we choose to be offline or not to produce in the weekends choosing profiles like u (2) , u (4) , u (6) . A different option, but with less probability is not to produce in peak hours, by choosing u (3) or u (5) .
Given the optimal profiles for the baseline model P and the alternative models Q ε * we can calculate the profits we make along the decision tree. To be precise, we denote by i τ ∈ N τ any node index at stage τ = 0, . . . , T . Since N 0 = 1, a possible path to follow forward in time up to stage t, will go through any sequence of nodes Fig. 9 On the left the profit tree by following the optimal profiles and on the right the distribution of the final profits obtained following every path (1, i 1 , . . . , i t−1 , i t ). The profit from i τ to i τ +1 is the profit made in stage τ and is written as h i τ ,i τ +1 τ . This profit is obtained with probability p τ (i τ , i τ +1 ). Therefore, the accumulated profit until stage t − 1 is (h 0,i 1 0 + · · · + h i t−1 ,i t t−1 ) with probability p 0 (1, i 1 ) · · · p t−1 (i t−1 , i t ), when we end in node i t ∈ N t . Figure 9 shows the accumulated profits following the tree of optimal profiles u * P as well as the distribution of the final profits.
If we include ambiguity, then the optimal profits change as well as their probabilities. Figure 10 shows the profit trees together with the final distribution with respect to the correspondent alternative model. We observe that for larger ε the alternative models put more weight at lower profits.

Conclusion
In this paper, we have shown how a realistic valuation of a power plant can be done by solving a multistage Markovian decision problem. The value is defined as the (discounted) expected net profit, that one can get from the operation of the plant, if an optimal production plan is implemented. In this valuation process, all relevant purchasing costs and selling prices are included in the model. The number of feasible production plans is finite and thus a discrete multistage optimization problem has to be solved. We use the classical backward algorithm for the Markovian control problem and a forward algorithm for determining an estimate of the achievable profit and its distribution. The novelty of the paper is twofold. First, we adopt a multiscale approach, where decisions are made on a coarser scale than costs are calculated. This allows us to keep the computational effort tractable. Second, we do not only consider the baseline model for the random factors, but rather a set of models (the ambiguity set) which are close to the baseline model. This allows to incorporate the fact that probability distributions for future costs and revenues are not known precisely. The more models, and especially the more unfavourable models are included in the ambiguity set, the smaller is the robust value of the plant. We demonstrate how the final value under model ambiguity depends on the degree of uncertainty about the correct price and cost model. Our distance model for the ambiguity set depends on the state of the system, taking into account that what is close for two price vectors depends also on the fact whether these prices are relevant for the state at hand. We also noticed that the optimal production strategy not only depends on the degree of ambiguity, but also gets more diversified for larger ambiguity, in contrast to some bang-bang solutions in unambiguous models.
The parameters of the multivariate lognormal distribution of the GBM process (7), i.e., the expectation vector μ(t) and the (components of the) variance-covariance matrix (t) of the associated multivariate normal distribution are given by: log(F e 0,t ) − 1 2 t 0 a 2 11 (s) + a 2 12 (s) + a 2 13 (s) + a 2 14 (s)) ds log(F Two versions of a scenario tree with different transition probabilities. The nested distance is 10. The Wasserstein distance between the first stage subtrees is 4, that between each of the two pairs of second stage subtrees is 2 A remark on the nested distance Definition 5.1 The nested distance dI(P,P) between two R m -valued filtered stochastic processes (P, F) and (P,F ) is defined as the optimal value of the following mass transportation problem, which optimizes over the set of all joint distributions that respect the given conditional marginals: inf π ω −ω π(dω, dω) Figure 11 illustrates that the nested distance (between two trees resulting from a variation of the transition probabilities) cannot be bounded by only considering the Wasserstein distance between the subtrees with matching node values. In the given example, the maximum Wasserstein distance in the second stage is 2, the first case Wassersein distance is 4. Thus, 4 + 2 · 2 = 8 would still be smaller than the nested distance between the two trees, which is 10.