Zaawansowana Analiza Szeregów Czasowych
Projekt Zaliczeniowy:
Analiza porównawcza ryzyka dla portfela zrównoważonego przy wykorzystaniu modeli klasy GARCH
Michał Skowronek & Adam Kordeczka
03.09.2017 r.
Wstęp
Celem niniejszego raportu jest przedstawienie analizy porównawczej ryzyka portfelowego wykonanej dla portfela zrównoważonego. Analiza zostala oparta o wyestymowanie warunkowej wariancji oraz wartości narażonej na ryzyko. Do tego celu wykorzystano modele klasy GARCH
Metodologia doboru spółek do portfela
Jednym z głównych wymagań projektu było zbudowanie portfela zrównoważonego. Idea jego konstrukcji polega na przyporządkowaniu jednakowych wag dla każdego z aktywów. W skład portfela weszły akcje następujących spółek: (notowane na NASDAQ i NYSE):
- akcje spółki Microsoft (MSFT.US, NASDAQ)
- akcje spółki Amazon (AMZN.US, NASDAQ)
- akcje spółki NVIDIA (NVDA.US, NASDAQ)
- akcje spółki Statoil (STO.US, NYSE)
- kurs złota do dolara amerykańskiego (XAUUSD, NYSE)
Pierwsze trzy spółki pochodzą z sektora IT. Zdecydowano się na ich wybór głównie z dwóch powodów. Pierwszy z nich wynikał z istotnego, długoterminowego trendu wzrostowego. Drugi to wewnętrzne przekonanie autorów projektu, że rozwijane przez owe spółki rozwiązania chmurowe (Microsoft Azure, Amazon Web Services) wzmocnią ich i tak już silną pozycje na rynku. Spółka Statoil to norweskie przedsiębiorstwo energetyczno-naftowe, którego udziałowcem większościowym jest Rząd Norwegii (67%). Nie bez znaczenia pozostaje fakt, że Norwegia posiada narodowy fundusz inwestycyjny, który co roku przynosi pokaźne stopy zwrotu. Taki udziałowiec nadaje spółce obiecujący status firmy stabilnej, nienarażonej na nieodpowiedzialne inwestycyjne decyzje zarządu. Ostatnie aktywo - kurs złota do dolara - to wybór czysto protekcjonalny. Panuje silne przekonanie, że złoto to zawsze bezpieczna inwestycja, zwłaszcza w czasie recesji. Wniosek ten znajduje swoje potwierdzenie przy porównaniu np. cen indeksu SP500 w ceny złota chociażby na początku XXI wieku.
Analiza opisowa
Dane zostały pobrane z portalu www.stooq.pl. Poniższe dwa wykresy przedstawiają zależności ceny zamknięcia kursów akcji od czasu dla wspólnego horyzontu czasowego, który został ustalony na okres od 01.03.2005 do 31.08.2017 (taki zakres był uwarunkowany chęciami uwzględnienia kryzysu z 2008 roku oraz późniejszych szybkich wzrostów dla spółek z sektora IT).
Na podstawie wartości wybranych aktywów skonstruowano portfel zrównoważony biorąc średnią arytmetyczną z logarytmicznych stóp zwrotu dla każdego z indeksów. Wykres zależności od czasu dla tak skonstruowanego portfela przedstawiono na rysunku poniżej:
Na wykresie widać wyraźnie kilka stylizowanych faktów charakterystycznych dla finansowych szeregów czasowych. Po pierwsze: grupowanie wariancji. Małe oraz duże zmiany logarytmicznych stóp zwrotu mają tendencje do występowania seriami. Najlepiej jest to widoczne podczas kryzysu z 2008 roku. Kolejny fakt stylizowany to efekt dźwigni. Zjawisko to określa tendencje do wzrostu wariancji w czasie gdy ceny aktywa spadają (tzw. asymetria warunkowej wariancji).
Sprawdźmy teraz, czy trzeci ze stylizowanych faktów jest obecny - zestawmy rozkład stóp zwrotu z portfela z rozkładem normalnym.
Powyższy wykres przedstawia histogram rozkładu logarytmicznych stóp zwrotu z naszego portfela, wraz z naniesioną teoretyczną krzywą rozkładu normalnego. Wyraźnie widzimy na nim, że wspomniany rozkład wykazuje się leptokurtycznością - świadczy o tym przede wszystkim znacznie podwyższony szczyt funkcji gęstości. Ponadto rozkład zdaje się mieć grube ogony. Sprawdźmy, czy będą one lepiej widoczne na wykresie Q-Q:
Postać wykresu Q-Q w zdecydowany sposób świadczy o tym, że badany rozkład ma grube ogony. Obliczmy jeszcze statystyki naszego rozkładu
## zwroty
## nobs 3071.000000
## NAs 0.000000
## Minimum -0.095275
## Maximum 0.099572
## 1. Quartile -0.006103
## 3. Quartile 0.008003
## Mean 0.000620
## Median 0.000880
## Sum 1.904628
## SE Mean 0.000258
## LCL Mean 0.000114
## UCL Mean 0.001126
## Variance 0.000205
## Stdev 0.014301
## Skewness -0.197807
## Kurtosis 5.911800Obliczone statystyki potwierdzają wcześniejsze wnioski - kurtoza i skośność znacznie odbiegają od wartości właściwych dla rozkładu normalnego. Na koniec przeprowadźmy test Jarque-Bera. Jego hipoteza zerowa mówi, że testowany rozkład jest normalny.
##
## Jarque Bera Test
##
## data: na.omit(portfel_mean$zwroty)
## X-squared = 4500.9, df = 2, p-value < 2.2e-16Widzimy, że hipoteza ta jest silnie odrzucana przez wynik testu. Ostatecznie stwierdzmay zatem, że rozkład stóp zwrotu z badanego portfela jest leptokurtotyczny, a co za tym idzie wszystkie trzy stylizowane fakty są obecne.
Sprawdźmy teraz autokorelację stóp zwrotu i ich kwadratów:
Na pierwszym wykresie widzimy, że wśród zwrotów z portfela wystepuje autokorelacja - opóźnienia 1, 2, 3, 5, 7, 11, 17 są istotne statystycznie. Z kolei drugi pokazuje silną autokorelację kwadratów logarytmicznych stóp zwrotu dla wszystkich opóźnień. Może to świadczyć o występowaniu efektów ARCH. Otrzymane wizualne rezultaty potwierdźmy testami statystycznymi:
## lag Autocorrelation D-W Statistic p-value
## 1 0.2902019 1.419573 0
## 2 0.2955578 1.408789 0
## 3 0.2468385 1.506168 0
## 4 0.3047129 1.390363 0
## 5 0.3441656 1.311414 0
## 6 0.3322022 1.335287 0
## 7 0.3183671 1.362885 0
## 8 0.1862127 1.627155 0
## 9 0.2829457 1.433619 0
## 10 0.2961132 1.407239 0
## Alternative hypothesis: rho[lag] != 0Za pomocą testu Durbina-Watsona sprawdzamy czy w naszym szeregu kwadratów stóp zwrotu występuje autokorelacja. Wartości statystyki testowej są znacznie poniżej dwóch, p-value testu dla opóźnień od pierwszego do dziesiątego dają podstawy do odrzucenia hipotezy zerowej. Między kwadratami stóp zwrotu występuje dodatnia autokorelacja. Przeprowadźmy teraz test L-M na obecność efektów ARCH:
##
## ARCH LM-test; Null hypothesis: no ARCH effects
##
## data: portfel_mean$zwroty
## Chi-squared = 587.64, df = 5, p-value < 2.2e-16Wartość p-value dla testu LM jest bliska zeru. Odrzucamy zatem hipotezę zerową, która mówi o braku efektów ARCH w testowanym szeregu na rzecz hipotezy alternatywnej, która je potwierdza.
Estymacja modelu GARCH
Zgodnie z zaleceniami do projektu podejmiemy próbę estymacji modelu GARCH(p,q) w celu późniejszego porówniania, z jego rozszerzeniami, wartośći narażonych na ryzyko. Z wykresu autokorelacji zwrotów z portfela wynika, że istotne statystycznie są pierwsze trzy opóźnienia. Spójrzmy jeszcze jak kształtują się autokorelacje cząstkowe:
Autokorelacje cząstkowe dla pierwszego i drugiego opóźnienia są istotne statystycznie. Zdecydowano się oszacować model GARCH(p,q) dla różnych wartości p oraz q. Przetestowano wszytskie kombinacje parametrów ze zbioru {1,2,3} (łącznie 9). Poniższa tabela zawiera zestawienie wartości kryteriów informacyjnych dla każdego z wyestymowanych modeli.
## p q AIC BIC
## 1 1 1 -5.950000 -5.944110
## 2 2 1 -5.949244 -5.941390
## 3 3 1 -5.948451 -5.938633
## 4 1 2 -5.949184 -5.941330
## 5 2 2 -5.948611 -5.938794
## 6 3 2 -5.948246 -5.936465
## 7 1 3 -5.948326 -5.938509
## 8 2 3 -5.947772 -5.935991
## 9 3 3 -5.949495 -5.935751Najniższa wartość obu kryteriów informacyjnych została osiągnięta dla modelu GARCH(1,1), który ostatecznie wykorzystamy do dalszych analiz. Sprawdźmy jego diagnostykę:
##
## Title:
## GARCH Modelling
##
## Call:
## garchFit(formula = ~garch(1, 1), data = portfel_mean_nona$zwroty,
## cond.dist = "norm", include.mean = FALSE, trace = FALSE)
##
## Mean and Variance Equation:
## data ~ garch(1, 1)
## <environment: 0x7fba5f5c4b80>
## [data = portfel_mean_nona$zwroty]
##
## Conditional Distribution:
## norm
##
## Coefficient(s):
## omega alpha1 beta1
## 3.3665e-06 6.8470e-02 9.1235e-01
##
## Std. Errors:
## based on Hessian
##
## Error Analysis:
## Estimate Std. Error t value Pr(>|t|)
## omega 3.366e-06 8.316e-07 4.048 5.16e-05 ***
## alpha1 6.847e-02 1.022e-02 6.699 2.10e-11 ***
## beta1 9.124e-01 1.329e-02 68.648 < 2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Log Likelihood:
## 9139.226 normalized: 2.975977
##
## Description:
## Sat Sep 9 19:43:56 2017 by user:
##
##
## Standardised Residuals Tests:
## Statistic p-Value
## Jarque-Bera Test R Chi^2 487.8113 0
## Shapiro-Wilk Test R W 0.983884 0
## Ljung-Box Test R Q(10) 14.75807 0.1411345
## Ljung-Box Test R Q(15) 15.93069 0.386673
## Ljung-Box Test R Q(20) 18.49046 0.5551333
## Ljung-Box Test R^2 Q(10) 13.70032 0.1871048
## Ljung-Box Test R^2 Q(15) 14.72576 0.4713455
## Ljung-Box Test R^2 Q(20) 18.49936 0.554549
## LM Arch Test R TR^2 13.40663 0.3401904
##
## Information Criterion Statistics:
## AIC BIC SIC HQIC
## -5.950000 -5.944110 -5.950002 -5.947884Wszystkie wartości wyestymowanych parametrów są istotne statystycznie (odrzucamy hipotezę zerową o zerowych wartościach współczynników). Współczynnik \(\omega\), oraz wyliczona na podstawie współczynników wariancja bezwarunkowa są dodatnie, nie występują żadne anomalie.
## omega alpha1 beta1
## 0.00000337 0.06847002 0.91235113Test Jarque-Bera stwierdza, że rozkład reszt nie jest normalny (praktycznie zerowe p-value, efekt większego prawdopodobieństwa,w stosunku do rozkładu normalnego, występowania nietypowych zmian kursu, charakterystyczne dla danych finansowych). Wyniki statystyk Ljunga-Boxa mają wysokie wartości p-value, nie mamy zatem podstaw do odrzucenia hipotezy zerowej o braku autokorelacji wśród reszt modelu. p-value dla testu ARCH jest wysokie co nie daje podstaw do odrzucenia hipotezy zerowej o braku występowania efektów ARCH. Nasz model zdaje się mieć dobrą specyfikacje. Poniższy wykres przedstawia wartości bezwzględne stóp zwrotu z nałożonymi wyestymowanymi w modelu warunkowymi wariancjami:
Na koniec oszacujmy jeszcze wariancję bezwarunkową zgodnie ze wzorem: \[ Var(u_t) = \frac{\omega}{1-(\alpha_1 + \beta_1)}\] Podstawiając oszacowane w modelu parametry otrzymujemy:
## omega
## 0.0001757144Model zdaje się dość dobrze uwzględniać stylizowany fakt grupowania wariancji (zwłaszcza w okolicach kryzysu z 2008 roku). GARCH(1,1) zdaje się zatem dobrze modelować zmienność naszego portfela, jednak należy wziąć pod uwagę dwie rzeczy. Po pierwsze, standaryzowane reszty z modelu nie podlegają założonemu tu przez nas rozkładowi normalnemu. Po drugie musimy pamiętać, że standardowe modele GARCH nie biorą pod uwagę jednego ze stylizowanych faktów, który zaobserwowaliśmy na etapie analizy opisowej - efektu dźwigni. By zaradzić tym dwóm problemom, zastosujmy odpowiednie rozszerzenia modelu GARCH(1,1).
Estymacja modelu E-GARCH
W rozdziale Analiza opisowa stwierdzono, że nasz szereg czasowy przejawia charakterystyczny efekt dźwigni. Jedym z rozszerzeń modelu GARCH, które stosunkowo dobrze uwzględnia wspomniany fakt stylizowany, jest model E-GARCH. Wykorzystuje się w nim asymetryczną krzywą wpływu informacji na wariancję warunkową (lewe ramie ma rosnąć szybciej niż prawe). Ponadto model ten nie wymaga nieujemności współczynników. Postanowiono wyestymować ten model w takim samym wariancie jak GARCH czyli (1,1). Dodatkowo, jako model średniej postanowiono dołączyć model ARFIMA(1,0,1) (pewna generalizacja modelu ARIMA). Wartości p = 1 oraz q = 1 wybrano na podstawie występowania istotnych statystycznie autokorelacji zwykłych oraz cząstkowych dla pierwszych opóźnień. Poniżej podsumowanie modelu:
##
## *---------------------------------*
## * GARCH Model Fit *
## *---------------------------------*
##
## Conditional Variance Dynamics
## -----------------------------------
## GARCH Model : eGARCH(1,1)
## Mean Model : ARFIMA(1,0,1)
## Distribution : norm
##
## Optimal Parameters
## ------------------------------------
## Estimate Std. Error t value Pr(>|t|)
## ar1 0.528316 0.034432 15.3438 0
## ma1 -0.530108 0.034302 -15.4543 0
## omega -0.137160 0.003822 -35.8827 0
## alpha1 -0.067939 0.005784 -11.7464 0
## beta1 0.983524 0.000461 2133.0579 0
## gamma1 0.113548 0.011531 9.8472 0
##
## Robust Standard Errors:
## Estimate Std. Error t value Pr(>|t|)
## ar1 0.528316 0.008001 66.0302 0e+00
## ma1 -0.530108 0.008015 -66.1357 0e+00
## omega -0.137160 0.008461 -16.2105 0e+00
## alpha1 -0.067939 0.008834 -7.6908 0e+00
## beta1 0.983524 0.000927 1060.5448 0e+00
## gamma1 0.113548 0.022699 5.0023 1e-06
##
## LogLikelihood : 9163.08
##
## Information Criteria
## ------------------------------------
##
## Akaike -5.9636
## Bayes -5.9518
## Shibata -5.9636
## Hannan-Quinn -5.9593
##
## Weighted Ljung-Box Test on Standardized Residuals
## ------------------------------------
## statistic p-value
## Lag[1] 0.5321 0.4657
## Lag[2*(p+q)+(p+q)-1][5] 1.4690 0.9984
## Lag[4*(p+q)+(p+q)-1][9] 4.1095 0.6652
## d.o.f=2
## H0 : No serial correlation
##
## Weighted Ljung-Box Test on Standardized Squared Residuals
## ------------------------------------
## statistic p-value
## Lag[1] 1.005 0.3161
## Lag[2*(p+q)+(p+q)-1][5] 3.059 0.3964
## Lag[4*(p+q)+(p+q)-1][9] 5.442 0.3675
## d.o.f=2
##
## Weighted ARCH LM Tests
## ------------------------------------
## Statistic Shape Scale P-Value
## ARCH Lag[3] 1.896 0.500 2.000 0.1686
## ARCH Lag[5] 3.770 1.440 1.667 0.1959
## ARCH Lag[7] 5.096 2.315 1.543 0.2150
##
## Nyblom stability test
## ------------------------------------
## Joint Statistic: 1.3262
## Individual Statistics:
## ar1 0.1213
## ma1 0.1213
## omega 0.6534
## alpha1 0.1955
## beta1 0.6133
## gamma1 0.1193
##
## Asymptotic Critical Values (10% 5% 1%)
## Joint Statistic: 1.49 1.68 2.12
## Individual Statistic: 0.35 0.47 0.75
##
## Sign Bias Test
## ------------------------------------
## t-value prob sig
## Sign Bias 0.4052 0.6854
## Negative Sign Bias 0.9561 0.3391
## Positive Sign Bias 0.4807 0.6307
## Joint Effect 1.3188 0.7247
##
##
## Adjusted Pearson Goodness-of-Fit Test:
## ------------------------------------
## group statistic p-value(g-1)
## 1 20 73.44 2.441e-08
## 2 30 91.21 2.367e-08
## 3 40 95.54 1.202e-06
## 4 50 102.67 1.132e-05
##
##
## Elapsed time : 0.3369091Warto zwrócić uwagę na fakt, że wszystkie współczynniki są istotnie statystycznie (również człony regresyjne i średniej ruchomej). Przyjrzyjmy się wyestymowanym parametrom:
## ar1 ma1 omega alpha1 beta1 gamma1
## 0.52831552 -0.53010758 -0.13715951 -0.06793861 0.98352355 0.11354754Model potwierdza zaobserwowany wcześniej efekt dźwigni - parametr alpha1 jest mniejszy od zera. W modelu jest on współczynnikiem proporcjonalności między \(log(\sigma^2_t)\) a \(\frac{\epsilon_t}{\sqrt{\sigma_{t-1}^2}}\) czyli zadaje zależność między logarytmem wariancji a szokiem w chwili t-1. Wartości testu Ljunga-Boxa sugerują brak autokorelacji reszt modelu, nie mamy również podstaw do odrzucenia hipotezy o braku efektów ARCH (test LM). Sprawdźmy jeszcze empiryczny rozkład reszt modelu:
##
## Jarque Bera Test
##
## data: model_eGarch@fit$residuals
## X-squared = 4490.4, df = 2, p-value < 2.2e-16Odrzucamy hipotezę zerową na rzecz hipotezy alternatywej, reszty modelu nie pochodzą z rozkładu normalnego. Znajduje to również potwierdzenie na poniższym empirycznym histogramu reszt:
Poniżej mamy wykres zależności wariancji warunkowej od opóźnionego szoku czyli krzywą wpływu informacji:
Jak widać, jej kształt jest zgodny z założeniami modelu, istnieje asymetria między lewym a prawym zboczem. Reakcja warunkowej wariancji portfela na napływające na rynek informacje jest zdecydowanie silniejsza dla tych złych niż dobrych (zmiana szoku po lewej stronie zera daje znacznie większy skok w wariancji warunkowej niż ten sam przyrost szoku po stronie prawej).
Estymacja modelu TGARCH
Zastosujmy teraz kolejne rozszerzenie modelu GARCH - tzw. threshold GARCH (TGARCH). Jest to kolejny z modeli, który posiada zdolność do identyfikowania efektu asymetrii wariancji warunkowej - o jego obecności mówić nam będzie istotność parametru \(\eta\). Tak jak poprzednio rozszerzenie będziemy stosować względem modelu GARCH(1,1), przy uwzględnieniu modelu średniej ARFIMA(1,1).
##
## *---------------------------------*
## * GARCH Model Fit *
## *---------------------------------*
##
## Conditional Variance Dynamics
## -----------------------------------
## GARCH Model : fGARCH(1,1)
## fGARCH Sub-Model : TGARCH
## Mean Model : ARFIMA(1,0,1)
## Distribution : norm
##
## Optimal Parameters
## ------------------------------------
## Estimate Std. Error t value Pr(>|t|)
## ar1 -0.761978 0.057199 -13.3216 0.0e+00
## ma1 0.756375 0.057801 13.0859 0.0e+00
## omega 0.000233 0.000060 3.8942 9.9e-05
## alpha1 0.059218 0.010751 5.5084 0.0e+00
## beta1 0.938462 0.011894 78.9006 0.0e+00
## eta11 0.641055 0.112368 5.7050 0.0e+00
##
## Robust Standard Errors:
## Estimate Std. Error t value Pr(>|t|)
## ar1 -0.761978 0.006394 -119.1732 0.000000
## ma1 0.756375 0.007178 105.3781 0.000000
## omega 0.000233 0.000116 2.0027 0.045205
## alpha1 0.059218 0.021919 2.7017 0.006899
## beta1 0.938462 0.024592 38.1613 0.000000
## eta11 0.641055 0.172305 3.7205 0.000199
##
## LogLikelihood : 9168.519
##
## Information Criteria
## ------------------------------------
##
## Akaike -5.9671
## Bayes -5.9553
## Shibata -5.9671
## Hannan-Quinn -5.9629
##
## Weighted Ljung-Box Test on Standardized Residuals
## ------------------------------------
## statistic p-value
## Lag[1] 0.2279 0.6331
## Lag[2*(p+q)+(p+q)-1][5] 1.1717 0.9999
## Lag[4*(p+q)+(p+q)-1][9] 3.6311 0.7740
## d.o.f=2
## H0 : No serial correlation
##
## Weighted Ljung-Box Test on Standardized Squared Residuals
## ------------------------------------
## statistic p-value
## Lag[1] 1.152 0.2831
## Lag[2*(p+q)+(p+q)-1][5] 3.376 0.3426
## Lag[4*(p+q)+(p+q)-1][9] 5.823 0.3191
## d.o.f=2
##
## Weighted ARCH LM Tests
## ------------------------------------
## Statistic Shape Scale P-Value
## ARCH Lag[3] 2.329 0.500 2.000 0.1270
## ARCH Lag[5] 4.024 1.440 1.667 0.1718
## ARCH Lag[7] 5.473 2.315 1.543 0.1809
##
## Nyblom stability test
## ------------------------------------
## Joint Statistic: 1.7583
## Individual Statistics:
## ar1 0.4986
## ma1 0.5051
## omega 0.7208
## alpha1 0.8587
## beta1 0.8480
## eta11 0.3494
##
## Asymptotic Critical Values (10% 5% 1%)
## Joint Statistic: 1.49 1.68 2.12
## Individual Statistic: 0.35 0.47 0.75
##
## Sign Bias Test
## ------------------------------------
## t-value prob sig
## Sign Bias 0.5586 0.5765
## Negative Sign Bias 0.9076 0.3641
## Positive Sign Bias 0.6490 0.5164
## Joint Effect 1.3265 0.7228
##
##
## Adjusted Pearson Goodness-of-Fit Test:
## ------------------------------------
## group statistic p-value(g-1)
## 1 20 71.64 4.892e-08
## 2 30 88.63 5.903e-08
## 3 40 92.05 3.534e-06
## 4 50 96.88 5.486e-05
##
##
## Elapsed time : 1.019946Wyświetlone podsumowanie świadczy o tym, że wyestymowany model TGARCH(1,1) jest dobrze dopasowany do danych. Po pierwsze, wszystkie parametry dopasowania są istotne statystycznie:
## Estimate Std. Error t value Pr(>|t|)
## ar1 -0.7619783516 5.719867e-02 -13.321611 0.000000e+00
## ma1 0.7563748991 5.780086e-02 13.085876 0.000000e+00
## omega 0.0002330889 5.985487e-05 3.894234 9.850938e-05
## alpha1 0.0592184599 1.075055e-02 5.508412 3.620853e-08
## beta1 0.9384617789 1.189423e-02 78.900601 0.000000e+00
## eta11 0.6410549656 1.123679e-01 5.704966 1.163664e-08Parametr \(\eta_{11}\) świadczący o asymetrii warunkowej wariancji jest istotny, co jest zgodne z wcześniejszymi rezultatami. Ponadto rezultaty testu Ljunga-Boxa świadczą o prawdopodobnym braku autokorelacji wśród standaryzowanych reszt, jak i wśród ich kwadratów. Sprawdźmy wykresy:
Z pierwszego wykresu wynika, że dla standaryzowanych reszt istotne statystycznie są 3 pojedyńcze opóźnienia, które nie pozwalają nam zidentyfikować żadnej konkretnej zależności. Z kolei dla kwadratów reszt istotna statystycznie jest jedynie autokorelacja dla 10 opóźnienia.
W podsumowaniu widzimy także wyniki testu LM na obecność efektów ARCH - hipoteza o ich występowaniu wśród standaryzowanych reszt jest odrzucana. Sprawdźmy jeszcze, czy reszty podlegają założonemu rozkładowi normalnemu: 
Widzimy, że reszty nie podlegają rozkładowi normalnemu - QQ plot wskazuje na wyraźne grube ogony. Udało się nam dopasować dwa modele uwzględniające efekt dźwigni, a tym samym na dwa różne sposoby potwierdzić obecność asymetrii warunkowej wariancji. Reszty z modelu nadal nie podlegają jednak zakładanemu rozkładowi normalnemu. Przetestujmy teraz kolejne rozszerzenie modelu GARCH(1,1) - model GARCH-t, który zakłada, że reszty podlegają rozkładowi t-studenta, który lepiej modeluje leptokurtyczność.
## Estymacja modelu GARCH-t
##
## *---------------------------------*
## * GARCH Model Fit *
## *---------------------------------*
##
## Conditional Variance Dynamics
## -----------------------------------
## GARCH Model : sGARCH(1,1)
## Mean Model : ARFIMA(0,0,0)
## Distribution : std
##
## Optimal Parameters
## ------------------------------------
## Estimate Std. Error t value Pr(>|t|)
## omega 0.000002 0.000001 1.7550 0.079252
## alpha1 0.072833 0.013315 5.4699 0.000000
## beta1 0.915632 0.015503 59.0634 0.000000
## shape 6.417288 0.674302 9.5169 0.000000
##
## Robust Standard Errors:
## Estimate Std. Error t value Pr(>|t|)
## omega 0.000002 0.000005 0.54239 0.58755
## alpha1 0.072833 0.052459 1.38839 0.16502
## beta1 0.915632 0.057452 15.93737 0.00000
## shape 6.417288 1.049277 6.11591 0.00000
##
## LogLikelihood : 9211.2
##
## Information Criteria
## ------------------------------------
##
## Akaike -5.9962
## Bayes -5.9884
## Shibata -5.9962
## Hannan-Quinn -5.9934
##
## Weighted Ljung-Box Test on Standardized Residuals
## ------------------------------------
## statistic p-value
## Lag[1] 0.2103 0.6465
## Lag[2*(p+q)+(p+q)-1][2] 0.2227 0.8404
## Lag[4*(p+q)+(p+q)-1][5] 1.3809 0.7692
## d.o.f=0
## H0 : No serial correlation
##
## Weighted Ljung-Box Test on Standardized Squared Residuals
## ------------------------------------
## statistic p-value
## Lag[1] 0.1004 0.7513
## Lag[2*(p+q)+(p+q)-1][5] 0.6013 0.9410
## Lag[4*(p+q)+(p+q)-1][9] 2.4433 0.8460
## d.o.f=2
##
## Weighted ARCH LM Tests
## ------------------------------------
## Statistic Shape Scale P-Value
## ARCH Lag[3] 0.4660 0.500 2.000 0.4948
## ARCH Lag[5] 0.8835 1.440 1.667 0.7678
## ARCH Lag[7] 2.0980 2.315 1.543 0.6966
##
## Nyblom stability test
## ------------------------------------
## Joint Statistic: 59.1102
## Individual Statistics:
## omega 4.7935
## alpha1 0.8342
## beta1 0.9662
## shape 0.5156
##
## Asymptotic Critical Values (10% 5% 1%)
## Joint Statistic: 1.07 1.24 1.6
## Individual Statistic: 0.35 0.47 0.75
##
## Sign Bias Test
## ------------------------------------
## t-value prob sig
## Sign Bias 0.1769 0.85963
## Negative Sign Bias 0.9611 0.33657
## Positive Sign Bias 1.7908 0.07343 *
## Joint Effect 6.8009 0.07852 *
##
##
## Adjusted Pearson Goodness-of-Fit Test:
## ------------------------------------
## group statistic p-value(g-1)
## 1 20 36.15 0.01011
## 2 30 42.89 0.04656
## 3 40 45.97 0.20595
## 4 50 52.40 0.34366
##
##
## Elapsed time : 0.2396719Dopasowany model GARCH-t(1,1) jest lepszy od pozostałych dotychczasowych modeli w sensie kryteriów informacyjnych. Tak jak dla pozostałych modeli odpowiednie testy odrzucają autokorelację standaryzowanych reszt i występowanie wśród nich efektów ARCH. Występuje jednak problem - parametr \(\omega\) jest na granicy istotności. Sprawdźmy jeszcze Q-Q plot, aby zweryfikować, czy nowy zakładany warunkowy rozkład reszt jest prawidłowy. 
Widzimy, że rozkład reszt z modelu GARCH-t(1,1) pozostaje w dobrej zgodności z założonym rozkładem. Na podstawie wyestymowanego modelu GARCH-t(1,1) możemy wnioskowac zatem, że dobrym pomysłem może być reestymacja wcześniej dopasowanych rozszerzeń (EGARCH i TGARCH), tym razem zakładając, że standaryzowane reszty podlegają rozkładowi t-Studenta.
Reestymacja modeli EGARCH i TGARCH
EGARCH
W pierwszym kroku ponownie wyestymujmy model EGARCH(1,1), uwzględniając model średniej ARFIMA(1,1).
##
## *---------------------------------*
## * GARCH Model Fit *
## *---------------------------------*
##
## Conditional Variance Dynamics
## -----------------------------------
## GARCH Model : eGARCH(1,1)
## Mean Model : ARFIMA(1,0,1)
## Distribution : std
##
## Optimal Parameters
## ------------------------------------
## Estimate Std. Error t value Pr(>|t|)
## ar1 0.46255 0.029516 15.6711 0
## ma1 -0.46307 0.029420 -15.7399 0
## omega -0.12713 0.007265 -17.4998 0
## alpha1 -0.07825 0.011115 -7.0402 0
## beta1 0.98512 0.000817 1205.7678 0
## gamma1 0.13845 0.019977 6.9306 0
## shape 6.80718 0.885747 7.6852 0
##
## Robust Standard Errors:
## Estimate Std. Error t value Pr(>|t|)
## ar1 0.46255 0.007851 58.9166 0
## ma1 -0.46307 0.007537 -61.4384 0
## omega -0.12713 0.006057 -20.9874 0
## alpha1 -0.07825 0.012570 -6.2254 0
## beta1 0.98512 0.000667 1475.8980 0
## gamma1 0.13845 0.025569 5.4149 0
## shape 6.80718 1.069792 6.3631 0
##
## LogLikelihood : 9232.051
##
## Information Criteria
## ------------------------------------
##
## Akaike -6.0078
## Bayes -5.9941
## Shibata -6.0079
## Hannan-Quinn -6.0029
##
## Weighted Ljung-Box Test on Standardized Residuals
## ------------------------------------
## statistic p-value
## Lag[1] 0.5095 0.4754
## Lag[2*(p+q)+(p+q)-1][5] 1.4215 0.9990
## Lag[4*(p+q)+(p+q)-1][9] 3.8853 0.7176
## d.o.f=2
## H0 : No serial correlation
##
## Weighted Ljung-Box Test on Standardized Squared Residuals
## ------------------------------------
## statistic p-value
## Lag[1] 0.261 0.6094
## Lag[2*(p+q)+(p+q)-1][5] 1.521 0.7348
## Lag[4*(p+q)+(p+q)-1][9] 3.253 0.7165
## d.o.f=2
##
## Weighted ARCH LM Tests
## ------------------------------------
## Statistic Shape Scale P-Value
## ARCH Lag[3] 0.8075 0.500 2.000 0.3689
## ARCH Lag[5] 1.8471 1.440 1.667 0.5060
## ARCH Lag[7] 2.6119 2.315 1.543 0.5896
##
## Nyblom stability test
## ------------------------------------
## Joint Statistic: 1.6101
## Individual Statistics:
## ar1 0.13981
## ma1 0.13978
## omega 0.83248
## alpha1 0.52778
## beta1 0.80384
## gamma1 0.03457
## shape 0.14262
##
## Asymptotic Critical Values (10% 5% 1%)
## Joint Statistic: 1.69 1.9 2.35
## Individual Statistic: 0.35 0.47 0.75
##
## Sign Bias Test
## ------------------------------------
## t-value prob sig
## Sign Bias 0.53729 0.5911
## Negative Sign Bias 0.06797 0.9458
## Positive Sign Bias 1.05894 0.2897
## Joint Effect 1.12658 0.7707
##
##
## Adjusted Pearson Goodness-of-Fit Test:
## ------------------------------------
## group statistic p-value(g-1)
## 1 20 35.45 0.01233
## 2 30 40.94 0.06969
## 3 40 53.34 0.06275
## 4 50 53.37 0.30988
##
##
## Elapsed time : 0.502692Widzimy, że otrzymano bardzo dobre dopasowanie - wszystkie uwzględnione parametry są istotne statystycznie, testy na autokorelację standaryzowanych reszt i ich kwadratów nie pozwalają przyjąć hipotezy zerowej o występowaniu autokorelacji, a test LM stwierdza, że efekty ARCH zostały wyjaśnione.
Brak widocznych zmian w ACF w stosunku do wcześniejszego modelu. Zobaczmy, czy krzywa wpływu informacji zmieniła się w stosunku do modelu EGARCH(1,1) przy założonym rozkładzie normalnym reszt.
Wykresy wykazują brak widocznych zmian.
TGARCH
Reestymujemy teraz model TGARCH(1,1) z modelem średniej ARIFMA(1,1), zakładając, że standaryzowane reszty mają rozkład t-Studenta.
##
## *---------------------------------*
## * GARCH Model Fit *
## *---------------------------------*
##
## Conditional Variance Dynamics
## -----------------------------------
## GARCH Model : fGARCH(1,1)
## fGARCH Sub-Model : TGARCH
## Mean Model : ARFIMA(1,0,1)
## Distribution : std
##
## Optimal Parameters
## ------------------------------------
## Estimate Std. Error t value Pr(>|t|)
## ar1 0.999339 0.000932 1.0722e+03 0.000000
## ma1 -0.996501 0.000005 -1.8935e+05 0.000000
## omega 0.000182 0.000068 2.6717e+00 0.007548
## alpha1 0.070350 0.010879 6.4668e+00 0.000000
## beta1 0.933483 0.012558 7.4331e+01 0.000000
## eta11 0.607644 0.016661 3.6471e+01 0.000000
## shape 6.730220 0.767256 8.7718e+00 0.000000
##
## Robust Standard Errors:
## Estimate Std. Error t value Pr(>|t|)
## ar1 0.999339 0.000919 1.0873e+03 0.000000
## ma1 -0.996501 0.000007 -1.5149e+05 0.000000
## omega 0.000182 0.000107 1.6896e+00 0.091107
## alpha1 0.070350 0.018748 3.7525e+00 0.000175
## beta1 0.933483 0.021515 4.3387e+01 0.000000
## eta11 0.607644 0.021237 2.8612e+01 0.000000
## shape 6.730220 0.787598 8.5452e+00 0.000000
##
## LogLikelihood : 9240.816
##
## Information Criteria
## ------------------------------------
##
## Akaike -6.0136
## Bayes -5.9998
## Shibata -6.0136
## Hannan-Quinn -6.0086
##
## Weighted Ljung-Box Test on Standardized Residuals
## ------------------------------------
## statistic p-value
## Lag[1] 0.6476 0.4210
## Lag[2*(p+q)+(p+q)-1][5] 1.8210 0.9827
## Lag[4*(p+q)+(p+q)-1][9] 4.5051 0.5700
## d.o.f=2
## H0 : No serial correlation
##
## Weighted Ljung-Box Test on Standardized Squared Residuals
## ------------------------------------
## statistic p-value
## Lag[1] 0.3524 0.5527
## Lag[2*(p+q)+(p+q)-1][5] 1.3136 0.7856
## Lag[4*(p+q)+(p+q)-1][9] 3.0215 0.7554
## d.o.f=2
##
## Weighted ARCH LM Tests
## ------------------------------------
## Statistic Shape Scale P-Value
## ARCH Lag[3] 0.9232 0.500 2.000 0.3366
## ARCH Lag[5] 1.6409 1.440 1.667 0.5561
## ARCH Lag[7] 2.5833 2.315 1.543 0.5954
##
## Nyblom stability test
## ------------------------------------
## Joint Statistic: 1.7309
## Individual Statistics:
## ar1 0.2707
## ma1 0.2603
## omega 1.0143
## alpha1 0.9364
## beta1 1.0749
## eta11 0.5927
## shape 0.5447
##
## Asymptotic Critical Values (10% 5% 1%)
## Joint Statistic: 1.69 1.9 2.35
## Individual Statistic: 0.35 0.47 0.75
##
## Sign Bias Test
## ------------------------------------
## t-value prob sig
## Sign Bias 0.3969 0.6914
## Negative Sign Bias 0.0406 0.9676
## Positive Sign Bias 0.9228 0.3562
## Joint Effect 0.8563 0.8360
##
##
## Adjusted Pearson Goodness-of-Fit Test:
## ------------------------------------
## group statistic p-value(g-1)
## 1 20 25.80 0.1360
## 2 30 38.40 0.1137
## 3 40 43.46 0.2870
## 4 50 52.10 0.3541
##
##
## Elapsed time : 1.232938Wszystkie otrzymane parametry są istotne statystycznie, wśród standaryzowanych reszt i ich kwadratów brak autokorelacji, test LM wskazuje na brak efektów ARCH. Porównując otrzymane wartości kryteriów informacyjnych z wartościami uzyskanymi w modelu TGARCH(1,1) z założonym normalnym rozkładem reszt, widoczna jest poprawa wszystkich ich wartości.
Porównajmy jeszcze krzywe wpływu informacji dla modeli o różnych zakładanych rozkładach standaryzowanych reszt. 
Otrzymane rezultaty dla obydwu reestymowanych modeli są bardzo obiecujące - zastosowanie rozszerzeń EGARCH i TGARCH pozwoliło uwzględnić asymetrię warunkowej wariancji, a połączenie tych dwóch modeli z założeniem o innym niż normalny rozkładzie standaryzowanych reszt pozwoliło polepszyc rezultaty dopasowania. W dalszym toku analizy do estymacji Value at Risk i analizy wrażliwości będziemy wykorzystywać więc 3 modele - EGARCH(1,1), TGARCH(1,1) (obydwa przy założeniu podlegania przez reszty rozkładowi t-Studenta i uwzględnieniu modelu średniej postaci ARIFMA(1,1)), oraz zwykły model GARCH(1,1).
Analiza VaR w okresach in sample oraz out of sample
Zakres notowań dla naszego portfela zaczyna się 02.03.2005 r. i kończy 31.08.2017 r. Jego konstrukcja była pomyślana z nastawieniem raczej na zysk długoterminowy niż na szybką sprzedaż. W związku z tym postanowiono wydzielić zbiory in sample oraz out of sample według następującej metodologii. Pierwszy ze zbiorów będzie obejmował okres od początku zakresu notowań do 31.08.2016 r. Chcemy w ten sposób uwzględnić kilka charakterystycznych momentów dla każdego z szeregów: krach finansowy z 2008 roku (duży klaster wariancji na wykresie stóp zwrotu z portfela ) - efekt dobrze widoczny dla Nvidii, Microsoftu oraz kursu ceny złota do dolara, słabiej natomiast dla Amazona oraz Statoila oraz początek szybkiego wzrostu kursu spółki Amazon oraz Nvidia na początku roku 2016. Drugi zbiór, out of sample, będzie rozciągał się na okres 01.09.2016 do 31.08.2017. W celu przeprowadzenia dalszych analiz musimy ponownie oszacować wymienione wczesniej modele: EGARCH(1,1), TGARCH(1,1) oraz GARCH(1,1) na próbce in sample. Ponadto przed estymacją wartości narażonej na ryzyko musimy zadbać o standaryzację naszych stóp zwrotu.
Poniższy wykres przedstawia oszacowane wartości narażone na ryzyko wraz z szeregiem stóp zwrotu dla każdego z wybranych wcześniej modeli.
Wykres sugeruje zgodne z oczekiwaniami zachowanie wartośći VaR w momentach największej zmienności stóp zwrotu (zwłaszcza rok 2008). Poniższa tabela zestawia odsetki oraz zliczenia przekroczeń wartości narażonej na ryzyko dla próbki in sample przy poziomie ufności 99%:
## EGRACH GARCH11 TGARCH
## Liczba dni przekraczających VaR 14.000000000 17.000000000 15.000000000
## Odsetek przekroczeń 0.004558776 0.005535656 0.004884402Wyniki pokazują, że najbradziej restrykcyjny okazał się model E-GARCH (0.45% przekroczeń VaR), a najmniej restrykcyjny model GARCH11 (0.55% przekroczeń).
W kolejnym kroku przetestujmy zachowanie VaR w próbce out of sample.
Wykres pokazuje, że modele EGARCH oraz TGARCH praktycznie identycznie estymują VaR. Obserwujemy rówież zastanawiające opóźnienie w dynamice VaR wszystkich modeli względem rzeczywistej zmienności stóp zwrotu na początku listopada. Tabela niżej zestawia zliczenia oraz częstości niedoszacowań wartości narażonej na ryzyko:
## EGRACH GARCH11 TGARCH
## Liczba dni przekraczających VaR 2.000000000 1.000000000 2.000000000
## Odsetek przekroczeń 0.008130081 0.004065041 0.008130081Na próbce out of sample modele wypadły bardzo podobnie: EGARCH oraz TGARCH mają po dwa przekroczenia wartości narażonej na ryzyko. Wyznaczone odsetki są jednak znacznie poniżej założonej wartości 1% niedoszacowanych przypadków. Modele wydają się być za bardzo restrykcyjne.
Analiza wrażliwości
Dla wszystkich 3 wyestymowanych wczesniej modeli w ostatnim kroku postanowiono przeprowadzić analizę wpływu wybranego do estymacji okna czasowego na zdolność modelu do tworzenia prognoz out of sample. Zacznijmy od wydłużenia okna estymacji i stworzenia krótszej niż wcześniejsza prognozy, tworząc sześciomiesięczny okres out of sample, równy połowie wykorzystanego wcześniej okresu.
Wykres zestawia przebieg zwrotów ze skonstruowanego portfela z obliczonym w półrocznym okresie out of sample VaR dla wybranych modeli. Wykresy pozostają zgodne z wcześniejszymi rezultatami - estymaty VaR dla modeli EGARCH(1,1) i TGARCH(1,1) są do siebie bardzo podobne (mają niemal identyczny kształt, jednak estymata z EGARCH(1,1) jest nieco bardziej konserwatywna). Z kolei model GARCH(1,1) w niektórych miejscach nieco lepiej reagował na zmiany wariancji stóp zwrotu, lecz ogólnie rzecz biorąc był bardziej restrykcyjny - dla modeli TGARCH i E-GARCH obserwujemy jedną sytuację przekroczenia przez zwroty zakładanej wartości VaR i dwie sytuacje, w których stopy zwrotu zbliżają się do wyestymowanej wartości VaR, zaś dla modelu GARCH we wszystkich tych 3 sytuacjach przebieg VaR pozostaje daleko od stóp zwrotu. Warto zwrócić uwagę, że wciąż widoczny jest efekt delikatnego opóźnienia wyliczonego VaR względem rzeczywistej wariancji. Obliczmy, jaką wartość procentową wszystkich punktów na wykresie stanowi jedno przekroczenie.
## EGARCH GARCH11 TGARCH
## Liczba dni przekraczających VaR 1.000000000 0 1.000000000
## Odsetek przekroczeń 0.007874016 0 0.007874016Widać, że jak wskazywała analiza wizualna, otrzymaliśmy po jednym przekroczeniu dla modeli EGARCH(1,1) i TGARCH(1,1), oraz 0 przekroczeń dla GARCH(1,1). Zatem przyjmując półroczny okres out of sample dla żadnego z modeli zatem nie przekraczamy założonego VaR wynoszącego 1%.
Pójdźmy teraz w drugą stronę i zamiast skracać wcześniej założony roczny okres out of sample o połowę, wydłużmy go o taką samą wartość, otrzymując okres póltoraroczny.
Wyrysowany wykres przedstawia podobne do wcześniejszych rezultaty. Wciąż widoczne jest małe (lecz wyraźne) przesunięcie VaR względem szoków, a wzajemna relacja między trzema modelami jest utrzymana - modele EGARCH i TGARCH zachowują się bardzo podobnie, przy czym EGARCH jest bardziej konserwatywny, zaś GARCH najbardziej przeszacowuje VaR. Sprawdźmy procentową liczbę przekroczeń:
## EGARCH GARCH11 TGARCH
## Liczba dni przekraczających VaR 3.000000000 2.00000000 3.000000000
## Odsetek przekroczeń 0.008042895 0.00536193 0.008042895Widzimy, że wydłużając okres uzyskaliśmy więcej przekroczeń. Procentowe wartości dla wszystkich modeli wzrosły, lecz nie przekroczyły założonej wartości VaR (1%). Możemy wnioskować, że wraz z wydłużaniem okresu prognozy out of sample będziemy stopniowo zbliżać się do tej wartości. Przeprowadźmy jeszcze jedną prognozę out of sample na znacznie dłuższy okres (~6 lat). Ze względu na trudności techniczne przy wydłużaniu okresu out of sample dla założonego modelu TGARCH(1,1) na potrzeby tej prognozy wykorzystano jedynie modele EGARCH i GARCH. Ze względu na bardzo duże podobieństwo zachowania modeli TGARCH i EGARCH dla dotychczas wykorzystanych okresów out of sample możemy jednak podejrzewać, że na dłuższej próbce taka zależność utrzymałaby się.
Analiza powyższego wykresu wskazuje, że w dłuższym horyzoncie czasowym prognozy out of sample warunkowej wariancji z modeli GARCH(1,1) i EGARCH(1,1) stają się bardziej podobne do siebie niż dla krótszych okresów. Widać także, że przebieg zwrotów znacznie częściej zbliża się do estymowanego VaR - jest to zgodne z podejrzeniem wysuniętym przy poprzedniej estymacji, mówiącym, że przy dłuższych horyzontach czasowych będziemy coraz bardziej zbliżać się do założonego poziomu VaR. Sprawdźmy przekroczenia:
## EGARCH GARCH11
## Liczba dni przekraczających VaR 6.000000000 6.000000000
## Odsetek przekroczeń 0.004070556 0.004070556Widzimy, że procentowo otrzymujemy mniej przekroczeń niż dla poprzedniego okresu out of sample, jednak tak jak już wspomniano, z analizy wizualnej wynika, że znacznie częściej występuje sytuacja, w której przebieg zwrotów zbliża się do założonego poziomu VaR, lecz go nie przekracza.
Podsumowanie
W niniejszym projekcie udało się skonstruowac portfel zrównoważony złożony z notowań 3 spółek technologicznych notowanych na NASDAQ, spółki z sektora petrochemicznego notowanego na NYSE, oraz kursu złota do dolara amerykańskiego. Na stopach zwrotu tego portfela wyestymowano kilka różnych modeli z rodziny GARCH. Najlepsze okazały się modele GARCH(1,1), oraz dwa jego rozszerzenia, EGARCH(1,1) i TGARCH(1,1). Obydwa rozszerzenia zawierały ARIFMA(1,1) jako model średniej. Dodatkowo w oparciu o estymację modelu GARCH-t(1,1), zdecydowano się do modeli EGARCH(1,1) i TGARCH(1,1) dodać założenie warunkowego rozkładu standaryzowanych reszt postaci t-Studenta. Portfel konstruowano z nastawieniem na uzyskanie długoterminowego zysku, w związku z czym wstępnie ustalono roczny okres out of sample, którym później manipulowano i porównywano rezultaty. Tak przeprowadzona analiza wrażliwości sugerowała, że wszystkie 3 modele były dość restrykcyjne, jednak przy dłuższych prognozach rzeczywista wartość narażona na ryzyko dążyła do założonej.