26.03

2015

Produkcja energii elektrycznej w Polsce — analizy i prognozy

Autor: adam

Zapraszam do lektury kolejnego wpisu poświęconego praktycznej analizie szeregów czasowych. Dzisiaj przedstawiamy przygotowane przez Hannę Loch studium przypadku na bazie danych dotyczących produkcji energii elektrycznej w Polsce (dla ostatnich kilku lat).

Przeprowadzona analiza pozwoliła odpowiedzieć na szereg ważnych pytań, np. jak kształtowało się zapotrzebowanie na energię elektryczną w ostatnich latach i czego możemy spodziewać się w kolejnych okresach. Szczegółowo porównano także skuteczność różnych metod prognozowania, z uwzględnieniem krótszego i dłuższego horyzontu. Oprócz bardziej standardowych metod (takich jak ARIMA, wygładzanie wykładnicze czy modele dekompozycji) do konstrukcji prognoz zostały zastosowane również sieci neuronowe, które dla tych danych okazały się jedną z najbardziej efektywnych metod prognozowania.

Na koniec zaznaczmy jeszcze, że poniższy artykuł (z założenia) nie jest tutorialem. który ma nauczyć od podstaw stosowania metod analizy i prognozowania szeregów czasowych. Osoby zainteresowane poznaniem podstaw metodologicznych, szczegółów dotyczących określonych metod analizy szeregów czasowych czy też konkretnych algorytmów dostępnych w pakiecie R odsyłamy do naszej książki: Analiza i prognozowanie szeregów czasowych. Praktyczne wprowadzenie na podstawie środowiska R, A. Zagdański, A. Suchwałko, PWN 2015.

1 Wprowadzenie

W poniższym opracowaniu przedstawione zostało studium przypadku (case study) w celu zaprezentowania wybranych metod analizy i prognozowania szeregów czasowych. Wybrane zostały metody i modele najbardziej odpowiednie dla danego szeregu; omówione zostaną m.in. proste metody dekompozycji, pozwalające na wyodrębnienie w szeregu regularnych składowych (takich jak trend długoterminowy i sezonowość), modele ARIMA, STL (Seasonal Decomposition of Time Series by Loess), ETS (ExponenTial Smoothing) oraz modelowanie przy użyciu sieci neuronowych.

Wykorzystane dane dotyczą miesięcznej całkowitej produkcji energii elektrycznej (brutto) w gigawatogodzinach w Polsce w okresie od stycznia 2008 do czerwca 2014. Są one dostępne na stronie Eurostatu. Wszystkie analizy zostały przeprowadzone w środowisku R.

1.1 Pytania biznesowe i analiza danych

Analiza danych dotyczących produkcji elektrycznej może nam powiedzieć sporo o tym, jak w przeszłości kształtowało się zapotrzebowanie na energię, jakie były długofalowe trendy oraz okresowe wahania. Są to niezwykle ważne informacje dla osób i przedsiębiorstw związanych z branżą energetyczną; zarówno producentów i dostawców energii, jak i na przykład dostawców węgla do elektrowni węglowych.

Mając już pewną wiedzę odnośnie historycznego przebiegu takiego szeregu czasowego, jakim jest produkcja energii, możemy się pokusić o próbę odpowiedzi na pytanie, jak będzie wyglądało przyszłe zapotrzebowanie na energię elektryczną. Ta informacja jest potrzebna między innymi wcześniej wymienionym producentom i dostawcom, ale też może być używana w wielu innych branżach, chociażby do zestawień ekonomiczno-gospodarczych (może być to pewnym wskaźnikiem rozwoju). Oczywiście analizy analogiczne do przedstawionych w tym opracowaniu można przeprowadzić dla danych związanych z produkcją innych surowców, na przykład gazu czy ropy naftowej.

Poniżej przedstawimy analizę danych wraz z prognozami, które mogłyby być zastosowane w praktyce.

1.2 Podstawowe własności danych

Analizę szeregu czasowego najlepiej jest zacząć od przyjrzenia się wykresowi danych, które będziemy analizować. Wykres naszych danych dotyczących produkcji prądu w Polsce został przedstawiony na rysunku 1.

Rysunek 1. Całkowita produkcja energii elektrycznej w Polsce w okresie 01.2008 - 07.2014

Rysunek 1. Całkowita produkcja energii elektrycznej w Polsce w okresie 01.2008 – 07.2014

Jak łatwo zauważyć na wspomnianym wykresie, nasze dane ulegają pewnym sezonowym zmianom — co roku w okresie zimowym produkcja jest zdecydowanie większa niż latem. Jest to z pewnością związane po pierwsze z krótszymi dniami i koniecznością dłuższego używania oświetlenia (co w pomieszczeniach biurowych czy przemysłowych skutkuje znaczącym zwiększeniem zużycia, a co za tym idzie, produkcji), a po drugie korzystaniem w tym czasie z elektrycznego ogrzewana bądź zasilania kotłów grzewczych.

Wzorzec sezonowy możemy dokładniej zidentyfikować patrząc na tak zwany wykres sezonowy, przedstawiony na rysunku 2. Na jego podstawie możemy jednoznacznie stwierdzić występowanie wyraźnej sezonowości w danych — ze wzrostem w zimowych miesiącach i spadkiem w letnich.

Rysunek 2. Całkowita produkcja energii elektrycznej - wykres sezonowy

Rysunek 2. Całkowita produkcja energii elektrycznej – wykres sezonowy

By zobaczyć, jak na przestrzeni lat zmieniała się produkcja energii elektrycznej w poszczególnych miesiącach, możemy się również przyjrzeć wykresowi miesięcznemu, przedstawionemu na rysunku 3. Czarne linie na tym wykresie przedstawiają średnie wartości produkcji prądu w danym miesiącu (średnia jest brana ze wszystkich lat). Jak widać, w przypadku większości miesięcy mieliśmy do czynienia z podobnymi zmianami — spadkiem w 2009 roku, a następnie wzrostem.

Rysunek 3. Całkowita produkcja energii elektrycznej - wykres miesięczny

Rysunek 3. Całkowita produkcja energii elektrycznej – wykres miesięczny

Patrząc również na wykres 1 nie jesteśmy w stanie jednoznacznie określić, czy mamy do czynienia z pewnym ogólnym trendem wzrostowym (lub spadkowym) czy niekoniecznie. W dalszej części dokładniej zbadamy obecność długoterminowego trendu w danych.

1.3 Inne własności danych

Dane związane z wolumenem produkcji energii elektrycznej są danymi miesięcznymi — jest to sumaryczne zużycie w danym miesiącu. Gdy spojrzymy na wykres zaprezentowany na rysunku 2, widzimy, że w każdym roku następował pewien spadek w lutym. Zwróćmy uwagę, że może być to spowodowane faktem, iż w lutym mamy mniejszą liczbę dni i przez to sumaryczna produkcja energii elektrycznej może być mniejsza. Należałoby zatem dostosować nasze dane do zmieniającej się długości miesięcy (przez podzielenie wartości przez liczbę dni w odpowiednim miesiącu, a następnie pomnożenie przez średnią liczbę dni). Zastosowanie takiej korekty efektów kalendarzowych (month-length adjustment) może uchronić nas przed błędnym zinterpretowaniem jako wahań sezonowych tych spadków lub wzrostów, które są związane wyłącznie z długością danych miesięcy.

Na rysunku 4 przedstawione są wykresy przed i po zastosowaniu korekty uwzględniającej różną długość miesięcy. Jak można zauważyć, krzywe na wykresie po dostosowaniu są gładsze, pozbawione skoków — na przykład nie ma gwałtownego spadku w lutym, który był wcześniej obecny. Spodziewamy się, że omawiana korekta ułatwi nam proces modelowania, dlatego w dalszej analizie będziemy wykorzystywali już tylko dane po jej zastosowaniu.

Rysunek 4. Całkowita produkcja energii elektrycznej przed i po zastosowaniu korekty kalendarzowej

Rysunek 4. Całkowita produkcja energii elektrycznej przed i po zastosowaniu korekty kalendarzowej

1.4 Identyfikacja regularnych wzorców w danych

Omawiając wzorzec wahań sezonowych, widoczny w naszych danych, wspomnieliśmy o wyodrębnieniu z szeregu jego regularnych składowych. Takie podejście, nazywane dekompozycją szeregu, zakłada, iż szereg czasowy składa się z trzech elementów: pewnego trendu (o którym też wspomnieliśmy — ogólnej tendencji, na przykład wzrostowej), sezonowości oraz losowych fluktuacji.

Model dekompozycji możemy zapisać symbolicznie: $$y_t = T_t + S_t + E_t \ \mbox{(dekompozycja addytywna)}$$ bądź $$y_t = T_t \times S_t \times E_t \ \mbox{(dekompozycja multiplikatywna)}$$ gdzie $T_t$ – trend, $S_t$ – sezonowość, $E_t$ – losowe fluktuacje.

Takie podejście jest rzeczywiście bardzo intuicyjne — prognoza na dany miesiąc jest sumą (bądź iloczynem) prognozy trendu, prognozy sezonowości i prognozy losowych fluktuacji (o których najczęściej zakłada się, iż mają średnią 0). W przypadku naszych danych dotyczących produkcji prądu, składnik trendu odpowiada długofalowym tendencjom — przeważnie zwiększaniu produkcji (na przykład ze względu na coraz większą liczbę używanych urządzeń elektrycznych), składnik sezonowy — zimowym wzrostom i letnim spadkom produkcji, natomiast fluktuacje losowe — wszystkim innym zdarzeniom (na przykład otwarciu kilku nowych zakładów produkcyjnych w Polsce). Jeżeli chodzi o wybór wariantu dekompozycji (tzn. dekompozycję addytywną bądź multiplikatywną) jest to zależne od charakteru danych. W praktyce, warto sprawdzić oba warianty i wybrać ten lepiej dopasowany.

Przechodząc już do wyodrębnienia poszczególnych składowych, zaczniemy od składowej trendu. Możemy próbować odgadnąć postać trendu z wykresu szeregu czasowego (rysunek 1), aczkolwiek ta metoda sprawdzi się tylko w wypadku prostych i bardzo wyraźnych trendów. Bardziej odpowiednim podejściem, które możemy wykorzystać do zidentyfikowania składowej trendu jest metoda średniej kroczącej (moving average), czyli pewnego rodzaju wygładzanie szeregu.

Średnią kroczącą w chwili $t$ określamy jako średnią z pewnej liczby sąsiednich obserwacji (liczba tych obserwacji to tzw. rząd średniej kroczącej). Dokładniej, średnią kroczącą rzędu $m$ możemy zapisać wzorem: $$\hat{T_t} = \frac{1}{m} \sum_{i=-k}^{k} y_{t+i},$$ gdzie $m=2k+1$.

Jak widać w powyższym wzorze, do średniej bierzemy $k$ poprzedzających i $k$ następujących obserwacji, dlatego nie możemy określić wartości średniej kroczącej dla wszystkich momentów czasowych (precyzyjniej, nie możemy jej określić dla kilku najbardziej skrajnych miesięcy). Oczywiście, jeżeli z jakiś przyczyn musielibyśmy uzyskać te wartości, można byłoby np. zastosować niesymetryczną bądź jednostronną średnią kroczącą.

Na rysunku 5 zaprezentowane zostało dopasowanie średniej kroczącej rzędu 5, 9 i 12. Jak widać na wspomnianym wykresie, średnia krocząca rzędu 12 najlepiej oddaje ogólne tendencje i nie ulega wpływom sezonowym. W przypadku sezonowych szeregów czasowych najczęściej rząd użytej średniej kroczącej jest równy okresowi — to właśnie zapobiega obecności wahań sezonowych w oszacowanym trendzie.

Rysunek 5. Dopasowanie trendu przy użyciu średniej kroczącej

Rysunek 5. Dopasowanie trendu przy użyciu średniej kroczącej

Ponieważ do uzyskania składowej trendu używamy dwustronnej średniej kroczącej, nie jesteśmy w stanie prognozować (ekstrapolować) przy jej użyciu trendu w dalszej perspektywie. Jednak otrzymana krzywa może nam podpowiedzieć jakiś prosty, parametryczny model trendu (na przykład linowy). Nawet w przypadku naszych danych możemy spróbować przyjąć liniowy model trendu (i zobaczyć, jak sprawdzą się prognozy wykonane przy takim założeniu) — czasami dobrze jest pozwolić na trochę gorsze, mniej dokładne prognozy, zyskując w zamian na prostocie modelu. Jest to szczególne ważne w przypadku, gdy prognozy musimy konstruować wielokrotnie i w krótkim czasie.

Po wyznaczeniu trendu możemy wyodrębnić składową sezonową. W tym celu wyznaczamy tzw. indeksy sezonowe — wartości przypisane każdemu z miesięcy, informujące o tym, o ile obserwacja z danego miesiąca różni się od średniej. By je otrzymać, należy z danych po usunięciu (wyeliminowaniu) trendu wyznaczyć średnią z wartości dla konkretnego miesiąca we wszystkich latach. Na przykład indeks sezonowy dla miesiąca sierpnia to średnia produkcja prądu (po usunięciu trendu) z tego miesiąca we wszystkich latach. Wykres indeksów sezonowych dla poszczególnych miesięcy przedstawia rysunek 6.

Rysunek 6. Wykres indeksów sezonowych w produkcji energii elektrycznej

Rysunek 6. Wykres indeksów sezonowych w produkcji energii elektrycznej

Po wyodrębnieniu tej składowej możemy zobaczyć, jak wyglądają same losowe fluktuacje (szereg z usuniętym trendem i sezonowością), które zostały przedstawione na rysunku 7.

Rysunek 7. Losowe fluktuacje - produkcja energii elektrycznej w Polsce

Rysunek 7. Losowe fluktuacje – produkcja energii elektrycznej w Polsce

Interpretując wykres losowych fluktuacji możemy podejrzewać, iż widoczne na nim najbardziej skrajne wartości — odpowiadające lutemu i marcowi 2012 — są wywołane jakimiś konkretnymi wydarzeniami. Zwróćmy jednak uwagę, że owe odstające obserwacje (anomalie) nie były szczególnie widoczne na wykresie szeregu czasowego (rysunek 1).

2 Konstrukcja prognoz

Jak już wcześniej wspominaliśmy, umiejętność prognozowania zapotrzebowania na energię elektryczną w kolejnych okresach może być bardzo istotna i często jest to główny powód analizy historycznych danych. Po dokładnej analizie własności szeregu czasowego możemy lepiej zrozumieć samo zjawisko i na tej podstawie skonstruować odpowiednie prognozy, oparte na dostosowanych do danych modelach.

Chcąc zbadać dokładność prognoz, które otrzymamy, wyodrębnimy dwa podzbiory danych, tzw. podzbiór uczący i testowy. Pierwszą część (podzbiór uczący) użyjemy do zbudowania odpowiednich modeli, na podstawie których dokonamy predykcji na następne okresy. Prognozy te następnie porównamy z wartościami z podzbioru testowego, czyli prawdziwymi danymi z okresu prognozowania.

W naszym przypadku podzbiór uczący będą tworzyły dane z okresu styczeń 2008 – czerwiec 2013, a podzbiór uczący z okresu lipiec 2013 – czerwiec 2014 (czyli będziemy prognozować na okres pełnego roku).

2.1 Ocena dokładności prognoz

Chcąc wybrać z otrzymanych prognoz najlepszą, czyli taką, która jest najbardziej zbliżona do danych rzeczywistych, musimy przyjąć pewne kryteria, oceniające popełniony błąd predykcji. W analizie zwykle korzysta się z kilku różnych kryteriów, my ograniczymy się do następujących dwóch:

  1. MAPE – mean absolute percentage error, czyli średni bezwzględny błąd procentowy. Wartość MAPE mówi nam, o ile procent średnio różniła się wartość prognozowana od rzeczywistej.
  2. MASE – mean absolute scaled error, czyli średni bezwzględny błąd skalowany. Jest to porównywanie otrzymanego błędu z błędem, który popełnilibyśmy stosując dla podzbioru uczącego naiwną metodę prognozowania.

Dla każdej z przygotowanych prognoz będziemy podawać wartości obu miar; ponadto porównamy błędy wyliczone dla horyzontu 6 i 12 miesięcy, by zaprezentować różnicę w dokładności predykcji w krótszym i dłuższym horyzoncie czasowym. Oczywiście za najlepszą prognozę będziemy uważać tę obarczoną najmniejszym błędem.

2.2 Przedziały predykcyjne

Przewidywane przez nas wartości to tzw. prognozy punktowe, czyli konkretna wartość, jaką prawdopodobnie przyjmie analizowany szereg czasowy w danym momencie. Oczywiście trudno jest przewidzieć na przykład zapotrzebowanie na energię elektryczną z dokładnością do 1kWh, dlatego często interesuje nas nie konkretna wartość, a pewien przedział wartości, w którym z pewnym (odpowiednio dużym) prawdopodobieństwem znajduje się ta rzeczywista. Te przedziały nazywamy przedziałami predykcyjnymi.

Warto zaznaczyć, że długości przedziałów predykcyjnych różnią się, w zależności od prawdopodobieństwa (poziomu ufności) z nimi związanego. Jeżeli bowiem chcemy mieć 95% pewność, że rzeczywista wartość znajdzie się w naszym przedziale, to naturalnie musi być on szerszy niż przedział, co do którego chcemy mieć tylko 70% pewność (musimy uwzględnić większe ryzyko). W praktyce najczęściej stosowane są dwa poziomy pewności: 95% i 80%. Takie też będą prezentowane w naszej pracy.

Przedziały predykcyjne pełnią również inną ważną funkcję — dzięki nim możliwa jest ocena ryzyka towarzyszącego prognozom punktowym, co sprawia, że są uzupełnieniem kryteriów dokładności prognoz, takich jak MAPE i MASE. Dodatkowo mając wyznaczone przedziały predykcyjne możemy przeanalizować różne scenariusze dotyczące przyszłych wartości i zaplanować różne strategie biznesowe.

2.3 Proste metody prognozowania

Zaczniemy od najprostszej metody prognozowania — tak zwanej metody naiwnej. Jest to bardzo prosty sposób predykcji: prognozowane dla kolejnych okresów wartości są równe ostatniej z zaobserwowanych. Zatem w naszym wypadku na każdy z miesięcy w okresie lipiec 2013 – czerwiec 2014 przewidywalibyśmy produkcję prądu równą tej z czerwca 2013, czyli 12291 gigawatogodzin. Prognoza tą metodą (wraz z uwzględnieniem 95% i 80% przedziałów ufności) jest przedstawiona na wykresie 8.

Rysunek 8. Prognoza przy użyciu metody naiwnej

Rysunek 8. Prognoza przy użyciu metody naiwnej

Jednak już na samym początku analizy położyliśmy duży nacisk na to, iż w naszych danych występują wyraźne wahania sezonowe. Nieuwzględnienie tej zmienności powoduje dużą niedokładność metody naiwnej (którą można zauważyć zarówno na wykresie 8, jak i w tabeli 1). Zauważmy, że w czerwcu wartość produkcji prądu jest jedną z najniższych w całym roku, zatem prognoza jest bardzo zaniżona w stosunku do rzeczywistych wartości.

Jedną z możliwości uzyskania lepszej prognozy jest użycie sezonowej odmiany metody naiwnej. W tym wypadku zamiast brać ostatnią zaobserwowaną wartość, do prognozowania używamy ostatnią wartość, ale dla danego miesiąca, np. prognozowana produkcja energii elektrycznej w lipcu 2013 będzie równa faktycznej produkcji w lipcu 2012. Prognozy skonstruowane na bazie tej metody zostały przedstawione na rysunku 9, zaś błędy predykcji umieszczone są w tabeli 1.

Rysunek 9. Prognoza przy użyciu sezonowej metody naiwnej

Rysunek 9. Prognoza przy użyciu sezonowej metody naiwnej

Jak widać na rysunku 9, prognozy otrzymane na bazie sezonowej metody naiwnej są o wiele bardziej zbliżone do rzeczywistych wartości. Również błędy predykcji są znacznie mniejsze w porównaniu do metody naiwnej.

Trzecią, dość prostą metodą, którą możemy próbować wykorzystać w tym przypadku (choć nie uwzględnia ona sezonowości, przez co może być gorsza niż omówiona przed chwilą sezonowej metody naiwnej), jest metoda oparta na średniej. W tym przypadku jako prognozę dla przyszłych miesięcy przyjmujemy wartość równą średniej z podzbioru uczącego. Bazując na wartości średniej unikamy sytuacji, gdy prognoza na cały rok jest równa najmniejszej wartości (tak, jak to było w przypadku metody naiwnej), ale nie uwzględniamy wahań sezonowych. Prognozy skonstruowane na bazie średniej są przedstawione na rysunku 10, zaś błędy MAPE i MASE zostały zaprezentowane w tabeli 1.

Rysunek 10. Prognoza przy użyciu średniej

Rysunek 10. Prognoza przy użyciu średniej

Tabela 1. Błędy predykcji dla prostych metod prognozowania
h = 6 | MAPE h = 6 | MASE h = 12 | MAPE h = 12 | MASE
Metoda naiwna 8.14 2.03 7.34 2.20
Sezonowa metoda naiwna 3.29 0.79 4.24 1.21
Metoda średniej 5.59 1.37 5.94 1.72

2.4 Konstrukcje prognoz na bazie dekompozycji

Jak już wcześniej wspomnieliśmy (przy okazji omawiania regularnych wzorców występujących w danych), możemy założyć, iż nasz szereg składa się z trzech składowych (trendu, sezonowości i losowych fluktuacji) i budować prognozy wykorzystując takie przedstawienie danych.

Na początku spróbujemy dopasować model, w którym jako trend będziemy przyjmować funkcję wielomianową (na przykład liniową czy kwadratową), natomiast sezonowość będziemy uwzględniać wykorzystując wspomniane indeksy sezonowe. Na wykresach 11, 1213 zaprezentowane zostały odpowiednio wyniki predykcji dla dekompozycji bez trendu (z samą sezonowością), z trendem będącym wielomianem 1-ego stopnia (liniowym) oraz trendem będącym wielomianem 3-ego stopnia (wielomian 2 -ego stopnia w tym przypadku daje bardzo zbliżone wyniki do wielomianu stopnia 1). W tabeli 2 przedstawione zostały wartości błędów MAPE i MASE dla tych prognoz.

Rysunek 11. Prognoza na podstawie dekompozycji - bez składowej trendu

Rysunek 11. Prognoza na podstawie dekompozycji – bez składowej trendu

Rysunek 12. Prognoza na podstawie dekompozycji z trendem liniowym

Rysunek 12. Prognoza na podstawie dekompozycji z trendem liniowym

Rysunek 13. Prognoza na podstawie dekompozycji z trendem wielomianowym 3. stopnia

Rysunek 13. Prognoza na podstawie dekompozycji z trendem wielomianowym 3. stopnia

Jak widać w tabeli 2, w długim horyzoncie czasowym (12 miesięcy) najlepszy wynik daje uwzględnienie tylko sezonowości (lepszy niż poprzednio omówione metody naiwne). Z kolei dekompozycja z trendem liniowym w okresie 6 miesięcy również daje bardzo dobre wyniki (MAPE równe 1.24); można to zauważyć na wykresie 12, gdzie prognoza punktowa w ciągu pierwszych kilku miesięcy niemal pokrywa się z wartościami rzeczywistymi.

Tabela 2. Błędy predykcji dla dekompozycji
h = 6 | MAPE h = 6 | MASE h = 12 | MAPE h = 12 | MASE
Dekompozycja bez trendu 4.25 1.01 3.02 0.87
Dekompozycja z trendem liniowym 1.24 0.30 4.05 1.16
Dekompozycja z trendem wielomianowym 4.02 0.97 3.74 1.05

Omawialiśmy teraz proste modele oparte na dekompozycji; poniżej przedstawiamy bardziej zaawansowany model oparty na tym podejściu, czyli STL.

2.5 STL

Pełna nazwa tej metody to Seasonal Decomposition of Time Series by Loess, czyli sezonowa dekompozycja szeregu czasowego, dokonywana przy pomocy algorytmu lokalnej regresji (Loess). Wyznaczanie poszczególnych składowych odbywa się w bardziej skomplikowany sposób (niż średnia krocząca i indeksy sezonowe) — dobieramy model na podstawie nieliniowych zależności pomiędzy obserwacjami, osobno dla trendu i sezonowości.

STL pozwala na dostosowanie wielu parametrów, między innymi liczby uwzględnianych w modelu obserwacji (może być ona różna dla trendu i sezonowości). Z jednej strony jest to dużą zaletą tej metody — algorytm STL może być skuteczną metodą dla danych z nieregularną (na przykład 5 miesięczną) sezonowością. Ponadto wzorzec wahań sezonowych może zmieniać się w czasie (nie zakładamy, że sezonowość jest stała w czasie, tak jak w przypadku indeksów sezonowych). Jest i druga strona medalu — względnie duża elastyczność algorytmu STL powoduje, że nie jest to metoda łatwa do zastosowania i wymaga od analityka odpowiedniego doświadczenia.

Rysunek 14 przedstawia prognozę wykonaną za pomocą algorytmu STL. Parametry zostały tak dobrane, by zminimalizować błąd na podzbiorze uczącym (choć warto pamiętać, że czasami nadmiernie dopasowany model daje złe prognozy). W tabeli 3 zostały zaprezentowane wartości MAPE i MASE dla obu horyzontów czasowych. Jak widać prognozy są lepsze niż w przypadku metod naiwnych oraz zbliżone do najlepszych prognoz skonstruowanych na bazie klasycznej dekompozycji, odpowiednio dla krótkiego i dłuższego horyzontu (tzn. dla horyzontu h=6 metoda STL daje podobną dokładność jak dekompozycja z trendem liniowym, a dla h=12 podobną jak dekompozycja bez trendu długoterminowego).

Rysunek 14. Prognoza na podstawie algorytmu STL

Rysunek 14. Prognoza na podstawie algorytmu STL

Tabela 3. Błędy predykcji dla metody STL
h = 6 | MAPE h = 6 | MASE h = 12 | MAPE h = 12 | MASE
STL 1.73 0.42 3.30 0.95

2.6 Wygładzanie wykładnicze

Innym algorytmem, który możemy wykorzystać do konstrukcji prognoz, jest wygładzanie wykładnicze. Polega ono na przedstawieniu każdej obserwacji jako ważonej średniej poprzednich obserwacji; dokładniej rzecz ujmując, wagi przy poszczególnych obserwacjach maleją wykładniczo — im starsza obserwacja, tym mniejsza jej waga.

Istnieją różne algorytmy wygładzania wykładniczego, dedykowane różnym typom danych. Dla danych bez trendu i sezonowości najbardziej odpowiedni jest Simple Exponential Smoothing (SES), dla danych tylko z trendem algorytm Holta, zaś dla danych z trendem i sezonowością algorytm Holta-Wintersa. Jak już ustaliliśmy wcześniej, nasze dane zawierają zarówno trend, jak i sezonowość, więc powinniśmy skorzystać z ostatniej metody; poniżej zaprezentujemy jednak również prognozy za pomocą SES, by pokazać, jak ważny jest wybór odpowiedniego wariantu algorytmu.

Na rysunku 15 zaprezentowana została prognoza za pomocą SES; widać, że jest ona zbliżona do prognozy metodą naiwną. Natomiast na rysunku 16 pokazana jest prognoza uzyskana za pomocą addytywnej wersji algorytmu Holta-Wintersa. W tabeli 4 zaprezentowane są błędy obu metod. Jak widać algorytm Holta-Wintersa daje znacznie lepsze wyniki (analogiczne do metody STL, choć lepsze dla krótszego horyzontu czasowego).

Rysunek 15. Prognoza na bazie algorytmu SES

Rysunek 15. Prognoza na bazie algorytmu SES

Rysunek 16. Prognoza za pomocą algorytmu Holta-Wintersa

Rysunek 16. Prognoza za pomocą algorytmu Holta-Wintersa

Tabela 5. Błędy predykcji dla wygładzania wykładniczego
h = 6 | MAPE h = 6 | MASE h = 12 | MAPE h = 12 | MASE
SES 8.14 2.03 7.34 2.20
Algorytm Holta-Wintersa 1.52 0.37 3.30 0.95

2.7 ETS

Uogólnieniem wygładzania wykładniczego są modele ETS, czyli Error Trend Seasonality (nazywane też ExponenTial Smoothing). Jest to przykład tzw. modelu przestrzeni stanów (state space) — opisujemy stan, w jakim znajduje się pewien układ. Model ETS składa się z równania różnicowego opisującego ewolucję oraz założenia o rozkładzie czynnika losowego. Rozważane wcześniej algorytmy SES i Holta-Wintersa są także jednymi z wariantów modeli ETS.

Na wykresie 17 przedstawione są prognozy, wyznaczone za pomocą modelu ETS, natomiast w tabeli 6 odpowiadające im błędy predykcji. Jak widać, są to jedne z najlepszych prognoz, jakie do tej pory skonstruowaliśmy.

Warto wspomnieć, że dobrana automatycznie do naszych danych postać modelu to ETS(M,N,A) — pierwsza litera oznacza multiplikatywne uwzględnienie losowych fluktuacji, druga brak składowych odpowiedzialnych za trend, natomiast trzecia wskazuje na addytywne uwzględnienie sezonowości. Uwzględniając fakt, iż spośród prostych metod dekompozycji najlepszy wynik dało dopasowanie samej sezonowości (bez trendu), możemy pokusić się o stwierdzenie, iż ekstrapolacja trendu na kolejne okresy (na podstawie modelu dopasowanego na zbiorze uczącym) nie daje zadowalających wyników.

Rysunek 17. Prognoza na bazie modelu ETS

Rysunek 17. Prognoza na bazie modelu ETS

Tabela 6. Błędy predykcji dla metody ETS
h = 6 | MAPE h = 6 | MASE h = 12 | MAPE h = 12 | MASE
ETS 1.53 0.37 3.17 0.91

2.8 Prognozowanie na bazie modelu ARIMA

Następnym modelem, który chcemy przedstawić, jest model ARIMA, czyli AutoRegressive Integrated Moving Average. Jego ideą jest przedstawienie wartości szeregu w danej chwili jako kombinacji liniowej wcześniejszych wartości szeregu (część autoregresyjna) oraz fluktuacji losowych w danym i poprzednich momentach czasowych (część związana ze średnią ruchomą). Ponadto, aby dopasować model konieczne może okazać się zróżnicowanie danych, czyli rozważanie różnic pomiędzy sąsiednimi obserwacjami zamiast oryginalnych obserwacji.

Wybór odpowiedniej postaci modelu, czyli rzędów modelu, nie jest łatwym zadaniem. W wielu przypadkach najlepiej zdać się na automatyczny wybór modelu przez pakiety statystyczne (które wykorzystują do tego pewne kryteria numeryczne). Niestety wadą tego rozwiązania jest to, iż bardzo często metody te znajdują jedynie "lokalnie najlepszy" model (tzn. najlepszy spośród kilku o zbliżonych parametrach), który nie jest najlepszym możliwym.

Podobnie jak w przypadku innych metod, także ARIMA ma swoją sezonową wersję — model sARIMA, który będzie bardziej odpowiedni do naszych danych. Na rysunkach 1819 zostały przedstawione prognozy wykonane za pomocą modelu sARIMA, odpowiednio dopasowanego automatycznie ($sARIMA(2,0,2)(1,0,1)_{[12]}$) i ,,ręcznie” przez analityka ($sARIMA(2,0,2)(0,1,1)_{[12]}$). W tabeli 7 zostały przedstawione błędy predykcji dla obu modeli. W obu przypadkach, przed dopasowaniem modelu została zastosowana transformacja logarytmiczna danych.

Rysunek 18. Prognoza za pomocą automatycznie dobranego modelu ARIMA

Rysunek 18. Prognoza za pomocą automatycznie dobranego modelu ARIMA

Rysunek 19. Prognoza za pomocą ręcznie dobranego modelu ARIMA

Rysunek 19. Prognoza za pomocą ręcznie dobranego modelu ARIMA

Tabela 7. Błędy predykcji dla modeli ARIMA
h = 6 | MAPE h = 6 | MASE h = 12 | MAPE h = 12 | MASE
sARIMA – dopasowanie automatyczne 4.75 1.15 4.18 1.20
sARIMA – dopasowanie ręczne 2.29 0.55 2.75 0.80

Jak widać we wspomnianej tabeli, model dopasowany ,,ręcznie” jest znacznie lepszy niż ten dopasowany automatycznie. Tak jak wspomnieliśmy, czasami metody automatyczne nie są najlepsze; więcej o dobieraniu optymalnych parametrów modelu ARIMA można przeczytać na przykład w książce Forecasting: principles and practice.

2.9 Prognozowanie za pomocą sieci neuronowych

Ostatnią metodą, którą chcemy przedstawić, jest prognozowanie za pomocą sieci neuronowych, które zyskały dużą popularność w ostatnich latach. Pozwalają one bardzo dobrze modelować nieliniowe zależności w danych, co w niektórych przypadkach jest niezbędne.

Sieć neuronowa składa się z warstw neuronów, czyli elementów przetwarzających w odpowiedni sposób sygnały wejściowe. W przypadku szeregów czasowych, każdy neuron wylicza pewną ważoną sumę danych wejściowych, czyli poprzednich wartości szeregu. Naszym celem jest uzyskanie takiego zestawu wag (dla każdego z neuronów i każdej obserwacji), dla których jako rezultat (wynik działania sieci) otrzymamy prognozowaną wartość szeregu. Więcej na ten temat można przeczytać na przykład w książce R.Hyndman’a.

Zaletą sieci neuronowej jest jej duża elastyczność (możliwość modelowania skomplikowanych, nieliniowych zależności) oraz brak założeń dotyczących danych (między innymi nie muszą one pochodzić z rozkładu gaussowskiego). W tym przypadku nie można jednak łatwo konstruować przedziałów predykcyjnych, odpowiadających prognozom punktowym. Z tego powodu pominiemy wyznaczanie tych przedziałów w naszej analizie.

Na rysunku 20 został przedstawiony wykres prognoz, wyznaczonych za pomocą sieci, natomiast w tabeli 8 wartości MAPE i MASE. Dopasowana sieć dla każdego momentu czasowego używa dwóch wcześniejszych obserwacji oraz obserwacji sprzed roku (w ten sposób uwzględnia sezonowość). W drugiej, ukrytej warstwie sieci znajdują się 4 neurony.

Rysunek 20. Prognoza za pomocą modelu sieci neuronowych

Rysunek 20. Prognoza za pomocą modelu sieci neuronowych

Tabela 8. Błędy predykcji dla sieci neuronowych
h = 6 | MAPE h = 6 | MASE h = 12 | MAPE h = 12 | MASE
Sieć neuronowa 1.43 0.34 1.98 0.54

2.10 Która metoda jest najlepsza?

W tabeli 9 zostały przedstawione błędy dla wszystkich zastosowanych powyżej metod konstrukcji prognoz.

Tabela 9. Błędy predykcji – porównanie wszystkich metod
h = 6 | MAPE h = 6 | MASE h = 12 | MAPE h = 12 | MASE
Metoda naiwna 8.14 2.03 7.34 2.20
Sezonowa metoda naiwna 3.29 0.79 4.24 1.21
Metoda średniej 5.59 1.37 5.94 1.72
Dekompozycja bez trendu 4.25 1.01 3.02 0.87
Dekompozycja z trendem liniowym 1.24 0.30 4.05 1.16
Dekompozycja z trendem wielomianowym 4.02 0.97 3.74 1.05
STL 1.73 0.42 3.30 0.95
SES 8.14 2.03 7.34 2.20
Algorytm Holta-Wintersa 1.52 0.37 3.30 0.95
ETS 1.53 0.37 3.17 0.91
sARIMA – dopasowanie automatyczne 4.75 1.15 4.18 1.20
sARIMA – dopasowanie ręczne 2.29 0.55 2.75 0.80
Sieć neuronowa 1.43 0.34 1.98 0.54

Jak widać w dłuższym horyzoncie czasowym (12 miesięcy) najlepsze efekty daje użycie sieci neuronowych; również gdy mówimy o prognozach na okres o połowę krótszy, ta metoda jest jedną z najlepszych (podobne wyniki daje prosta dekompozycja z trendem liniowym). Ważne jednak by nie ulegać złudzeniu, iż prognozowanie za pomocą sieci zawsze daje tak dobre wyniki, niezależnie od danych i horyzontu czasowego.

Zdecydowanie natomiast należy odrzucić proste metody, takie jak metoda naiwna i metoda średniej, które nie uwzględniają sezonowości i przez to nie są dobrze dopasowane do naszych danych. Tak naprawdę każda z bardziej zaawansowanych metod, takich jak ARIMA, wygładzanie wykładnicze, STL czy ETS, o ile wybrany został jej odpowiedni wariant, dała dość dobre wyniki (błąd procentowy MAPE poniżej 3,5%).

Dla większości prognoz (poza tymi na bazie sieci neuronowych) możemy porównać również skonstruowane przedziały predykcyjne. Porównanie, przedstawione w tabeli 10, zostało przygotowane w następujący sposób: spośród wyznaczonych przedziałów predykcyjnych dla $h=1,2,…12$ została wzięta średnia szerokość przedziału (zarówno 80% i 95%), a następnie porównana ze średnim przedziałem dla metody naiwnej (na przykład średnia długość przedziału dla modelu ETS stanowi około 25% długości przedziału dla metody naiwnej).

Tabela 10. Porównanie średnich szerokości przedziałów predykcyjnych względem metody naiwnej
Metoda naiwna
Metoda naiwna 100 %
Sezonowa metoda naiwna 35.9 %
Metoda średniej 76.84 %
Dekompozycja bez trendu 38.17 %
Dekompozycja z trendem liniowym 27.61 %
Dekompozycja z trendem wielomianowym 32.51 %
STL 24.38 %
SES 99.75 %
Algorytm Holta-Wintersa 24.62 %
ETS 24.64 %
sARIMA – dopasowanie automatyczne 31.79 %
sARIMA – dopasowanie ręczne 30.07 %

Jak widać w powyższej tabeli, wszystkie wcześniej wspomniane zaawansowane metody (ARIMA, wygładzanie wykładnicze, STL, ETS) mają zbliżone długości przedziałów predykcyjnych; jest to jeszcze jeden powód, dla którego możemy je zakwalifikować do grupy metod, które okazały się skuteczne dla naszych danych.

Spróbuj ponownie