Long-Term Generation Scheduling for Cascade Hydropower Plants Considering Price Correlation between Multiple Markets

: The large-scale cascade hydropower plants in southwestern China now challenge a multi-market environment in the new round of electricity market reform. They not only have to supply the load for the local provincial market, but also need to deliver electricity to the central and eastern load centers in external markets, which makes the generation scheduling much more complicated, with a correlated uncertain market environment. Considering the uncertainty of prices and correlation between multiple markets, this paper has proposed a novel optimization model of long-term generation scheduling for cascade hydropower plants in multiple markets to seek for the maximization of overall beneﬁts. The Copula function is introduced to describe the correlation of stochastic prices between multiple markets. The price scenarios that obey the Copula ﬁtting function are then generated and further reduced by using a scenario reduction strategy that combines hierarchical clustering and inconsistent values. The proposed model is applied to perform the long-term generation scheduling for the Wu River cascade hydropower plants and achieves an increase of 106.93 million yuan of annual income compared with the conventional scheduling model, without considering price scenarios, showing better performance in e ﬀ ectiveness and robustness in multiple markets.


Introduction
On 15 March 2015, the Communist Party of China's (CPC's) Central Committee and State Council proposed opinions on further deepening the reforms of the electric power system, marking the beginning of a new round of electricity market reform in China [1], after a series of reforms in 1985, 1996, and 2002 [2]. This new round of reform aims to break the monopoly, introduce competition, and achieve optimum allocation of power resources in a broader scope. Over the past four years, most provinces and regions have been actively trying to get involved in the reform and have continuously improved the market operation rules. Although the market is still constrained by many factors, the hydropower plants in some provinces of southwestern China have successfully participated in the yearly and monthly electricity markets [3,4]. It is foreseeable that a competitive and stable electricity market environment will be formed in the future.
With extremely rich hydropower resources, many huge cascaded hydropower plants have been constructed in southwestern China. As the main sending end for the West-to-East Electricity Transmission (WEET) Project, these hydropower plants face a special multi-market environment.
The remainder of this paper is organized as follows: In Section 2, the MCP scenarios for multiple markets are described. Section 3 is devoted to the detailed presentation of the proposed optimization model. In Section 4, the solution method for the proposed model is introduced. In Section 5, the case study and optimization results are carefully discussed. Finally, Section 6 outlines the main conclusions.

Description of MCP Scenarios for Multiple Markets
The market electricity price derives from a complex process with uncertain factors [23,24], and is generally obtained through prediction methods in the stage of reporting transactions [25]. However, the prediction price is subject to many factors, such as the accuracy of the input data and the fitting parameters, and there may be a certain prediction error, objectively. Directly using the deterministic prediction value will bring participants many more market risks. In this paper, the MCP is firstly described as a variable of its prediction value and random prediction error. Then, considering the correlation between the multiple markets, a joint distribution function of the MCP for the multi-market environment is constructed and further discretized into several independent and identical distribution scenarios. Finally, these scenarios are reduced based on the balance of the result's accuracy and the computation workload.

Building MCP Scenarios Based on the Copula Theory
For any time period of t, the predicted MCP of market i is denoted byλ i,t , and the prediction error is ε i,t . The actual MCP λ i,t can be described as a random variable ofλ i,t and ε i,t . That is: where, i = 1, 2, · · · , I indicates the serial number of the market, and the total number of markets, I; λ max i,t and λ min i,t , indicate the upper and lower bounds for λ i,t , respectively. Due to the complex correlation between supply and demand of multiple markets, the MCPs of these markets can fluctuate greatly and are related to each other, thus the MCP of each market is not appropriate for separate description. In this paper, the Copula theory in econometrics is introduced to describe the correlation of MCPs in multiple markets, and a method for generating combined MCP scenarios (i.e., MCP scenarios for multi-markets) based on the Copula function is proposed. The Copula theory was firstly proposed by Sklar in 1959 to describe complex nonlinear correlations between random variables, which is now widely used in the financial sector as it has no restrictions on the marginal distribution of variables [26]. The Copula function can put the marginal distributions together to form a joint distribution, and its strict definition is given later by Nelsen [27]. That is, there is a function C that makes: H(x 1 , x 2 , · · · , x n ) = C(H 1 (x 1 ), H 2 (x 2 ), · · · , H n (x n )) (2) where, H is the joint distribution function of the random vectors x 1 , x 2 , · · · , x n , and its marginal distribution is H 1 (x 1 ), H 2 (x 2 ), · · · , H n (x n ), respectively. It has been proved that if all these marginal distributions are continuous, the Copula function C is unique. Defining H −1 i (x i ), i ∈ {1, 2, · · · , n} as the inverse function of H i (x i ), the Copula function can be expressed as: C(u 1 , u 2 , · · · , u n ) = H H −1 1 (u 1 ), H −1 2 (u 2 ), · · · , H −1 n (u n ) where, u i (i = 1, 2, · · · , n) is a variable obeying uniform distribution in [0, 1]. Obviously, the Copula function C(u 1 , u 2 , · · · , u n ) is a joint distribution function constructed by n marginal distributions of uniform distribution in [0, 1]. There are five commonly used Copula functions, including Gaussian-, t-, Gumbel-, Clayton-and Frank-Copula functions. Their expressions, characteristics, and density distribution plots are listed in Table 1.

Type Expression Characteristics Density Distribution Plot
Gaussian Distributed symmetrically, without reflecting the tail correlation   , ; Distributed symmetrically, without reflecting the tail correlation As mentioned above, the cascade hydropower plants in southwestern China mainly participate in two electricity markets, i.e., the provincial and external markets. Based on the Copula theory, the combined MCP scenarios were generated using the following steps: (1) Construct MCP marginal cumulative distribution functions for the provincial and external markets. As the MCP , i t λ is related to its prediction value , λ i t and prediction error , σ . The probability density function (PDF) is: Distributed symmetrically, reflecting the tail correlation   , ; Distributed symmetrically, without reflecting the tail correlation As mentioned above, the cascade hydropower plants in southwestern China mainly participate in two electricity markets, i.e., the provincial and external markets. Based on the Copula theory, the combined MCP scenarios were generated using the following steps: (1) Construct MCP marginal cumulative distribution functions for the provincial and external markets. As the MCP , i t λ is related to its prediction value , λ i t and prediction error , σ . The probability density function (PDF) is: Distributed asymmetrically, the upper tail related and the lower tail independent Distributed symmetrically, without reflecting the tail correlation As mentioned above, the cascade hydropower plants in southwestern China mainly participate in two electricity markets, i.e., the provincial and external markets. Based on the Copula theory, the combined MCP scenarios were generated using the following steps: (1) Construct MCP marginal cumulative distribution functions for the provincial and external markets. As the MCP , i t λ is related to its prediction value , λ i t and prediction error , σ . The probability density function (PDF) is: Distributed asymmetrically, the upper tail independent and the lower tail related Distributed symmetrically, without reflecting the tail correlation As mentioned above, the cascade hydropower plants in southwestern China mainly participate in two electricity markets, i.e., the provincial and external markets. Based on the Copula theory, the combined MCP scenarios were generated using the following steps: (1) Construct MCP marginal cumulative distribution functions for the provincial and external markets. As the MCP , i t λ is related to its prediction value , λ i t and prediction error , ε i t , the distribution function for , i t ε needs to be firstly determined. At the early stage of reform, the monthly MCP data were accumulated for a short time and the distribution of prediction error was not clear. Therefore, according to the central limit theorem, suppose that ( )   , ; Distributed symmetrically, without reflecting the tail correlation As mentioned above, the cascade hydropower plants in southwestern China mainly participate in two electricity markets, i.e., the provincial and external markets. Based on the Copula theory, the combined MCP scenarios were generated using the following steps: (1) Construct MCP marginal cumulative distribution functions for the provincial and external markets. As the MCP , i t λ is related to its prediction value , λ i t and prediction error , ε i t , the distribution function for , i t ε needs to be firstly determined. At the early stage of reform, the monthly MCP data were accumulated for a short time and the distribution of prediction error was not clear. Therefore, according to the central limit theorem, suppose that namely, the MCP of market i is subject to the Gaussian distribution with the mean value ,i t λ and variance , i t σ . The probability density function (PDF) is: As mentioned above, the cascade hydropower plants in southwestern China mainly participate in two electricity markets, i.e., the provincial and external markets. Based on the Copula theory, the combined MCP scenarios were generated using the following steps: (1) Construct MCP marginal cumulative distribution functions for the provincial and external markets. As the MCP λ i,t is related to its prediction valueλ i,t and prediction error ε i,t , the distribution function for ε i,t needs to be firstly determined. At the early stage of reform, the monthly MCP data were accumulated for a short time and the distribution of prediction error was not clear. Therefore, according to the central limit theorem, suppose that the MCP of market i is subject to the Gaussian distribution with the mean valueλ i,t and variance σ i,t . The probability density function (PDF) is: Then, the marginal cumulative distribution function H is obtained by the integral of h(λ i,t ). That is: (2) Identify parameters for the Copula function. According to the characteristics of marginal distribution functions, appropriate Copula functions were firstly selected, and the unknown parameter α was estimated by the maximum likelihood estimate method, and the estimated value is: is the density function of the Copula function C(u 1 , u 2 ; α).
(3) Check the fitting performance of the Copula function. The squared Euclidean distance d 2 is used and its expression is: where, n is the number of samples, C(u 1 , u 2 ; α) is the Copula function, and G(u 1 , u 2 ; α) is the empirical distribution function. The smaller the value is, the better the fitting performance is.
(4) Generate the combined MCP scenarios. S × I dimensional data samples were generated by means of uniform discretization. S and I indicate, respectively, the number of samples and random variables. For each sample u s 1 , u s 2 , its combined MCP scenario vector ξ s = λ s 1,t , λ s 2,t could be obtained by the inverse operation of the marginal distribution function, i.e., λ i,t = H −1 i (u i ), and the probability of the scenario is P r (ξ s ) = 1/S.

Reducing MCP Scenarios
Generally, the results will be more accurate with more simulated scenarios, but the computation scale and time will increase accordingly. To guarantee the result accuracy and simultaneously strike a balance with the computing workload, the number of scenarios needed to be properly determined. Clustering is a commonly used scenario reduction method, but the number of reserved clusters is difficult to determine. In this paper, the hierarchical clustering method combined with inconsistent values is introduced to reduce the MCP scenarios, and the silhouette value was used to measure the rationality of clustering results. The process of reducing scenarios is detailed below: (1) Take each MCP scenario as a cluster and calculate the distance between every two clusters. Then, the absolute distance is adopted.
(2) Merge the nearest two clusters into a new cluster and calculate the distance between the new cluster and other clusters.
(4) Calculate the inconsistent values after each merging and express them as {θ k , k = 1, 2, · · · , S − 1}, and where, θ k is the inconsistent value after the kth merging. S is the number of total original scenarios, d k is the distance between the two clusters of the kth merging, d k is the mean value of the distances calculated in the kth merging, and σ k is the standard deviation of the distances calculated in the kth merging. The inconsistent value was calculated by the inconsistent function in MATLAB. Each inconsistent value was a measure of separation between the two merged clusters, compared to the separation between subclusters merged within those clusters. If only two original MCP scenarios are merged in the kth merging, the inconsistent value is 0. Defining ∆θ k = θ k − θ k−1 , the larger the ∆θ k was, the worse clustering result of the kth merging was, and the better the clustering result of the k − 1th merging was. Generally, the larger the number of reserved clusters, the more the reserved scenarios remain consistent with the original MCP scenarios, and the computation time will increase accordingly. The final number of reserved clusters could be determined by reference to the variation of inconsistent values {θ k , k = 1, 2, · · · , S − 1}.
(5) The silhouette value proposed by Rousseeuw [28] is often used to evaluate the result of the cluster analysis, which can be calculated by Equation (9): where, β i is the silhouette value of scenario i, i = 1, 2, · · · , S, r is the cluster index that scenario i belongs to, a i is the average distance between scenario i and other scenarios in the cluster r, The occurrence probability of scenario ξ s is derived from Equation (11): where, Ω s is the set of scenarios involved in cluster s, ξ s,m is the scenario m involved in cluster s, λ s,m 1,t is the electricity price of market i in period t under scenario m, λ s i,t is the electricity price of market i in period t under representative scenario ξ s , m is the scenario index involved in cluster s, 1 ≤ m ≤ M, t is the time period, 1 ≤ t ≤ T, i is the market index, 1 ≤ i ≤ N, and N is the total number of sub markets for the participation of hydropower cascades.

Optimization Model
In this paper, the cascade hydropower plants are assumed as one generator in the electricity market and cleared by the MCP. The annual market is considered and each time period represents one month. The operation objective was to maximize the market benefits. Considering the different MCP scenarios in the future period, the objective function can be defined as Equation (12). where, e s i,t is the supply energy for market i in the time period t under the representative scenario ξ s , S t is the number of reserved scenarios in period t. The hydropower generation cost mainly involved the constant cost with little influence on the optimization results, and so it was not considered in the optimization model. For the cascade hydropower plants in southwestern China, the provincial and external markets are considered, i= 1 for the provincial market and i= 2 for the external market.
To get the proposed model more applicable to the scheduling in practice, the energy ratio for different electricity markets and the minimum supply energy for the provincial market, during the whole scheduling period, were also included in the constraints. The constraints involved in the model were as follows, in which the hydraulic constraints also can be seen in [8,9,13].
(1) Water balance equation: where, V m,t and V m,t+1 are the storage of reservoir m at the beginning and end of the time period t, respectively, m 3 ; Q m,t is the total inflow of reservoir m in the period t, q m,t , d m,t , and I m,t are the power release, non-power release, and natural inflow of reservoir m in the time period t, respectively, m 3 /s, m = 1, 2, · · · , M is the reservoir index, m − 1 is the index of the upstream reservoir of reservoir m, and ∆T t is the number of seconds in the period t.
(2) Initial and end reservoir level constraints: where, Z m,1 and Z m,T are the initial and end reservoir levels of planning horizon, respectively, which are specified before the optimization process, m. In the case study, Z m,1 and Z m,T were specified based on the historical records for practical operation.
where, V m,t and V m,t are the minimum and maximum storage for reservoir m at the beginning of the time period t, respectively, m 3 . (6) Market-to-market energy ratio constraint: where, θ is the energy ratio of cascade hydropower plants involved in the provincial and external markets, which is agreed by the power producers in the river basin, the provincial development and reform commission, and other parties, e i,t is the supply energy of cascade hydropower plants for market i in the time period t, MWh, and E 1 and E 2 are the supply energy of cascade hydropower plants in the provincial and external markets, respectively, during the whole optimization horizon, MWh. (7) Minimum supply energy for provincial market: where, e min 1,t is the minimum supply energy of the cascade hydropower plants in the provincial market in period t, MWh.
(8) Relationship between fore-bay water level and reservoir storage: where, Z m,t is the fore-bay water level for reservoir m in period t, m V m,t is the storage for reservoir m in the period t, m 3 , and the coefficients (a i,m , i = 1, 2, 3, 4) are fitting parameters for reservoir m.
(9) Relationship between tailrace level and total release: where, ZT m,t is the tailrace level for reservoir m in the time period t, m, and the coefficients (b i,m , i = 1, 2, 3, 4) are fitting parameters for reservoir m.

Solution Method
The solution method typically depends on the characteristics of the optimization model. The commonly used solution methods for cascaded hydropower plants fall into three categories: Mathematical programming, dynamic programming and its improved algorithms, and artificial intelligence algorithms. The aforesaid optimization problem, including Equation (19), Equation (21), Equation (22), and other nonlinear constraints, is a typical nonlinear problem. Especially, the proposed model cannot meet the Markovian character, because the market-to-market energy ratio constraint, i.e., constraint (6), is included. Consequently, the dynamic programming and similar algorithms are not applicable. Moreover, the intelligent algorithms have problems, such as long computing time, easily trapped into local optimization, and unstable calculation results. Therefore, nonlinear programming (NLP) in mathematical programming is chosen in this paper.
In this paper, the solution procedure involves two main steps. At step 1, the combined MCP scenarios are generated and reduced by using MATLAB software; at step 2, based on the reduced scenarios, the NLP is implemented by the Global Solver, available in LINGO 17.0 (http://www.lindo. com/), to solve the proposed model. LINGO is a special software for solving mathematical programming problems, such as linear programming, non-linear programming, quadratic programming, and integer programming problems. It has the advantages of easy operation, fast computation, and stable calculation results. The basic computational principle and detailed computational procedure of the Global Solver for solving NLP problems have been given in literatures [29,30]. It mainly includes three steps: (1) With the branch-and-bound algorithm, the entire feasible solution spaces are segmented into series sub-problems; (2) a slack problem for each sub-problem can be obtained by convexity and linearization transformation and solved by the linear programming method; and (3) the global optimal solution is obtained by traversing all the sub-problems (i.e., the whole feasible space). Figure 1 shows the detailed calculating process. To set the parameters and calibrate the Global Solver, read the manual provided in LINDO.
scenarios are generated and reduced by using MATLAB software; at step 2, based on the reduced scenarios, the NLP is implemented by the Global Solver, available in LINGO 17.0 (http://www.lindo.com/), to solve the proposed model. LINGO is a special software for solving mathematical programming problems, such as linear programming, non-linear programming, quadratic programming, and integer programming problems. It has the advantages of easy operation, fast computation, and stable calculation results. The basic computational principle and detailed computational procedure of the Global Solver for solving NLP problems have been given in literatures [29,30]. It mainly includes three steps: (1) With the branch-and-bound algorithm, the entire feasible solution spaces are segmented into series sub-problems; (2) a slack problem for each subproblem can be obtained by convexity and linearization transformation and solved by the linear programming method; and (3) the global optimal solution is obtained by traversing all the subproblems (i.e., the whole feasible space). Figure 1 shows the detailed calculating process. To set the parameters and calibrate the Global Solver, read the manual provided in LINDO.

Study System
The proposed model was applied to the cascade hydropower plants of the Wu River in the Guizhou Province, China. Seven hydropower plants have been built on the mainstream of the Wu River, with a total installed capacity of 8345 MW. Figure 2 shows the distribution chart of the Wu River cascade, with reservoir characteristics in Table 2. Four hydropower plants have long-term regulation abilities in the studied cascade, and they are Hongjiadu, Dongfeng, Wujiangdu, and Goupitan. The storage capacities of the other three reservoirs were much smaller, and they were

Study System
The proposed model was applied to the cascade hydropower plants of the Wu River in the Guizhou Province, China. Seven hydropower plants have been built on the mainstream of the Wu River, with a total installed capacity of 8345 MW. Figure 2 shows the distribution chart of the Wu River cascade, with reservoir characteristics in Table 2. Four hydropower plants have long-term regulation abilities in the studied cascade, and they are Hongjiadu, Dongfeng, Wujiangdu, and Goupitan. The storage capacities of the other three reservoirs were much smaller, and they were addressed as runoff hydropower plants in the long-term operation. Thus, the cascade power decisions were distributed through the operation of these four large reservoirs. In this paper, the long-term provincial and external electricity markets are considered, assuming that the energy ratio between the markets is 3:1. The estimated electricity MCP and the standard deviation of estimated error are treated as known and given in Table 3. The monthly local inflow data in 2016 were used for simulation and shown in Figure 3. In the later text, market 1 represents the provincial electricity market and market 2 represents the external electricity market. All the basic data of the study system was collected from the dispatching center of the Wu River basin. long-term provincial and external electricity markets are considered, assuming that the energy ratio between the markets is 3:1. The estimated electricity MCP and the standard deviation of estimated error are treated as known and given in Table 3. The monthly local inflow data in 2016 were used for simulation and shown in Figure 3. In the later text, market 1 represents the provincial electricity market and market 2 represents the external electricity market. All the basic data of the study system was collected from the dispatching center of the Wu River basin.     addressed as runoff hydropower plants in the long-term operation. Thus, the cascade power decisions were distributed through the operation of these four large reservoirs. In this paper, the long-term provincial and external electricity markets are considered, assuming that the energy ratio between the markets is 3:1. The estimated electricity MCP and the standard deviation of estimated error are treated as known and given in Table 3. The monthly local inflow data in 2016 were used for simulation and shown in Figure 3. In the later text, market 1 represents the provincial electricity market and market 2 represents the external electricity market. All the basic data of the study system was collected from the dispatching center of the Wu River basin.

Scenarios Analysis
The Copula function should be adopted to fit the MCPs in the provincial and external markets. In this paper, it is assumed that the prediction error is the Gaussian distribution with the mean value of 0. Therefore, the probability density distribution was symmetric for the MCPs in the provincial and external markets. Three symmetrically distributed Copula functions in Table 1, including the bivariate Gaussian-Copula function, bivariate t-Copula function and bivariate Frank-Copula function were candidate functions, but the final closed function should be subject to the fitting performance.
There were 1000 groups of samples generated discretely. Based on the marginal distribution functions of MCPs in the provincial and external markets, the parameter estimation result of different Copula functions is shown in Table 4. The fitting performance of different Copula functions is tested by Equation (7), and the result is shown in Table 5. It can be seen that the values of the Gaussian-Copula function were the smallest, which indicates the better fitting performance of the Gaussian-Copula function in each month compared with the t-and Frank-Copula functions. Therefore, the Gaussian-Copula function was used to describe the correlation between the MCPs in the provincial and external markets. Taking April as an example, the Gaussian-Copula density and distribution function are shown in Figure 4.

Scenarios Analysis
The Copula function should be adopted to fit the MCPs in the provincial and external markets. In this paper, it is assumed that the prediction error is the Gaussian distribution with the mean value of 0. Therefore, the probability density distribution was symmetric for the MCPs in the provincial and external markets. Three symmetrically distributed Copula functions in Table 1, including the bivariate Gaussian-Copula function, bivariate t-Copula function and bivariate Frank-Copula function were candidate functions, but the final closed function should be subject to the fitting performance.
There were 1000 groups of samples generated discretely. Based on the marginal distribution functions of MCPs in the provincial and external markets, the parameter estimation result of different Copula functions is shown in Table 4. The fitting performance of different Copula functions is tested by Equation (7), and the result is shown in Table 5. It can be seen that the values of the Gaussian-Copula function were the smallest, which indicates the better fitting performance of the Gaussian-Copula function in each month compared with the t-and Frank-Copula functions. Therefore, the Gaussian-Copula function was used to describe the correlation between the MCPs in the provincial and external markets. Taking April as an example, the Gaussian-Copula density and distribution function are shown in Figure 4.   It is well known that the more discrete scenarios are, the more precise the description of the original distribution will be. However, the typicality of each scenario will be poor, and the workload for computation will be more. Upon the test, there are 1000 sets of scenarios generated by monthly discretization. After scenario reduction, the minimum number of scenarios from January to December is 47, and the maximum number is 64. A diagram of silhouette values is shown in Figure 5. It was found that the silhouette values were positive for each original scenario after clustering, and more than 95% silhouette values were greater than 0.2 for the original scenarios. This indicated that the classification number determined by the inconsistent values was appropriate, and the scenario set after clustering could be used to describe the original scenario set.  It is well known that the more discrete scenarios are, the more precise the description of the original distribution will be. However, the typicality of each scenario will be poor, and the workload for computation will be more. Upon the test, there are 1000 sets of scenarios generated by monthly discretization. After scenario reduction, the minimum number of scenarios from January to December is 47, and the maximum number is 64. A diagram of silhouette values is shown in Figure  5. It was found that the silhouette values were positive for each original scenario after clustering, and more than 95% silhouette values were greater than 0.2 for the original scenarios. This indicated that the classification number determined by the inconsistent values was appropriate, and the scenario set after clustering could be used to describe the original scenario set.

Effectiveness Analysis
The effectiveness of the proposed model was verified by the comparison with the conventional model. In the conventional model, the electricity price factor was not included, but only to achieve the maximum expected power generation in the scheduling horizon, maintain an adequate security of supply, and arrange the cascade plants to deliver electricity outwards in a certain proportion in case of spilled water in the flood seasons. In 2016, the natural inflow concentrated from June to August, accounted for 50% of the total inflow in the year. There was a risk of spilled water, and so the cascade hydropower plants were arranged to deliver electricity outwards, i.e., participate the external electricity market. Other constraints, such as the final water level, maximum and minimum power output, combined MCP scenarios, and so on, were the same for the two models. Figure 6 shows the water level of the cascade hydropower plants with regulation capacity, and the power output process were also obtained. It was found that the water level change trend of each

Effectiveness Analysis
The effectiveness of the proposed model was verified by the comparison with the conventional model. In the conventional model, the electricity price factor was not included, but only to achieve the maximum expected power generation in the scheduling horizon, maintain an adequate security of supply, and arrange the cascade plants to deliver electricity outwards in a certain proportion in case of spilled water in the flood seasons. In 2016, the natural inflow concentrated from June to August, accounted for 50% of the total inflow in the year. There was a risk of spilled water, and so the cascade hydropower plants were arranged to deliver electricity outwards, i.e., participate the external electricity market. Other constraints, such as the final water level, maximum and minimum power output, combined MCP scenarios, and so on, were the same for the two models. Figure 6 shows the water level of the cascade hydropower plants with regulation capacity, and the power output process were also obtained. It was found that the water level change trend of each hydropower plant was consistent for both the proposed model and conventional model. Generally, before the flood seasons, the water level of each hydropower plant fell, vacating the reservoir. During the flood seasons, the water was stored, the water level rose, and the generation water head increased. After the flood seasons, the water level reduced back to the control water level. These trends were in line with the cascade hydropower operation rules.  Figure 7 shows the supply energy for the two different markets obtained by the two comparison models under the assumptions that the same MCP scenarios were adopted. According to the analysis, on the basis of ensuring the power supply in the province market during the scheduling horizon, the proposed model responded better to the changing of market electricity prices and arranged the electricity delivery plan rationally. The power generation in the dispatching periods was 23.98 billion kWh, and the power generation revenue was 7.11 billion yuan. They were 24.07 billion kWh and 7.00 billion yuan, respectively, for the conventional model. It was found that, based on the same energy storage at the end of scheduling period, the proposed model resulted in an increase of 106.93 million yuan for the total income, although the power generation was less, compared with the conventional model. This proved that the electricity price factors would greatly influence the expected benefits in the market of the cascade hydropower plants. In terms of power generation in different months, the power generation of the proposed model during the dry seasons was much more than the conventional model, and the water level drops of each hydropower plant before the flood seasons were greater. This was because the power market was of a high proportion of hydropower and price in dry seasons was much higher. The hydropower plants with regulation ability preferred to use a self-regulating storage capacity to achieve the redistribution of water, so as to respond to price changes and improve their benefits. For example, in the proposed model, the water level drop at the Hogjiadu plant (i.e., the first stage of the cascade plants) was serious before the flood seasons, and it was close to the lowest operation water level. When the flood seasons came, the Hongjiadu plant was concentrated to store water and almost no power was generated. The water level of the downstream plants, including Dongfeng, Wujiangdu, and Goupitan, changed lightly between the flood seasons and dry seasons, which, to some extent, benefitted from inflow compensation of the Hongjiadu plant. Therefore, these downstream hydropower plants could keep a high water head for high generation efficiency. Figure 7 shows the supply energy for the two different markets obtained by the two comparison models under the assumptions that the same MCP scenarios were adopted. According to the analysis, on the basis of ensuring the power supply in the province market during the scheduling horizon, the proposed model responded better to the changing of market electricity prices and arranged the electricity delivery plan rationally. The power generation in the dispatching periods was 23.98 billion kWh, and the power generation revenue was 7.11 billion yuan. They were 24.07 billion kWh and 7.00 billion yuan, respectively, for the conventional model. It was found that, based on the same energy storage at the end of scheduling period, the proposed model resulted in an increase of 106.93 million yuan for the total income, although the power generation was less, compared with the conventional model. This proved that the electricity price factors would greatly influence the expected benefits in the market of the cascade hydropower plants. Figure 6. The optimization process of regulating reservoirs (M1 represents the proposed model, M2 represents the conventional model). Figure 7 shows the supply energy for the two different markets obtained by the two comparison models under the assumptions that the same MCP scenarios were adopted. According to the analysis, on the basis of ensuring the power supply in the province market during the scheduling horizon, the proposed model responded better to the changing of market electricity prices and arranged the electricity delivery plan rationally. The power generation in the dispatching periods was 23.98 billion kWh, and the power generation revenue was 7.11 billion yuan. They were 24.07 billion kWh and 7.00 billion yuan, respectively, for the conventional model. It was found that, based on the same energy storage at the end of scheduling period, the proposed model resulted in an increase of 106.93 million yuan for the total income, although the power generation was less, compared with the conventional model. This proved that the electricity price factors would greatly influence the expected benefits in the market of the cascade hydropower plants.

Robustness Analysis
The robustness of the proposed model was verified by the comparison with different price scenarios. Five typical MCP scenarios were selected for the comparison. Scenario 1 embodied the situation for good prediction, and actual electricity price was higher than the prediction value. Scenario 2 embodied the situation for good prediction and the actual electricity price was lower than the prediction value. Scenario 3 embodied the situation for perfect prediction and the actual electricity price was exactly the same as the prediction value. Scenario 4 embodied the situation for poor prediction and the actual electricity price was higher than the prediction value. Scenario 5 embodied the situation for poor prediction and the actual electricity price was lower than the prediction value. For different scenarios, the indicators for evaluating the prediction result, including the mean absolute percent error (MAPE), the mean average error (MAE) and the relative correlation coefficient (R), are listed in Table 6.
The comparison results are shown in Table 7. It was found that the benefits of other MCP scenarios were in the range of Scenario 1 and Scenario 5 for both models. Compared with the conventional model, the proposed model achieved better power generation benefits in different MCP scenarios. Even in the adverse scenario with a serious deviation in the price prediction, the proposed model still had relatively higher benefits. For instance, in Scenarios 4 and 5, the benefits were 7.135 billion yuan and 6.961 billion yuan for the proposed model. However, they were 6.959 billion yuan

Robustness Analysis
The robustness of the proposed model was verified by the comparison with different price scenarios. Five typical MCP scenarios were selected for the comparison. Scenario 1 embodied the situation for good prediction, and actual electricity price was higher than the prediction value. Scenario 2 embodied the situation for good prediction and the actual electricity price was lower than the prediction value. Scenario 3 embodied the situation for perfect prediction and the actual electricity price was exactly the same as the prediction value. Scenario 4 embodied the situation for poor prediction and the actual electricity price was higher than the prediction value. Scenario 5 embodied the situation for poor prediction and the actual electricity price was lower than the prediction value. For different scenarios, the indicators for evaluating the prediction result, including the mean absolute percent error (MAPE), the mean average error (MAE) and the relative correlation coefficient (R), are listed in Table 6. The comparison results are shown in Table 7. It was found that the benefits of other MCP scenarios were in the range of Scenario 1 and Scenario 5 for both models. Compared with the conventional model, the proposed model achieved better power generation benefits in different MCP scenarios. Even in the adverse scenario with a serious deviation in the price prediction, the proposed model still had relatively higher benefits. For instance, in Scenarios 4 and 5, the benefits were 7.135 billion yuan and 6.961 billion yuan for the proposed model. However, they were 6.959 billion yuan and 6.827 billion yuan for the conventional model, decreased by 176.022 million yuan and 132.952 million yuan, respectively. This indicated that the proposed model could better adapt to different price scenarios in the future, and could enhance the profitability and robustness for the generation scheduling.

Conclusions
The large-scale cascade hydropower plants in southwestern China challenge a new multi-market environment. Considering the price uncertainty and the correlations between multiple markets, this paper has proposed a novel optimization model of long-term generation scheduling for cascade hydropower plants to seek for the overall benefits. The main contributions and conclusions were as follows: (1) The correlation of MCP between multiple markets was described. The prediction price was treated as a random variable that consists of its prediction value and a random prediction error. Then, the MCP scenarios for the hybrid market environment were generated by using the Copula theory and further reduced by the hierarchical clustering method combined with inconsistent values.
(2) A novel long-term generation scheduling model was constructed based on the stochastic programming theory and the reduced MCP scenarios, in which the market-to-market supply energy ratio constraint and the minimum constraint of supply energy in the provincial market were creatively taken into consideration, so as to meet the actual engineering needs.
(3) The scenarios analysis showed that the Gaussian-Copula function has better fitting performance in describing the correlation between the MCPs of the provincial and external markets than other functions. The proposed reduction strategy for the MCP scenarios was effective, after the reduction of scenarios. The minimum number of scenarios from January to December was 47, and the maximum number was 64.
(4) The effectiveness analysis showed that the proposed model could better respond to the changing of MCP signals, and so the arrangement of power generation in different months and markets was more reasonable. This resulted in an increase of 106.93 million yuan for the total income, although the power generation was less, compared with the conventional model. (5) The robustness analysis showed that the proposed model could better adapt to different price scenarios in the future, even in the adverse scenario, with a serious deviation in the price prediction. This enhanced profitability and robustness for the generation scheduling.
It should be mentioned that the market reform in China is still at an early stage and the available samples are in a relatively small quantity. This paper only analyzed qualitatively the price correlation between the provincial and external markets but did not give the specific influence factors and their different roles. With continuous accumulation of data in the future, the interaction characteristics between the markets, consumers, and producers can be analyzed in detail, so as to better describe the price correlation between multiple markets.