# Electricity Load and Price Forecasting Using a Hybrid Method Based Bidirectional Long Short-Term Memory with Attention Mechanism Model

Authors: William Gomez; Fu-Kwun Wang; Zemenu Endalamaw Amogne
Journal: International Journal of Energy Research
 (2023)

---

## Abstract

The accuracy of forecasting short- or medium-term electricity loads and prices is critical in the energy market. The amplitude and duration of abnormally high prices and load spikes can be detrimental to retailers and production systems. Therefore, predicting these spikes to effectively manage risk is critical. In this paper, a novel hybrid method that combines ensemble empirical mode decomposition (EEMD) algorithm and a bidirectional long short-term memory with attention mechanism (BiLSTM-AM) model is proposed to predict electricity loads and prices. A simple approach is proposed to determine the number of intrinsic mode functions (IMFs) that decompose raw data using EEMD to avoid overdecomposition, irrelevant components, and high computational cost. Each selected mode is then modeled with BiLSTM-AM to obtain a predicted sequence. These sequences are summed and then reverted to obtain the final predicted value. The proposed method is validated using two datasets (PMJ and Australian Energy Market Operator) with different time intervals to demonstrate the generality and robustness of the forecasts, especially in temporal valley or peak forecasting. The results show that the proposed method outperforms other methods in prediction accuracy and spike-capturing ability, with EEMD reducing the mean absolute percentage error (MAPE) by 53%, 54%, and 60%, respectively. In the three forecast periods, the average MAPE andR̲2 are 0.097 and 0.92, respectively. Furthermore, we use Kolmogorov-Smirnov predictive accuracy (KSPA) test and model confidence set (MCS) test to validate the superiority of the proposed model. The results demonstrate its suitability, reliability, and performance in short- and medium-term forecasting.


---

## Body

## 1. Introduction

The energy market is one of the most important markets in the world, and its price fluctuations affect the global economy. In the energy industry, reliable power system forecasts allow market participants, especially production sectors, to base their scheduling and production plans on forecast information. Investors are interested in short-term and long-term economic price trends. In addition to helping reduce volatility in speculative markets, energy forecasts can also help predict major financial and economic events, such as the 2009 financial crisis or energy shortages from COVID-19, to ensure better industrial and economic progress.The ability to forecast energy market loads and prices can enable consumers to adjust their consumption or supply while still achieving their economic and environmental goals. The production, transmission, distribution, and markets of electrical energy require accurate price and load forecasts from suppliers, investors, financial institutions, and others. Electricity operators can develop effective market plans to increase the benefits of energy management through accurate energy price forecasts [1]. Therefore, developing predictive models that can provide accurate prices is a major priority. Different types of literature recommend various forecasting methods, ranging from traditional statistical time series to modern data analytic models. However, the performance of a single technique-based predictive model is insufficient due to inherent limitations. Energy prices are difficult to predict due to rapid peaks and troughs caused by changes in supply and demand due to intraday system constraints [2].Therefore, a generalized forecasting model will help in grid monitoring, optimal dispatch, and energy production and conversion control. Not much works has attempted to model electricity load and price with the same model to examine its generality, flexibility, and reliability. Different time-range datasets were also tested to test its ability to handle different variational frequency data and better capture valley-peak. We applied a simple technique to determine the number IMFs to avoid over-decomposition and irrelevant components and reduce computation time.High demand for electric vehicles and microgrids will lead to increased energy distribution and frequent market surges. So, electricity load and price forecasting is becoming more challenging. Much of the recent literature looks more at models that can improve accuracy and capture abnormal price situations.Over the years, machine learning (ML) models have made great contributions to high accuracy forecasting of electricity load and price forecasting [3]. The nonlinear learning capabilities and fault tolerance of these models make them excellent decision-making tools for rapidly changing energy markets [4]. Lu et al. [5] combined three machine learning methods, including least absolute shrinkage and selection operator (LASSO) for features extraction, principal component analysis (PCA) for dimensionality reduction, and Bayesian ridge regression (BRR) for predictive load. Khan et al. [6] applied a hybrid model combining multilayer perceptron, support vector regression (SVR), and CatBoost. Bhinge et al. [7] developed a nonparametric model based on a Gaussian process regression (GPR) approach to predict energy consumption with a relative error of less than 6% on test data. A gradient boosting regression trees (GBRT) method was developed for 15-minute electric load forecasting in time domain industry energy [8]. Zhao et al. [9] used a kernel extreme learning machine (KELM) model to improve predictive performance. Tschora et al. [10] applied some machine learning models called ×AI to study the ability of different machine learning techniques to accurately predict electricity prices.Recently, deep learning models have received more attention due to their ability to handle energy time-series data’s nonlinear and chaotic properties. Several methods, including long short-term memory (LSTM) [11], hybrid methods based on autoencoders and variational mode decomposition (VMD) [12], reinforcement learning [13], support vector machine (SVM) [14], multiobjective optimization [15], and combining convolution, gated recurrent unit (GRU), and LSTM [16] have been applied to short-term prediction.Modeling data using attention mechanisms has achieved great success in many areas of deep learning [17], recognition of speech [18], and energy time series data [19]. Recurrent neural networks with attention mechanisms are more accurate than traditional neural networks. Combining multiple models to form these methods shows that they perform better than a single model. Decomposition methods include frequency scale resolution using wavelet transform analysis [20] and empirical mode decomposition (EMD) [21] for nonlinear and nonstationary data. The task of combining models through decomposition algorithms such as EMD is the process of finding hyperparameters and the number of recommended intrinsic mode functions (IMFs). This process involves an optimization analysis of the trade-off between accuracy and computation time. Determining spikes and peaks in the energy market is necessary for a wide range of decisions, such as risk management, transmission congestion, and outage monitoring, so we apply an ensemble empirical mode decomposition (EEMD) method.Highly unpredictable price spikes can occur due to many complex factors affecting intraday grid conditions. We present a novel method based on ensemble empirical mode decomposition and bidirectional long short-term memory with attention mechanism (BiLSTM-AM) to the electricity load and price forecasting problem. BiLSTM-AM can capture nonlinear and temporal trends in energy prices and demand. The BiLSTM-AM model is used to predict normal economic situations and record price and load spikes above a predetermined threshold. The proposed model is validated using two electricity loads and one electricity price data. The contributions of this paper are summarized as follows:(1)
A Bayesian optimization (BO) algorithm based on random forest regressors is used to tune the hyperparameters of the proposed hybrid EEMD and BiLSTM-AM model(2)
The proposed model can be used to monitor electricity load and price(3)
Provide a simple and effective method to determine the number of IMFs, avoiding over-decomposition, irrelevant components, and high calculation costsThe rest of this article is organized as follows. The literature related to the energy market forecasting models are reviewed in Section2. Section 3 describes the proposed model. Section 4 provides the analysis results. Finally, we make conclusions and further research directions.


## 2. Related Work

Price and load spikes frequently occur in the energy market that can significantly impact economic activity. Therefore, predicting normal price and load fluctuations as well as spike situations, helps gauge the risks to the economy associated with changes in energy markets. In energy markets, capturing spikes is vital as they are a core component of market attention. The BiLSTM model, coupled with an attention mechanism, allows us to model load and price spikes and provide accurate predictions. A number of methods in the recent literature have been used to address price or load spikes. Singhal and Swarup [22] applied an artificial neuron network (ANN) to capture normal, small, and large price spikes to predict market-clearing prices (MCPs) in day-ahead energy markets. Clements et al. [23] investigated the regional relationship of load spikes in Australia. Christensen et al. [24] used autoregressive conditional hazard to handle spikes by treating prices as discrete-time point processes.Several methods, such as similarity-day approach [13], SVM [14], and ANN and fuzzy logic methods [25], have progressed in load forecasting to help overcome the limitations of traditional forecasting models. Kiartzis et al. [26] proposed an expert system-based peak prediction model. Polson et al. [2] presented a combination of deep learning and extreme value theory (DL-EVT) to deal with sharp peaks, troughs, and sudden price changes in energy markets due to the fluctuations in demand and supply caused by intraday constraints. Still, this approach fails to capture the ups and downs in the levels. Radovanovic et al. [1] developed a holistic approach to recovering many energy market structures and forecasting node prices based solely on publicly available data, particularly historical prices, system loads, and a grid-wide mix of generation types. Huang et al. [21] presented an adaptive signal time-frequency processing method called EMD, which is particularly suitable for analyzing and processing nonlinear and nonstationary signals. Wang et al. [29] used the time-frequency spectrum obtained from the EEMD to overcome the mode mixing problem and more accurately reflect time series conditions. Wu et al. [30] presented an improved EEMD with LSTM for crude oil price forecasting. Zheng et al. [31] proposed empirical mode decomposition and long short-term memory (EMD-LSTM) with the Xgboost algorithm for short-term load prediction. Regarding the hybrid methods, Meng et al. [19] proposed a load prediction method based on EMD and long short-term memory with an attention mechanism. Huang [32] proposed a hybrid deep neural network model for short-term electricity price prediction. Bedi and Toshniwal [33] adopted a deep learning model based on EEMD to predict crude oil prices. Liu et al. [34] used EMD with ANN for wind speed prediction. Table 1 summarizes some previous related work on energy market forecasting.Table 1
Previous related work on electricity price and load.
Data typeReferenceModelsData sourcePrice[1]ALG−M^+ARIMA,ALG−M^+dayagoSouthwest power Pool (SPP) market[2]F-ARIMAPJM[23]SEMPPAustralian energy market operator (AEMO)[27]FFNNElectricity markets of mainland Spain and California[28]RNNPJM[32]SEPNetNew York, USA (ENGIE group database)Load[2]DL-EVTPJM[13]SD-RL-BPNNKorea[14]SVMEastern Province in Saudi Arabia[15]MOFA-GRNNAustralia[30]SD-EMD-LSTMISO New EnglandNote: F-ARIMA: Fourier-autoregressive integrated moving average; SEMPP: multivariate self-exciting point process model; FFNN: feedforward neural network; RNN: recurrent neural network; SEPNet: combing variational mode decomposition (VMD), a convolutional neural network (CNN), and gated recurrent unit (GRU); DL-EVT: deep learning with extreme value theory; SD-RL-BPNN: similar days reinforcement learning back-propagation neural network; SVM: support vector machine; MOFA-GRNN: modified generalized regression neural network with a multiobjective firefly algorithm; SD-EMD-LSTM: similar days-empirical mode decomposition-long short-term memory.


### 2.1. LSTM and BiLSTM Models

The availability of data and computing power has made deep learning an integral part of time series forecasting. In machine learning, an LSTM model is a special type of recurrent neural network (RNN) developed by Hochreiter and Schmidhuber [35]. It quickly gained popularity, particularly for solving time series prediction issues in many research areas. These networks were introduced to eliminate long-term dependency and vanishing gradient problems. In regular RNNs, the problem mainly occurs when previous information is connected with new information. As a result, LSTM units are used instead of hidden layers to better process timing-related inputs. LSTM units consist of forget gate ft, input gate it , cell state C~t, a memory unit Ct, output gate Ot, and hidden layer ht.
(1)ft=σWfhht−1+Wfxxt+bf,(2)it=i1t∗i2t,where i1t=σWihht−1+Wixxt+bi and i2t=tanhWghht−1+Wgxxt+bg.
(3)C~t=tanhWθhht−1+Wθxxt+bθ,(4)Ct=ft∗Ct−1+it∗C~t,(5)Ot=σWohht−1+Woxxt+bo,(6)ht=tanhCt∗Ot.These equations contain an offset term (b) and a weight coefficient matrix (W), a hyperbolic tangent activation function (tanh), and a sigmoid activation function (σ).The concept of BiLSTM can be defined as making a neural network that has information in both directions, backward (from future to past) and forwards (from past to future) [36]. Compared to regular LSTM where input can only flow in one direction, BiLSTMs can let input flow in both directions [37]. The bidirectional recurrent neural networks (BRNNs) concept is simple to grasp; backward LSTMs are calculated similarly to forward LSTMs, except that the time direction is reversed. The forward layer, the backward layer and the hidden layer’s final output are
(7)h⟶=fw⟶1xt+w⟶2ht−1,h⟵=fw⟵1xt+w⟵2ht−1,yi=gwo1∗h⟶+wo2∗h⟵+bo,where h⟶ is the forward output, h⟵ is the backward output, and yi is the final output of LSTM. It is the concatenation of the two hidden states h⟶andh⟵ that determines the final output yt of the BiLSTM.In this study, the BiLSTM-AM model is used to capture spike situations in the energy market, the periodicity in the data, and past and future effects on current data.


### 2.2. Decomposition Algorithms

EMD is an adaptive time-space analytic approach for nonstationary and nonlinear series processing. EMD executes operations that partition a series into IMFs without leaving the time domain. Huang et al. [21] proposed the EMD method, which divides the complex original signal into a sequence of IMFs with varying amplitudes and a residual. In this case, the series is decomposed into a finite number of IMFs. Since numerous economic and social factors cause a change in energy market price and demand, they show various characteristics of the original signal at different time scales. The term IMF refers to a monocomponent function with one instantaneous frequency that must satisfy the two conditions as follows:
(i)
The local extreme value points and zero-crossing points must be equal, or at least one difference, during the whole period(ii)
There is a zero-mean value for the envelopes defined by the upper and lower envelope at every point in the time seriesThe EMD decomposition steps are as follows:(1)
All maxima and minima in the time seriesxt are calculated using cubic spline interpolation(2)
Do a local mean value calculation for the upper and lower envelopesmt=Ut+Lt/2.(3)
mt is subtracted from the original time series xt to obtain the intermediate signal ht=xt−mt.(4)
Input seriesht are repeated until the mean of them approaches zero, and the ith IMF is computed, denoted Cit, satisfying the two conditions (i) and (ii) mentioned above. ht=Cit=IMF1; If not, consider ht as the new signal to perform (1) to (4). If hkt does not meet the definition of IMF, repeat k times, then extract IMF1: hkt=C1t=IMF1.(5)
A new signalRt=xt−IMF1can be obtained by subtracting IMF1 from original signal xt. the above steps 1 to 4 are repeated to get IMF1.(6)
Repeat steps 1-5 by usingRt as the input series until the termination condition is met, where Rt satisfies monotonicity; then, EMD decomposition of the original signal endsUsing the number ofIMFs=n, the original signal xt is reconstructed as follows:
(8)xt=∑i=1nIMFit+Rt,where Rtis the last residue.To overcome the mode mixing problem of the EMD method, Wu and Huang [38] suggested a noise-assisted EMD algorithm called ensemble empirical mode decomposition (EEMD). EEMD is an improvement on the traditional EMD algorithm. It introduces some noise into the process so that the IMFs produced on the other end should be more accurate. Unlike the standard EMD method, this method allows a better scale separation. With EEMD, white noise is added to the original signal according to the EMD algorithm. The signal is evenly distributed at the extreme point interval throughout the band, reducing mode mixing effects.EEMD algorithm [39] involves the following steps:
(1)
Add white noiseεit to the original signal xt with a standard normal distribution N0,1. Yit=xt+εit, for i=1,2,.,.,.,N, where N is the number of ensembles(2)
The EMD algorithm decomposed the transformed signalYit into a set of IMFs and a residual
(9)Ynt=∑m=1M−1IMFmit+rMit,where M−1 is the total number of IMFs produced. IMFmiis the mth IMF and rMi the residual at ith trial(3)
Steps (1) and (2) are repeated for N trials. A different white noise seriesεit is added to the original signal in each trial(4)
The last IMF of the EEMD algorithmIMFmave is computed by averaging the total m IMF in relation to N trials: IMFmavet=1/N∑i=1NIMFmit, therefore IMF1=1/N∑i=1NIMF1it and Residual r1t=xt−IMF1The results produced by EEMD depend on the number of ensemble N and the added noiseβ which should satisfy the following relation [18].
(10)α=β√N,where αis the final standard deviation of the residual.


## 2.1. LSTM and BiLSTM Models

The availability of data and computing power has made deep learning an integral part of time series forecasting. In machine learning, an LSTM model is a special type of recurrent neural network (RNN) developed by Hochreiter and Schmidhuber [35]. It quickly gained popularity, particularly for solving time series prediction issues in many research areas. These networks were introduced to eliminate long-term dependency and vanishing gradient problems. In regular RNNs, the problem mainly occurs when previous information is connected with new information. As a result, LSTM units are used instead of hidden layers to better process timing-related inputs. LSTM units consist of forget gate ft, input gate it , cell state C~t, a memory unit Ct, output gate Ot, and hidden layer ht.
(1)ft=σWfhht−1+Wfxxt+bf,(2)it=i1t∗i2t,where i1t=σWihht−1+Wixxt+bi and i2t=tanhWghht−1+Wgxxt+bg.
(3)C~t=tanhWθhht−1+Wθxxt+bθ,(4)Ct=ft∗Ct−1+it∗C~t,(5)Ot=σWohht−1+Woxxt+bo,(6)ht=tanhCt∗Ot.These equations contain an offset term (b) and a weight coefficient matrix (W), a hyperbolic tangent activation function (tanh), and a sigmoid activation function (σ).The concept of BiLSTM can be defined as making a neural network that has information in both directions, backward (from future to past) and forwards (from past to future) [36]. Compared to regular LSTM where input can only flow in one direction, BiLSTMs can let input flow in both directions [37]. The bidirectional recurrent neural networks (BRNNs) concept is simple to grasp; backward LSTMs are calculated similarly to forward LSTMs, except that the time direction is reversed. The forward layer, the backward layer and the hidden layer’s final output are
(7)h⟶=fw⟶1xt+w⟶2ht−1,h⟵=fw⟵1xt+w⟵2ht−1,yi=gwo1∗h⟶+wo2∗h⟵+bo,where h⟶ is the forward output, h⟵ is the backward output, and yi is the final output of LSTM. It is the concatenation of the two hidden states h⟶andh⟵ that determines the final output yt of the BiLSTM.In this study, the BiLSTM-AM model is used to capture spike situations in the energy market, the periodicity in the data, and past and future effects on current data.


## 2.2. Decomposition Algorithms

EMD is an adaptive time-space analytic approach for nonstationary and nonlinear series processing. EMD executes operations that partition a series into IMFs without leaving the time domain. Huang et al. [21] proposed the EMD method, which divides the complex original signal into a sequence of IMFs with varying amplitudes and a residual. In this case, the series is decomposed into a finite number of IMFs. Since numerous economic and social factors cause a change in energy market price and demand, they show various characteristics of the original signal at different time scales. The term IMF refers to a monocomponent function with one instantaneous frequency that must satisfy the two conditions as follows:
(i)
The local extreme value points and zero-crossing points must be equal, or at least one difference, during the whole period(ii)
There is a zero-mean value for the envelopes defined by the upper and lower envelope at every point in the time seriesThe EMD decomposition steps are as follows:(1)
All maxima and minima in the time seriesxt are calculated using cubic spline interpolation(2)
Do a local mean value calculation for the upper and lower envelopesmt=Ut+Lt/2.(3)
mt is subtracted from the original time series xt to obtain the intermediate signal ht=xt−mt.(4)
Input seriesht are repeated until the mean of them approaches zero, and the ith IMF is computed, denoted Cit, satisfying the two conditions (i) and (ii) mentioned above. ht=Cit=IMF1; If not, consider ht as the new signal to perform (1) to (4). If hkt does not meet the definition of IMF, repeat k times, then extract IMF1: hkt=C1t=IMF1.(5)
A new signalRt=xt−IMF1can be obtained by subtracting IMF1 from original signal xt. the above steps 1 to 4 are repeated to get IMF1.(6)
Repeat steps 1-5 by usingRt as the input series until the termination condition is met, where Rt satisfies monotonicity; then, EMD decomposition of the original signal endsUsing the number ofIMFs=n, the original signal xt is reconstructed as follows:
(8)xt=∑i=1nIMFit+Rt,where Rtis the last residue.To overcome the mode mixing problem of the EMD method, Wu and Huang [38] suggested a noise-assisted EMD algorithm called ensemble empirical mode decomposition (EEMD). EEMD is an improvement on the traditional EMD algorithm. It introduces some noise into the process so that the IMFs produced on the other end should be more accurate. Unlike the standard EMD method, this method allows a better scale separation. With EEMD, white noise is added to the original signal according to the EMD algorithm. The signal is evenly distributed at the extreme point interval throughout the band, reducing mode mixing effects.EEMD algorithm [39] involves the following steps:
(1)
Add white noiseεit to the original signal xt with a standard normal distribution N0,1. Yit=xt+εit, for i=1,2,.,.,.,N, where N is the number of ensembles(2)
The EMD algorithm decomposed the transformed signalYit into a set of IMFs and a residual
(9)Ynt=∑m=1M−1IMFmit+rMit,where M−1 is the total number of IMFs produced. IMFmiis the mth IMF and rMi the residual at ith trial(3)
Steps (1) and (2) are repeated for N trials. A different white noise seriesεit is added to the original signal in each trial(4)
The last IMF of the EEMD algorithmIMFmave is computed by averaging the total m IMF in relation to N trials: IMFmavet=1/N∑i=1NIMFmit, therefore IMF1=1/N∑i=1NIMF1it and Residual r1t=xt−IMF1The results produced by EEMD depend on the number of ensemble N and the added noiseβ which should satisfy the following relation [18].
(10)α=β√N,where αis the final standard deviation of the residual.


## 3. Proposed Method

We propose a novel deep learning framework based on integrated ensemble empirical mode decomposition and bidirectional long short-term memory with attention to forecasting energy market data. The purpose of combining the attention mechanism with BiLSTM is to capture the influential part of the input sequence during the prediction process. With the introduction of the attention mechanism, BiLSTM becomes more accurate by assigning different probability weights to inputs to highlight more important factors. The normalization term in the RNN keeps all summed inputs of the rescaling layer constant, making the hidden-to-hidden dynamics more stable [17–18].The many-to-one attention mechanism is calculated as follows:(11)Attentionweight:αfb=expScorehf,hb∑i=1nexpScorehf,hb′.Attentionvector:πt=σCtht=vtanhWht+Wht+b,where the
(12)Contextvector:C=∑i=1nαfbhi,(13)Scorehf,hb=hfTWhbLuong′smultiplicativestylevαTvtanhW1hf+W2hbBahdanan′sadditivestyle.The hidden layer state value of the final output is calculated asHt′=GC,ht,xt.As shown in Figure1, the proposed BiLSTM-AM model consists of five layers: input vector, forward LSTM hidden layer, backward LSTM hidden layer, attention layer, fully connected layer, and output layer. The input sequence data is passed to a forward LSTM hidden layer and a backward LSTM hidden layer, which combine to output a processed vector. The attention layer computes a weight vector based on the LSTM model and then combines the weight vector with the shallow output to compute the predictions of the fully connected layer. The structural architecture of the proposed BiLSTM-AM is shown in Figure 1.Figure 1
Architecture of BiLSTM with attention mechanism.We transform the data to the range between [0,1] using the min-max normalization method. Due to the noisy nature of time series data, to reduce the noise impact of the collected data for better forecasting, we normalize and then decompose using EEMD method prior to modeling. The normalization formula is(14)Xn=X−XminXmax−Xmin,where the normalized data is Xn, the overall sample minimum value is Xmin, and the overall sample maximum value is Xmax. After making predictions using EEMD-BiLSTM-AM, the predicted values from the sum of IMFs predicted values are restored to the original values using the following equation;
(15)Yt=Xt′Xmax−Xmin+Xmin,where Xt′ is the forecast value from the proposed model and Yt is the final reconstructed prediction value.Figure2 provides an illustration of the proposed method. The proposed method consists of four parts:
(1)
The data is split into training and test sets, and preprocessed the training data using Eq. (14). Then, we apply the EEMD to decompose the data(2)
Since the number of mode components is difficult to determine, much research has been done on how to determine the number of components to decompose from the original data. A simple approach based on the model evaluation method is proposed to determine the number of mode components. The iteration process applies the LSTM model to compare the root mean square error (RMSE) of the mode components (MCn) and MCn+1 by Eq. (16). When the RMSEn+1 of MCn+1 is greater than the RMSEn of MCn, the iteration process is stopped and the optimal number of IMFs is obtained
(16)RMSE=1N∑t=1Nyt−ft2,where N is the total number of predicted values, yt represents the actual value, and ft is the predicted value. Algorithm 1 gives the optimal IMF number selection process.(3)
We model each sequence by using the BiLSTM-AM model. The predicted subseries are sum and we use Equation (15) for reverse normalization to obtain the predicted value(4)
The trained BiLSTM-AM model with the sliding window approach is used to obtain the final predicted value for the test dataAlgorithm 1: The optimal IMF number selection procedure.
Input: Dataset X, m (# of modes by default)1: Divide the dataset into a training datasetXtr and a test dataset Xte based on the prediction period P.2:Fxi⟵IMFsx and rx.3: Assign: Number of modesn=1,⋯,m. Where m is the number of IMFs by default.4:for n =1,…m do5:  decomposeFx using EEMD algorithm6:  trainIMFsx and rx using LSTM7:  calculate the RMSE at each number of IMFs8:ifRMSEn+1≤RMSEnthen9:  repeat steps 5 and 610:else11:   the decomposition will stop, andn will be the optimum number of IMFs12:endif13:endforOutput: Optimum number of IMFsNote:m=12 for PJM price and load, m=13 for NSW load.Figure 2
Flow chart of the proposed method.


## 4. Datasets

We use two datasets of hourly and half-hourly electricity loads and prices for model validation. The hourly dataset is from PJM [39]. The price of electricity in the real-time (RT) market reflects the actual cost of supplying electricity to a location, given the actual constraints of the grid. PJM uses hourly RT energy market to be locational marginal pricing (LMP). Wholesale electricity prices are based on LMP, loads, generation patterns, and physical constraints on the transmission system, so they reflect the value of electricity in different regions. The raw data for daily-ahead (DA) and RT price from January 1, 2020 to June 30, 2021. Figure 3(a) shows a plot of load data in megawatt-hour (MWh). Figure 3(b) shows the RT price diagram.Figure 3
PJM electricity load and price graphs.
(a)(b)The second dataset is from the Australian Energy Market Operator (AEMO) [40], which is aggregated half-hourly demand and price data. Australia’s AEMO manages electricity and gas networks in five states. We use the load data from New South Wales (NSW) for our predictive analysis. Data range from January 2020 to August 2021. Figure 4(a) shows the load plot for NSW half-hourly data, and Figure 4(b) compares the time-series relationship between price and load for the first ten days from January 1 to January 10, 2021. Both the NSW and PJM interconnection load diagrams depict relatively homogeneous and periodic characteristics, while the price graphs show non-homogeneous and rarely periodic characteristics.Figure 4
NSW electricity load and price graphs.
(a)(b)Table2 shows the summary statistics of the three datasets. The skewness values of all variables are greater than 0, indicating a right-skewed effect.Table 2
Descriptive statistical analysis.
DataSizeMeanSDSkewnessKurtosisMinimumMaximumPJM load2059130992.726233.560.850.8619255.0056391.00PJM price1312222.6014.749.37198.660.04537.79NSW load290077773.721264.870.803.595221.1313634.63


## 5. Analysis Results

Using our proposed method in Figure2, the optimal number of IMFs for PJM load under 2, 10, and 30-day forecast periods are 3, 4, and 4, respectively (see Table 3). For the 2-day forecast of PJM load, the RMSEs of different numbers of IMFs, such as 1, 2, 3, and 4, are 323.44, 290.51, 254.01, and 267.87, respectively. It shows that the number of IMFs=3 for the PJM load provides the smallest RMSE. Therefore, the optimal number of IMFs for PJM load under 2-day prediction using the LSTM model is 3. The optimal number of IMFs for PJM price under 2, 10, and 30-day forecast periods are 5, 3, and 3, respectively.Table 3
Experimental results of LSTM model with different numbers of IMFs.
(a)PJM load2-day10-day30-dayIMF #1 +12+13+14+11+12+13+14+15+11+12+13+14+15+1RMSE323.44290.51254.01267.87323.01276.46258.27253.62264.19341.61321.42266.55258.37258.39Time10.3613.2519.1625.059.5712.4518.5624.3531.399.5212.0616.1323.5728.13(b)PJM price2-day10-day30-dayIMF #1+12 +13+14+15+16+11+12+13+14+11+12+13+14+1RMSE4.654.104.094.064.034.224.384.264.074.274.484.254.144.17Time7.4711.5415.1421.1627.3735.327.5411.4514.3220.195.579.6413.2918.11(c)NSW load2 -day10-day30-dayIMF #1+12+13+14+15+11+12+13+14+15+11+12+13+14+15+1RMSE65.9644.7844.5844.2345.8857.4744.7444.4744.2245.3659.2545.7544.0243.9644.56Time19.3129.2138.0547.0654.0718.2432.1949.2854.3559.2917.0928.5137.4548.0453.46For a fair comparison, all models are hyperparameter tuned to achieve their best predicted values. Hyperparameters of the proposed model and other models used for comparison were tuned using a Bayesian optimization (BO) algorithm based on a random forest regressor as a surrogate model. Table4 describes the optimal hyperparameters for two days predictive analysis. N+1 means the number of IMFs=N, and 1 is the residual.Table 4
Optimal hyperparameters of the proposed model for NSW 2 days prediction.
ParametersIMF1IMF2IMF3IMF4ResidualNum_layers21131Neuros[128,32]256430[200,64,16]300Dropout0.0100.100.04OptimizerAdamAdamAdamAdamAdamLearning_rate0.00010.0010.002010.0007420.005For predictive analysis, the data is first normalized by Equation (14) and then split into training and test data. A multistep prediction approach with a sliding window is used for the prediction. For example, we take a 10-day forecast: PJM price (Jan./1/2020–Jun./20/2021) as training data and (Jun./21/2021–Jun./30/2021) as test data; and NSW load (Jan./1/2020 at 0 : 30 am–20/8/2021 at 7 : 00 am) as training data and (Aug./20/2021 at 7 : 30 am–Aug./30/2021 at 23 : 30 pm) as test data. The same segmentation procedure as above was applied for other prediction periods, with the length of the prediction period as the test set and the rest as the training set.All DL models use 10% of the training set as a validation set to evaluate the model’s fitness while tuning its hyperparameters based on the training data set. The $100/MWh threshold is generally considered to be the threshold for exceeding an extreme price event in Australia, so a spike is defined asPt≥$100MWh [23].The performance of our proposed model and benchmark models is evaluated by using RMSE: mean absolute percentage error (MAPE),R2: mean absolute error (MAE) and mean square error (MSE), where MAPE, R2, MAE, and MSE are defined as
(17)MAPE=1N∑t=1Nyt−ftytx100%,R2=1−∑t=1Nyi−ft2∑t=1Nyi−y¯2,MAE=∑t=1Nyt−ytN,MSE=∑t=1Nyt−yt2N,where N is the total number of predicted values, yt represents the actual electricity price or load, ft is the predicted value,y¯ is the mean of all actual values, and R2 is the coefficient of determination. Different evaluation measures tend to produce different results in predicting effects. Therefore, we also apply the model confidence set (MCS) test [41] and Kolmogorov-Smirnov predictive accuracy (KSPA) test [42] to compare the predictive ability of each model.To validate the proposed model and EMD-BiLSTM-AM, the predicted subsequences are summed using an ensemble summation strategy, and all predicted values are renormalized to the original data scale. The spike threshold for PJM price was set to $100 per MWh, the NSW load peak was set to 11,000 MWh, and the PJM load peak was set to 39,000 MWh, as shown in Figures3 and 4.We compared the proposed model with nine other models in terms of prediction accuracy and ability to capture peaks. The comparison models include four machine learning models (support vector regression (SVR), extreme learning machine (ELM), gradient boosting regression tree (GBRT), and back-propagation neuron network (BPNN)), three deep learning models (LSTM, BiLSTM, and BiLSTM-AM), one classical autoregressive integrated moving average (ARIMA) model, and a hybrid model (EMD-BiLSTM-AM). The training lengths for the three forecast periods (2, 10, and 30 days) for PJM prices are 13,074 hours, 12,882 hours, and 12,402 hours, respectively.Python is used to implement all analyses of the model. TensorFlow is installed in the Spyder 3.7.3 environment. Bayesian optimization (BO) based on random forest regression was used to tune the hyperparameters of all other models except the ELM model with default parameters implemented using the scikit-learn neural network package and automatic ARIMA. The proposed model and other DL models used for comparison apply a “ReLU” activation function, the loss function is the mean square error (MSE), Adam as the optimizer, alookback=6, epochs=100, and batchsize=4 for the load dataset, 3 batch size for the price dataset. All experiments are performed on a desktop with processors: Intel(R) Core(TM) i7-10750H CPU @ 2.60GHz and NVIDIA GeForce GTX 1660 Ti GPU; OS: Windows x64; Memory: 16.00 GiB.Tables5 and 6 show the prediction accuracy of the eight models under different forecast periods in terms of MAPE, RMSE, and R2 for PJM and NSW loads. The proposed model has a MAPE of less than 1% in the three forecast periods. The proposed model has the lowest RMSE for all different forecast periods. Using the R2 metric, the proposed model can achieve 99.6% to 99.8% accuracy in the three forecast periods. The results show that the proposed model achieves the best results on benchmarks for all the prediction periods, followed by other deep learning models. All models can capture peaks in load data, while some cannot capture price spikes. Regardless of the time interval difference, EEMD-BiLSTM-AM accurately predicts the actual electricity load while better capturing peaks. Therefore, the proposed model can predict the electrical load at different time intervals.Table 5
Prediction results of different models for PJM load.
Models2-day10-day30-dayMAPERMSER2MAPERMSER2MAPERMSER2SVR2.3131178.1310.8801.9661086.5010.9481.732874.9270.945ARIMA4.2421427.7820.7142.3071071.0990.9491.967843.1880.947GBRT2.8291027.6200.9081.801960.1280.9591.901968.4800.959ELM5.0321570.2960.7863.8611418.4940.9113.1001140.2300.907BPNN3.1711241.5580.8662.4321099.2500.9472.207913.2890.940LSTM1.9061062.3020.9021.860972.8640.9591.714785.7090.956BiLSTM1.547959.5960.9201.492903.3850.9651.699725.9100.963BiLSTM-AM2.2651144.2770.8961.693913.0580.9641.621742.5780.961EMD-BiLSTM-AM1.136548.8600.9741.666620.9150.9831.290484.7150.981EEMD-BiLSTM-AM0.775232.9890.9960.729258.4730.9970.636228.9450.996Table 6
Prediction results of different models for NSW load.
Models2-day10-day30-dayMAPERMSER2MAPERMSER2MAPERMSER2SVR1.438136.0660.9831.317150.1740.9761.476141.7560.985ARIMA1.386142.6450.9821.305149.5700.9801.326153.0250.982GBRT1.189123.3850.9871.287136.9880.9811.155130.5630.987ELM2.208217.0350.9611.489162.0970.9691.937214.6870.966BPNN1.541155.5480.9801.579174.9860.9601.475159.3550.981LSTM1.356135.8170.9851.276143.2690.9861.287148.4250.984BiLSTM1.279132.4970.9861.248141.1620.9871.377147.7410.984BiLSTM-AM1.108133.7660.9891.312144.340.9861.368148.160.984EMD-BiLSTM-AM0.73168.7820.9940.80593.9980.9941.088139.4210.986EEMD-BiLSTM-AM0.48048.2160.9980.60262.4640.9970.55356.1810.998On the MAPE metric for NSW load data, applying the EEMD method to BiLSTM-AM can improve the prediction results by 53%, 54%, and 60% for 2-day, 10-day and 30-day forecasts, compared to using BiLSTM. Figure5 shows the 2-day test set load forecast for NSW with good accuracy for the spike and peak capture levels for the proposed model. Machine learning models and classical statistical model ARIMA perform poorly in valley-peak prediction as compared to the hybrid models. SVR and ELM also miss the peak time predictions. This suggests that hybrid deep learning models combined with decomposition methods can help improve predictive performance, especially on high-frequency data, as shown in Figure 5 for the 30-minute interval data for NSW.Figure 5
Prediction results for 2 days for NSW electricity load.Table7 shows each model’s peak time and amplitude analysis in Figure 5(b). The table shows that the proposed model has the closest amplitude to the actual peak and captures the exact time of the peak. All other models captured the exact time index of the peak except SVR, ARIMA, and ELM but failed to capture the amplitude. The proposed model underestimates the load peak by only 1.98 MW, which shows the prediction power of the proposed model in capturing the peak and spike.Table 7
Peak points and time difference analysis of 2-day prediction for NSW load.
ModelsActual peak loadPredicted peakLoad differenceActual peak timePredicted peak timeTime differenceSVR10,151.6710,529.31-377.648/29/21 19 : 008/29/21 19 : 3030 minARIMA10,151.6710,447.81-296.1398/29/21 19 : 008/29/21 20 : 0060 minGBRT10,151.6710190.85-39.17788/29/21 19 : 008/29/21 19 : 00—ELM10,151.6710,448.04-296.378/29/21 19 : 008/29/21 19 : 3030 minBP10,151.6710,501.46-349.798/29/21 19 : 008/29/21 19 : 00—LSTM10,151.6710,393.57-241.908/29/21 19 : 008/29/21 19 : 00—BiLSTM10,151.6710,370.99-219.328/29/21 19 : 008/29/21 19 : 00—BiLSTM-AM10,151.6710,341.87-190.208/29/21 19 : 008/29/21 19 : 00—EMD-BiLSTM-AM10,151.6710,071.3780.308/29/21 19 : 008/29/21 19 : 00—EEMD-BiLSTM-AM10,151.6710,149.691.988/29/21 19 : 008/29/21 19 : 00—As shown in Figure6, the spike threshold level is set at $100/MWh, and the 10-day forecast from June 21st-30th, 2021 shows a one-time spike situation as on June 29th, 2021. ELM, EMD-BiLSTM-AM, SVR, and the proposed model can capture the spikes. LSTM also captures the spikes, but not in amplitude compared to others, while the rest fail to capture the spikes. Compared with the other models in Table 8, the proposed model outperforms the comparative models in all different forecast periods and better captures the peak time and amplitude, as shown in Figure 6. Although the EMD-BiLSTM-AM model captures the same peak level as the proposed model, it still has lower prediction accuracy than the proposed model.Figure 6
Prediction for PJM electricity price.Table 8
The results of different models for PJM price under different forecast periods.
Models2-day10-day30-dayMAPERMSER2MAPERMSER2MAPERMSER2SVR0.10210.9360.5720.11711.5590.6000.14311.7930.451ARIMA0.22511.6980.4570.19012.1270.5600.22014.5380.330GBRT0.20317.8020.2410.12112.0020.5680.13211.3680.489ELM0.18011.8990.4940.17211.9370.5740.20311.6610.463BPNN0.12010.0060.6420.12410.6700.6360.14910.5570.560LSTM0.15411.2800.5450.16410.9300.6430.14510.4190.571BiLSTM0.13211.3930.5360.10810.9930.6390.14710.6680.551BiLSTM-AM0.13912.0820.4780.11810.6200.6630.15610.4060.572EMD-BiLSTM-AM0.1787.8590.7790.1226.8530.8890.1297.2030.842EEMD-BiLSTM-AM0.0895.1340.9060.1014.8900.9290.1004.3920.924On average, the MAPE andR2 of the proposed model for the three test periods were 0.097 and 0.920, respectively. EEMD-BiLSTM-AM outperforms other models in all performance metrics used. We can see in Table 8 that using EEMD-BiLSTM-AM improves the MAPE of the prediction results by 50%, 17%, and 37% compared to EMD with BiLSTM-AM in 2-day, 10-day, and 30-day predictions, respectively.From the above analysis, we conclude that the proposed method can predict data with different fluctuation frequencies, time intervals, and range, i.e. load and price data types.Figures5 and 6 compare the prediction accuracy of nine benchmark models and our proposed model on different data types and different time horizons. EEMD-BiLSTM-AM shows higher forecasting performance and less peak-to-valley amplitude.Table9 provides the computation times in Python, including hyperparameter tuning, model training, and prediction. Hybrid deep learning models take longer to compute than deep learning and machine learning models. However, the hybrid deep learning models outperform other models in terms of prediction performance and valley-peak capture.Table 9
Model computation time in seconds.
ModelsPJM priceNSW load2 days10 days30 days2 days10 daysSVR8.0169.84511.92230.70131.908ARIMA76.13570.54268.287260.041215.041GBRT24.75025.00024.62522.8281323.4375ELM4.3133.9537.03112.94514.870BPNN9.70113.65914.09135.9806737.89465LSTM479.035505.580593.8041304.9211263.869BiLSTM574.893608.733693.3002594.8162232.46BiLSTM-AM406.360451.976527.8751984.8281864.757EMD-BiLSTM-AM2495.8092709.8362907.4486272.2675803.947EEMD-BiLSTM-AM2593.9082860.8993009.2946666.8256298.189The MCS test is used to select the best model with a given level of confidence [42]. The KSPA test is a nonparametric test based on the Kolmogorov-Smirnov (KS) test [43]. This test was designed to determine whether the forecast error distributions were statistically significantly different and to distinguish between more predictive and less predictive models.We provide statistical test analysis results of 2 days NSW load and 30 days PJM price forecast data. The KSPA error distribution and empirical cumulative distribution function are shown in Figures7 and 8 and, respectively. The two-sided (p value) of the KSPA test for all compared models was less than 0.01 for the 2 days NSW load data. This confirms the statistical difference between the proposed and competing models and shows that the proposed model outperforms the others. One-sided KSPA tests are also performed to show the low errors of the proposed model and its maximum prediction accuracy. The calculated statistical errors of the predictions show that the proposed model outperforms the comparison models. The random error deviation is better explained by the proposed EEMD-BiLSTM-AM model, which shows the smallest error and the best prediction performance. Our proposed model better captures random bias with minimum error. The p value results of MCS test are shown in Table 10. We set a threshold of 0.25 for the p value to indicate that the model survived under the loss function, which is underlined. A larger p value means a more accurate model. According to the above analysis, the prediction ability of the EEMD-BiLSTM-AM model is the best.Figure 7
Absolute error distribution plots for 2 days NSW forecast data.
(a)
Statistical and machine learning models(b)
Deep learning and hybrid deep learning modelsFigure 8
Empirical cumulative distribution function of error plots.
(a)
2 days NSW load forecast data(b)
30 days PJM price forecast dataTable 10
P-values of MCS test using different models for NSW load and PJM price.
Models2 days NSW load10 days PJM priceMAEMAPERMSEMSEMAEMAPERMSEMSESVR0.00000.00000.00000.00000.00000.51280.28340.0000ELM0.00000.00000.00000.00000.00000.00000.00000.0000GBRT0.00000.00000.00000.00000.00000.00000.00000.0000BPNN0.00000.00000.00000.00000.00000.00000.00000.0000ARIMA0.00000.00000.00000.00000.00000.00000.00000.0000LSTM0.00000.00000.00000.00000.00000.00000.00000.0000BiLSTM0.00000.00000.00000.00000.00000.00000.00000.0000BiLSTM-AM0.00000.00000.00000.00000.00000.00000.00000.0000EMD-BiLSTM-AM0.00001.00000.00000.00000.00001.00000.00000.0000EEMD-BiLSTM-AM1.00001.00001.00001.00001.0001.00001.00001.0000Note: A higherp value means a more accurate model.


## 6. Conclusions

Energy market participants have always been nervous about the economic consequences of soaring electricity prices. Many forecasting models in the energy market focus on predicting normal price and load movements without dealing with extreme load or price situations in the market. In this paper, we combine EEMD and BiLSTM-AM to predict electricity loads and prices while capturing spikes and peak points. The EEMD algorithm decomposes the data to reduce the effect of noise. We then applied the BiLSTM-AM model to obtain the occurrence and normal fluctuation values of price and load peaks from the decomposed data. Using three data with different time horizons and different volatility patterns, we examine the predictive performance, efficiency, and consistency of the proposed model by forecasting the short to medium term. The proposed method is compared with night benchmark methods, including four machine learning models (SVR, ELM, GBRT, and BPNN), three deep learning models (LSTM, BiLSTM, and BiLSTM-AM), a classical ARIMA model, and a hybrid method (EMD-BiLSTM-AM). The results show that the proposed method outperforms other models in terms of prediction accuracy and spike capture ability, with EEMD reducing the mean absolute percentage error (MAPE) by 53%, 54%, and 60%, respectively. Over the three forecast periods, the average MAPE andR2 are 0.097 and 0.92, respectively. Furthermore, the analytical results of MCS and KSPA tests show that the EEMD-BiLSTM-AM model has the best predictive ability.The analysis results show that the price data exhibits large extreme volatility and uncertainty in peak behavior. Its heterogeneity and many seasonality make modeling more challenging compared to loads.



---

*Source: 3815063-2023-02-03.xml*