Bootstrap prediction regions for daily curves of electricity demand and price using functional data

Bo otstrap prediction regions for daily curv es of electricit y demand and price using functional data Reb eca P el´ aez a , Germ´ an Aneiros a , Juan Vilar a a R ese ar ch Gr oup MODES, Dep artment of Mathematics, F aculty of Computer Scienc e, Universidade da Coru˜ na, A Coru˜ na, 15071, Sp ain Abstract The aim of this pap er is to compute one-da y-ahead prediction regions for daily curv es of electricit y demand and price. Three mo del-based pro cedures to construct general prediction regions are proposed, all of them using b o otstrap algorithms. The first proposed metho d considers an y L p norm for functional data to measure the distance b et w een curves, the second one is designed to take differen t v ariabilities along the curve in to accoun t, and the third one takes adv antage of the notion of depth of a functional data. The regression mo del with functional resp onse on whic h our prop osed prediction regions are based is rather general: it allo ws to include both endogenous and exogenous functional v ariables, as w ell as exogenous scalar v ariables; in addition, the effect of suc h v ariables on the resp onse one is modeled in a parametric, nonparametric or semi-parametric w ay . A comparativ e study is carried out to analyse the p erformance of these prediction regions for the electricit y market of mainland Spain, in year 2012. This work extends and complements the metho ds and results in Aneiros et al. (2016) (fo cused on curve prediction) and Vilar et al. (2018) (fo cused on prediction in terv als), which use the same database as here. K eywor ds: Bo otstrap, Electricit y mark ets, Load and price, F unctional time series, Prediction regions, Regression 2000 MSC: , 62F40, 62G08, 62M10, 91B84 1. In tro duction Prediction is one of the main aims of time series analysis. In fact, ha ving at hand accurate future v alues of some v ariable of in terest enables managers or p olicy mak ers to mak e prop erly informed decisions. In the case of electricity mark ets, giv en that electricity cannot b e stored, it is v ery imp ortan t for pro ducers to forecast future v alues of electricity demand to an ticipate to future demands and to av oid o verproduction. In the same wa y , kno wing prop er future v alues of electricity price allo ws agen ts and companies in v olv ed in the electricit y mark et to design appropriate strategies, contributing to an ticipate scenarios to mak e timely decisions. In the last decades, a lot of research w as fo cused on electricit y demand and price fore- casting. Most of suc h researc h considered the case where all the in volv ed v ariables are scalar (real-v alued v ariables), fo cusing mainly on hourly and daily demand (or price) forecasting at different forecast horizons (one da y ahead for hourly forecasts and several da ys ahead for daily forecasts). F or time-series metho dology related to electricity demand and price fore- casting, see for instance the reviews by Sugan thi and Samuel (2012), W eron (2014), Hong Pr eprint submitte d to IJEPES January 23, 2024 and F an (2016) and Now otarski and W eron (2018), as well as W eron (2006) for a monograph. See also Kay a et al. (2023) for a recent study based on deep learning approaches (artificial neural netw orks). No wada ys, with the dev elopment of mo dern tec hnology together with the high storage capacit y , it is usual to deal with functional data (function-v alued v ariables); that is, data ob- serv ed on a contin uum, often time but not only . Some examples of functional data are curves (daily electricity demand or price curv es, . . . ), images (magnetic resonance imaging, . . . ), etc. F rom the 1990s, F unctional Data Analysis (FDA) has emerged as one of the most chal- lenging tasks in Statistics, and many inference techniques (estimation, testing of hypotheses, v ariance analysis, prediction, classification, etc) w ere extended from multidimensional data (finite dimensional) to functional data (infinite dimensional). See, for instance, Ramsa y and Silv erman (2005) and F errat y and Vieu (2006) for some early monographs on parametric and nonparametric modelling of functional data, respectively , and the more recen t monographs b y Zhang (2014) and Kok oszka and Reimherr (2017) for v ariance analysis and an introduc- tion to sev eral topics in FDA, resp ectiv ely; see also Aneiros et al. (2022) for a comp endium of recent adv ances in FDA. The case of electricity demand and price forecasting from mo dels including functional data was also considered in the literature (which is m uch more limited than that of the scalar case). A pioneer w ork w as that of Vilar et al. (2012), where hourly electricit y demand and price were predicted from b oth functional nonparametric and semi-functional partially linear regression mo dels with scalar resp onse (hourly electricity demand or price) and some functional explanatory v ariable (lagged daily curv es of electricicit y demand or price); in the case of the semi-functional partially linear mo del, dummy v ariables w ere also included to indicate the type of day . See F erraty and Vieu 2004 and Aneiros-P´ erez and Vieu 2008 for some first theoretical results on the t wo referred regression mo dels, resp ectiv ely . The article by Vilar et al. (2012) was extended in Ra ˜ na et al. (2018), where additive mo dels w ere considered; in addition, some exogenous scalar co v ariates were included in the mo dels. Liebl (2013) prop osed a new statistical p ersp ectiv e for mo delling and forecasting electricity sp ot prices that accoun ts for the merit order mo del, where a functional factor mo del is used to parameterise the series of daily price-demand functions. Recently , Lisi and Shah (2020) in vestigated the forecasting p erformance of several mo dels for the 1-da y-ahead prediction of demand and price on four electricit y mark ets: nonparametric pro cedures are used to estimate the deterministic comp onen t, while the functional approach takes part in the mo delling of the residual sto c hastic comp onen t. All the articles referred to in the previous paragraph are fo cused on prediction of a scalar resp onse v ariable (in general, hourly or daily electricity demand or price) using information from some functional explanatory v ariable. An to c h et al. (2010) was one of the first pap ers dev oted to prediction of curv es of electricity consumption. More sp ecifically , they predicted the forthcoming w eekend –and weekda ys– consumption curves from functional linear regres- sion mo dels, where b oth the resp onse and explanatory v ariables w ere of functional nature. Aneiros et al. (2016) extended the pro cedures in Vilar et al. (2012) in tw o wa ys: on the one hand, some exogenous scalar co v ariates w ere added to the mo dels (temp erature when the aim w as demand prediction, and forecasted daily demand and wind pow er production for the case of price prediction); on the other hand, the resp onse v ariable was functional: daily electricity demand and price curv es. Among the conclusions obtained in Aneiros et al. (2016), it is 2 w orth highlighting: ‘to tak e information from exogenous cov ariates impro ves the forecasts, hourly predictions obtained from the discretisacion of the predicted curves are b etter than the predictions given b y the corresp onding functional mo dels with scalar resp onse, and the computational time to obtain suc h hourly predictions from discretisation is m uch low er than the needed time to implemen t the corresp onding scalar procedures’ . Some other related w orks dealing with functional resp onse (more sp ecifically , fo cused on prediction of in terest curv es in the electric p o wer mark et) are Aneiros et al. (2013), P aparo ditis and Sapatinas (2013), Chen and Li (2017), Portela-Gonz´ alez et al. (2018), El ´ ıas et al. (2022) and Barrien tos-Mar ´ ın et al. (2023). This pap er is motiv ated b y the conclusions in Aneiros et al. (2016) referred in the previous paragraph. The aim of this pap er is twice: on the one hand, to prop ose metho dology to construct prediction regions for a future curv e in a setting of functional time series; on the other hand, to put in practice such regions for the case of daily electricit y demand and price curves. A prediction region, at a predetermined confidence level 1 − α , for a random v ariable ζ N +1 is a random region in which ζ N +1 will tak e v alues with probability 1 − α . In the case that ζ N +1 is real v alued, the term ‘prediction in terv al’ is used instead of ‘prediction region’ . At this momen t, it is w orth b eing noted the main difference b et ween prediction in terv als and prediction regions. T o clarify this, let us consider the example w e will discuss in Section 4, where a curve to predict will b e the daily curv e of electricit y demand, sa y { ζ N +1 ( t ); t ∈ (0 , 24] } . If R α = { ( LR N +1 ( t ) , U R N +1 ( t )); t ∈ (0 , 24] } is a prediction region for ζ N +1 at confidence lev el 1 − α , then it verifies that P ( LR N +1 ( t ) < ζ N +1 ( t ) < U R N +1 ( t ); ∀ t ∈ (0 , 24]) = 1 − α. Nev ertheless, if one computes a prediction interv al for ζ N +1 ( t ), sa y ( LI N +1 ( t ) , U I N +1 ( t )), for eac h t ∈ (0 , 24], what is verified is that P ( LI N +1 ( t ) < ζ N +1 ( t ) < U I N +1 ( t )) = 1 − α, ∀ t ∈ (0 , 24] , but not P ( LI N +1 ( t ) < ζ N +1 ( t ) < U I N +1 ( t ); ∀ t ∈ (0 , 24]) = 1 − α. In particular, when one computes b oth R α and I α = { ( LI N +1 ( t ) , U I N +1 ( t )); t ∈ (0 , 24] } from observ ed data, one should exp ect that the area of R α is greater than the one of I α . The topic of prediction interv als for electricit y demand and price using functional data w as dealt, for instance, in Vilar et al. (2018) and Ra ˜ na et al. (2018). In particular, they constructed mo del-based prediction in terv als using b ootstrap pro cedures, and suc h interv als w ere applied to hourly electricity demand and price. The case of prediction regions in a con- text of functional time series was studied, from a theoretical p oint of view, in Zhu and P olitis (2017); as in the case of Vilar et al. (2018) and Ra ˜ na et al. (2018) for prediction interv als, the prediction regions prop osed in Zhu and Politis (2017) were based on b oth a regression mo del and b ootstrap tec hniques. In this pap er, we prop ose three mo del-based approaches to construct prediction regions us ing b o otstrap algorithms: the first approach extends the prop osal in Zh u and Politis (2017) from a functional autoregressive nonparametric mo del to a fairly general functional mo del (allowing to include b oth endogenous and exogenous func- tional v ariables, as well as exogenous scalar v ariables); the second approach is designed to tak e different v ariabilities along the curve in to account; finally , the third approach is based 3 on the notion of depth of a functional data. Then, these three approaches will b e used to compute and compare prediction regions for future daily curv es of electricity demand and price in the Spanish Electricity Mark et. The rest of this work is organised as follows. The three mo del-based prop osals to con- struct prediction regions are presented in Section 2. More sp ecifically: our prediction regions dep ends on a fairly general regression mo del, which is stated in Section 2.1; then, Section 2.2 is devoted to motiv ate b oth the construction of each prediction region and the corresp onding b ootstrap algorithms; b ecause of the regression function is unkno wn, Section 2.3 presen ts t w o particular regression functions and their corresp onding estimators; some prop osals to c ho ose the tuning parameters related to our pro cedures are given in Section 2.4. Our metho dology to construct prediction regions is applied to data related to the electricit y mark et of mainland Spain, in y ears 2011 and 2012, and they are presen ted in Section 3. The obtained results are rep orted and compared in Section 4. Finally , Section 5 pro vides some conclusions. 2. Prediction regions for functional data This section con tains the main methodological contribution of this paper: three prop osals to construct prediction regions for functional time series. Such prop osals are b oth mo del based and b o otstrap based. Firstly we establish the considered general mo del and then we construct the prop osed b o otstrap prediction regions. Some particular mo dels and a guide to select the tuning parameters are also included. 2.1. Gener al setting Let { ζ ( t ) } t ∈ R b e a real v alued con tinuous time seasonal sto c hastic pro cess with seasonal length τ . Then, assuming that suc h pro cess is observed on the in terv al I = ( a, b ] ⊂ R with b = a + N τ , one can cut { ζ ( t ) } t ∈ I in to N succesiv e random curves, { ζ i } N i =1 , where ζ i ( t ) = ζ ( a + ( i − 1) τ + t ) , with t ∈ (0 , τ ] . And finally , it seems natural to tak e adv an tage of the dep endence b et ween the random curv es, ζ i , in order to predict ζ N +1 using information given from the discrete time sto c hastic pro cess, { ζ i } N i =1 . Note that, mathematically , ζ i is a functional data that take v alues in the space of real functions defined in (0 , τ ]. Suc h space will b e denoted by F . It is a separable Hilb ert space, so it is endo wed with an inner pro duct ⟨· , ·⟩ and the corresp onding norm, ∥·∥ , defined b y ∥ ζ ∥ 2 = ⟨ ζ , ζ ⟩ . In addition to ζ i ( i = 1 , . . . , N ), some exogenous v ariables (say x i and ξ i ) could b e relev ant in predicting ζ N +1 . In this pap er, w e will obtain model-based prediction regions. More sp ecifically , w e will assume that the follo wing regression mo del holds: ζ i = r ( ζ i − 1 , x i , ξ i ) + ε i , i ∈ S ⊂ I = { 2 , . . . , N } , (1) where x i and ξ i are p -dimensional and q -dimensional v ectors of exogenous scalar and func- tional v ariables, resp ectiv ely . The unkno wn op erator r ( · ) denotes the regression function while { ε i } are the random functional errors with zero mean and indep enden t and identically distributed. F rom no w on, n S will denote the cardinal of S , while the vector of explanatory v ariables in mo del (1), ( ζ i − 1 , x i , ξ i ), will b e denoted by χ i . 4 Note that mo del (1) is a fairly general mo del, allo wing to include the effect of an en- dogenous v ariable and exogenous v ariables of scalar and/or functional nature. Note also that, in order to make clearer our prop osal, we ha v e considered only one endogenous v ariable (autoregression of order one), although the pro cedure could b e generalised to higher orders. 2.2. Building the pr e diction r e gions A prediction region for ζ N +1 at confidence level 1 − α is a random subset R α ⊆ F , such that P ( ζ N +1 ∈ R α ) = 1 − α. F rom this definition, and taking in to accoun t mo del (1), it follows that our prediction region dep ends on the op erator r ( · ), which is unkno wn. F rom now on, b r h ( · ) will denote a nonpara- metric estimator of r ( · ), with h b eing a smo othing parameter. Section 2.3 will present the regression functions, r ( · ), and the corresp onding estimators, b r h ( · ), explored in this study . In this pap er, three metho ds for computing prediction regions for ζ N +1 from information giv en by the sample { ( χ i , ζ i ) , i ∈ S } are prop osed. The three metho ds are based on an algorithm for b ootstrap resampling of the such sample. The algorithm is shown b elo w: A lgorithm for Bo otstr ap R esampling (BR). 1. Compute b r b ( χ i ) , i ∈ S , that is, the prediction of ζ i | χ i , where b > 0 is a smoothing parameter. 2. Compute the cen tered residuals, b ϵ i,b = b ε i,b − b ¯ ε b , i ∈ S, (2) where b ε i,b = ζ i − b r b ( χ i ) and b ¯ ε b = n − 1 S P i ∈ S b ε i,b . 3. Dra w n S i.i.d. random v ariables from the empirical distribution of { b ϵ i,b } i ∈ S . Such v ariables are denoted b y b ε ∗ i , i ∈ S (b o otstrap pseudo-residuals). 4. Compute ζ ∗ i = b r b ( χ i ) + b ε ∗ i , i ∈ S . 5. Obtain { ( χ i , ζ ∗ i ) , i ∈ S } , a b o otstrap resample of the original sample { ( χ i , ζ i ) , i ∈ S } . The following Sections 2.2.1, 2.2.2 and 2.2.3 present our three prop osed metho ds to con- struct prediction regions. 2.2.1. L p -metho d for pr e diction r e gions In this section, we assume that the observ ed curv e ζ N +1 | χ N +1 b elongs to the normed func- tional space L p ( F ) ≡ ( F , ∥ · ∥ p ), where the norm ∥ · ∥ p is defined using a natural generalisation of the p -norm for finite-dimensional vector spaces. The pap er b y Zh u and P olitis (2017) w as one of the first articles in the statistical literature where mo del-based prediction regions w ere prop osed in a setting of functional time series. Suc h prediction regions w ere based on functional autorregression mo dels of order 1. The metho d we prop ose here extends the prop osal in Zhu and P olitis (2017) to the case of the general mo del in (1). 5 If we consider a v alue ρ α > 0 v erifying P     ζ N +1 | χ N +1 − b r h ( χ N +1 )    p < ρ α  = 1 − α, then, R ρ α = n φ ∈ L p ( F ) :    φ − b r h ( χ N +1 )    p < ρ α o is a prediction region for ζ N +1 | χ N +1 computed at the 1 − α confidence lev el (note that R ρ α is the ball around b r h ( χ N +1 ) of radius ρ α , whic h is denoted b y B ( b r h ( χ N +1 ) , ρ α )). But, in practice, R ρ α is infeasible b ecause ρ α dep ends on the distribution of ∥ ζ N +1 | χ N +1 − b r h ( χ N +1 ) ∥ p , whic h is unknown. W e prop ose to appro ximate the distribution of ∥ ζ N +1 | χ N +1 − b r h ( χ N +1 ) ∥ p b y means of a b ootstrap pro cedure, which takes adv entage of the equality (see (1)) ζ N +1 | χ N +1 − b r h ( χ N +1 ) = r ( χ N +1 ) − b r h ( χ N +1 ) + ε N +1 | χ N +1 . (3) In the follo wing, a b o otstrap approximation for ∥ r ( χ N +1 ) − b r h ( χ N +1 ) + ε N +1 | χ N +1 ∥ p is given b y ∥ b r b ( χ N +1 ) − b r ∗ h ( χ N +1 ) + ε ∗ N +1 ∥ p . Then, our prop osed b ootstrap prediction region is R ρ ∗ α , where ρ ∗ α is such that P ∗     b r b ( χ N +1 ) − b r ∗ h ( χ N +1 ) + ε ∗ N +1    p < ρ ∗ α  = 1 − α. The algorithm to obtain the prediction region for ζ N +1 at confidence level 1 − α is explained b elo w. The Monte Carlo metho d is used to appro ximate the radius ρ ∗ α , so that the prediction region, R ρ ∗ α , has a confidence level appro ximately equal to 1 − α . A lgorithm of the L p metho d for pr e diction r e gions. 1. Compute the predictor for ζ N +1 giv en b y b ζ N +1 = b r h ( χ N +1 ) with bandwidth h . 2. Using the algorithm BR, obtain the b o otstrap resample ζ ∗ i = b r b ( χ i ) + b ε ∗ i , i ∈ S . (Note that ζ ∗ i dep ends on b ; for the sak e of clarity , w e hav e written ζ ∗ i = ζ ∗ b,i .) 3. Compute the b o otstrap predictor of ζ ∗ N +1 from the b o otstrap resample { ( χ i , ζ ∗ i ) , i ∈ S } as follows: b ζ ∗ N +1 = b r ∗ h ( χ N +1 ) . 4. Rep eating B times Steps 2-3, obtain the B b o otstrap predictors n b ζ ∗ ,j N +1 o B j =1 = n b r ∗ ,j h ( χ N +1 ) o B j =1 . 5. Dra w B i.i.d random v ariables n ε ∗ ,j N +1 o B j =1 from the empirical distribution of the cen tered residuals { b ϵ i,b } i ∈ S (see (2)). 6. F or j = 1 , . . . , B , obtain the future b o otstrap observ ation ζ ∗ ,j N +1 = b r b ( χ N +1 ) + ε ∗ ,j N +1 . 6 7. F or j = 1 , . . . , B , compute the b ootstrap error: E ∗ j = ζ ∗ ,j N +1 − b ζ ∗ ,j N +1 = b r b ( χ N +1 ) − b r ∗ ,j h ( χ N +1 ) + ε ∗ ,j N +1 . 8. F or j = 1 , . . . , B , obtain ρ ∗ j = ∥ E ∗ j ∥ p . 9. Sort the v alues ρ ∗ 1 , . . . , ρ ∗ B , obtaining ρ ∗ (1) , . . . , ρ ∗ ( B ) and select ρ ∗ α = ρ ∗ ([ B (1 − α )]) . 10. The (1 − α ) 100% prediction region for ζ N +1 | χ N +1 consists of all φ suc h that ∥ φ − b r h ( χ N +1 ) ∥ p ≤ ρ ∗ α . Remark 1. It is worth b eing note d that, given the p opularity of the L p norms, in this se ction we have c onsider e d that kind of norms; so the metho d was name d L p -metho d. A ctual ly, everything pr esente d in this se ction is valid if inste ad of the norm ∥ · ∥ p one uses any other norm in the functional sp ac e F . 2.2.2. λ -metho d for pr e diction r e gions The λ -metho d for prediction regions is based on finding the v alue of λ α > 0 suc h that P     ζ N +1 ( t ) | χ N +1 − b r h ( χ N +1 )( t )    < λ α σ ( t ) , ∀ t ∈ (0 , τ ]  = 1 − α, where σ 2 ( t ) = V ar  b r h ( χ N +1 )( t )  ∀ t ∈ (0 , τ ]. Th us, the theoretical prediction region at confidence level 1 − α is defined by R λ α ,σ =  φ ∈ F : φ ( t ) ∈  b r h ( χ N +1 )( t ) − λ α σ ( t ) , b r h ( χ N +1 )( t ) + λ α σ ( t )  , ∀ t ∈ (0 , τ ]  . As in the case of our first prop osed prediction region, R ρ α (see Section 2.2.1), R λ α ,σ is in- feasible: λ α dep ends on b oth the unkno wn parameter σ ( t ) and the unknown distribution of | ζ N +1 ( t ) | χ N +1 − b r h ( χ N +1 )( t )    . W e prop ose a b o otstrap algorithm to approximate λ α and σ ( t ) b y means of λ ∗ α and σ ∗ ( t ). More specifically , according to (3), the v alues of λ ∗ α and σ ∗ ( t ) m ust satisfy that p ( λ ∗ α ) = P ∗     b r b ( χ N +1 )( t ) − b r ∗ h ( χ N +1 )( t ) + ε ∗ N +1 ( t )    < λ ∗ α σ ∗ ( t ) , ∀ t ∈ (0 , τ ]  = 1 − α, (4) where b is some auxiliary smo othing parameter. The algorithm to obtain the bo otstrap prediction region for ζ N +1 | χ N +1 at confidence lev el 1 − α is explained b elo w. The Mon te Carlo metho d is used to approximate σ ∗ ( t ), and an iterativ e metho d is used to approximate the v alue of λ ∗ α so that the prop osed prediction region, R λ ∗ α ,σ ∗ , has a confidence level appro ximately equal to 1 − α . 7 A lgorithm of the λ -metho d for pr e diction r e gions. 1. Compute the predictor for ζ N +1 giv en b y b ζ N +1 = b r h ( χ N +1 ) with bandwidth h . 2. Using the algorithm BR, obtain the b ootstrap resample ζ ∗ i = b r b ( χ i ) + b ε ∗ i , with i ∈ S . 3. Compute the b o otstrap predictor of ζ ∗ N +1 from the b o otstrap resample { ( χ i , ζ ∗ i ) , i ∈ S } as follows: b ζ ∗ N +1 = b r ∗ h ( χ N +1 ) . 4. Rep eating B times Steps 2-3, obtain the B b o otstrap predictors n b ζ ∗ ,j N +1 o B j =1 = n b r ∗ ,j h ( χ N +1 ) o B j =1 . 5. Dra w B i.i.d random v ariables { ε ∗ ,j N +1 } B j =1 from the empirical distribution of the cen tered residuals { b ϵ i,b } i ∈ S (see (2)). 6. Appro ximate the standard deviation of b r h ( χ N +1 )( t ) by σ ∗ ( t ) ≃   1 B B X j =1 b r ∗ ,j h ( χ N +1 )( t ) − 1 B B X k =1 b r ∗ ,k h ( χ N +1 )( t ) ! 2   1 / 2 . 7. Use an iterativ e metho d to obtain an approximation of the v alue λ ∗ α defined in (4). 8. The (1 − α ) 100% prediction region for ζ N +1 | χ N +1 consists of all φ suc h that | φ ( t ) − b r h ( χ N +1 ) ( t ) | ≤ λ ∗ α σ ∗ ( t ) , for all t ∈ (0 , τ ]. Iter ative metho d to appr oximate λ ∗ α . The iterative metho d to approximate the v alue of λ ∗ α > 0 so that the b o otstrap prediction region has a confidence level appro ximately equal to 1 − α is explained b elo w. This algorithm, which follows ideas of Cao et al. (2010), allo ws to quickly and efficiently appro ximate the parameter λ ∗ α . Let n b ζ ∗ ,j N +1 o B j =1 = n b r ∗ ,j hb ( χ N +1 ) o B j =1 b e the b ootstrap predictors of ζ N +1 and { ε ∗ ,j N +1 } B j =1 the b ootstrap mo del errors. Define the Mon te Carlo appro ximation of p ( λ ) in (4), for an y λ > 0, as follows: p ( λ ) ≃ 1 B B X j =1 I     b r b ( χ N +1 )( t ) − b r ∗ ,j h ( χ N +1 )( t ) + ε ∗ ,j N +1 ( t )    < λσ ∗ ( t ) , ∀ t ∈ (0 , τ ]  . (5) Let λ L , λ H > 0 be suc h that p ( λ L ) ≤ 1 − α ≤ p ( λ H ) and let η > 0 be a tolerance, for example, η = 10 − 4 . 1. Obtain λ M = λ L + λ H 2 and compute Mon te Carlo appro ximations of p ( λ L ), p ( λ M ) and p ( λ H ) according to (5). 2. If p ( λ M ) = 1 − α or p ( λ H ) − p ( λ L ) < η , then λ ∗ α = λ M . Otherwise, (a) If 1 − α < p ( λ M ), then λ H = λ M and return to Step 1. (b) If p ( λ M ) < 1 − α , then λ L = λ M and return to Step 1. 8 2.2.3. Depth-b ase d metho d for pr e diction r e gions El ´ ıas et al. (2022) prop osed depth-based metho ds for functional time series forecasting. In this section, some of such ideas are used to construct prediction regions. In an y case the metho d prop osed in this pap er differs from the one prop osed in El ´ ıas et al. (2022) since ours inv olv es b o otstrap pro cedures in the computation of prediction regions. The prop osed algorithm is defined b elo w. A lgorithm of the depth-b ase d metho d for pr e diction r e gions. 1. Using the algorithm BR, obtain the b ootstrap resample ζ ∗ i = b r b ( χ i ) + b ε ∗ i , with i ∈ S . 2. Compute the b o otstrap predictor of ζ ∗ N +1 from the b o otstrap resample { ( χ i , ζ ∗ i ) , i ∈ S } as follows: b ζ ∗ N +1 = b r ∗ h ( χ N +1 ) . 3. Rep eating B times Steps 1-2, obtain the B b o otstrap predictors n b ζ ∗ ,j N +1 o B j =1 = n b r ∗ ,j h ( χ N +1 ) o B j =1 . 4. Dra w B i.i.d random v ariables n ε ∗ ,j N +1 o B j =1 from the empirical distribution of the cen tered residuals { b ϵ i,b } i ∈ S (see (2)). 5. F or j = 1 , . . . , B , obtain the future b o otstrap observ ation ζ ∗ ,j N +1 = b r ∗ ,j h ( χ N +1 ) + ε ∗ ,j N +1 . 6. Obtain the C = [(1 − α ) B ] deepest b ootstrap future observ ations { ζ ∗ , ( j ) N +1 } C j =1 where ζ ∗ , ( j ) N +1 is the j -th deep est curve. 7. Compute the lo wer and upp er limits of the prediction region: L N +1 ( t ) = min { ζ ∗ , ( j ) N +1 ( t ) : j = 1 , . . . , C } and U N +1 ( t ) = max { ζ ∗ , ( j ) N +1 ( t ) : j = 1 , . . . , C } , resp ectiv ely . 2.3. Pr e diction metho ds for functional time series Recalling the general formulation of the regression mo del given in (1), different particular mo dels can b e assumed for the unknown regression function r ( · ) and differen t estimation metho ds can b e considered. In this pap er, t wo regression mo dels with functional resp onse will b e used: the F unctional Non Parametric (FNP) mo del and the Semi-F unctional Partial Linear (SFPL) mo del. Both mo dels ha ve been used in Aneiros et al. (2013) for the functional prediction of residual demand curves in the electricity mark et and in Aneiros et al. (2016) for predicting daily electricit y demand and price curves. F or the case of the FNP mo del, Nadara ya-W atson type estimators were used; for the case of the SFPL mo del, a combina- tion of least squares t yp e estimators and Nadaray a-W atson type estimators were considered. Those mo dels and estimators will b e also used in this pap er, and they are presented in the t wo follo wing sections. 9 2.3.1. The functional non p ar ametric mo del Consider that the vector of explanatory v ariables in the general mo del (1) is reduced to χ i = ζ i − 1 (that is, no exogenous v ariables are considered). Then, mo del (1) b ecomes the functional autoregressive model of order 1 ζ i = r ( ζ i − 1 ) + ε i , i ∈ S. (6) Giv en that we only assume general conditions on r ( · ) (for instance, smo othness), (6) is, in particular, a FNP mo del. Therefore, the op erator r ( · ) should b e estimated from some nonparametric approach. W e will use the nonparametric Nadaray a-W atson type estimator, whic h is defined as follows: b r F N P h ( ζ ) = X i ∈ S w h ( ζ , ζ i − 1 ) ζ i , (7) where w h ( · , · ) are the Nadara ya-W atson weigh ts giv en b y w h ( ζ , ζ i − 1 ) = K ( d ( ζ , ζ i − 1 ) /h ) P j ∈ S K ( d ( ζ , ζ j − 1 ) /h ) , with K : [0 , ∞ ) → [0 , ∞ ) b eing a k ernel function and h > 0 a smo othing parameter. In addition, d ( · , · ) denotes the semi-metric used to measure the pro ximity b et ween t wo curves in F . Then, a one-step-ahead prediction for the curve ζ N +1 can b e obtained b y means of b ζ N +1 = b r F N P h ( ζ N ). Studies on the asymptotic prop erties of the estimator (7) can b e seen in F erraty et al. (2012) and Zhu and P olitis (2017), among others. In Aneiros et al. (2016) the mo del (6) and the estimator (7) are used to predict daily curv es of electricity demand and price. The case of scalar resp onse, that is, considering ζ i ( t ) as resp onse v ariable in (6) instead of ζ i , was dealt in Masry (2005) and Delsol (2009) among others. Finally , the case of indep endence b et w een the resp onse and explanatory v ariables in general FNP mo dels w as considered in F erraty et al. (2011). 2.3.2. The semi-functional p artial line ar mo del Consider now that χ i = ( ζ i − 1 , x i ) and r ( χ i ) = x ⊤ i β + m ( ζ i − 1 ) (8) in the general mo del (1). This means, in addition to an endogenous functional v ariable, exogenous escalar ones are included in the mo del. Then, we ha v e that ζ i = x ⊤ i β + m ( ζ i − 1 ) + ε i , i ∈ S. (9) As in the previous section, w e assume that m ( · ) is an unknown smo oth op erator and β = ( β 1 , . . . , β p ) ⊤ is a v ector of unkno wn functional parameters. Model (9) is, in particular, a SFPL mo del that extends the previous FNP one (6) by adding linear effects of some exogenous scalar v ariables in the regression function. Thus, the predictions obtained from mo del (9) could improv e the predictions obtained from mo del (6). Mo del (9) was prop osed in Aneiros et al. (2013) to predict curves of residual demand in electricit y mark ets. In order to put in practice suc h prediction, they also proposed estimators 10 for β and m ( · ) based on a com bination of least squares and Nadara y a-W atson t yp e estimates. Sp ecifically , such estimators are given b y b β h = ( f X T h f X h ) − 1 f X T h e ζ h (10) and c m S F P L h ( ζ ) = X i ∈ S w h ( ζ , ζ i − 1 )  ζ i − x T i b β h  , (11) resp ectiv ely , where f X h = ( I − W h ) X (12) and e ζ h = ( I − W h ) ζ , (13) with W h = ( w h ( ζ i , ζ j )) i +1 ,j +1 ∈ S , X = ( x i ) ⊤ i ∈ S = ( x ij ) i ∈ S 1 ≤ j ≤ p and ζ = ( ζ i ) i ∈ S . No w, from (8), (10) and (11), we ha v e that b r S F P L h ( x , ζ ) = x ⊤ b β h + c m S F P L h ( ζ ) is a natural estimator for r ( x , ζ ) defined in (8) in the SFPL mo del given in (9). Finally , the one-step-ahead forecast of the curve ζ N +1 is given b y b ζ N +1 = b r S F P L h ( x N +1 , ζ N ) . It is worth b eing noted that the estimators (10) and (11) were studied from a theoret- ical p oin t of view in Aneiros-P´ erez and Vieu (2008) in the case of scalar resp onse (that is, considering ζ i ( t ) as resp onse v ariable in (9) instead of ζ i ). 2.4. T uning p ar ameters As usual when one deals with nonparametric estimation, our prop osals dep end on several tuning parameters that must be chosen in some appropriate wa y . This section is devoted to indicate how suc h parameters are chosen in our application to real data (see Section 4). In the b o otstrap algorithms introduced in Section 2.2, an auxiliary bandwidth, b , is needed, as usual when w orking with b o otstrap pro cedures in nonparametric regression. The bandwidth b is used to construct the residuals to b e resampled and w as theoretically pro ven to b e larger than the bandwidth h used to smo oth the b o otstrap sample (see F erraty et al. (2010) and Ra ˜ na et al. (2016) for further details). In this pap er, the choice of the bandwidth b follo w ed the ideas of Vilar et al. (2018) for prediction interv als. See also F erraty et al. (2011) and Aneiros et al. (2016) for more details. The L p metho d dep ends on the chosen norm ∥ · ∥ p and the depth-based metho d dep ends on the depth measuremen t. Regarding the norm to b e used, the usual norms ∥ · ∥ 1 and ∥ · ∥ 2 allo w us to mathematically define the prediction region and to chec k whether or not a giv en curv e b elongs to this region. The disadv antage of these norms is that they do not allow a graphical representation of the region. Cho osing the norm ∥ · ∥ ∞ could provide a useful graphical represen tation of the resulting prediction regions. F or this reason, in Section 4 the norm ∥ · ∥ ∞ will b e used. The metho d for obtaining prediction regions based on L ∞ has the dra wback of not taking in to accoun t the v olatility of the functional curv es. But this method 11 will b e compared with the other tw o prop osals that do consider this volatilit y . Concerning the choice of depth measuremen t, in Section 4, the depth-based metho d uses the Random T ukey depth prop osed b y Cuesta-Alb ertos and Nieto-Reyes (2008), which is not very time consuming, since it is p ossible to obtain similar results with it to those obtained with more in volv ed depths by taking only a few one-dimensional pro jections. Finally , the forecasting pro cedures exp osed in Section 2.3 dep end on v arious tuning pa- rameters that m ust b e selected from the data. Sp ecifically , the FNP and SFPL predictors dep end on the bandwidth h , the semi-metric d ( · , · ) and the kernel function K ( · ). T o c ho ose h , d ( · , · ) and K ( · ), the suggestions giv en in Aneiros et al. (2013) were follo w ed. F or the bandwidth, h , the k -nearest-neigh b ours metho d w as considered. The n umber of neighbours, k , was selected b y means of the local cross-v alidation metho d. The semimetric, d ( · , · ), is c hosen based on a seminorm. Since our data are rough, we consider a seminorm based on functional principal comp onen t analysis. Finally , w e use the Epanec hniko v kernel defined b y K ( u ) = 0 . 75(1 − u 2 ) I (0 , 1) , although it is known that the c hoice of kernel has little influence on the results. 3. The data Our goal is to compute prediction regions of the daily curves of electricity demand and price, one day-ahead, corresp onding to mainland Spain, year 2012. Some information will b e taken from past curves (endogeneus functional v ariables). In addition, exogeneus scalar co v ariates could be used. Sp ecifically , our database con tains information related to years 2011 and 2012: hourly observ ations of electricity demand and price, maximum daily temp erature and daily wind p o wer pro duction. This data set w as previously analysed in v arious articles from differen t mo dels and with differen t goals: Aneiros et al. (2016) obtained p oin twise predictions for hourly electricit y demand and price from FNP (6) and SFPL (9) models, b y c hanging the functional resp onse ζ i b y the scalar one ζ i ( j ) , j = 1 , . . . , 24; based on the same mo dels, Vilar et al. (2018) completed suc h study by computing prediction interv als for future hourly electricit y demand and price; Ra ˜ na et al. (2018) also rep orted b oth p oin twise prediction and prediction in terv als for hourly electricit y demand and price, but based on additive mo dels instead of FNP and SFPL ones; Aneiros et al. (2016), in addition to the aforementioned p oint wise (scalar) prediction, dealt the case of prediction of daily curves of electricity demand and prices from mo dels (6) and (9) (now, maintaning the functional resp onse ζ i ). T o our kno wledge, there are no analyses of this data set that present prediction regions for the future curves of daily electricity demand and price. The motiv ation for c ho osing these data to illustrate the prop osals of this work is t wice: on the one hand, to b e able to compare them with those prop osed in Vilar et al. (2018) for prediction inte rv als; on the other hand, to complete and finalise the serie of pap ers related to the analysis of these data presented in the previous paragraph. Eac h daily functional data comes from the 24 hourly observ ations of demand (measured in MWh, Megaw att-hour) or price (measured in Cen t/k Wh), in each da y . Note that smo othing tec hniques are used to conv ert the 24 hourly data in a functional observ ation (a curv e). Scalar co v ariates are related to daily demand, temp erature and wind p ow er pro duction. 12 Electricit y demand and price curves show daily and w eekly seasonalit y and are affected by the w eek end calendar. Figure 1 sho ws the daily electricit y demand and price curv es along the y ear 2012. Figures 2 and 3 sho w the differences in electricit y demand and price, resp ectiv ely , according to the da y t yp e: w eekda y (Monda y to F rida y), Saturda y and Sunda y . A more comprehensiv e analysis of this effect can b e found in Vilar et al. (2018). Due to this feature detected in the curv es, the analyses carried out in this work will b e done separately for eac h t yp e of da y . 5 10 15 20 10000 20000 30000 Electricity demand Hour MWh 5 10 15 20 0 20 40 60 80 Electricity price Hour Cent/kWh Figure 1: Electricit y demand (left) and price (right) daily curves along the year 2012. 5 10 15 20 10000 15000 20000 25000 30000 Hour MWh 5 10 15 20 10000 15000 20000 25000 30000 Hour MWh 5 10 15 20 10000 15000 20000 25000 30000 Hour MWh Figure 2: Electricit y demand daily curves along the year 2012 distinguishing weekda ys (left), Saturdays (middle) and Sunda ys (righ t). 13 5 10 15 20 0 20 40 60 80 Hour Cent/kWh 5 10 15 20 0 20 40 60 80 Hour Cent/kWh 5 10 15 20 0 20 40 60 80 Hour Cent/kWh Figure 3: Electricity price daily curv es along the y ear 2012 distinguishing w eekdays (left), Saturdays (middle) and Sunda ys (righ t). Both demand and price daily curv es presen t some outliers. T o iden tify these at ypical observ ations, metho ds prop osed in Ra ˜ na et al. (2015) and Vilar et al. (2016) are used. In Ra ˜ na et al. (2015) a test based on functional depths is p erformed to determine if a curve is an outlier or not. In Vilar et al. (2016), the authors consider robust principal comp onen t analysis to lo ok for outliers in the pro jection of the first comp onen t. The identified outliers are replaced in the sample b y a weigh ted moving a verage of the surrounded days. See Vilar et al. (2018) for more details. There exist exogeneus v ariables that could improv e the forecasts and, therefore, the pre- diction regions of the daily curves of electricit y demand and price. See T aylor and Buizza (2003), T aylor et al. (2006) and Hyde and Ho dnett (2015), among others, for a discussion ab out the impact of meteorological factors on the electricit y demand and the effect of the electricit y demand and wind p o w er pro duction on the electricit y price. In Aneiros et al. (2016) and Vilar et al. (2018), a detailed analysis of these data and the influence of certain co v ariates on them is carried out. In general, there is a high demand for electric heating in cold w eather. Th us, the tem- p erature influences the electricit y demand and the maximum daily temp erature (in Celsius degrees, ◦ C ) in Spain is considered for the study . Ho wev er, this v ariable has a nonlinear effect on demand; demand as a function of temp erature exhibits a U-shap ed pattern. Since the prediction metho ds to b e used in this study (see Section 2.3) include the information of the exogenous v ariables in a linear wa y , some transformation of the maximum temp erature that verifies this assumption is needed. In this pap er, we consider the HDD (Heating Degree Da ys) and CDD (Co oling Degree Da ys) v ariables, which exert a linear effect on the electricit y demand. They are defined as follo ws: HDD( t ) = max { 20 − T ( t ) , 0 } and CDD( t ) = max { T ( t ) − 24 , 0 } . where T ( t ) is the maxim um daily temperature in da y t . See Cancelo et al. (2008) and Aneiros et al. (2016) for more details on the justification b ehind the construction of these v ariables. 14 F or the building of prediction regions of the daily electricit y price curves, b oth daily forecasted demand and wind p o wer pro duction will b e considered as scalar co v ariates. The System Op erator in the Spanish Electricity Mark et, REE, monitors the demand and its generator structure and pro vides the amount of demand cov ered b y wind p o wer during each ten minutes p erio d of the year. Then, the hourly or daily demand co vered by wind p o w er can b e obtained. I n this case, b oth cov ariates do ha ve a linear effect on the electricity price. It is imp ortan t to note that unobserved cov ariates (temp erature and wind pro duction for the day to b e forecast) are going to b e included in the regression mo dels. On the one hand, in order to apply the SFPL pro cedure, it is necessary to hav e accurate forecasts of the v alues of these co v ariates for the following day . On the other hand, sophisticated meteorological mo dels are kno wn to provide very go od forecasts of temp erature and wind p o wer. How ev er, neither these mo dels nor these forecasts are publicly a v ailable. Since the predictiv e p o wer of the SFPL mo del could b e greatly diminished or masked if not v ery go o d forecasts of these co v ariates are incorp orated in the mo del, we chose to include the ideal forecasts giv en by the v alues themselves. 4. Application to the electricit y mark et In this section, prediction regions for the daily curves of electricit y demand and price are computed using the three b o otstrap metho ds presented in Section 2.2: L p -metho d, λ - metho d and depth-based metho d. Although the metho ds can b e extended to any forecasting metho d, in this section w e will use those asso ciated with the FNP and SFPL mo dels presen ted in Section 2.3 (see (6) and (9), resp ectiv ely). Electricit y demand and price are b oth considered to b e contin uous time sto c hastic pro- cesses with seasonal length τ = 24, and the notation { ζ ( t ) } t ∈ R is used to refer an y of them (units for t are hours). The aim of this section is to put in practice the prediction regions prop osed in Sections 2.2.1, 2.2.2 and 2.2.3. More sp ecifically , to compute prediction regions for each daily curve, ζ N +1 , of electricity demand or price in y ear 2012 from information given b y the previous 365 days. This set of curves corresp onding to the 365 curv es preceding ζ N +1 is determined by the follo wing set of index: I = { N − 364 , N − 363 , . . . , N − 1 , N } . How ev er, the dynamic of the curves of demand or price dep ends on the type of day to b e predicted: w eekdays, Saturda ys or Sundays (see Figures 2 and 3). This fact suggests to compute three differen t regression mo dels according to the type of day: to forecast a curve corresp onding to either Sunday or Saturday , information from the previous curve (i.e. from the curve ob- serv ed on previous Saturda y or F riday , resp ectiv ely) will b e used; mean while, if one wishes to forecast a curv e corresponding to a weekda y , information will b e tak en from the curve observ ed on the previous weekda y (note that the previous weekda y to a Monday is a F riday). Therefore, the following notation is in tro duced to distinguish the v arious scenarios: • S W = { 1 , . . . , #( S 0 W ) } , where S 0 W = { i ∈ I : i is a weekda y } is the set of indices corresp onding to weekda ys in the 365 curv es preceding ζ N +1 . In addition we denote, ζ i = ζ ( i ) , for i ∈ S W , where ( i ) is the element in S 0 W whose rank is i , • S Sat = { i ∈ I : i is a Saturda y } is the set of indices corresp onding to Saturda ys in the 365 curves preceding ζ N +1 , 15 • S Sun = { i ∈ I : i is a Sunday } is the set of indices corresp onding to Sundays in the 365 curves preceding ζ N +1 . Note that the notation ab o ve aims to b e able to use, in eac h of the three scenarios, our notation corresp onding to the algorithms for prediction regions proposed in Section 2. In the case that ζ N +1 corresp onds to Saturda y or Sunday , it is suffice to consider S = S Sat or S = S Sun , respectively . If ζ N +1 corresp onds to a w eekday , S = S W should b e considered and, in addition, curves { ζ i } should b e used instead of { ζ i } . Note also that the num ber of curv es a v ailable to predict ζ N +1 dep ends on the type of day , since the sets S W , S Sat and S Sun do not ha ve the same cardinal. Sp ecifically , #( S W ) = 261, #( S Sat ) = 52 and #( S Sun ) = 53. Therefore, the sample size to estimate the daily demand or price curv e is muc h larger if one w ants to estimate a da y of the w eekday t yp e than if it is Saturday or Sunday . In this study , as mentioned in Section 2.4, L p metho d uses the norm L ∞ and the depth- based metho d uses the random T ukey depth. The random T uk ey depth is av ailable at the fda.usc package from the Comprehensiv e R Archiv e Net work (see F ebrero-Bande and Oviedo de la F uente (2012)). The n umber of resamples for the b o otstrap pro cedures is B = 500. 4.1. Performanc e me asur ements A prediction region p erforms w ell if its co v erage is close to the nominal one, (1 − α ) 100%, and has a small area or av erage width. The following v alues measure the p erformance of a set of J prediction regions (one prediction region for eac h curv e ζ N j +1 , j = 1 , . . . , J ) and allo w for the comparison of results. Denoting L N j +1 the lo wer limit and U N j +1 the upp er limit of the prediction region for ζ N j +1 , with j = 1 , . . . , J , the functional cov erage is the p ercen tage of b ootstrap regions that con tain to the corresp onding daily curv es of electricity demand or price, and it is defined as follo ws: F C ov = 100 × 1 J J X j =1 I  ζ N j +1 ( t ) ∈  L N j +1 ( t ) , U N j +1 ( t )  , ∀ t ∈ (0 , 24]  . The mean p oint wise cov erage is the mean of the p ercentage of time grid v alues for whic h the prediction regions contain the corresp onding daily curves of electricity demand or price. It is giv en b y P C ov = 100 × 1 J J X j =1 1 24 24 X t =1 I  ζ N j +1 ( t ) ∈  L N j +1 ( t ) , U N j +1 ( t )  ! . A verage width of the b o otstrap prediction region is defined by AW idth = 1 J J X j =1   1 24 24 X t =1  U N j +1 ( t ) − L N j +1 ( t )    . Winkler or interv al score (see Winkler (1972) and Gneiting et al. (2007)) is also used to compare the b eha viour of the metho ds. F or classical prediction in terv als, it is defined as the length of the interv al plus a p enalty if the theoretical v alue is outside the interv al. Thus, it com bines width and cov erage. F or v alues that fall within the in terv al, the Winkler score is 16 simply the length of the in terv al. So, lo w scores are asso ciated with narro w interv als. When the theoretical v alue falls outside the interv al, the p enalt y is prop ortional to how far the observ ation is from the interv al. Given a prediction interv al, (L N +1 ( t ) , U N +1 ( t )), for ζ N +1 ( t ) ( t is fixed), the form ula of the Winkler score is as follows: WS N +1 ( t ) = U N +1 ( t ) − L N +1 ( t ) + 2 α (L N +1 ( t ) − ζ N +1 ( t )) I  ζ N +1 ( t ) < L N +1 ( t )  + 2 α ( ζ N +1 ( t ) − U N +1 ( t )) I  ζ N +1 ( t ) > U N +1 ( t )  . Since w e are w orking with prediction regions, the functional Winkler score is prop osed as a criteria for the comparison. It is an extension of the classic Winkler score to the functional con text, defined as follo ws: FWS N +1 = δ (L N +1 , U N +1 ) + 2 α min n δ (L N +1 , ζ N +1 ) , δ (U N +1 , ζ N +1 ) o I  ∃ t ∈ (0 , τ ] : ζ N +1 ( t ) < L N +1 ( t ) or ζ N +1 ( t ) > U N +1 ( t )  , where δ ( · , · ) is a semi-metric. Then, when one has several prediction regions, one considers the mean of their functional Winkler scores: FWS = 1 J J X j =1 FWS N j +1 . In our study , the semi-metric to b e considered is δ ( ξ , χ ) = ∥ ξ − χ ∥ 1 = Z (0 ,τ ] | ξ ( t ) − χ ( t ) | dt. (14) Section 4.2 con tains the results for the electricity demand curves. Section 4.3 contains an analogous study for the electricity price curves. 4.2. Pr e diction r e gions for the ele ctricity demand daily curves Previous w orks ha v e proposed metho ds to obtain prediction in terv als for hourly electricity demand and price. In Vilar et al. (2018), prediction interv als are obtained for eac h hour, one-da y ahead, of the electricity demand and price, ζ N +1 ( t ) with t ∈ { 1 , . . . , 24 } , in the year 2012. They prov ed the go o d p erformance of the prediction in terv als in the problem of the next-da y forecasting of electricity demand and price. T able 1 contains a comparison b etw een the metho ds prop osed in this pap er and the use of the prediction interv als from Vilar et al. (2018) as prediction regions (after smo othing the former) for the electricity demand curves in the year 2012. The functional non parametric forecasting metho d is considered for this analysis. The results in T able 1 show that the functional co verage of the prediction regions obtained b y smo othing prediction interv als (PI) is far from the nominal cov erage: for a nominal co verage of 95%, the metho d based on PI pro vides functional cov erages of around 70%, down ev en to 54%. It is also observ ed that the functional Winkler score of the metho d based on p oin t wise prediction interv als is higher than in the three metho ds for functional prediction regions. This is to b e exp ected since computing a p oin twise prediction interv al for each of 17 the 24 p oin ts considered leads to narro w in terv als. But globally , they constitute a narro w region with a high probability that some p oin t of the real curve is not co vered. Therefore, their use as a metho d for computing functional prediction regions is discarded. W eekdays Saturdays Sunda ys Y ear L ∞ -metho d F Co v 93.8 90.4 84.9 92.0 PCo v 98.8 97.5 97.4 98.4 A Width 8360.0 11594.6 10294.8 9099.7 FWS 13433.3 15399.1 16896.6 14214.1 λ -metho d F Cov 95.4 90.4 81.1 92.6 PCo v 99.2 98.2 96.5 98.7 A Width 7684.8 11272.5 9414.3 8444.9 FWS 10684.2 14957.7 16970.9 12201.7 Depth metho d F Cov 92.3 82.7 69.8 87.7 PCo v 98.6 97.0 93.4 97.6 A Width 6932.8 9128.5 7780.8 7367.5 FWS 13014.0 20035.7 27992.8 16180.7 PI F Cov 70.1 53.8 64.2 66.9 PCo v 93.7 88.2 90.4 92.4 A Width 5339.4 4384.5 5389.3 5211.0 FWS 21005.0 24717.9 25680.9 22209.7 T able 1: F unctional cov erage (in %), point wise co verage, width and functional Winkler score of the prediction regions for the daily electricity demand based on the FNP mo del using the L ∞ -metho d, the λ -metho d, the depth-based metho d and the prediction regions obtained by smo othing the prediction interv als (PI) with α = 0 . 05 for eac h t yp e of da y in 2012. T able 2 shows the cov erage, width and Winkler score of the prediction regions obtained with the prop osed methods and using the FNP forecasting mo del for the 2012 electricity demand curv es at 80% confidence level. The condition that the en tire daily curve is con tained in the prediction region is a demanding one, so w orking with high confidence levels leads to v ery wide and, in many cases, uninformativ e prediction regions. Therefore, it is reasonable to work with α = 0 . 20 in this context. Results for 95% w ere already shown in T able 1. The co v erage p ercen tages of the three metho ds are reasonable for weekda ys, but they are low on Saturdays and Sundays, esp ecially on the latter. A t the same time, the forecast regions obtained for Saturda ys and Sundays are wider than for weekda ys. This p erformance is motiv ated b y the sample size of I W , I S at and I S un , that is, the num b er of preceding curves to b e used in the prediction for each group. When predicting demand curves on Saturda ys and Sundays, the sample size of I S at and I S un , resp ectiv ely , is m uch smaller than the sample size of I W for weekda ys. In addition, Sunda y’s daily electricit y demand curv es are highly v ariable and mak e the problem of forecasting and construction of prediction regions difficult. In the comparison of the three methods it is observ ed that λ -metho d pro vides the smallest FWS. F urthermore, this metho d has the adv an tage of taking into account the volatilit y of the demand curv e. The other t wo metho ds presen t similar results, sligh tly b etter for the L ∞ -metho d whose drawbac k is providing confidence regions of constant width throughout the day . 18 W eekdays Saturdays Sunda ys Y ear L ∞ -metho d F Co v 77.8 80.8 66.0 76.5 A Width 5640.9 7570.4 7231.6 6145.4 FWS 8937.3 10006.1 12278.0 9572.9 λ -metho d F Cov 73.9 76.9 66.0 73.2 A Width 5285.8 7084.4 6708.9 5747.4 FWS 8798.3 10234.9 11586.0 9406.1 Depth metho d F Cov 76.2 76.9 52.8 72.9 A Width 5851.9 7955.4 5750.0 6136.0 FWS 9212.4 11140.7 12276.3 9930.1 T able 2: F unctional cov erage (in %), width and functional Winkler score of the prediction regions for the daily electricity demand based on the FNP model using the L ∞ -metho d, the λ -method and the depth-based metho d at 80% confidence lev el for eac h t yp e of da y in 2012. Figure 4 shows the forecast and prediction region of electricity demand on a w eekday using the λ -metho d at the 80% confidence lev el and the L ∞ -metho d at the 95% confidence lev el based on the FNP mo del. The prediction region obtained at 80% is more informative, as it has a smaller width. The λ -metho d also takes into accoun t the volatilit y of electricit y demand throughout the da y . On the contrary , the L ∞ -metho d pro vides a region of larger width b ecause the 95% confidence lev el and also constan t for each hour. 5 10 15 20 5000 15000 25000 35000 Electric demand, Spain − 2012 Hour MWh 5 10 15 20 5000 15000 25000 35000 Electric demand, Spain − 2012 Hour MWh Figure 4: Daily curv es of electricit y demand corresponding to weekda ys (grey solid lines), curve of electricity demand for 2nd April 2012 (red solid line), prediction curve (black solid line) and prediction region (blac k dashed lines) b y means of λ -metho d at 80% of confidence (left panel) and L ∞ -metho d at 95% of confidence (righ t panel) based on the FNP forecasting model. The metho ds for obtaining prediction regions for the electricit y demand daily curves based on the SFPL regression mo del are discussed b elow. In this w ork, the vector of scalar co v ariates included in the SFPL mo del giv en in (9) to forecast the daily curv es of electricity demand is x = ( x 1 , x 2 ) T = (HDD , CDD) T , where HDD and CDD are the Heating Degree Da ys and Co oling Degree Da ys, resp ectiv ely , in tro duced in Section 3. 19 T able 3 shows the functional co v erage, width and Winkler score of the prediction regions obtained with the prop osed metho ds, and using the SFPL forecasting mo del, for the 2012 electricit y demand curves at confidence levels 80% and 95%. The results obtained with the SFPL regression mo del are very similar to those obtained with the FNP mo del. Comparing T ables 1, 2 and 3, it is concluded that the functional co verage of the 80% prediction regions based on the SFPL regression mo del is low er than that of the regions based on the FNP , but so is the region width, which leads to lo wer FWS. When working with a nominal co verage of 95% the effect is the opp osite. Regarding the metho ds for obtaining prediction regions, the w orst b eha viour with low er co verages and larger widths is observed for the depth-based metho d. The L ∞ and λ metho ds pro vide similar results, sligh tly b etter for the L ∞ -metho d. α = 0 . 05 W eekdays Saturdays Sunda ys Y ear L ∞ -metho d F Cov 90.8 88.5 79.2 89.0 A Width 8469.0 10318.0 8700.9 8765.3 FWS 16084.8 19116.7 23472.6 17585.4 λ -metho d F Cov 90.8 86.5 69.8 87.1 A Width 7793.2 9845.6 7609.6 8058.2 FWS 14724.5 19648.2 27294.4 17244.3 Depth metho d ACo v erage 87.0 76.9 62.2 82.0 A Width 7526.2 8039.3 6967.7 7518.2 FWS 16398.0 20954.7 29681.6 18969.0 α = 0 . 20 W eekdays Saturdays Sunda ys Y ear L ∞ -metho d F Cov 74.3 80.8 56.6 72.7 A Width 5284.8 6824.9 6171.3 5632.0 FWS 8800.6 9213.6 11820.8 9296.6 λ -metho d F Cov 71.3 75.0 47.2 68.3 A Width 5063.3 6477.5 5382.0 5310.4 FWS 8740.9 9438.6 11394.9 9224.4 Depth metho d FCo v 72.8 73.1 47.2 69.1 A Width 5402.4 6953.5 5935.1 5699.9 FWS 8914.3 10473.1 13074.3 9738.2 T able 3: F unctional cov erage (in %), width and functional Winkler score of the prediction regions for the daily electricity demand based on the SFPL mo del using the L ∞ -metho d, the λ -metho d and the depth-based metho d for each kind of da y in 2012 with α = 0 . 20 and α = 0 . 05. Figure 5 shows the forecast and prediction region of electricity demand on a w eekday using the λ -metho d at the 80% confidence lev el and the L ∞ -metho d at the 95% confidence lev el based on the SFPL mo del. The prediction region obtained at 80% is more informativ e, as it has a smaller width. The L ∞ -metho d pro vides a region of larger width b ecause the 95% confidence level. F urthermore, the λ -metho d takes in to account the volatilit y of electricity demand throughout the day , mean while L ∞ -metho d pro vides a region of larger width. Figure 6 sho ws the forecast and prediction region of electricit y demand on a Saturda y and a Sunday using the λ -metho d based on the SFPL mo del at the 80% confidence level. The shortcomings of forecasting on suc h days are obvious. With less data av ailable and higher 20 v olatility , as illustrated b y the cloud of curves in Figure 6, these b o otstrap metho ds p erform w orse and the resulting prediction regions do not alw ays contain the true daily demand curve. 5 10 15 20 5000 15000 25000 35000 Electric demand, Spain − 2012 Hour MWh 5 10 15 20 5000 15000 25000 35000 Electric demand, Spain − 2012 Hour MWh Figure 5: Daily curv es of electricit y demand corresponding to weekda ys (grey solid lines), curve of electricity demand for 2nd April 2012 (red solid line), prediction curve (black solid line) and prediction region (blac k dashed lines) b y means of λ -metho d at 80% of confidence (left panel) and L ∞ -metho d at 95% of confidence (righ t panel) based on the SFPL forecasting model. 5 10 15 20 5000 15000 25000 35000 Electric demand, Spain − 2012 Hour MWh 5 10 15 20 5000 15000 25000 35000 Electric demand, Spain − 2012 Hour MWh Figure 6: Daily curves of electricit y demand corresp onding to Saturda ys (grey solid lines), curve of electricit y demand for 7th Jan uary 2012 (red solid line), prediction curve (black solid line) and prediction region (black dashed lines) b y means of λ -method based on the SFPL forecasting model at 80% of confidence (left panel) and daily curves of electricit y demand corresp onding to Sundays (grey solid lines), curv e of electricity demand for 8th Jan uary 2012 (red solid line), prediction curv e (blac k solid line) and prediction region (black dashed lines) b y means of λ -method based on the SFPL forecasting model at 80% of confidence (right panel). 4.3. Pr e diction r e gions for the ele ctricity pric e daily curves In this section, an analogous study is carried out for the electricity price data. Prediction regions for the electricity price daily curves are computed using the L ∞ -metho d, λ -metho d and depth-based metho d and the forecasting mo dels FNP and SFPL. 21 T able 4 sho ws the functional cov erage, width and functional Winkler score of the predic- tion regions obtained with the prop osed methods, and using the FNP forecasting mo del, for the 2012 electricit y price curves. This table compares the different metho ds for obtaining prediction regions based on the FNP , also taking in to accoun t the fixed confidence level. At 95% of confidence, the L ∞ -metho d and λ -metho d hav e acceptable functional co verage. The width of the regions obtained with the λ -metho d is smaller, whic h results in a low er FWS. Therefore, the λ -metho d seems to b e the b est alternativ e in this scenario. The depth-metho d giv es v ery small prediction region widths, th us decreasing their co verage. Therefore, it is the one that pro vides the w orst FWS. The smaller sample size and the higher v ariabilit y of the curv es cause the results for Sunda ys to b e w orse. Fixed the confidence lev el at 80%, the three metho ds giv e similar results in this scenario. α = 0 . 05 W eekdays Saturdays Sunda ys Y ear L ∞ -metho d FCo v 95.4 90.4 86.8 93.4 A Width 59.8 54.4 63.9 59.6 FWS 92.4 105.1 145.9 102.0 λ -metho d F Cov 93.5 86.5 90.6 92.1 A Width 47.7 47.5 58.4 49.2 FWS 81.2 114.7 110.5 90.2 Depth metho d F Cov 85.1 76.9 75.5 82.5 A Width 41.3 38.5 48.5 41.9 FWS 104.3 133.2 172.9 118.3 α = 0 . 20 W eekdays Saturdays Sunda ys Y ear L ∞ -metho d FCo v 78.2 76.9 64.2 76.0 A Width 33.6 33.7 38.9 34.4 FWS 54.8 55.5 74.6 57.8 λ -metho d F Cov 68.2 73.1 60.4 67.8 A Width 29.9 31.2 38.4 31.4 FWS 55.4 54.7 75.5 58.2 Depth metho d F Cov 70.2 67.3 62.2 68.6 A Width 30.7 31.6 40.9 32.3 FWS 54.0 60.3 82.1 59.0 T able 4: F unctional cov erage (in %), width and functional Winkler score of the prediction regions for the daily electricity price based on the FNP mo del using the L ∞ -metho d, the λ -metho d and the depth-based metho d with α = 0 . 05 for eac h kind of day in 2012. Figure 7 sho ws the forecast and prediction region of electricity price on a weekda y using the λ -metho d at the 80% confidence level and the L ∞ -metho d at the 95% confidence lev el based on the FNP mo del. The volatilit y of the daily electricit y price curves is considerable, as can b e seen in Figure 7. This mak es the p erformance of the forecasting mo dels and the metho ds for computing prediction regions worse. At a confidence level of 95%, the metho ds yield prediction regions of large width and, therefore, not v ery informative. This is emphasised in the case of the L ∞ -metho d when providing prediction regions of constant width. 22 5 10 15 20 0 20 40 60 80 100 Electric price, Spain − 2012 Hour Cents/Kwh 5 10 15 20 0 20 40 60 80 100 Electric price, Spain − 2012 Hour Cents/Kwh Figure 7: Daily curves of electricity price corresp onding to w eekdays (grey solid lines), curve of electricit y demand for 2nd April 2012 (red solid line), prediction curve (black solid line) and prediction region (blac k dashed lines) b y means of λ -metho d at 80% of confidence (left panel) and L ∞ -metho d at 95% of confidence (righ t panel) based on the FNP forecasting model. The metho ds for obtaining prediction regions for the electricit y price daily curves based on the SFPL regression mo del are discussed b elo w. The v ector of scalar cov ariates included in the SFPL model to forecast the daily curv es of electricit y price is x = ( x 1 , x 2 ) T = (D , PP) T where D is the forecasted daily demand and PP is the wind p o w er pro duction. Other co v ariates that ha v e a linear relationship with the resp onse v ariable can be incorp orated in to the model. T able 5 shows the functional co v erage, width and Winkler score of the prediction regions obtained with the prop osed metho ds, and using the SFPL forecasting mo del, for the 2012 electricit y demand curves at confidence level 80% and 95% confidence levels. The main conclusions drawn ab o ve for the prediction region metho ds based on the FNP mo del can b e extrap olated to the metho ds based on the SFPL model, although the differences b et w een the L ∞ -metho d and the λ -metho d based on FNP mo del are attenuated when using SFPL mo del. The results for Saturdays and Sundays are still w orse than for weekda ys. A ccording to the FWS v alues obtained, it seems that the results are b etter with the SFPL mo del than with the FNP mo del, esp ecially at 80% confidence. The p erformance of the depth metho d based on the SFPL mo del is worse than the other t wo metho ds, although it is less w orse than when using the FNP mo del. 23 α = 0 . 05 W eekdays Saturdays Sunda ys Y ear L ∞ -metho d Co v erage 93.5 84.6 69.8 88.8 Width 49.5 40.6 41.9 47.1 FWS 88.8 101.8 174.4 103.1 λ -metho d Co verage 88.1 84.6 60.4 83.6 Width 40.3 34.2 37.5 39.1 FWS 93.2 83.7 190.7 105.9 Depth metho d Cov erage 78.9 61.5 45.3 71.6 Width 34.2 29.5 34.5 33.6 FWS 105.9 138.8 232.2 128.9 α = 0 . 20 W eekdays Saturdays Sunda ys Y ear L ∞ -metho d Co v erage 73.2 63.5 45.3 67.8 Width 27.9 26.5 29.2 27.9 FWS 49.1 52.6 68.7 52.5 λ -metho d Co verage 63.6 59.6 39.6 59.6 Width 25.6 23.2 26.9 25.4 FWS 50.6 45.3 69.3 52.6 Depth metho d Cov erage 59.0 50.000 39.6 54.9 Width 25.8 24.5 28.5 26.0 FWS 54.3 54.3 74.9 57.3 T able 5: F unctional cov erage (in %), width and functional Winkler score of the prediction regions for the daily electricit y price based on the SFPL model using the L ∞ -metho d, the λ -metho d and the P-based method with α = 0 . 05 for eac h kind of day in 2012. Figure 8 shows the forecast and prediction region of electricity price on a Saturday and a Sunda y using the λ -method based on the SFPL model at the 80% confidence level. As illustrated by Figure 6, the greater volatilit y of the daily electricity price curv es is increased on Saturdays and Sunda ys and the b o otstrap metho ds p erform w orse. 24 5 10 15 20 0 20 40 60 80 100 Electric price, Spain − 2012 Hour Cents/Kwh 5 10 15 20 0 20 40 60 80 100 Electric price, Spain − 2012 Hour Cents/Kwh Figure 8: Daily curves of electricity price corresponding to Saturda ys (grey solid lines), curv e of electricit y demand for 7th Jan uary 2012 (red solid line), prediction curve (black solid line) and prediction region (black dashed lines) b y means of λ -method based on the SFPL forecasting model at 80% of confidence (left panel) and daily curves of electricit y demand corresp onding to Sundays (grey solid lines), curv e of electricity demand for 8th Jan uary 2012 (red solid line), prediction curv e (blac k solid line) and prediction region (black dashed lines) b y means of λ -method based on the SFPL forecasting model at 80% of confidence (right panel). 4.4. Computation time An imp ortant asp ect of these techniques to tak e in to accoun t is their computational cost. T able 6 contains the computational time required to obtain the prediction region for the electricity price or demand curv e o ver a day with each of the three b ootstrap metho ds, L ∞ -metho d, λ -metho d and depth-based metho d, based on the t wo considered regression mo dels, FNP and SFPL. The pro cedures for deriving the functional prediction regions do not differ significantly in terms of computational time. The time differences in practice lie in the regression mo del chosen. Although the SFPL mo del allows more information from the data to b e included in the estimation through exogeneous scalar cov ariates and, therefore, pro vides more accurate results in some of the analysed scenarios, it has the disadv antage of b eing computationally more exp ensive than the FNP mo del. Specifically , obtaining the one-da y-ahead prediction of a daily electricity price or demand curve using the FNP mo del requires less time than using the SFPL mo del. Thus, the metho ds for obtaining FNP mo del- based prediction regions are also faster than the SFPL mo del-based metho ds. This results in large differences b etw een b o otstrap metho ds based on one or the other regression mo del, as shown in T able 6. L ∞ -metho d λ -metho d depth metho d FNP 1.8 1.8 2.0 SFPL 145.2 132.0 160.5 T able 6: Computation time (in seconds) of L ∞ -metho d, λ -method and depth-based metho d for the compu- tation of the prediction region for one d ay using FNP and SFPL models. 25 5. Conclusions The ob jectiv e of this pap er is to solv e the functional regression problem where the resp onse is the daily curv e of electricit y demand or price, in order to obtain one-da y-ahead prediction regions for these curves. T o construct the prediction regions, three b ootstrap procedures, L ∞ - metho d, λ -metho d and depth-based metho d, are proposed. They are based on assuming some functional resp onse forecasting mo del and b o otstraping the residuals. T o this exten t, these metho ds are very general, since any functional resp onse regression mo del could b e assumed and an y estimator could b e used to fit the assumed mo del. F urthermore, they are able to capture t wo sources of v ariability , the one due to the estimation of the regression function and the one due to the mo del error. In this work, these three general metho ds are particularised to a functional autoregressive mo del that is nonparametrically estimated by the Nadaray a- W atson functional estimator (FNP) and a partial linear mo del with functional resp onse that includes an autoregressiv e part and linearly introduces scalar cov ariates (SFPL). The three proposed tec hniques hav e been used to compute prediction regions for the daily curv es of b oth electricity demand and price in the Spanish market. Predictions hav e b een obtained for all days of the y ear 2012 using information from the previous 365 days. Due to the differences b et ween the types of da y , weekda y , Saturda y and Sunda y , it w as necessary to fit three differen t mo dels to predict the electricity demand or price for each da y . The database used in this pap er has b een analysed in other works on the electricit y mark et. In particular, in Vilar et al. (2018), they solved the problem of computing p oint wise prediction interv als and their results are compared with those of this pap er. The prediction regions of the daily curv es obtained from the results of the p oin t wise (hourly) prediction in terv als are very narrow and they hav e a very lo w cov erage, clearly m uch lo wer than the one assumed. Therefore, that technique is discarded and direct pro cedures based on mo dels with functional resp onse, as the three prop osals of this pap er, are needed. The results obtained in the study allow us to conclude that, in general, the L ∞ -metho d and the λ -metho d are reasonable options for computing prediction regions for the electricity demand and price curv es in Spain. Both ha ve fairly similar b ehaviour, with minor differences. It is imp ortan t to highlight that the λ -metho d takes in to account the volatilit y of the curves, giving rise to prediction regions of non-constant width, an adv antage that this metho d has o ver the L ∞ -metho d. The depth-based metho d p erformed consistently worse than the other t wo throughout the study . In general, these observ ations remain consisten t when using b oth the FNP regression mo del and the SFPL mo del. If anything, it is noticed that the use of the SFPL mo del attenuates the differences b et ween the metho ds for computing prediction regions. The prop osed pro cedures for deriving the functional prediction regions do not differ sig- nifican tly in terms of computational time. The time differences in practice lie in the chosen regression mo del. Although the SFPL mo del allo ws more information from the data to b e included in the estimation through exogeneous scalar co v ariates and, therefore, pro vides more accurate results in some of the analysed scenarios, it has the disadv an tage of b eing computa- tionally m uch more exp ensiv e than the FNP mo del. It should also b e noted that regression mo dels with functional resp onse where an estimation of the daily curve is directly obtained are faster than regression mo dels with scalar resp onse that ha ve to be fitted 24 times, once for eac h hour or p oin t in the time grid considered. This is another adv an tage of the techniques 26 presen ted here, in addition to the improv emen ts in co verage men tioned ab o v e. A c knowledgemen ts This researc h has b een supp orted by MICINN Gran t PID2020-113578RB-100, by the Xun ta de Galicia (Grup o de Referencia Comp etitiv a ED431C-2020-14 and Centro Singular de Inv estigaci´ on de Galicia ED431G 2019/01), all of them through the ERDF. References Aneiros, G., Horov´ a, I., Huvsk ov´ a, M., Vieu, P ., 2022. On functional data analysis and related topics. Journal of Multiv ariate Analysis 189, 104861. Aneiros, G., Vilar, J.M., Cao, R., Mu˜ noz San-Ro que, A., 2013. F unctional prediction for the residual demand in electricit y spot mark ets. IEEE T ransactions on Po w er Systems 28, 4201–4208. Aneiros, G., Vilar, J.M., Ra ˜ na, P ., 2016. Short-term forecast of daily curves of electricity demand and price. Electrical Po w er and Energy Systems 80, 96–108. Aneiros-P ´ erez, G., Vieu, P ., 2008. Nonparametric time series prediction: A semi-functional partial linear mo deling. Journal of Multiv ariate Analysis 99, 834–857. An to c h, J., Prc hal, L., Rosa, M.R.D., Sarda, P ., 2010. Electricity consumption prediction with functional linear regression using spline estimators. Journal of Applied Statistics 37, 2027–2041. Barrien tos-Mar ´ ın, J., Marulanda, L.M., Duque, F.V., 2023. Analyzing electricit y demand in Colom bia: A functional time series approach. In ternational Journal of Energy Economics and Policy 13, 75–84. Cancelo, J.R., Espasa, A., Grafe, R., 2008. F orecasting the electricit y load from one day to one week ahead for the spanish system op erator. International Journal of F orecasting 24, 588–602. Cao, R., F rancisco-F ern´ andez, M., Quinto, E., 2010. A random effect multiplicativ e het- eroscedastic mo del for bacterial growth. BMC Bioinformatics 11. Chen, Y., Li, B., 2017. An adaptive functional autoregressive forecast mo del to predict electricit y price curves. Journal of Business & Economic Statistics 35, 371–388. Cuesta-Alb ertos, J., Nieto-Reyes, A., 2008. The random tukey depth. Computational Statis- tics and Data Analysis 52, 4979–4988. Delsol, L., 2009. Adv ances on asymptotic normality in non-parametric functional time series analysis. Statistics 43, 13–33. El ´ ıas, A., Jim ´ enez, R., Shang, H.L., 2022. On pro jection metho ds for functional time series forecasting. Journal of Multiv ariate Analysis 189, 104890. 27 F ebrero-Bande, M., Oviedo de la F uente, M., 2012. Statistical computing in functional data analysis: The r pac kage fda.usc. Journal of Statistical Softw are 51, 1–28. F erraty , F., Laksaci, A., T adj, A., Vieu, P ., 2011. Kernel regression with functional resp onse. Electronic Journal of Statistics 5, 159–171. F erraty , F., Laksaci, A., T adj, A., Vieu, P ., 2012. Nonparametric regression for functional resp onse and functional regressor under dep endance. Comptes Rendus Mathematique 350, 717–720. F erraty , F., V an Keilegom, I., Vieu, P ., 2010. On the v alidity of the b o otstrap in non- parametric functional regression. Scandinavian Journal of Statistics 37, 286–306. F erraty , F., Vieu, P ., 2004. Nonparametric mo dels for functional data, with application in regression, time-series prediction and curve discrimination. Journal of Nonparametric Statistics 16, 111–125. F erraty , F., Vieu, P ., 2006. Nonparametric F unctional Data Analysis. Theory and Practice. Springer, New Y ork. Gneiting, T., Balab daoui, F., Raftery , A.E., 2007. Probabilistic F orecasts, Calibration and Sharpness. Journal of the Roy al Statistical So ciet y Series B: Statistical Metho dology 69, 243–268. Hong, T., F an, S., 2016. Probabilistic electric load forecasting: a tutorial review. In ternational Journal of F orecasting 32, 914–938. Hyde, O., Ho dnett, P ., 2015. Mo deling the effect of weather in short-term load forecasting. Mathematical Engineering in Industry 6, 155–169. Ka ya, M., Karan, M.B., T elatar, E., 2023. Electricity price estimation using deep learning approac hes: An empirical study on T urkish markets in normal and Covid-19 p erio ds. Exp ert Systems With Applications 224, 120026. K okoszka, P ., Reimherr, M., 2017. In tro duction to F unctional Data Analysis. Chapman and Hall, London. Liebl, D., 2013. Mo deling and forecasting electricit y sp ot prices: A functional data p ersp ec- tiv e. The Annals of Applied Statistics 7, 1562–1592. Lisi, F., Shah, I., 2020. F orecasting next-day electricity demand and prices based on func- tional mo dels. Energy Systems 11, 947–979. Masry , E., 2005. Nonparametric regression estimation for dependent functional data: asymp- totic normality . Sto c hastic Pro cesses and their Applications 115, 155–177. No wotarski, J., W eron, R., 2018. Recent adv ances in electricit y price forecasting: A review of probabilistic forecasting. Renew able and Sustainable Energy Reviews 81, 1548–1568. 28 P aparo ditis, E., Sapatinas, T., 2013. Short-term load forecasting: The similar shap e func- tional time-series predictor. IEEE T ransactions on Po w er Systems 28, 3818–3825. P ortela-Gonz´ alez, J., Ro que, A.M.S., P ´ erez, E.A., 2018. F orecasting functional time series with a new Hilb ertian ARMAX mo del: Application to electricity price forecasting. IEEE T ransactions on Po w er Systems 33, 545–556. Ra ˜ na, P ., Aneiros, G., Vilar, J.M., 2015. Detection of outliers in functional time series. En vironmetrics 26, 178–191. Ramsa y , J., Silv erman, B., 2005. F unctional Data Analysis. Springer, New Y ork. Ra ˜ na, P ., Aneiros, G., Vilar, J., Vieu, P ., 2016. Bo otstrap confidence interv als in functional nonparametric regression under dep endence. Electronic Journal of Statistics 10, 1973–1999. Ra ˜ na, P ., Vilar, J., Aneiros, G., 2018. On the use of functional additive mo dels for electricit y demand and price prediction. IEEE Acces 6, 9603–9613. Sugan thi, L., Sam uel, A., 2012. Energy mo dels for demand forecasting - A revie w. Renewable and Sustainable Energy Reviews 16, 1223–1240. T aylor, J.W., Buizza, R., 2003. Using weather ensamble predictions in electricit y demand forecasting. In ternational Journal of F orecasting 19, 57–70. T aylor, J.W., de Menezes, L., McSharry , P .E., 2006. A comparison of univ ariate metho ds for forecasting electricity demand up to a day ahead. International Journal of F orecasting 22, 1–16. Vilar, J.M., Cao, R., Aneiros, G., 2012. F orecasting next-day electricity demand and price using nonparametric functional metho ds. Electrical P ow er and Energy Systems 39, 48–55. Vilar, J.M., Ra ˜ na, P ., Aneiros, G., 2016. Using robust fp ca to identify outliers in functional time series, with applications to the electricity mark et. SOR T 40, 321–348. Vilar, J.M., Ra ˜ na, P ., Aneiros, G ., 2018. Prediction in terv als for electricit y demand and price using functional data. Electrical Po w er and Energy Systems 96, 457–492. W eron, R., 2006. Mo deling and F orecasting Electricit y Loads and Prices: A Statistical Approac h. John Wiley and Sons, New Y ork. W eron, R., 2014. Electricity price forecasting: A review of the state-of-the-art with a lo ok in to the future. International Journal of F orecasting 30, 1030–1081. Winkler, R.L., 1972. A decision-theoretic approach to in terv al estimation. Journal of the American Statistical Asso ciation 67, 187–191. Zhang, J.T., 2014. Analysis of V ariance for F unctional Data. Chapman and Hall, London. Zh u, T., Politis, D., 2017. Kernel estimates of nonparametric functional autoregression mo dels and their b o otstrap approximation. Electronic Journal of Statistics 11, 2876–2906. 29