Wybrane pakiety do analizy finansowych szeregów czasowych w R

Jakiś czas temu pisałem o zestawie dokumentów przygotowanych przez Krzyśka Trajkowskiego. Zastanawiałem się później, w jaki sposób na dostępność materiałów wpływa format np. pdf czy html. Można się zastanawiać, można sprawdzić. Gdy Krzysiek przygotował kolejny dokument skonwertowałem go do html’a i … dodałem do listy gościnnych wpisów.
Konwersja została wykonana z użyciem konwertera pandoc i była ‘prawie’ bezbolesna, tzn te same źródła, które były użyte do kompilacji do pdf’a, po drobnych zmianach posłużyły do kompilacji do html’a. Jedna z niewielu zmian to zastąpienie środowiska lstlisting tagiem <pre lang=”rsplus”>

Tak więc, dziś będzie gościnny wpis Krzysztofa Trajkowskiego, dotyczący wybranych narzędzi dostępnych w R pozwalających na analizę pewnych szeregów czasowych w finansach.
Zapraszamy!

Wersja pdf znajduje się tutaj, wersja HTML poniżej.

Wybrane pakiety do analizy finansowych szeregów czasowych w R

Krzysztof Trajkowski

Pozyskiwanie danych

Do pozyskania aktualnych danych giełdowych (obiekt OHLC) wykorzystamy funkcję getSymbols{quantmod}. Poniżej przykład wczytania danych giełdowych spółki Adobe Systems Inc. z datą końcową 2013-03-31.

Aby korzystać z szeregu funkcji do przetwarzania danych dostępnych w pakiecie quantmod (np. pozyskiwanie informacji o szeregu, agregacja danych) zalecane jest aby szereg czasowy był obiektem klasy xts{xts} lub zoo{zoo}.

Poniżej kilka przykładów funkcji do pozyskiwania informacji o szeregu ADBE.

Możemy również wyodrębnić dane z interesującego nas okresu.

Zamiast opcji days możemy stosować także inne warianty np. years, quarters, months i weeks.

Pozostałe dostępne opcje to: to.yearly, to.quarterly, to.monthly i to.daily.

Na podstawie wcześniej wczytanych danych tzn. ADBE dotyczących spółki Adobe obliczymy logarytmiczne stopy zwrotu z wykorzystaniem cen zamknięcia.

Teraz wczytamy kilka obiektów OHLC dla głównych indeksów giełdowych: niemiecki DAX, francuski CAC40, brytyjski FTSE, japoński NIKKEY oraz amerykański NASDAQ.

Oprócz danych dotyczących cen akcji lub indeksów giełdowych możemy także wczytywać dane dotyczące np. kursów walutowych, cen metali kolorowych i surowców oraz dane ekonomiczne różnych krajów.

Osoby zajmujące się analizą techniczną z pewnością zinteresuje możliwość wykonania np. wykresów świecowych. Warto podkreślić, że dla dużego zakresu danych taki wykres upodabnia się do wykresu liniowego.

Wykres liniowy dla spółki Adobe od 03-01-2007 do 28-03-2013.’ [TAF1]

Wykres świecowy dla spółki Adobe od 03-01-2007 do 28-03-2013.’ [TAF2]

Świece są bardziej widoczne gdy wykonamy wykres dla krótszego okresu np. dla stycznia 2013 roku.

Wykres świecowy dla spółki Adobe - marzec 2013.’ [TAF3]

Jak można łatwo zauważyć wykres świecowy dostarcza nam o wiele więcej informacji niż wykres liniowy. Jest to szczególnie widoczne gdy obserwujemy dane o małym zakresie czasowym (rys. [TAF3]). Każda świeca zawiera informacje zawarte w szeregu czasowym typu OHLC tj. o cenie otwarcia, maksymalnej, minimalnej i zamknięcia. Z kolei jej kolor jest uzależniony od tego czy cena otwarcia była większa od ceny zamknięcia – kolor pomarańczowy czy też nie – kolor zielony.

Gdy w kodzie pominiemy opcję TA=NULL dodatkowo pojawi się wykres słupkowy obrazujący wielkość obrotów w danym okresie czasu.

Wykres świecowy dla spółki Adobe - marzec 2013.’ [TAF4]

Zwróćmy także uwagę na pojawiające się wartości liczbowe w górnym lewym rogu każdego okienka wykresu (rys. [TAF4]). Na wykresie świecowym oraz słupkowym jest podana informacja odpowiednio o cenie zamknięcia i wielkości obrotów w ostatnim okresie tj. dla 28-03-2013.

Ciekawą opcją jest również możliwość dodawania wszystkich funkcji dostępnych w pakiecie TTR. Na rysunku [TAF5] przykład wykresu świecowego oraz wstęga Bollingera. Natomiast w dolnej części wykresy MACAD i słupkowy wielkości obrotów.

Wykresy dla spółki Adobe - rok 2012.’ [TAF5]

Budowa szeregu czasowego

W tym podrozdziale pokażemy jak narysować wykres szeregu czasowego – pakiet ggplot2 oraz jak zbudować szereg czasowy – pakiet xts. Zastosowanie pakietu xts do budowy szeregu czasowego umożliwi modyfikację danych za pomocą wielu funkcji dostępnych w bibliotece quantmod. Dane które będziemy wykorzystywać dotyczą cen akcji spółki Telekomunikacja Polska s.a. i zostały pozyskane ze strony www.stooq.pl w formie pliku .csv.

Dzienne ceny zamknięcia dla spółki TPsa od 18-11-1998 do 31-01-2013.’ [cana1]

Aby uzyskać logarytmiczne stopy zwrotu wykorzystamy pakiet xts oraz quantmod.

Dzienne stopy zwrotu spółki TPsa od 18-11-1998 do 28-03-2013.’ [return1]

Właściwości finansowych szeregów czasowych

Rozkład prawdopodobieństwa

Każdy rozkład stopy zwrotu (np. z akcji) charakteryzuje się pewnymi cechami które mają wpływ na kształt rozkładu prawdopodobieństwa zmiennej. Do tych cech można zaliczyć:

  • skośność — miara asymetrii gęstości prawdopodobieństwa w stosunku do rozkładu normalnego,

  • leptokurtoza — wyższy szczyt gęstości prawdopodobieństwa w porównaniu z rozkładem normalnym,

  • grube ogony — gęstość prawdopodobieństwa na końcach rozkładu jest większa niż w rozkładzie normalnym.

Sprawdźmy zatem czy rzeczywiście zmienna rtp1 nie ma rozkładu normalnego.

Na podstawie powyższych obliczeń dochodzimy do wniosku, że rozkład stóp zwrotu nie ma rozkładu normalnego. Wniosek ten potwierdzają także dwa wykresy: kwantylowy oraz gęstości.

Wykres kwantylowy oraz gęstości - dzienne stopy zwrotu.’ [QD1]

Połączymy teraz nasze dane w większe jednostki czasowe (agregacja szeregów czasowych) i sprawdzimy jak to przekształcenie szeregu wpłynie na rozkład prawdopodobieństwa. W pierwszej kolejności dokonamy agregacji danych dziennych do tygodniowych, następnie do miesięcznych i kwartalnych.

Na podstawie powyższych obliczeń możemy wysnuć wniosek, że im większa agregacja zmiennej (stopa zwrotu) tym bardziej rozkład prawdopodobieństwa upodabnia się do rozkładu normalnego.

W dalszej części tego opracowania będziemy analizować dzienne stopy zwrotu z akcji spółki Telekomunikacja Polska s.a. Zatem interesująca może się okazać odpowiedź na pytanie: jakiemu rozkładowi prawdopodobieństwa może podlegać zmienna rtp1.

Na podstawie testu Andersona-Darlinga możemy stwierdzić, że badana stopa zwrotu ma rozkład skośny t-Studenta. Również wykresy (rys. [QD2]) potwierdzają ten wniosek.

Wykres kwantylowy oraz gęstości - dzienne stopy zwrotu.’ [QD2]

Niejednorodność wariancji

Efekt skupiania zmienności (volatility clustering) to kolejna własność stóp zwrotu. Do oceny czy w badanym szeregu występuje zjawisko grupowania wariancji wykorzystamy test Engle’a.

Zatem na poziomie istotności α = 0, 05 odrzucamy hipotezę zerową która zakłada jednorodność wariancji. Potwierdzeniem efektu ARCH może być również autokorelacja kwadratów stóp zwrotu.

Autokorelacja kwadratów stóp zwrotu spółki TPsa.’ [AUT2]

Z problemem grupowania wariancji wiąże się kolejna własność szeregów finansowych tzw. efekt dźwigni (leverage effects). Na rysunku [CS1] można zaobserwować silniejszą reakcję zmienności na duże spadki cen, niż na wzrosty.

Ceny i stópy zwrotu spółki TPsa.’ [CS1]

Efekt długiej pamięci

Na podstawie wartości współczynnika d (ułamkowy pierwiastek jednostkowy) lub H (wykładnik Hursta) możemy ocenić czy w badanym szeregu czasowym występuje efekt długiej pamięci (long memory processes). Warto podkreślić, że między omawianymi parametrami występuje pewien związek który można opisać równaniem:


d=H-0,5

Poniżej przedstawimy własności szeregów czasowych w zależności od wartości parametru d oraz H:

  • szereg antypersystentny (krótka pamięć):
    d ∈ ( − 0, 5; 0) lub H ∈ (0; 0, 5)

  • szereg losowy (biały szum):
    d = 0 lub H = 0, 5

  • szereg persystentny (długa pamięć):
    d ∈ (0; 0, 5) lub H ∈ (0, 5; 1)

Estymacja współczynnika integracji ułamkowej d.

Estymacja wykładnika Hursta H.

Ponieważ oszacowana wartość niecałkowitego pierwiastka jednostkowego jest równa  =  − 0. 00618375 oraz wartość wykładnia Hursta wynosi Ĥ = 0. 4975533 to możemy stwierdzić, że w badanej stopie zwrotu nie występuje efekt długiej pamięci. Warto podkreślić, że w sytuacji odwrotnej tzn. gdy w badanym szeregu wystąpi efekt długiej pamięci można zastosować model z warunkową wartością oczekiwaną ARFIMA z wykorzystaniem np. funkcji fracdiff{FGN}. Dodajmy jeszcze, że w przypadku krótkiej pamięci funkcja autokorelacji ACF opada szybko tzn. w tempie wykładniczym, natomiast gdy występuje długa pamięć opadanie funkcji ACF jest wolniejsze tzn. hiperboliczne.

Funkcja autokorelacji stóp zwrotu.

Autokorelacja stóp zwrotu spółki TPsa.’ [AUT1]

Jak się okazuje (rys. [AUT1]) szereg finansowy jest ciągiem nieskorelowanych zmiennych losowych. Zjawisko autokorelacji w badanym szeregu nie występuje.

Należy także wspomnieć o bardzo bogatej kolekcji metod estymacji współczynników d oraz H dostępnych w pakiecie fArma.

Modele ARCH/GARCH

Do estymacji standardowego modelu GARCH możemy użyć funkcji garch{tseries}. Dla modelu ARCH(1) [Engle, 1982] składnia jest następująca:

Natomiast dla modelu GARCH(1,1) [Bollerslev, 1986]:

Zwróćmy uwagę na to, że funkcja garch{tseries} zawsze zakłada rozkład normalny reszt. Jak się okazało to założenie zostało naruszone – test Jarque-Bera. W związku z tym należy oszacować model, który uwzględni rozkład inny niż normalny np. rozkład skośny t-Studenta.

Modele APARCH

Gdy w szeregu czasowym występuje asymetria lub efekt dźwigni dobrym rozwiązaniem jest zastosowanie modelu GARCH lub APARCH z odpowiednim rozkładem innowacji. Modele tej klasy możemy estymować za pomocą funkcji garchFit{fGarch}. Opcja cond.dist określa jaki rozkład prawdopodobieństwa będzie wykorzystany w modelu: std, sstd, ged, sged, snig, QMLE, snorm. Rozkład norm czyli rozkład normalny, jest ustawiony jako domyślny. Jest też możliwość uwzględnienia modelu ARMA np. ARMA(2,1)-APARCH(1,1). W takim przypadku formuła określająca model będzie miała postać ~arma(2,1)+aparch(1,1).

  • APARCH – Asymmetric Power Autoregressive Conditional Heteroscedasticity:

Dodajmy jeszcze, że include.delta=T oraz leverage=T są ustawieniami domyślnymi. Oznacza to, że parametry δ oraz γ w modelu APARCH będą oszacowane. W przeciwnym wypadku tzn. dla include.delta=F oraz leverage=F parametry δ oraz γ nie będą szacowane. Wtedy domyślnymi wartościami będą delta=2 i gamma=0. Oczywiście tym parametrom można przypisywać dowolne wartości. Wynika z tego, że w modelu APARCH jest zagnieżdżonych kilka innych modeli z rodziny GARCH. Poniżej wymienimy niektóre z nich:

  • ARCH — Autoregressive Conditional Heteroscedasticity:

  • GARCH — Generalized Autoregressive Conditional Heteroscedasticity:

  • TS-GARCH

  • GJR-GARCH

  • T-GARCH — Threshold GARCH:

  • N-ARCH — Nonlinear ARCH:

Poniżej sprawdzimy czy model sstd-N-ARCH(1,1) będzie dobrze dopasowany do badanej zmienej rtp1. Zauważmy, że include.delta=T oraz gamma=0 (dla opcji leverage=F) są wartościami domyślnymi i nie musimy ich uwzględniać w kodzie.

Dokonajmy jeszcze estymacji drugiego modelu tzn. sstd-GARCH(1,1), ale z wykorzystaniem formuły ~garch zamiast ~aparch:

Ponieważ logarytm wiarygodności jest większy niż w modelu sstd-N-ARCH to do dalszych analiz wykorzystamy model sstd-GARCH. Poniżej kilka podstawowych testów diagnostycznych dla modelu sstd-GARCH(1,1).

W omawianym modelu fit2 założyliśmy rozkład skośny t-Studenta (na powyższym wydruku jest badana zgodność reszt tylko z rozkładem normalnym) więc sprawdzimy to założenie przy pomocy testu Andersona-Darlinga.

Dodatkowo mamy do wyboru kilkanaście wykresów diagnostycznych.

Poniżej przedstawimy wykresy o numerach 9 i 13 czyli wykres standaryzowanych reszt oraz wykres kwantylowy dla standaryzowanych reszt.

Wykres reszt oraz kwantylowy.’ [ga1]

Inne modele

Dużo większe możliwości od omawianych powyżej bibliotek posiada pakiet rugarch. Oprócz wyżej wymienionych metod mamy możliwość estymacji modeli które uwzględniają np. efekt długiej pamięci. Poniżej kilka przykładów estymacji modeli z rodziny GARCH.

  • I-GARCH -

  • EWMA

  • E-GARCH

  • GJR-GARCH

  • T-ARCH

  • F-GARCH

Ponieważ w analizowanym szeregu rtp1 nie występuje efekt długiej pamięci, więc oszacujemy model sstd-GARCH(1,1).

Jak można zauważyć oprócz wyników estymacji danego modelu zostały podane także wyniki kilku testów diagnostycznych. Zaznaczmy, że wynik dowolnego testu oraz kilka innych informacji możemy uzyskać wpisując odpowiednią komendę.

Podobnie jak w pakiecie fGarch mamy możliwość wygenerować kilkanaście wykresów diagnostycznych dla modelu oszacowanego z wykorzystaniem biblioteki rugarch.

Wykresy diagnostyczne.’ [ga1]

Na koniec zaznaczmy, że nie sposób przedstawić wszystkie modele z rodziny GARCH, cały czas pojawiają się nowe metody. Ciekawą propozycją są modele beta-t-EGARCH dostępne w bibliotece betategarch. Z kolei zwolenników metod bayesowskich z pewnością zainteresuje biblioteka bayesGARCH. Wiele interesujących informacji można znaleźć na bologu http://www.unstarched.net/blog/. Są tam prezentowane ciekawe przykłady z wykorzystaniem różnych pakietów środowiska R np. biblioteki rugarch. Bardzo dobry opis możliwości tego pakietu jest dostępny na stronie http://cran.r-project.org/web/packages/rugarch/vignettes/Introduction_to_the_rugarch_package.pdf. Osoby zainteresowane modelami wielowymiarowymi GARCH (M-GARCH) powinny zwrócić uwagę np. na pakiet gogarch uogólnione modele ortogonalne GARCH (GO-GARCH) oraz ccgarch modele korelacji warunkowej GARCH (CC-GARCH) np. DCC-GARCH (dynamic conditional correlations GARCH) i CCC-GARCH (constant conditional correlation GARCH).

Krzysztof Trajkowski.

2 myśli na temat “Wybrane pakiety do analizy finansowych szeregów czasowych w R”

Dodaj komentarz

Twój adres e-mail nie zostanie opublikowany. Pola, których wypełnienie jest wymagane, są oznaczone symbolem *

Możesz użyć następujących tagów oraz atrybutów HTML-a: <a href="" title=""> <abbr title=""> <acronym title=""> <b> <blockquote cite=""> <cite> <code class="" title="" data-url=""> <del datetime=""> <em> <i> <q cite=""> <strike> <strong> <pre class="" title="" data-url=""> <span class="" title="" data-url="">