Prediction of Extreme Conditional Quantiles of Electricity Demand: An Application Using South African Data

Abstract

It is important to predict extreme electricity demand in power utilities as the uncertainties in the future of electricity demand distribution have to be taken into consideration to achieve the desired goals. The study focused on the prediction of extremely high conditional quantiles (between 0.95 and 0.9999) and extremely low quantiles (between 0.001 and 0.05) of electricity demand using South African data. The paper discusses a comparative analysis of the additive quantile regression model with an extremal mixture model and a nonlinear quantile regression model. The estimated quantiles at each level were then combined using the median approach. The comparisons were carried out using daily peak electricity demand data ranging from January 1997 to May 2014. Proper scoring rules were used to compare the three models, and the model with the smallest score was preferred. The results could be useful to system operators including decision-makers in power utility companies by giving insights and guidance for future electricity demand patterns. The prediction of extremely high quantiles of daily peak electricity demand could help system operators know the possible largest demand that will enable them to supply adequate electricity to consumers and shift demand to off-peak periods. The prediction of extreme conditional quantiles of daily peak electricity demand in the context of South Africa using additive quantile regression, nonlinear quantile regression, and extremal mixture models has not been performed previously to the best of our knowledge.

1. Introduction

In an attempt to diversify the energy mix, the South African government has developed a Renewable Energy Independent Power Producer Procurement Program (REIPPPP). The REIPPPP has proven to be a very successful program especially in bringing renewable energy projects to commercial operations over the past five years [1]. The country also came up with a new plan, namely the National Development Plan (NDP). The NDP is a plan for infrastructure development from 2013 to 2030. Moreover, the NDP is also an excellent guiding economic plan that set the GDP growth target per annum for the country to be able to meet its economic, social and political objectives [2].

In identifying the national goals relevant to establishing renewable energy policy objectives, South Africa needs to identify the key goals for the nation and how the electricity sector fits in among its priorities. The country also needs to take urgent actions to ensure the sustainability of renewable energy and energy efficiency by 2030. According to [3], renewable and energy efficiency have a positive impact on electricity demand during peak hours.

It is important to accurately predict extremely high quantiles of electricity demand. Uncertainties related to electricity demand have to be taken into account when predicting electricity demand. The superiority of extreme conditional quantile models over the least squares (conditional mean) model in the extremes of the conditional distribution is well established in the literature. Hence, accurate prediction of electricity demand at extreme levels could produce more useful information to decision-makers on the sustainability of renewable energy, energy generation and energy purchases. Furthermore, the accurate prediction of extreme electricity demand distributions would have a significant impact on mitigating load shedding and overloading and allow energy-efficient storage.

The main contribution of this paper is the use of additive quantile regression (AQR), extremal mixture (EM) and nonlinear quantile regression (NLQR) models in estimating the extremely high and extremely low quantiles of electricity demand using South African data. Such a study has not been carried out elsewhere to the best of our knowledge.

The highlights of the study are summarised as follows:

• The study carried out a comparative analysis of EM, AQR and NLQR models in predicting extremely high and low daily peak electricity demand;
• The identification of how electricity demand will change in the distribution networks in five to fifteen years going forward;
• The prediction of extremely high quantiles of DPED could help system operators know the possible largest demand that will enable them to supply adequate electricity to consumers;
• Knowing the possible largest demand of electricity at a given point in time can help system operators shift demand to off-peak periods.

The rest of the paper is structured as follows: Section 2 reviews several approaches proposed in the literature with some previous studies. Section 3 focuses on the semiparametric extremal mixture models, additive quantile regression model, nonlinear quantile regression model, and a brief explanation of the scoring rules. Section 4 presents the data used in the study including the exploratory data analysis. Section 5 gives a discussion of the results. Finally, the study is then concluded in Section 6.

2. Literature Review

Many statistical applications focus on either the lower or upper quantiles of the distributions. Consequently, the theory of extreme value techniques has been extensively used in many disciplines for the past several decades. The extremal quantile regression (EQR) model can also be utilised to model extreme values as an alternative to classical extreme value theory (EVT) techniques. EVT is widely used as a tool for risk management. Several estimators have been proposed for the estimation of extreme conditional quantiles, which include [4,5,6,7,8], among others. The literature has extensively dealt with EQR, which consists of one-stage and two-stage approaches. A study by [4] used a two-step procedure based on two case studies in estimating extreme conditional quantiles. The procedure was evaluated using small sample simulation for both heavy-tailed and right-bounded distributions. The limitation of quantile regression (QR) in estimating extreme conditional quantiles is that a one-step extreme conditional quantile procedure based on QR underestimates these conditional quantiles. The introduction of a two-step extreme conditional quantile procedure suggested a way of overcoming this limitation [4]. However, the results showed that the two-step procedure did not prove to be useful for practical purposes [4], and hence, it was not considered in this paper.

Although many approaches have been discussed in the literature on the estimation of extreme quantiles, which includes those of [9,10,11,12], among others, the comparison of the performances of additive quantile regression (AQR), extremal mixture (EM) and nonlinear quantile regression (NLQR) models at extremely high and low quantiles using South African electricity demand data has not been performed to the best of our knowledge. The present study took an alternative approach from that of [9,10,11,12,13] in that it compared extreme conditional quantiles between 0.95 and 0.9999.

A study by [14] built on EVT to develop EQR models. It studied the asymptotic theory of extremal quantile regression. In this context, the study obtained the large sample properties of extremal (extreme order and intermediate order) quantile regression for the class of linear quantile regression models with conditional tails of the response variable. This was restricted to the domain of minimum attraction and closed under the tail equivalence across conditioning values. Based on the results from [14], the present paper built on this procedure to estimate extremal conditional quantiles. This was done using the scoring rules for evaluating the predictive accuracy of the developed models. Scoring rules were considered in this study, as they are significantly important in: (1) parameter estimation, (2) evaluating the predictive performance of extreme conditional quantiles and AQR models; and (3) obtaining the probabilities of rare events.

The use of scoring rules in forecasting evaluation with an emphasis on extreme events has been proposed in the recent literature. The paper by [15] compared the predictive performance of the forecasting models using the continuous ranked probability score (CRPS) and logarithmic score (LogS). The LogS was computed under the assumption of a normal predictive distribution with mean and variance considered in the experiment, while CRPS was computed using numerical methods. For evaluating probability forecasts with scoring rules, [16] emphasised the scoring rules as a comparative forecast evaluation of probability models that cover a wide range of situations in statistical applied works. The authors used several examples including economics examples, where the predictive distribution was given as a simulated sample to compute the CRPS and LogS by applying Bayesian forecasts of the United States (U.S.) gross domestic product (GDP) growth rates. However, the authors did not include the Dawid–Sebastiani score (DSS) in their comparative analysis.

To treat the asymptotic theory of extreme conditional quantile estimators, models such as parametric, semiparametric and nonparametric quantile regression can be used. For instance, [17] proposed kernel smoothing for extremal quantile regression. The paper used nonparametric regression quantiles obtained by inverting a kernel estimator of the conditional distribution. The main idea of the paper was to extend the asymptotics of the extreme conditional quantile estimator in a nonparametric regression model. Nonparametric models have also received more considerable attention in papers such as [4] and [18]. The extreme conditional quantiles estimator’s approach is one way of modelling extremes by fitting the generalised extreme value distribution (GEVD) or generalised Pareto distribution (GPD) where the location ($\mu$), shape ($\xi$) and scale ($\sigma$) parameters depend on either parametric or nonparametric covariates. According to [7], this approach captures the covariate effects at different tails of the response distribution. However, the limitation of extreme conditional quantile estimators is instability with heavy-tailed distributions due to scarce data in the tails of the distributions [7]. Researchers, [19,20] have developed semiparametric and parametric models by establishing the asymptotic normality of a simple estimation procedure, which combines QRs in the tails or EQRs with the minimum distance, as well as reviewing the technique and proposing some extensions based on the point process of high-level exceedances, respectively. The approach suggested by [21] uses a stagewise pack of a fit testing procedure to propose the excess distribution function. They introduced the procedure for selecting the sufficiently high threshold, which is the $\tau$ quantile. Their simulation results showed that the proposed adaptive procedure nearly captures the best choice for testing Pareto degrees of freedom. It was also shown that the choice based on the detection of the lack of fit point introduces a significant bias.

The study undertaken by [22] considered modelling the extremal behaviour of stock market data to compare the tail distribution of the generalised logistic distribution (GLD), GEVD and GPD. The empirical results from the study showed that the three distributions were asymptotically equivalent in their tails. Hence, the GLD plays the central role as it possesses important characteristics in classical extreme value theory for distribution properties, stability and convergence. The study by [18] proposed a standard Pareto model using a simulation study. It used the mean-squared relative error (MSRE) of adaptive quantile estimators in 0.99, 0.999 and 0.9999 quantiles. The findings suggested that the developed model provides excellent approximation in the Fréchet domain of attraction. This current paper differs from the previous one by [18], as it used three scoring rules in evaluating the extremely low and high quantiles of daily peak electricity demand (DPED). The GPD is the distribution that can be used to model the under demand and over demand prediction of peak electricity demand providing a basis for risk assessment and quantification with forecasting uncertainty [23]. However, for the GPD, only the observations above the threshold are considered, and therefore, a large amount of information is not wasted.

This present study was motivated by [24], who used partially linear additive quantile regression with an application to the unit commitment problem to forecast short-term electricity demand during peak hours. The South African data from January 2009 to June 2012 were considered for their analysis. Their results showed that the hourly peak load demand occurred between 18:00 and 20:00. The study revealed that electricity demand forecasting during this period is important to system operators as they have to use power plants to balance the high demand with what can be supplied. Even though the study is useful to system operators in power utility companies in South Africa, the study did not consider extreme peak load forecasting. The extreme peaks could be useful in assessing the power reliability in South Africa, as indicated by [25]. The study by [26] used the peak-over-threshold (POT) model for modelling extreme peak electricity demand during a heatwave. The authors used South African data from 2000 to 2013 for their analysis. Their results showed that electricity demand increases as a result of heat during days of extremely high temperature, especially during October, November and February. The proposed model in their study could help the South African government know how much electricity is needed during a heatwave period. However, the study failed to extend the analysis of extreme peak electricity demand to seasons, for example spring and autumn. In South Africa, the impact of the COVID-19 pandemic is still emerging. However, it seems that the energy sector will also continue to be significantly impacted. The following methods were used to analyse the impacts of the COVID-19 pandemic on electricity demand forecasting: regressive and neural network models [27], hybrid multi-objective optimizer-based models [28], bidirectional long short-term memory (Bi-LSTM) models [29], autoregressive integrated moving average with exogenous (ARIMAX) and artificial neural network (ANN) models [30], among others.

The different methods applied have limitations as they failed to prove that they are useful for practical purposes especially when modelling extreme peaks in electricity demand. This study was based on three techniques, namely AQR, EM and NLQR, which were compared using the historical electricity demand data. These models are very important in predicting the extremely high and low peaks of electricity demand as they yield accurate predictions. The accurate prediction of the extreme electricity demand distribution is important to decision-makers in the electricity sector and should be monitored regularly.

A summary of some previous studies on the modelling and forecasting of extremely high quantiles is given in

Table 1. Summary of previous studies on the modelling and estimation of extremely high quantiles.

Authors	Data	Models	Main Findings
Muller et al. [31]	Marseilles hourly rainfall data from 1882–2003	GPD and SHYPRE hourly rainfall stochastic models.	Results show that both methods give similar results and have similar uncertainties. Moreover, the sensitivity of the GPD to the shape parameter is quite high.
Gardes and Girard [32]	French hourly rainfall data from 1993–2000	Nearest neighbour model.	The results show that the nearest neighbour Hill estimator gives the same weight to all the largest observations.
Sigauke et al. [9]	Eskom aggregated DPD data from 2000–2011	GSP distribution and GPD models.	The Q-Q plot of the GSP distribution incorporates most extreme observations in the tail slightly better than the GPD.
Cai and Reeve [33]	Venice sea-level data from 1931–1981	Semiparametric, QR, and parametric quantile function models.	The performances of the parametric and the semiparametric approaches are very similar at the lower quantile levels. However, the performance of a quantile function modelling approach may vary from dataset to dataset at high quantile levels.
Chavez-Demoulin et al. [34]	USB data from 27 June 2002 to 18 May 2010	NPOT and classical POT models.	The results of NPOT confirmed a rather precise and adapted estimation of high quantile-based risk measures for financial time series.
Diriba et al. [10]	Port Elizabeth weather station data from 1949–2013	GPD model.	The GPD model for the minimum daily winter temperature shows no improvement in the parameter estimates’ precision.
Gijbels et al. [35]	Hurricane data from 1971–2017	Semiparametric and nonparametric models.	The results show that the semiparametric model provides the smallest estimated prediction error compared to the nonparametric model.
Taylor [36]	Hourly Nord Pool market prices data from 2013–2018	AR-GJR-GARCH and AR models.	The results show that the AR-GJR-GARCH model performs better than the AR model for both wider and narrower quantile intervals.

.

3. Methodology

An overview of the semiparametric extremal mixture, AQR and NLQR models with their implementations based on extreme quantiles is discussed in this section. The additive quantile regression models discussed here are becoming increasingly popular in many applications as they are known to be robust and flexible.

3.1. Semi-Parametric Extremal Mixture Models

Let $X_{t_1},\dots ,X_{t_n}$ denote the DPED, where $t_i,i=1,\dots ,n$ is a sequence of times ($0\le t_1\le \dots \le t_n\le T_{\max }$). Suppose the random variable $X_{t_i}$ has a distribution function $F_{t_i}$. We seek to estimate extreme quantiles, i.e., $F_t^{-1}\left(\tau \right)$ for $0.950\le \tau \le 0.9999$. Now, if $F_t$ is in the domain of attraction of the Fréchet distribution, then the excess distribution function given in Equation (1) can be estimated by a Pareto distribution in Equation (2) [37]:

$$F_{t,\tau }\left(x\right)=1-\frac{1-F_t\left(x\right)}{1-F_t\left(\tau \right)},x\in [\tau ,\infty )$$

(1)

$$G_{\tau ,\theta }\left(x\right)=1-{\left(\frac{x}{\tau }\right)}^{-\frac{1}{\theta }},x\in [\tau ,\infty ),$$

(2)

where $\theta >0$ and $\tau \ge x_0$ ($\tau$ is the unknown threshold). Consider a semiparametric mixture model (bulk model and tail model) given by [18]:

$$F_{t,\tau ,\theta }\left(x\right)=\left\{\begin{array}{cc}{F}_t\left(x\right)\hfill & \mathrm{if}x\in [x_0,\tau ]\hfill \\ 1-(1-F_t\left(\tau \right))(1-G_{\tau ,\theta }\left(x\right))\hfill & \mathrm{if}x>\tau ,\hfill \end{array}\right.$$

(3)

where $\tau \ge x_0$ represents the threshold. For any $p\in (0,1)$, the extreme quantile of $X_t$ is given by:

$$q_p(t,h)=\left\{\begin{array}{cc}{F}_{t,h}^{-1}\left(p\right)\hfill & \mathrm{if}p<p_{\tau }\hfill \\ \tau {\left(\frac{1-p_{\tau }}{1-p}\right)}^{{\theta }_{t,h,\tau }}\hfill & \mathrm{otherwise},\hfill \end{array}\right.$$

(4)

where:

$$F_{t,h}\left(x\right)=\frac{1}{{\sum }_{i=1}^{n}K\left(\frac{t_i-t}{h}\right)}\times \sum _{i=1}^{n}K\left(\frac{t_i-t}{h}\right){\mathbf{I}}_{X_{t_i}\le x},$$

(5)

where $K(.)$ is a kernel function, h is the bandwidth and $\mathbf{I}$ is an indicator function.

Threshold Selection

In peaks-over-threshold models, the threshold is normally estimated first before fitting the desired model to the exceedances [38]. A sufficiently high threshold is vital to guarantee the stability of the parameters. If the threshold is incorrectly chosen at some value larger than $\tau$, the number of observations on which the distribution is fit becomes smaller, which leads to unstable parameter estimates. However, if the threshold is too high, it produces fewer excesses to estimate the scale and shape parameters, resulting in a higher variance [39]. Classically, a quantile-based approach to find an appropriate threshold is used. The study provides a threshold selection method given in Equation (3) that effectively detects whether a sample follows a certain distribution F. However, it is noted that other threshold selection methods might be equally valid, depending on the circumstances.

3.2. Additive Quantile Regression Model

The study used the quantile generalised additive model (quantGAM) based on the work of [40] and extended by [41], defined as:

$$y_{t,\tau }=\sum _{j=1}^p{s}_{j,\tau }\left(x_{tj}\right)+{\epsilon }_{t,\tau };\tau \in (0,1).$$

(6)

The smoothing function, s, is written as:

$$s_{j\tau }\left(x\right)=\sum _{k=1}^q{\beta }_{kj}{b}_{kj}\left(x_{tj}\right),$$

(7)

where ${\beta }_{kj}$ denotes the jth parameter and $b_{kj}\left(x\right)$ represents the jth basis function with the dimension of the basis being denoted by q. The parameter estimates of Equation (6) are obtained by minimising the function given in Equation (8) as:

$$q_{Y|X}\left(\tau \right)=\sum _{t=1}^n{\rho }_{\tau }\left(y_{t,\tau }-\sum _{j=1}^p{s}_{j,\tau }\left(x_{tj}\right)\right),$$

(8)

where $q_{Y|X}\left(\tau \right)$ is the extreme conditional quantile function of $\tau$ and ${\rho }_{\tau }\left(u\right)=u[\tau -\mathbf{I}(u<0)]$ is a check function. In this study, we are interested in estimating extreme conditional quantiles, i.e., for $\tau \in (0,1)$.

3.3. Nonlinear Quantile Regression

NLQR is an extension of linear quantile regression in which the models are nonlinear in their parameters, whereas the linear-in-parameters quantile regression model is given by:

$$q_{Y|X}\left(\tau \right)=X^T\beta \left(\tau \right),$$

(9)

The NLQR model (nonlinear in parameters) is given in Equation (10) and was discussed in detail in [42]:

$$q_{Y|X}\left(\tau \right)=g(X,{\beta }_0\left(\tau \right))$$

(10)

and the corresponding estimator is given by:

$$\widehat{\beta }\left(\tau \right)=\underset{b\in B}{argmin}\sum _{i=1}^n{\rho }_{\tau }\left(y_i-g(x_i,b)\right),$$

(11)

where $g(x_i,b)$ is a function with unknown parameters, $B\in {\mathbb{R}}^p$ and $\widehat{\beta }\left(\tau \right)$ is the unknown regression coefficient for the $\tau$th quantile. The NLQR model is formed by replacing the linear quantile regression model in Equation (8) with the quantile curve in Equation (10). The study by [43] proposed the interior point algorithm approach for computing NLQR estimates. The proposed estimation is used to solve Equation (11).

3.4. Combination of Estimated Extreme Quantiles

The median (Md) method is fairly easy to use and is not sensitive to outliers [44]. The Md uses the position rather than the specific value of each data entry. If the extreme values of a dataset change, the Md usually does not change and is not sensitive to outliers. Moreover, the Md does not change because it is only dependent on the middle observation’s value. It is given as:

$$L_{Md}=median(L_1,...,L_K)U_{Md}=median(U_1,...,U_K),$$

(12)

where $L_{Md}$ and $U_{Md}$ are the lower and upper median values; $L_1,\dots ,L_K$ and $U_1,\dots ,U_K$ are the sequences of the lower and upper values, respectively.

Table 2. Comparison of the models.

Models	Strengths	Weaknesses
M1 (AQR)	1. A hybrid model that combines GAMS with QR. 2. Estimation is distribution free. 3. Robust to outliers in the response variable.	1. Requires a smoothing function of the covariates. 2. Parameters are harder to estimate. 3. Does not give any details about the size of the high level of possible exceedances.
M2 (EM)	1. Semiparametric extremal mixture model. 2. Based on one covariate, which is $t=1,\dots ,n$.	1. Has limitations on accuracy and stability. 2. Very sensitive to numbers and the location of the measured points.
M3 (NLQR)	1. Inference is performed based on large sample approximation. 2. Robust to outliers in the response variable.	1. Requires a smoothing parameter. 2. Outliers only have an influence on quantile curves close to them, i.e., they affect extreme quantiles

presents a summary of the strengths and weaknesses of the proposed, additive quantile regression, extremal mixture and nonlinear quantile regression models.

3.5. Scoring Rules for Quantiles

The scoring rules were used to evaluate the probability forecasts by assigning a real number to the predictive distribution and observations. The score values were used to compare the models. The scoring rule $S(y,p)$ is defined as a special case of a loss function $L(y_t,q_{\tau })$, measuring the negative worth of behaviour when the variable Y turns out to be $y_t$, where $q_{\tau }$ and $y_t$ denote the quantile forecast and actual value of the DPED [45]. Moreover, the proper scoring rules are preferable if the calibration and sharpness are measured simultaneously. This study used the scoring rules designed for the comparative forecasting evaluation, namely: continuous ranked probability score (CRPS), logarithmic score (LogS), Dawid–Sebastiani score (DSS), pinball loss (PL) and interval width (IW).

3.5.1. Continuous Ranked Probability Score

The CRPS has recently attracted much attention in forecasting performance as it quantifies and takes into consideration both sharpness and calibration [46]. The CRPS is defined by:

$$\mathrm{CRPS}(F,y)={\int }_{-\infty }^{\infty }{\{F\left(x\right)-\mathbf{I}(y\le x)\}}^{2}dx,$$

(13)

where y is the actual observation, $\mathbf{I}(.)$ is an indicator function equal to one for $y\le x$ and $F\left(x\right)$ is the cumulative distribution function (CDF).

3.5.2. Logarithmic Score

The LogS is the only proper local score that assumes the regularity conditions of $f(y,\theta )$ for almost all the actual observations, where $\theta$ is a parameter [15]. The LogS is given by:

$$\mathrm{LogS}(F,y)=-\mathrm{Log}\left(f\right(y\left)\right),$$

(14)

where $f\left(y\right)$ is the probability density function (PDF) and F is a strictly proper scoring rule relative to the probability distribution.

3.5.3. Dawid–Sebastiani Score

The DSS depends on the first and second moments of the forecast. The DSS is given by [46]:

$$\mathrm{DSS}(F,y)=\frac{{(y-\mu F)}^2}{{\sigma }^{2}F}+2\mathrm{Log}\left(\sigma F\right),$$

(15)

where F is the predictive distribution of y with first and second moments $\mu F$ and ${\sigma }^{2}F$, respectively, ${\left(y-\mu F\right)}^2$ is a squared error score and $\mu F$ and ${\sigma }^{2}F$ denote the mean and variance of the predictive distribution F, respectively.

3.5.4. Pinball Loss Function

The pinball loss (PL) function is relatively easy to use and is given as:

$$PL\left(q_{\tau ,t}\right)=\left\{\begin{array}{cc}2(1-\tau )|y_t-q_{\tau ,t}|,\hfill & \mathrm{if}{y}_t<q_{\tau ,t},\hfill \\ 2\tau |y_t-q_{\tau ,t}|,\hfill & \mathrm{if}{y}_t\ge q_{\tau ,t},\hfill \end{array}\right.$$

(16)

where $q_{\tau ,t}$ is the quantile forecast at time t and $y_t$ is the observed value of the DPED at time t. The interpretation of $PL\left(q_{\tau ,t}\right)$ is made easier by the inclusion of the multiplier number 2 [47]. When $\tau =0.5, {PL}_{0.5,t}=|y_t-q_{\tau ,t}|$, which is the same as the absolute error. Hence, $PL{(}_{\tau ,t})$ is generally interpreted as an absolute error.

3.5.5. Estimated Intervals’ Widths

The interval width (IW) is the difference between the estimated upper and lower quantile values. It is given in Equation (17) as:

$${\mathrm{IW}}_t=q_{\tau ,t}-q_{(1-\tau ,t)},$$

(17)

where $q_{\tau ,t}$ and $q_{(1-\tau ,t)}$ are the upper and lower quantiles, respectively. The coverage probability (CP) is the probability procedure for constructing random regions to produce an interval covering the true value. The CP is used to evaluate the performance of estimated intervals. It is considered as a property of the interval-producing procedure, which is independent of the particular sample to which such a procedure is applied. It also evaluates the reliability of the estimated interval widths and is given in Equation (18) by:

$$\mathrm{CP}=\frac{1}{h}\sum _{t=1}^h{I}_t,I_t=\left\{\begin{array}{c}1if{y}_tϵ(L_t,U_t)\hfill \\ 0\mathrm{otherwise}\hfill \end{array}\right.,$$

(18)

where $L_t$ and $U_t$ are the lower and upper specification limits, respectively.

4. Empirical Results

4.1. Data and Software

The data used in this paper were from Eskom, South Africa’s power utility company. The data considered were Net Energy Sent Out (NESO). The NESO, measured in megawatts, is the rate at which electrical energy is delivered to customers. The NESO considered in this paper is the daily peak electricity demand (DPED), which is the maximum hourly demand in a 24 h period. Aggregated DPED data are used for the agricultural, commercial, domestic and industrial sectors of South Africa and were for the period January 1997 to May 2014, giving us a total of 6299 observations. The R statistical package Version 4.04 [48] was used in this study.

4.2. Exploratory Data Analysis

Table 3. Descriptive statistics of DPED.

Var	Min	Q1	Mean	Median	Q3	Max	Skew	Kurt
DPED	17,605	25,706	28,688	29,149	31,596	37,158	−0.232	2.287

presents the summary statistics of the DPED. The minimum and maximum DPED values during the sampling period were 17 605MW and 37 158MW, respectively. The skewness value was −0.232, showing that the distribution of the DPED was skewed to the left. The density (Panel (c)) and box (Panel (d)) plots given in Figure 1 confirm that the distribution of the DPED was negatively skewed.

Figure 1 (a) shows a plot of the DPED. The DPED data were not normally distributed as shown by Panel (b) of Figure 1.

A plot of the DPED superimposed with a nonlinear trend variable is given in Figure 2. A penalised cubic regression spline was used to obtain the nonlinear trend, which was used as a covariate for all the models, except for the extremal mixture model, where a linear trend was used as a covariate.

The penalised cubic regression smoothing spline model is given in Equation (19) as:

$$\pi \left(t\right)=\sum _{i=1}^n{\left(y_i-f\left(t_i\right)\right)}^2+\lambda \int {\left(f^{″}\left(t\right)\right)}^{2}dt,$$

(19)

where $y_i$ denotes the DPED and $\lambda$ is a smoothing parameter. The optimal estimate of the smoothing parameter $\lambda$ was found based on the generalised cross-validation (GCV). Considering the AQR model, we have:

$$q_{Y|X}\left(\tau \right)=\sum _{t=1}^n{\rho }_{\tau }\left(y_{t,\tau }-s_{\tau }\left(\mathrm{noltrend}\right)\right),\tau \in (0,1).$$

(20)

In this study, we considered four quantile values, which were 0.95, 0.99, 0.999 and 0.9999, respectively.

4.3. Results

The probabilistic accuracy measures (scoring rules) were calculated as follows: we fit a parametric distribution to the forecasts and then estimated its parameters. The study considered a comparative analysis of three models, i.e., the AQR (M1), EM (M2) and NLQR (M3) models at 0.95, 0.99 and 0.9999 quantile forecasts. To evaluate the forecast accuracy measures of the AQR, EM and NLQR models, we used the CRPS, LogS, DSS, interval width and the pinball losses, as described in Section 3.5. We used the empirical scores of the CRPS, Logs, DSS, PL and interval width to evaluate the forecast performances. The lower the values of the scoring rules, the better the prediction performance was. The scoring rules were then computed based on the fitted parametric distribution.

A comparison of the models at both extremely high and extremely low quantiles is given in

Table 4. Model comparisons (extremely high quantiles).

95.0$\mathrm{th}$ percentiles (0.95 quantile)
Models	CRPS	LogS	DSS	PL
M1 (AQR)	2144.546	9.6088	17.4418	165.9795
M2 (EM)	2069.789	9.5560	17.3701	209.2725
M3 (NLQR)	2155.875	9.6116	17.4495	161.9122
M4 (Median)	2155.875	9.6116	17.4495	165.5629
99.0$\mathrm{th}$ percentiles (0.99 quantile)
Models	CRPS	LogS	DSS	PL
M1 (AQR)	2125.457	9.5947	17.4256	43.2765
M2 (EM)	2131.232	inf	17.4315	56.8979
M3 (NLQR)	2163.629	inf	17.4598	42.6107
M4 (Median)	2163.629	inf	174598	43.2509
99.9$\mathrm{th}$ percentiles (0.999 quantile)
Models	CRPS	LogS	DSS	PL
M1 (AQR)	2172.784	inf	17.4759	5.529
M2 (EM)	2426.084	inf	17.7509	7.9267
M3 (NLQR)	2190.201	inf	17.4958	5.3599
M4 (Median)	2190.201	inf	17.4958	5.4869
99.99$\mathrm{th}$ percentiles (0.9999 quantile)
Models	CRPS	LogS	DSS	PL
M1 (AQR)	2168.997	inf	17.4747	0.6116
M2 (EM)	2945.86	inf	18.3882	1.0272
M3 (NLQR)	2202.221	inf	17.5149	0.6053
M4 (Median)	2202.221	inf	17.5149	0.6274
Interval widths for the 0.01 and 0.99 quantiles (CP = 0.98)
Models	Ave IW	Cov Prob	Below 0.01 quantile	Above the 0.99 quantile
M1 (AQR)	5364	0.9833	59	46
M2 (EM)	6807	0.9886	52	20
M3 (NLQR)	5312	0.9806	65	57
M4 (Median)	5385	0.9843	56	43

and

Table 5. Model comparisons (extremely low quantiles).

5.0$\mathrm{th}$ percentiles (0.05 quantile)
Models	CRPS	LogS	DSS	PL
M1 (AQR)	3826.125	10.3962	19.1919	226.0744
M2 (EM)	4225.651	10.5635	19.6090	287.7971
M3 (NLQR)	3789.174	10.3804	19.1481	220.9624
M4 (Median)	3789.174	10.3804	19.1481	224.6259
1.0$\mathrm{th}$ percentiles (0.01 quantile)
Models	CRPS	LogS	DSS	PL
M1 (AQR)	4529.206	10.6899	20.0098	64.0034
M2 (EM)	5089.073	10.9302	20.6903	79.2499
M3 (NLQR)	4523.984	10.6905	20.0264	63.6298
M4 (Median)	4523.984	10.6905	20.0264	64.45815
0.1$\mathrm{th}$ percentiles (0.001 quantile)
Models	CRPS	LogS	DSS	PL
M1 (AQR)	5050.386	10.9086	20.6212	7.793396
M2 (EM)	5257.614	11.0022	20.9055	8.342102
M3 (NLQR)	5033.644	10.9037	20.6112	7.756549
M4 (Median)	5033.644	10.9037	20.6112	7.781235

, respectively. At the 0.95 quantile, the EM model had the lowest CRPS, LogS and DSS, making it the best fitting model at this quantile level. The AQR model was the best fitting model at the 0.99, 0.999 and 0.9999 quantiles based on the CRPS and DSS evaluation metrics. Based on the prediction interval width coverage probability of 0.98 (see Table 4), all the models had valid coverage probabilities since they were all greater than 0.98. However, the EM model provided the largest coverage probability, which was 0.9886, and had the least number of observations below the 0.01 and 0.99 quantiles, respectively. As for the extremely low quantiles, the NLQR model had the smallest CRPS, LogS and DSS at the 0.05 quantile. At the 0.01 quantile, the AQR model had the smallest LogS and DSS values. The NLQR model was the best fitting model at the 0.001 quantile based on the CRPS, LogS and DSS evaluation metrics.

Figure 3, Figure 4 and Figure 5, respectively, show the DPED superimposed with the 0.95, 0.99 and 0.01 quantiles from the M1, M2 and M3 models, respectively. All figures reflect the peak identification in the electricity demand based on the 0.95, 0.99 and 0.999 quantiles of the DPED data. Furthermore, Figure 3, Figure 4 and Figure 5 show the black solid line depicting the DPED from January 1997 to May 2014, the red dotted line the M1 model distribution, the blue dotted line the M2 model and the green dotted line the M3 model distribution, respectively. Model M2 gave the best results at the 0.95 quantile.

Figure 6, Figure 7 and Figure 8, respectively, show three box plots of the 0.99, 0.999 and 0.9999 quantiles with the AQR, EM and NLQR models. The box plots for the EM model have shorter left tails and longer right tails compared to the other two models, AQR and NLQR, respectively. The medians of the AQR and NLQR models are all at the same level; however, their box plots show very different distributions for the DPED data.

Figure 9 shows the forecast distributions for the DPED data from 1997 to 2021. It also summarises the forecast distributions drawn at each specific date interval. For every five years, the histogram represents a sample from the model’s forecast distribution.

The Murphy diagrams in Figure 10, Figure 11 and Figure 12 show empirical scores and differences in scores for the M1 and M2, M1 and M3 and M2 and M3 models, respectively. A negative difference means that the regime-switching forecast is preferable. It must be noted that Murphy diagrams (Figure 10, Figure 11 and Figure 12) might lead to inconclusive situations in which neither of the three forecast methods dominates the other, and as a result, it would be unhelpful in decision-making. The p-values of the three models for three different forecasting horizons were analysed. Based on the information given in Section 3.4, the forecast accuracies of the 0.95, 0.99, 0.999 and 0.9999 quantiles are summarised. Since all p-values were less than 0.05, we, therefore, rejected the null hypothesis that there was no significant difference in the predictive abilities amongst the models. The forecasting performance of all the models was significant. Hence, the forecasting accuracy of M1 was better than that of the M2 and M3 models in extremely high quantiles, and the M3 model was better than the M1 and M2 models in extremely low quantiles.

5. Discussion

The current study was motivated by the need to accurately predict the largest possible demand for electricity at any given point in time. This can help system operators in shifting demand to off-peak periods including the scheduling of generating units during peak periods.

The study carried out a comparative analysis of the AQR, EM and NLQR models in predicting extremely high and extremely low quantiles of daily peak electricity demand. Other approaches, for example [4,13,32] among others, underestimated conditional quantiles or failed to estimate extreme conditional quantiles from their proposed models.

The additive quantile regression models are becoming increasingly popular in many applications as they are known to be flexible and robust. The AQR model does not require a predetermined function fit; however, it determines the best fit from the DPED data under extremely high quantiles. Table 4 shows the AQR model as a powerful method as suggested by the low evaluation metrics. Table 5 shows NLQR to be the best fitting model for extremely low quantiles.

The predictive performance of the models was evaluated based on three evaluation metrics, the CRPS, LogS and DS, respectively. The models were compared based on these scoring rules. The comparative forecast evaluation is generally applauded for giving critical features with out-of-sample forecasts in modelling comparison. The AQR model showed the smallest values of the scoring rules in all three extremely high quantiles, except for one quantile. All sets of scores’ (CRPS, LogS and DSS) values in Table 4 suggest that the M1 model was the best model, while Table 5 suggests that M3 was the best fitting model. Hence, the model M1 provided the highest predictive accuracy at the 0.95 and 0.99 quantiles, respectively, as given in both Figure 3 and Figure 4.

The accurate predictions in extreme conditional quantiles of electricity demand in this study are necessary for planning power systems and assessing investment projects in South Africa. Based on the Murphy diagrams, M1 (AQR model) had the highest predictive ability compared to the extremal mixture and nonlinear quantile regression models. The question of model comparisons in both extremely high and low quantiles in this study arises from the need to:

• Address the uncertainties that seem to be ignored in practice (for example, the uncertainty in the process that is generating the occurrence of the extreme events);
• Quantify the uncertainties in the estimated parameters of the distribution;
• Predict extremely high quantiles of daily peak electricity demand. This helps system operators know the possible largest demand, which will enable them to supply adequate electricity to consumers and shift load to off-peak periods.

A limitation of this study is that the sample data were for the period 1997–2014. It would have been good to cover the period of the COVID-19 pandemic. However, at the time of the study, the data covering the pandemic period were not available to the researchers. This would have helped in seeing the impact of the pandemic on electricity demand. Future research, when data are made available, should include this period. Other future research directions would include the inclusion of more covariates such as meteorological and economic including calendar variables so that we may capture their effect on the prediction of extremely high and extremely low electricity demand.

6. Conclusions

The paper presented a comparative analysis of the prediction of extremely high and extremely low daily peak electricity demand. An extremal mixture model in which a kernel density was fit to the bulk model and a Pareto distribution fit to the tail model was compared with an additive quantile regression and a nonlinear quantile regression models. The results indicated that the additive quantile regression model produces the most accurate predictions at both extremely high and low quantile levels. Furthermore, the additive quantile regression and nonlinear quantile regression models fit well based on the Murphy diagrams for the comparison of the 0.99 quantile compared to the other models. The results from this study could be useful to decision-makers in power utility companies in the prediction of extremely high and low electricity demand, thereby assisting them in the management of the risk of overprediction and underprediction.