Value at risk approach to producer's best response in an electricity market with uncertain demand

We deal with several sources of uncertainty in electricity markets. The independent system operator (ISO) maximizes the social welfare using chance constraints to hedge against discrepancies between the estimated and real electricity demand. We find an explicit solution to the ISO problem and use it to tackle the problem of a producer. In our model, production, as well as the income of a producer, are determined based on the estimated electricity demand predicted by the ISO, which is unknown to producers. Thus, each producer is hedging against the uncertainty of the prediction of the demand using the value-at-risk approach. To illustrate our results, a numerical study of a producer's best response given a historical distribution of both estimated and real electricity demand is provided.


Introduction
In this paper, we deal with the problem which arises in the deregulated electricity markets.Since nowadays the electricity cannot be effectively stored, see [1] for a deep discussion, one of the goals is to balance the aggregated demand and the supply from several producers over a short time period.In particular, we focus on the day-ahead market where the producers offer electricity deliveries a day before the real demand is observed.The producers provide the bidding curves stating the price for a particular electricity quantity that they are able to deliver.Since, in the real world, the aggregated demand is not precisely known and thus is uncertain, we incorporate it into our models as a random variable leading to stochastic optimization problems.An independent system operator (ISO) collects the bidding curves (bids) from several producers and according to an estimated distribution of the aggregated demand, it computes the production quantities to be dispatched.Under our settings, the ISO goal is to maximize social welfare while satisfying the demand with a high probability.In our problems, we do not react on past outcomes of the random variables, which would correspond to the 'wait-and-see' approach, but we hedge against future unknown outcomes, so our decisions must be made 'here and know' and thus are not random but deterministic.
We focus on modelling the various ways the ISO and the individual producers are facing the uncertainty of (future) electricity demand in the day-ahead market.We have two main points in this respect.First, the way the ISO and the producers manage the uncertain demand should reflect their different roles.The ISO has to balance the supply and the demand with high reliability, whereas a producer is more focused on the reliability of the profit.To this end, we model the ISO problem with a chance constraint for the demand satisfaction, and we introduce a chance-constrained problem (CCP) of a producer where the upper quantile of the random profit is maximized.
Second, we argue that the uncertainty of the demand has to be modelled using distinct random variables.Indeed, in real markets, producers are bidding earlier than the ISO is clearing the market.Although this time difference may be relatively small, it can nevertheless force producers to use, e.g. a less precise weather forecast to compute their optimal bids.Even if this was not the case, it is reasonable to assume that different tools and heuristics used by competing market participants to model this uncertainty yield different results.To simplify the notation, we use one random variable for the ISO and another one for all the producers; note that one may easily generalize the model by using an individual random variable for each producer.Such random variables may be dependent, but it is not necessary to specify this dependence, because these variables are used in independent models in our approach.
To be able to provide a numerical study, we had to estimate the probability distribution of the future electricity demand from the point of view of the ISO as well as producers.We were unable to provide realistic estimates, since both the data and methods used are publicly unavailable.Thus, we used publicly available data limited to point estimates and real observations of the demand available for several weeks, and fit the parameters of the lognormal distribution.We are aware that this is a simplifying approach and in practice we would need a more sophisticated method for the demand forecasting such as time series analysis, cf.[2,3], which can take into account seasonality effects and external factors affecting the demand, such as weather.Note that papers [4,5] dealt with the problem of energy pricing considering limited production capacity, unit commitment and fixed costs, however, not taking into account the uncertainty in the electricity demand.
General risky design equilibrium problems with stochastic elements were investigated in [6,7].Recently, [8] proposed a new stochastic-programming market-clearing mechanism to optimize pre-dispatch quantities given the probability distribution of the random demand and the costs of real-time deviations.
The previously proposed stochastic real-time clearing formulation in which generation capacity, demand and transmission line capacity are considered as random, has been extended in [9] by employing the social surplus function, which induces penalties between day-ahead and real-time quantities.
To deal with uncertain (random) demand, we employ the chance-constrained formulations of the problems of the ISO and each producer.CCPs, a standard tool of stochastic optimization, cf.[10,11], are usually used to get optimal solutions which are highly reliable with respect to stochastic parts of the optimization problems under uncertainty.Recent progress in this area includes sequential algorithm based on an exact penalty [12], optimality conditions and regularization [13,14], or new quantile cuts for MINLP reformulations [15].First-and second-order differentiability results under elliptically symmetric distributions, which can be directly employed in standard NLP solvers, have been derived in [16].Paper [17] compared four stochastic programming approaches (including chance constraints) to a large-scale unit-commitment problem which deals with computing the most cost-efficient production schedule while meeting customer load under the operational constraints.A chance-constrained economic dispatch model was presented by [18].The model integrates energy storage and high renewable penetration to satisfy renewable portfolio requirements.The Mordukhovich subdifferential of probabilistic/robust (probust) functions was characterized in [19].In our case, the stochastic problem of each producer is related to the Value at Risk problem, which was elaborated by several previous works, see, e.g.[20][21][22][23].However, due to the specific structure of the producer problems, we are able to derive and solve a nonlinear programming equivalent using the demand distribution function with decision dependent arguments, which is, according to our best knowledge, the first attempt in the area of CCP.
From the modelling point of view, the above-described problem leads to a multi-leader-common-follower problem, where the producers are considered as the leaders and the ISO is viewed as a common follower.Showing the existence of solutions of such a bi-level problem is typically difficult.Even if convexity is assumed at both levels, it is no more satisfied once the upper-level pay-off functions are composed with the solution map of the lower-level problem (the ISO).For a specific setting, all equilibria may be found analytically, see, e.g.[24,25], assuming, however, that the demand is deterministic.Alternatively, one may model electricity markets as supply function equilibrium (SFE), a concept introduced in [26] that naturally generalizes market models of both Cournot and Bertrand.Modelling the competition of producers by Nash equilibrium, the profit-maximizing support functions (i.e.inverses of bid functions) are smooth functions described by differential equations.This general model of market competition has been well adapted to the particular situation of electricity markets, see, e.g.[27,28] and the references therein.Note, however, that such an approach is orthogonal to the value-at-risk approach used here.Indeed, supply functions determined by maximization of (expected) profit are in no relation to optimal bid functions of producers using value-at-risk.As far as we know, an analytical solution to such a market is yet to be determined; our work may be considered as the first step in this direction.
To better take into account the uncertainty in the electricity demand, we simplified the model of a pay-as-clear day-ahead market in several aspects.Next, we will discuss our main assumptions: (1) To focus more on competition amongst producers, we don't model individual consumers.They may be probably included in a more detailed model by following, e.g.[29].(2) We consider one specific time period day-ahead (a quarter of an hour), which is independent of other time periods.Thus, the model would have to be considerably adapted to incorporate other market participants, e.g. market speculators selling and purchasing contracts for different periods in the day.(3) Transmission constraints are not taken into account; in other words, model is formulated at a single node of transmission network.Such constraints substantially complicate the analysis of the problem, see, e.g. the discussion in [28] where radial transmission network with local demand shocks is analysed, and existence and uniqueness of SFE in two-node networks is shown.(4) The aggregated electricity demand is considered to be in-elastic.The model may be generalized in this respect by following the direction of [30], thus modelling the linear elasticity of the demand.(5) We assume that the producers and the ISO are able to estimate the continuous probability distribution of the demand in each step of the modelling.(6) We model the production costs using convex quadratic functions, which is quite common as a reasonable simplification in the analysis of equilibria in electricity markets, see, e.g.[24,29,31,32].Such approximation captures well, at least qualitatively, the increasing marginal costs of electricity production.(7) We limit producers to bid functions that are convex and quadratic, following again, e.g.[24,29,31,32], thus obtaining approximation that is convenient for further mathematical analysis.In real markets, however, producers' bids are typically piecewise-linear functions with limited production capacity.
Some of the above limitations may possibly be overcome in the future, whereas others seem to be inevitable to facilitate the analysis below.In particular, even if one found an analytic solution to the ISO problem with capacity constraints (cf.Theorem 3.1), a statement analogous to Theorem 4.5 would still depend on the assumption of quadratic cost and bid functions, see the points 6 and 7 above.We, however, argue that such a simplified model still carries the main features of an electricity market, and focus on the application of the probability constraints in the electricity market modelling; it is an open question to which extent one could impose probability constraints using a more realistic model.
The paper is organized as follows.In Section 2, the basic notation, assumptions and market settings are introduced.Section 3 is focused on the optimal dispatch problem of the ISO.Section 4 deals with the chance-constrained profit maximization problem of a producer.A numerical study using the real data from the French electricity market is proposed in Section 5. Section 6 concludes the paper.

Notation and problem setting
First, we summarize the basic hypothesis that are considered in this work: we consider a pay-as-clear electricity market with N > 1 producers; we only consider producers, that is, the demand of consumers is aggregated.Finally, the transmission network is not taken into account, thus also thermal losses and 'local demand' are omitted.
By δ > 0, we denote the (aggregated) electricity demand, N = {1, . . ., N} is the set of producers, and q i ≥ 0 represents the non-negative production quantity of the ith producer.Considering q ∈ R N + , we use to denote the vector (q 1 , . . ., q i−1 , q i+1 , . . ., q N ), and the same convention is used also for other vectors hereinafter.For i ∈ N , we use a i , b i ≥ 0 to denote the coefficients of the ith producer's bid a i q i + b i q 2 i and A i ≥ 0, B i > 0 to denote the coefficients of the true production cost function A i q i + B i q 2 i .We use R ++ = R + \ {0}.For a one-dimensional random variable X on a probability space ( , F, P), we denote its distribution by μ which is defined as for all Borel measurable subsets A ⊆ R.This distribution induces the distribution function the inverse of which is the quantile function F −1 X defined by We say that a measurable real function f X is a density of X, if for all Borel measurable subsets A ⊆ R or, equivalently, if

Problem of ISO
Each producer provides the ISO with a quadratic bid a i q i + b i q 2 i .The ISO thus have knowledge of the bid vectors a ++ , however, the ISO is not aware of the true production cost parameters A i , B i , i = 1, . . ., N. Thus, knowing only the bid vectors, the ISO computes the production quantity to be dispatched to the producers q = (q 1 , . . ., q N ) ∈ R N + to maximize the so-called social welfare, see, e.g.[29,33].Assuming, moreover, that also the aggregated demand δ > 0 is given, the problem ISO(a, b, δ) reads ISO(a, b, δ) min Note that the market clearing price λ(a, b, δ) corresponds to the Lagrange multiplier of the demand satisfaction constraint in ISO(a, b, δ).The following result is fundamental for this work.
++ , the market clearing price λ(a, b, δ) and optimal production quantities q(a, b, δ) in problem ISO(a, b, δ) are the solutions of the system of equations and This statement is already known, see, for example, a more general setting of [25, Theorem 2.1].For the sake of completeness, we nevertheless include a concise proof.

Proof:
The Karush-Kuhn-Tucker (KKT) conditions for the problem ISO(a, b, δ) are stated as with the first two conditions considered for all i ∈ N .Since δ > 0, at any feasible point of ISO(a, b, δ) there has to be some i ∈ N such that q i > 0. Then one may easily verify that Linear Independence Constraint Qualification (LICQ) is satisfied everywhere, and so KKT conditions (3) are necessary and sufficient for the solution to ISO(a, b, δ).
To solve (3), we first show that λ ≥ 0. Let λ < 0 for a contradiction.Assumption δ > 0 implies that q j > 0 for some j ∈ N , and so ν j = 0.Then, a j + 2b j q j = λ < 0 contradicts a j , b j , q j ≥ 0. Next, let us observe that {i ∈ N : Indeed, for all i ∈ N such that ν i = 0 we have On the other hand, ν i > 0 implies q i = 0 and thus also λ = a i − ν i < a i .Now, we may verify formula (2); considering first i ∈ N such that a i > λ, we have ν i > 0 by (4), and so q i = 0. Similarly, a i ≤ λ implies ν i = 0, and combined with (5) leads to q i = λ−a i 2b i .Next, by substituting (2) into the last equation of (3) we directly verify also (1).
Finally, it remains to show that the implicit definition of the real function 1) is correct.Indeed, the left-hand side of (1) is continuous in λ, and also strictly monotone in λ using λ > 0 implied by assumption δ > 0.
The fact that λ(a, b, δ) is well defined by (1) may also be seen from the following remark.

Remark 3.2:
Consider the setting of Theorem 3.1.If we moreover assume, without loss of generality, that a i ≤ a j for i < j, we may restate (1) equally as For details see [25,Remark 3].In this article, however, we will not assume any ordering of producers.Now, we turn our attention to the ISO problem with stochastic demand represented by a positive random variable D ISO on the probability space ( , F, P).
In the remained of the paper, we impose the following assumption on the distribution of the random demand: (CD) Let the distribution of random demand be absolutely continuous, i.e. it has a probability density function with respect to the Borel measure.
Since we deal with one particular part of the day (a quarter of an hour), we are not using any time index.We assume that the main goal of the ISO is to establish equilibrium between supply and demand with great reliability to avoid high costs associated with the supply failure.Such a problem may be formulated as a chance-constrained problem where a probability p ∈ ]0, 1[ is prescribed to satisfy the demand: As the individual chance constraint above has a structure of the so-called separable (random) right-hand side, cf.[11], one can easily reformulate it using the quantile function of D ISO , thus obtaining a deterministic constraint.Consequentially, the explicit solution of ISO(a, b, δ) stated in Theorem 3.1 remains valid even for SD-ISO(a, b) provided δ is replaced by . This important observation will be used in Section 5 to solve the stochastic ISO problem to get the optimal dispatch quantities and clearing price.

Problem of producer
In this section, we illuminate the point of view of a particular producer i ∈ N supposing that its true production cost function is given by A i q i + B i q 2 i with known A i ≥ 0 and B i > 0. Producer i ∈ N then aims to maximize his profit function with respect to his decision variables a i ≥ 0, b i > 0, with the remaining bid coefficients (a −i , b −i ) fixed.Furthermore, we assume that the electricity demand δ is not known when the producer's bid is submitted to the ISO.Instead, we consider stochastic demand is given by a positive random variable D on the probability space ( , F, P) with a probability density function f (δ).We stress that the producers and the ISO can represent demand with random variables having different distributions since the ISO can use more recent information to predict the demand for the considered time period in the next day.Now, producer i can solve the following chance-constrained problem, where the profit m i that can be reached with a probability p i ∈ ]0, 1[ with respect to the random demand D is maximized: Note that this formulation is related to the value at risk (VaR), where the analogous chance constraint is imposed on random losses resulting from investments on financial markets, cf.[20,23].Alternatively, one may consider losses above the quantile leading to the measure known as conditional value at risk (CVaR), see [21,22].

Existence of solution
In this subsection, we show that problem Theorem 4.1: Before proving the above theorem (on p. 10), we show several auxiliary lemmas.First, we extend the domain of the profit function π i to include points having b i = 0, and establish continuity of thus adjusted function; we shall denote such function still by π i .Should the optimal bid function of producer i be such that b i = 0, it has to be interpreted as 'limiting' bid function, cf.remarks following equation ( 15) in [25], since such b i is not considered in the problem of the ISO(a, b, δ).

Lemma 4.2:
Thus adjusted function π i is well-defined.Moreover, given the fixed values of a −i , b −i and δ, function

Proof of Lemma 4.2:
We order producers as in Remark 3.2, i.e. for all k, l ∈ N such that k < l it holds a k ≤ a l (slightly abusing the notation by denoting the producer in question still by i).Next, define function q i for any x ∈ R by q , and observe that it is bounded by δ.Moreover, denoting λ = λ(a 1 , . . ., a i−1 , b 1 , . . ., b i−1 , δ) and using (1) we may verify that q i (x) > 0 for all x < λ, and q i (x) < 0 for all x > λ.First, we will show that lim For any ãi ≥ 0 and b i > 0 one may rewrite Equation (1) equivalently as By evaluating both sides of (12) for ãi converging to a i and b i converging to 0 from above, we will prove (11) by contradiction.Let us denote the left-hand side of (11) by L; for L < min{ λ, a i } the left-hand side of ( 12) is converging to 0, whereas the right-hand side is, eventually, positive, leading to a contradiction.For L > λ, the right-hand side of ( 12) is, eventually, negative, which is a contradiction with non-negativity of the left-hand side.Finally, for L > a i , the left-hand side of ( 12) is, eventually, above δ, thus we get a contradiction with boundedness of the right-hand side.Further, we will show that lim For a i > λ, Equation ( 11) implies that, eventually, (8).Thus formula ( 13) is justified.Otherwise, i.e. for a i ≤ λ, the clearing price converges to a i due to (11).Then . Next, by employing identity valid for any ãi and b i > 0 one establishes convergence of q i (ã i ) to q i (a i ), and so (14) by the means of ( 8).
We will complete the proof by showing that function π i , now defined for all 2), ( 6) and (8).Continuity at points (a i , 0) with respect to sequences b i → 0+ is due to the definition, see (10).Thus, one has only to verify continuity of π i at points (a i , 0), i.e. of the function given by ( 13) and ( 14), with respect to a i .To this end observe first that q i ( λ) = 0, thus there is no jump when passing between ( 13) and (14).Finally, formula (14) is continuous in a i due to continuity of q i (a i ).
In the rest of this subsection, we will denote by λ(a −i , b −i , δ) a value of a function defined on Let us henceforth denote by m i (a −i , b −i , p i ) the supremum of the objective function in P i (a −i , b −i , p i ).We will establish lower and upper bounds on To show the first inequality, observe that λ(a, b, δ) ≤ a i implies π i (a, b, δ) = 0 using (2); the latter stems from (15).Now, finding a i high enough such that To show that m i (a −i , b −i , p i ) is finite, we will now find an upper bound on π i (a, b, δ) being valid for all δ > 0 such that π i (a, b, δ) > 0. Necessarily, one then has A i < λ(a, b, δ) due to (8), thus also A i < λ(a −i , b −i , δ) using (15).Assuming b i > 0 for the moment, we observe from (8) . The first inequality is due to (15); to show the latter it suffices to maximize the middle formula w.r.t.q i ∈ R.Moreover, such an upper bound on π i holds also for b i = 0 using continuity of π i as stated in Lemma 4.2.Further, denoting α = j =i a j 2b j and β = j =i β for all δ, and since Thus, given some p i ∈ ]0, 1[ , for a i , b i ∈ R + and m i > 0 to be a feasible point of P i (a −i , b −i , p i ), it has to hold In other terms, β of D; we thus obtained an upper bound on m i .

Now, we define probability function ρ
omitting fixed parameters a −i and b −i in the notation.Note that X i (m i ) is an intersection of the set of feasible points of Proof: To show the statement, we may alternatively prove that for m i > 0 setvalued mapping and S(y) ⊃ S(x) for all x > y > 0.
Observe first that the monotonicity of S(m i ) is directly due to monotonicity of ρ i (a i , b i , m i ) in m i , see (14).
Next, consider a realization δ > 0 of random demand D and observe that π i (a, b, δ) ≥ m i implies q i (a, b, δ) > 0 due to (2), and so λ(a, b, δ) > a i .For the case of b i > 0, one may use Equations ( 2) and ( 8) to obtain We see that ρ i (a i , b i , m i ) ≥ p i implies P λ(a, b, D) ≥ max {a i , √ 2m i b i } ≥ p i .Considering the quantile function F −1 D of D, the last formula may be rewritten as where we included also the previously avoided values of b i = 0.

Proof of Theorem 4.1:
has a solution.This solution is, however, also a solution of P i (a −i , b −i , p i ), cf.(17).

Problem reformulation
In this subsection, we reformulate problem P i (a −i , b −i , p i ) to facilitate the numerical experiment in Section 5. First, we analyse what values δ > 0 of the demand yield π i (a, b, δ) ≥ m i .To this end, we define functions for all possible values of a i , b i ∈ R + and m i > 0 in the following way: Then, we may reformulate the probability function ρ i .

Theorem 4.5:
By comparing ( 18) with ( 16), one may deduce that functions we necessarily have q i (a, b, δ) > 0 and so using (2) also and λ(a, b, δ) > a i .Then, denoting α = a i − A i , β = 2b i − B i , and and substituting ( 20) into (8), we may reformulate inequality (19) equivalently as λ(a, b, δ) ∈ (a i , b i , m i ).Evaluating functions inf{ (a i , b i , m i )} and sup{ (a i , b i , m i )}, with the discriminant of the left-hand side of the inequality in (21) being thus showing (18).Further, observe that all above arguments are valid also for b i = 0 due to continuity of profit π i with respect to b i → 0+, see Lemma 4.2.
Next, we define functions yielding what value of the aggregated demand corresponds to clearing prices λ 1,2 (a i , b i , m i ), see Equation (1), and conclude this section with a corollary playing a key role in the numerical experiment.

and denote by F D the distribution function of the demand D. Then either the optimal profit in P
Proof: Using Theorem 4.5, it suffices to observe that the probability function can be expressed as

Numerical study
In the numerical study, we apply the above-derived results to real data from the French electricity market.We derive the optimal bidding curves for five artificial producers and report the corresponding optimal dispatch quantities, which are provided by the ISO.We employ the real data to estimate the distribution of the random demand.Note that these estimates are different for the producers and the ISO as discussed below.Our goal is also to illustrate that the reformulations obtained in the previous section lead to problems which can be solved by standard software tool such as Matlab.
The sequence of steps is from the perspective of producer i ∈ N as follows: estimate the future demand distribution D; use it to solve producer's optimization problem, thus obtaining bidding parameters a i , b i to be submitted to the ISO; after the clearing process of the ISO using D ISO , the producer obtains the dispatch order and the payment.The sequence of steps from the perspective of the ISO: generate a forecast of the demand distribution D ISO ; obtain the bids a, b from the producers; use these data to clear the market day-ahead; announce the dispatch and pay the producers according to the clearing price.In the case that the demand realization D ISO = δ does not match the planned supply, the difference is then compensated in the intraday market.The situation is outlined in Figure 1.
We will now introduce a naive approach to estimating parameters of D and D ISO .Table 1 contains the point estimates and the real data observed between January 3 and February 28, 2017.These days correspond to Tuesdays, Wednesdays, and Thursdays; we wanted to avoid Mondays and Fridays when the demand can differ considerably.Note that our model can be used to deal with any particular time period for which we are able to construct reliable estimate of the random demand distributions.The table contains the point estimates Dt of the demand used by producers, the point estimates DISO t used by the ISO, the real consumptions on the day D ISO t , and the clearing prices.We assume that the point estimates are i.i.d.realizations of the demand forecasts for the next day.We are aware that this is a highly simplifying assumption and in practice we would need more sophisticated approach as we have discussed in the Introduction.Based on  the observations, we have estimated the parameters of the lognormal distribution, cf.Table 2.Note that our approach is not limited to this particular probability distribution, but we can use any distribution with positive support, e.g.Gamma or inverse-Gaussian.As already discussed in the Introduction, we further construct the producers' estimate of D based on point estimates Dt that are published by the ISO.Such a naive method is necessary to illustrate our approach since one may not use private data and/or models of producers.For a producer, the parameters of the lognormal distribution of D are estimated using the pairs Dt , DISO t .The expected value corresponds to the sample mean of the producer demand estimate D, i.e.The values of the parameters μ, σ 2 for a producer are then obtained by the following arithmetic operations valid for the lognormal distribution using the desirable expected value and variance Analogous approach is used to estimate the parameters of D ISO with the pairs of observations DISO t , D ISO t in the place of Dt , DISO t , see Table 2. Considering five producers, we will solve problem (9) for each producer given the initial values of the bidding coefficients a i , b i and the production cost coefficients A i , B i , see Table 3.Note that producer 1 is considered as the largest one with the smallest linear cost coefficient and the highest quadratic one, whereas producer 5 is the smallest with corresponding cost curve.We will employ the reformulated form of the problem of producer, see (23), assuming p i = 0.9 for all i.
As we mentioned above, we are not able to get the equilibrium over all producers strategies.Instead, we are going to compare three approaches which demonstrate different information propagation/availability of the other producers strategies.We consider the following approaches: (1) compute âi , bi given (a −i , b −i ) for all i = 1, . . ., N , i.e. all producers perform optimization independently and know only initial/unoptimized bidding coefficients of other producers, which means they have only basic information about other producers strategies, (2) compute only â3 , b3 given (a −3 , b −3 ), i.e.only one (randomly selected) middle-sized producer optimizes its profit and the others use basic unoptimized coefficients as the inputs for the ISO, (3) compute âi , bi given â1 , . . .âi−1 , a i+1 , . . ., a N , b1 , . . .bi−1 , b i+1 , . . ., b N for all i = 1, . . ., N , i.e. producer i is given the optimal bids of producers 1, . . ., i − 1, which corresponds to partially informed producers.Note that the last producer N knows the optimal strategies of other producers, i.e. it is perfectly informed.
The problems are solved by optimization solvers available in Matlab.The implementation of the producer problem (23) includes the following steps: (0) Writing the lower and upper bounds δ 1,2 (a, b, m i ) as Matlab functions.
Defining the nonlinear constraint and objective function as anonymous functions using @() declaration.Setting the solver option to SQP algorithm.(1) Fixing the parameters (a −i , b −i ) according to the selected approach (1)-( 3).
(2) Selecting the starting points for the decision vector (a i , b i , m i ).
(3) Setting the lower and upper limits for the decision vector.(4) Solving the problem for each producer using fmincon solver.
The obtained optimal bids âi , bi , i = 1, . . ., N are then reported to the ISO, which produces the dispatch orders qi and the clearing price λ by solving SD-ISO(a, b) problem (7).This optimization problem is solved by quadprog solver.The clearing price λ is obtained as the value of the Lagrange multiplier corresponding to the constraint, which is part of the quadprog output.Table 3 contains the numerical results.We can observe that although the delivered bids are all different in all cases, the clearing price is stable.If we compare the dispatched quantities, we can observe that the producer 1 is expected to deliver the most ( m1 = 21.36),although not as much more than the producer 2 in the first approach ( m2 = 19.47).The difference is then greater in the other approaches.As expected, the smallest producer will deliver the least in all approaches.
We also investigated the development of optimal values of parameters for producer 3 with respect to the changes of the probabilistic level p i ∈ [0.5, 0.99], see Figure 2. Note that a −3 and b −3 are hold fixed and that the optimal values from the previous iteration are used as the starting points for the update.We can observe that the behaviour is quite stable for m3 , whereas â3 , b3 rapidly change for high probabilistic level p 3 where the coefficient for the linear term â3 rapidly increases whereas b3 decreases.We can deduce that the optimization tries to stabilize the bidding function by reducing the influence of the quadratic term and by strengthening the linear term.
Sensitivity analysis of optimal value m 3 was also performed with respect to its production cost parameters A 3 , B 3 , see Figure 3.As the production costs of producer 3 increase, the profit quantile decreases, which is certainly in line  with expectations.When the coefficient for the quadratic term B 3 increases, the decrease of the profit quantile m3 is more significant than that associated with growth of the linear term A 2 .Perhaps less expected is the behaviour of the profit quantile relative to the expected bidding function of another producer.The development of the profit quantile m3 with respect to the producer's bid coefficients a 2 , b 2 can be found in Figure 4.As the coefficients of the bidding function of producer 2 increase, so does the quantile of producer 3's profit m3 .This can be interpreted as that if producer 2 is more expensive, producer 3 can benefit from it and achieve higher profits with a high probability by optimizing its bidding coefficients.When we increase the coefficient for the quadratic term b 2 , the growth of the profit quantile m3 is more significant than for the linear one a 2 .

Conclusions
In this paper, we have investigated two closely connected problems appearing on deregulated electricity markets which are subject to uncertainty.We have focused on the stochastic demand and employed the chance-constrained formulations for the problems of the ISO and producers.We have shown that due to the structure of the ISO problem, it is possible to use an earlier result and derive an explicit solution for the production quantities.For each producer, we have formulated a value at risk problem with the maximization of profit which can be reached with a certain level of probability.Then, we have derived an explicit reformulation of the probability function, which enables to solve the problem using a nonlinear programming solver.In the numerical study, we have illustrated our approach using real data from the French market.We have solved the producers as well as the ISO problems and compared several approaches to sharing the information about producers' bidding functions showing their impact on the optimal production quantities.
represents the Mean Square Prediction Error (MSPE) which is estimated as a sum of the sample variance of D Var( DISO , D) = Var( D) + MSE( DISO , D).

Figure 3 .
Figure 3. Optimal solutions of problem (9) for producer 3 -sensitivity with respect to the production cost parameters A 3 , B 3 .

Figure 4 .
Figure 4. Optimal solutions of problem (9) for producer 3 -sensitivity with respect to the producer's bid coefficients a 2 , b 2 .
then a maximizer of P i (a −i , b −i , p i ) exists due to Lemma 4.3, thus we may further assume m i > 0. Next, for any 0 < mi < m i a non-empty set X i ( mi ) is closed since ρ i in (17) is upper-semicontinuous.This stems from (16) using continuity of π i in (a i , b i ) due to Lemma 4.2.Moreover, X i ( mi ) is also bounded using Lemma 4.4, and so problem max

Table 1 .
Input data a : point estimates of the demand at various stages ( Dt used by producers, DISO

Table 2 .
Parameters estimates of the lognormal distribution.

Table 3 .
Starting values of the coefficients (a i , b i ) and results of the producer (â i , bi , mi ) and the ISO problems (q i , λ) for considered approaches.