Modeling hourly (or even more frequent) electricity load is a well-studied problem and many approaches are suggested. Many of these approaches can be handled using the UCM procedure. At this time I am using this forum for illustrating one such approach discussed in Case Study 3 in Chapter 9 of "Time Series Modelling with Unobserved Components," by
Pelagatti, M. M. (2015). Boca Raton, FL: CRC Press. The electricity load data set used in this example is taken from
http://www.ucmbook.info/ , which the website associated with the book.
I am also including a link to one of my workshops on UCM and SSM procedures. This is a link to the slides of a workshop at
"The 36th International Symposium on Forecasting Santander, Spain | 19-22 June 2016":
https://forecasters.org/wp-content/uploads/gravity_forms/7-621289a708af3e7af65a7cd487aee6eb/2016/07/...
See the example "Hourly Electricity Load in Italy" in the Illustrations section for additional explanations.
Note that each PROC UCM call can take up to 15 minutes to finish (processes 9 years of hourly data).
/* ------------------Program-----------------------------------------------*/
libname wk 'C:\udrive\ucmBook\AX2016\programs'; * directory location of the elettroload data set;
data test(drop=h1-h24);
set wk.elettroload;
array hr{24} h1-h24;
do hour=1 to 24;
eload = hr[hour];
dthour = intnx('dthour', dhms(date, 0, 0, 0), hour-1);
output;
end;
format dthour datetime.;
run;
proc sort data=test;
by hour date;
run;
data test(drop=j twoPi yf days year);
set test;
array cosf{16} cos1-cos16;
array sinf{16} sin1-sin16;
year = year(date);
yf = mdy(1, 1, year);
days = date - yf; *day within that year;
twoPi = 2*constant('pi')/365.25;
do j=1 to 16;
cosf[j] = cos(twoPi*j*days);
sinf[j] = sin(twoPi*j*days);
end;
run;
/* Analysis as per Case Study 3 in Chapt 9 of this book. */
proc ucm data=test;
*where hour=1;
id date interval=day;
by hour;
irregular;
level;
cycle period=7 rho=1 noest=(period rho);
cycle period=3.5 rho=1 noest=(period rho);
cycle period=2.3333 rho=1 noest=(period rho);
randomreg cos1-cos16 sin1-sin16;
model eload = dec24 dec25 dec26 jan1 jan6 aug15 easterSat
easterSun easterMon easterTue holidays holySat
holySun bridgeDay endYear;
estimate back=14 plot=panel;
forecast back=14 lead=14 outfor=loadfor;
run;
*An alternate formulation that uses a trigonometric season component with
first 16 harmonics to model the day-of-the-year pattern;
title "Model that uses parsimonious trigonometric season";
proc ucm data=test;
*where hour=1;
id date interval=day;
by hour;
irregular;
level;
cycle period=7 rho=1 noest=(period rho);
cycle period=3.5 rho=1 noest=(period rho);
cycle period=2.3333 rho=1 noest=(period rho);
season length=365 type=trig keeph=1 to 16 by 1;
model eload = dec24 dec25 dec26 jan1 jan6 aug15 easterSat
easterSun easterMon easterTue holidays holySat
holySun bridgeDay endYear;
estimate back=14 plot=panel;
forecast back=14 lead=14 outfor=loadfor1;
run;
proc sgplot data=loadfor1;
title "Day of the year pattern for hour=1 in 2006";
where year(date) = 2006 && hour=1;
series x=date y=s_season;
run;
title;
*Yet another formulation that uses spline season to model
the day-of-the-year pattern;
title "Model that uses spline season";
proc ucm data=test;
*where hour=1;
id date interval=day;
by hour;
irregular;
level;
cycle period=7 rho=1 noest=(period rho);
cycle period=3.5 rho=1 noest=(period rho);
cycle period=2.3333 rho=1 noest=(period rho);
splineseason length=365 degree=3 knots=4 15 30 45 60 75 90 105 120 135 150 165
180 195 210 225 240 255 270 285 300 315
330 345 360;
model eload = dec24 dec25 dec26 jan1 jan6 aug15 easterSat
easterSun easterMon easterTue holidays holySat
holySun bridgeDay endYear;
estimate back=14 plot=panel;
forecast back=14 lead=14 outfor=ucmfor2;
run;
proc sgplot data=ucmfor2;
title "Day-of-the-year pattern for hour=1 in 2006";
where year(date) = 2006 && hour=1;
series x=date y=s_splseas;
run;
