#### Abstract

Electricity price risk assessment (EPRA) is essential for spot market analysis and operation. The statistical moments (i.e., the mean and standard deviation) of the price need to be assessed to support market risk control. This paper proposes a data-driven approach for EPRA based on the Gaussian process (GP) framework. Compared with the deep learning algorithms, GP has two merits: (1) the scale of training sample required is small and (2) the time-consuming hyperparameter tuning process is avoided. However, the direct application of GP for EPRA is not tractable due to the complicated discrete relationship between the system operating status and the electricity price. To deal with that, a data-driven EPRA framework is developed that contains a GP surrogate model for the direct current optimal power flow (DC-OPF) problem and a hybrid model-data-based hybrid electricity price calculation method. To guarantee the accuracy of EPRA, an adaptability criterion and a second verification process based on the Karush–Kuhn–Tucker (KKT) condition are developed to distinguish the samples with GP learning errors. Numerical results carried out on IEEE benchmark systems demonstrate that the proposed method can achieve exactly the same EPRA results as Monte Carlo (MC) simulation, which significantly improved the computational efficiency.

#### 1. Introduction

To reduce pollution and greenhouse gas emissions, a high share of renewable energy integration has become one of the basic characteristics of the smart grid [1–3]. With the development of renewable energy and the adoption of locational marginal pricing (LMP) methodology, the spot market is full of uncertainties, such as load deviation and renewable variation [4]. The abovementioned uncertainties cause the electricity price to ﬂuctuate violently, bringing significant operational and planning risks for electricity market participants.

Risk assessment can provide power system operators with prior knowledge and theoretical basis to ensure safe and stable power system operation [5–7]. In the spot market, electricity price risk assessment (EPRA) is crucial for independent system operators (ISOs) and market participants. However, it is more volatile and challenging to predict the fluctuations in electricity prices than the uncertainties of power production and consumption [8]. The reliable and efficient assessment method is the basis for spot market operation and risk control. Current studies focus on the risk caused by the electricity price fluctuation for the risk assessment in electricity markets. Reference [9] analyses the optimal electricity procurement problem for large consumers considering the electricity price fluctuation. Reference [10] proposes a value-at-risk (VaR) and conditional VaR (CVaR) assessment for electricity price risk based on historical data. Reference [4] uses the Monte Carlo simulation method in electricity price risk management. Reference [11] analyses the price risk of power portfolios in multimarkets based on the well-established mean-variance model.

For EPRA, the expectation and standard deviation of LMP need to be assessed to support market risk control of independent system operator*s* (ISOs) [12]. Generally, LMP can be obtained based on the direct current optimal power flow (DC-OPF) model, which is derived from the Lagrangian multipliers of the power balance constraint and transmission constraints, including an energy component and a congestion component [13]. Probabilistic optimal power flow (POPF) is able to comprehensively consider various uncertainties in the spot market and thus has become an effective tool to estimate LMP in the deregulated market [14, 15].

To solve the POPF problem, two main calculation approaches have been developed, namely, model-based and data-based approaches. The model-based methods can be roughly divided into analytical methods and simulation methods. Typical analytical methods, such as the point estimation method, construct representative samples according to the probability density function (PDF) of the uncertainty variables [16, 17]. EPRA results can be obtained according to the OPF solutions of representative samples, which is computationally efficient, but complicated mathematical derivations and strict assumptions are required. Typical simulation methods such as Monte Carlo (MC) simulations obtain EPRA results by using massive random generated samples that are carried out on the OPF model, which is reliable but computationally demanding [18, 19]. Recently, data-based machine learning methods have been widely applied in power system [20, 21], showing a promising way to achieve EPRA with high precision and fast computational speed. For POPF problem, the core idea of the data-based approach is to construct a data-driven surrogate model that treats the OPF problem as a functional mapping between the system operating status and the OPF solutions, thus greatly improving the computational efficiency of the POPF problem. In [22, 23], a deep neural network (DNN) approach for solving OPF problems was developed based on historical data and offline simulations. Reference [24] proposed a data-driven machine learning framework for the OPF problem considering the characteristics of the physical model. However, these data-driven approaches have several technical challenges between the POPF problem and EPRA. Several challenges need to be addressed. First, the discrete features of LMP are hard to be learned by the existing data-driven methods. Second, the data-driven methods, such as DNN-based approaches, usually require massive training samples, which may not align with the current spot market practice. Third, the inherent learning error of the data-driven methods is inevitable and may yield an unreliable EPRA result. To overcome the challenges mentioned above, our work combines the advantages of the two aforementioned approaches to develop a data-driven assisted electricity price risk assessment method based on the Gaussian process (GP) and the physical model of DC-OPF.

Compared with traditional DNN-based methods [25, 26], the GP is a novel machine learning technology that requires smaller training samples and fewer hyperparameters for learning [27, 28], making the GP align well with current industry practice. The GP is used extensively as a nonparametric regression tool in various scenarios, e.g., active learning [29], multitask learning [30, 31], manifold learning [32], and optimization [33]. However, the learning error of GP is also inevitable. Further advanced technology is required to accurately learn the features of LMP.

The objective of EPRA is to obtain the statistical moments of the LMP according to various uncertainties of the system operating status. Data-driven methods can build a surrogate model with cheap computation cost to replace the time-consuming LMP calculation process. Note that an efficient data-driven EPRA algorithm needs not only high precision and fast computing speed but also good generalization capability with limited training sample, which makes the unique characteristics of GP (e.g., fast training, less intervention, and small sample requirement) an ideal candidate. However, unlike the POPF problem, the relationship between the input (the system operating status) and the output (the LMP) is rather complex due to the discontinuous property of LMP. Hence, direct learning LMP using data-driven methods is intractable, which will be shown in the simulation results. Fortunately, the physical model of DC-OPF is known, and this motivates us to develop a new framework to achieve the LMP assessment based on the POPF results by including its physical characteristics.

To this end, a data-driven EPRA approach is proposed based on both physical models and historical data. Compared with existing methods, the proposed data-driven method combines the advantages of the model-based and data-based approaches to achieve a more efficient EPRA without accuracy loss. Specifically, we embed the GP surrogate model for DC-OPF into the model-based EPRA process to improve the computational efficiency of the traditional model-based method. By providing the strict judging criteria (adaptability criterion and a second verification) to determine the inaccurate samples obtained by the proposed data-driven method, the accuracy of EPRA is guaranteed. Note that the proposed approach is a general method for EPRA, even in a specific scenario with limited samples. It has the following advantages: (1) the accuracy of EPRA is maintained. The EPRA results obtained through our approach are exactly the same as those of the MC method. (2) The training sample size for learning the LMP has significantly reduced thanks to the GP. (3) The efficiency of EPRA is improved because a large proportion of the time-consuming POPF process is replaced by direct GP mapping.

The main contributions of this paper are summarized as follows:(1)A data-driven framework is proposed to reduce the scale and accelerate the computational speed of the EPRA problem. Specifically, to avoid directly learning the LMP, a GP surrogate model for the DC-OPF problem is developed to identify key information for LMP calculation (e.g., the marginal generators and congested transmission lines). Then, a model-data hybrid EPRA method is proposed by solving a set of linear equations. The proposed method can significantly improve the efficiency of the EPRA without compromising its accuracy.(2)Under this framework, a model-based adaptability criterion and a second verification for EPRA are developed to determine inaccurate samples. Before using the sample with marginal generators and congested transmission lines identified by the GP to calculate the LMP, physical model information is used to distinguish the samples with learning errors. Hence, the accuracy of EPRA is guaranteed.

The rest of the paper is organized as follows: the data-driven EPRA framework is developed in Section 2. Section 3 presents the proposed GP surrogate model for DC-OPF. Numerical results are analyzed in Section 4, and finally, Section 5 concludes the paper.

#### 2. The Data-Driven Framework for EPRA

In the spot market, the LMP arises from an economic dispatch. Speciﬁcally, the system operator solves a DC-OPF problem for the optimal economic generation that meets the variational load and renewable energy while satisﬁes the generation and transmission constraints [34]. In fact, there is a linear relationship between the LMP and the Lagrangian multiplier of the DC-OPF model. The relationship relies on the marginal generator and congested transmission line, which can be obtained through the DC-OPF solutions. Hence, the key idea of the proposed data-driven framework is to build a GP surrogate model for the DC-OPF problem to identify the marginal generator and congested transmission line, thus improving the computational efficiency of EPRA. The physical characteristics of the DC-OPF model are considered to ensure accuracy. In this section, the LMP formulation is first studied using a general DC-OPF formulation and its Karush–Kuhn–Tucke*r* (KKT) condition. Then, the data-driven approach is proposed for EPRA.

Note that this paper is focused on the LMP risk arised from the uncertainty of load and renewable energy. Within the proposed scope, we assume the topology of power grid is invariable. The uncertainties from equipment fault are ignored.

##### 2.1. A General DC-OPF Formulation

A general DC-OPF model for economic dispatch is formulated as follows:

Objective function:

Constraints:where *c*_{i} is the generator cost for production; *PTDF*_{li} represents the power transfer distribution factor of Bus *i* to Line *l*; is the transmission power flow of Line *l*, which is denoted as *PF*_{l}; and *P*_{i} and *D*_{i} are the generator output and demand quantity, respectively. Note that renewables are treated as negative loads in this paper, which are included in *D*_{i}.

The linear objective function of the DC-OPF model is designed to minimize the operating costs associated with supplying real power to meet the demand requirement. Equation (2) is the system power balance equation, and is the corresponding Lagrangian multiplier. The constraints in (3) and (4) limit the transmission line power flow, and and are the Lagrangian multipliers of the upper and lower transmission limit constraints, respectively. The constraints in (5) are the operational limits for the real generator power, and , are the Lagrangian multipliers of the upper and lower limits of the generator output constraints, respectively.

##### 2.2. Deduction for the LMP Formulation

To understand the internal relationship between the LMP and DC-OPF problem, the KKT condition is used to analyze the properties of LMP.

###### 2.2.1. The LMP Formulation

According to the KKT condition, we derive the relationships among the LMP, the Lagrangian multiplier of the power balance , and the dual multiplier of the transmission line limits and . Note that in the following analysis, the saddle point used by the KKT condition corresponds to the global optimum of the OPF model.

To obtain the LMP for the EPRA, the Lagrangian function of the DC-OPF models (1)–(5) is denoted by *LF*, as follows:

Then, the LMP for the load at Bus *i* is derived from the Lagrangian function in (6) as

From (7), we note that the LMP is obtained by the marginal generator through and the congested transmission line through and . The relationship between the marginal generators and congested transmission lines is discussed in the next section.

###### 2.2.2. The Relationship between Marginal Generators and Congested Transmission Lines

Based on the KKT condition, the following equation for generator *i* can be obtained:

For the marginal generators, , and equation (8) can be expressed aswhere and *I*_{MG} and *N*_{MG} are the set and the number of marginal generators, respectively. Then, the matrix form of the equations above is analyzed. The matrix form of (9) is

There are *N*_{MG} equations and (*L *+* *1) unknown variables in (10), which cannot yield a unique solution. To solve these equations, (*L *+* *1* *−* N*_{MG}) variables should be determined in advance. Hence, the information of the transmission line constraints is introduced. Note that for transmission lines without congestion, .

Hence, the number of congested transmission lines *N*_{CL} is

Then, the equations in (10) can be rewritten as

Based on the above discussion, we know that if the marginal generators and the congested transmission lines are known, the LMP can be calculated by solving a set of linear equations in (12). Hence, the identification of marginal generators and congested transmission lines is the key problem for EPRA.

##### 2.3. The Proposed Framework

In this section, a trustworthy data-driven framework is proposed for EPRA. Based on historical data, a GP surrogate model is developed for DC-OPF to identify the marginal generators and congested transmission lines, which will be introduced in the next section. After identification, the LMP can be obtained immediately without solving the time-consuming DC-OPF problem. Unfortunately, inherited learning error in data-driven methods, including the proposed GP surrogate model, is unavoidable, making the identification of marginal generators and congested transmission lines unreliable. To overcome this challenge, an adaptability criterion is developed based on the KKT condition of the physical model discussed in Section 2.2. The proposed framework for EPRA is illustrated in Figure 1. The three steps are described as follows: *Step 1: Marginal Generators and Congested Transmission Line Identification.* The power system operating data (e.g., *D*_{i}) are collected from historical data or by running a Monte Carlo simulation, and they are input into the trained GP surrogate model for DC-OPF. For each sample, the DC-OPF outputs (e.g., the generator output and transmission power flow) are immediately obtained. Then, the marginal generators and congested transmission lines can be identified with high precision and fast computation. *Step 2: Distinguishing the Error Samples.* According to the discussion presented in Section 2.2, for the DC-OPF problem, we note that the number of marginal generators *N*_{MG} is *n *+* *1 if the number of congested transmission lines *N*_{CL} is *n*. This natural property motivates us to adapt it to the proposed framework to distinguish the samples with learning errors. Hence, for each sample with marginal generators and congested transmission lines identified by the GP surrogate model, the proposed adaptability criterion is as follows: It should be noted that the proposed adaptability criterion is not strictly complete. There are a very few samples that meet the adaptability criterion but have learning errors. For these samples, the LMP error of all the buses can diverge significantly from the real value because the marginal cost of generation is changed. Hence, a second verification process is proposed based on historical data, as follows: where mean(*LMP*_{i}) is the mean value of LMP at Bus *i*, which is obtained from historical data. In this work, we set *p* as 50%. *Step 3: LMP Calculation and EPRA.* Based on equation (13), if the sample meets the adaptability criterion, the LMP can be calculated by solving a set of linear equations; otherwise, we should run DC-OPF to obtain the LMP. Then, with all the LMP data in hand, the EPRA can be completed by performing statistical analysis on the LMP data. Note that the EPRA results are obtained by the Monte Carlo (MC) simulation, which generates massive random samples that are carried out on the OPF model. The proposed method provides an effective tool to replace the time-consuming OPF calculation process to reduce the computationally demanding of the MC simulation. Hence, the convergence of the proposed method is the same as the traditional MC.

#### 3. GP Surrogate Model for DC-OPF

This section ﬁrst brieﬂy introduces the GP. Then, the GP surrogate model for DC-OPF is proposed to identify the marginal generators and the congested transmission lines. The basic idea of the proposed GP surrogate model is to use the GP to replace the time-consuming optimization process of DC-OPF, as illustrated in Figure 2. The complicated DC-OPF features can be represented by the mapping relationship *f*: *D*_{i} ⟶ *P*_{i}, *PF*_{l}, which is the learning target of the GP.

##### 3.1. A Brief Introduction to the GP Regression Method

The GP is generally used to solve hard regression and classification problems. It is attractive because of its flexible nonparametric nature and computational simplicity. In nonparametric statistics, the regularity of a relationship can be postulated without requiring the dataset to be focused on an easily describable class. This efficient property allows the GP to predict the functional behavior inside and outside of the input domain with a small sample size [35].

The GP is introduced for regression in this paper. We denote the regression function by *f*(·), which is the output of the GP surrogate model. Its corresponding input vector of *p* dimensions is denoted as **x**. For a GP regression problem, a finite collection of training sample inputs ** x** is denoted as [

**x**

_{1},

**x**

_{2},…,

**x**

_{n}]. Accordingly, the corresponding output

*f*(

**x**) can be denoted as [

*f*(

**x**

_{1}),

*f*(

**x**

_{2}),…,

*f*(

**x**

_{n})]. According to [36], the model output

*f*(

**x**) is expected to follow a joint multivariate normal probability distribution, as follows:where

*m*(·) represents the mean function and

*C*(·,·) is a kernel function representing the covariance function. Then, (15) can be rewritten aswhere

*X*is an

*n*×

*p*matrix denoted by [

**x**

_{1},

**x**

_{2},…,

**x**

_{n}]

^{T}. Then, considering that there is independent, identically distributed noise in the model output

**f**(

**X**), the realizations

**Y**can be formulated aswhere

*σ*

^{2}and

**I**

_{n}are the variances of the noise and an

*n*-dimensional identity matrix, respectively. Note that the noise is always accounted for in practical implementation in the GP output.

To infer the GP regression output with noise from the abovementioned sample set (**X**, **Y**), a Bayesian inference framework is introduced. It is well known that a Bayesian posterior distribution of the model output can be inferred from a Bayesian prior distribution of *y*(**x**)|**x** and the likelihoods obtained from the realizations. For a test sample input set **X**_{t}, the Bayesian prior distribution of **Y**_{t} can be expressed as

Combined with the training sample set (**X**, **Y**), the joint distribution of **Y** and **Y**_{t}|**X** can be formulated as follows:where **C**_{11} = **C**(**X**, **X**) + *σ*^{2}**I**_{n},**C**_{12} = **C**(**X**, **X**_{t}), **C**_{21} = **C**(**X**_{t}, **X**), and **C**_{22} = **C**(**X**_{t}, **X**_{t}) + *σ*^{2}**I**_{nt}. Then, based on the rules of the conditional Gaussian distribution, the Bayesian posterior distribution of **Y**_{t} can be inferred from the training sample set (**Y**, **X**) and test sample input set **X**_{t}. It follows a Gaussian distribution *N*(*μ*(**X**_{t}), Σ(**X**_{t})). The expected value of **Y**_{t} can be expressed as follows:

Thus far, the general form of GP regression has been derived. Equation (20) can be used as a surrogate model for a complicated DC-OPF model with a low computational cost. Further details about the GP can be found in [35, 36].

##### 3.2. Proposed GP Surrogate Model for DC-OPF

The basic DC-OPF model introduced in Section 2 can be further expressed as the following linear programming (LP) problem:where **y** is the output set, including the generation output and transmission line power flow, which are the key variables for determining the marginal generators and congested transmission lines, and c(·) is the associated cost function. **A** and **B** are the corresponding matrices concerning vectors ** x** and

**, respectively. With a random**

*y***x**, the DC-OPF in (21) can be cast as a POPF problem.

In the proposed approach, a GP surrogate model is developed for the DC-OPF problem with a lower computational cost to improve the effectiveness of POPF. To this end, the steps for the proposed GP surrogate method are illustrated below.

###### 3.2.1. Training Sample Generation

To construct the GP surrogate model described in (20) for the DC-OPF problem, the training sample set **D** = (**X**, **Y**) can be obtained from historical operational ISO data or by running an MC simulation where the uncertainty vector *ω* is sampled and the resulting DC-OPF output is calculated for a large number of samples. In particular, the Latin hypercube sampling method is implemented due to the small sample requirement of the GP. Here, each row of **X** is an *I-*dimensional uncertain input vector, including the load demands *D*_{i} of all buses. The output matrix **Y** contains columns of *I*_{G} + *I*_{L} output variables as the generator output **P** and transmission line power flow **PF** corresponding to each input vector **x**.

###### 3.2.2. GP Surrogate Model Construction

With the training dataset **D** = [**X**, **Y**], we choose the squared exponential (SE) covariance kernel function for our regression problem, i.e.,

The hyperparameters *ξ** **= ** *(*τ*, *l*) can be estimated by the Gaussian maximum likelihood estimator (MLE) method, which is optimal under the Gaussian assumption and is easy to implement [31]. Equation (17) with hyperparameters can be rewritten as

Based on the MLE, we obtain

The marginal log-likelihood can be expressed as

After utilizing a gradient-based optimizer, the hyperparameters are obtained while the GP surrogate model for DC-OPF is fully constructed.

###### 3.2.3. The Key Information Identification for EPRA

For the output of the proposed GP surrogate model (e.g., **P** and **PF**), the learning error is not avoided. To identify the key information for EPRA (e.g., the marginal generator and congested transmission line) and address the learning error effect, a relaxing factor *ε* is introduced. Then, the marginal generators and congested transmission lines can be identified according to the following equations:

Marginal generator *i*:

Congested transmission line *l*:

After obtaining the marginal generators and congested transmission lines of each sample, based on the framework proposed in Section 2, the EPRA can be achieved in short order without accuracy loss.

#### 4. Numerical Results

In this section, the proposed method is tested on the IEEE 30-bus system to illustrate its effectiveness, while the IEEE 118-bus system is used to demonstrate its scalability to the larger system. All simulations are performed on a PC equipped with an AMD Ryzen 5 3600X 6-Core Processor CPU @ 3.80 GHz with 16 GB RAM. The algorithm is implemented in MATLAB.

The following methods are compared:(i)M0: Monte Carlo simulation (benchmark)(ii)M1: a neural network method based on SAE [37](iii)M2: a neural network method based on SDAE [38](iv)M3: a stacked extreme learning machine [24](v)M4: the Gaussian process [36](vi)M5: the proposed method

The hyperparameter settings of each data-driven method are shown in Table 1, which are obtained according to the artificial experience and reference [39]. The learning error of M1∼M5, which is deﬁned as the average difference between the results obtained with M1∼M5 and M0, is evaluated as follows:where and *p* are the output of the data-driven method and its dimensions; is the output of Monte Carlo simulation as benchmark; *K* represents the number of testing samples; and is the MAPE index.

##### 4.1. Evaluation of the Proposed Approach on the IEEE 30-Bus System

###### 4.1.1. Learning Performances of DC-OPF and LMP

We first show that the EPRA problem of LMP assessment is more complicated than the OPF problem, which cannot be learned directly by data-driven methods. The learning performances of the LMP and DC-OPF outputs learned by M1–M4 are compared on IEEE 30 in Table 2. For M1–M3, the number of training samples is set as 10000. For M4 and M5, the number of training is set as 200. The number of testing samples is set as 10000 for all methods.

The results show that directly learning the LMP is intractable for data-driven methods because of its discontinuous property. Additionally, the DC-OPF problem is more comfortable to learn, making the proposed method based on the learning output of DC-OPF reasonable.

###### 4.1.2. Effectiveness of the Proposed Method

To demonstrate the benefits achieved by the proposed approach, we compare the LMP errors of M1–M5 in the IEEE 30-bus system, as shown in Table 3.

Several conclusions can be drawn, as follows:(1)Among all the methods for EPRA, the proposed method (M5) achieves the best accuracy, which is exactly the same as that of the benchmark method. In the proposed framework, the learning error of the GP is filtered out by the identification process and model-based adaptability criterion.(2)Compared with M0, the testing time decreases by 59.34%. This shows that the computational efficiency of LMP is significantly improved by the proposed method without accuracy loss.(3)Comparing M4 with M1–M3, the results show that the GP can achieve a similar accuracy with a small sample size. However, the learning errors of the data-driven methods are unavoidable, even with a large number of training samples.(4)For M1–M5, because of the superior quality of the GP (e.g., its small sample requirement), the training sample size of the proposed method is much smaller than those of M1–M3, which aligns well with the current industry practice.

The errors of the EPRA results are compared, as shown in Figure 3. For the proposed method, the EPRA results for the mean and standard deviation are very accurate. However, for M1∼M4, the estimation results for the mean are relatively accurate, while the standard deviations are far from the real value obtained by M0, so they cannot achieve an accurate EPRA. In particular, as seen in Figure 3, the error is abnormally large in some specific buses (e.g., Bus 8 and Bus 24). This is because the LMPs at these buses fluctuate much more than other buses, making it challenging to be directly learned by the data-driven methods, resulting in false information being passed to ISOs and market participants, leading to severe market risks.

**(a)**

**(b)**

##### 4.2. Results on the IEEE 118-Bus System

The IEEE 118-bus system is used to demonstrate the scalability of the proposed method. The learning error and EPRA results are given in Table 4 and Figure 4, respectively. The test sample number is set to be the same as for the IEEE 30-bus system.

**(a)**

**(b)**

The results show that the proposed method can guarantee EPRA accuracy even in a large case while improving the computational efficiency. It also shows that the statistical moments (i.e., the mean and standard deviation) diverge far from the real value at Buses 72∼110, which demonstrates that EPRA cannot be achieved by data-driven direct learning.

#### 5. Conclusions and Future Work

This paper proposes a data-driven framework for EPRA. Specifically, a GP surrogate model is developed to identify the marginal generators and congested transmission lines of DC-OPF. This paves the way for improving the efficiency of EPRA. Based on the KKT condition, an adaptability criterion is proposed to identify samples with learning errors. The simulation results show that the proposed method increases EPRA accuracy. Comparisons with recent data-driven methods show that the proposed approach can greatly improve the computational efficiency of EPRA without compromising its accuracy. It is also shown that direct learning for LMP may not be tractable due to the problem complexity. Future work will consider more uncertain boundary conditions for the power system.

As shown in this paper, although the inherent learning error in data-driven methods is unavoidable, helpful reference information can be analyzed to increase the computational speed of EPRA. Hence, the idea of improving the efficiency of spot market risk analysis through data-driven techniques is worthy of further study. This paper focuses on the EPRA in the spot market, which is commonly used in the Guangdong Electric Power Trading Center of China to manage LMP risk. Our study shows that the GP is capable of learning the complex relationship between a set of key information (e.g., the marginal generator and congested transmission line) and system operating conditions with minimal training samples and a fast training speed. Therefore, for more complicated problems considering bid price, generating capacity, and unit commitment, the GP is also a promising tool for extracting useful information from historical data. However, effectively utilizing GP techniques regarding the specific properties of these problems is worthy of further exploration. In addition, the physical models of spot market operation are known by the ISO. Therefore, improving the learning performance in combination with power domain expertise should be investigated further. To further speed up the calculation, it is necessary to improve the learning accuracy considering the characteristics of the physical model and to make the proposed method applicable to more scenarios according to the selection of samples in future work.

#### Data Availability

The data that support the findings of this study are available on request from the corresponding author. The data are not publicly available due to privacy or ethical restrictions.

#### Conflicts of Interest

The authors declare that they have no conflicts of interest.

#### Acknowledgments

This project was supported by the Science and Technology Project of Guangdong Electric Power Trading Center Co., Ltd. (GDKJXM20185365).