Wprowadzenie

Niniejsza praca stanowi raport z przeprowadzonej analizy porównawczej ryzyka, oszacowanego za pomocą modeli klasy GARCH. Szeregi czasowe charakteryzują się wieloma cechami (m.in. efekt grubych ogonów, leptokurtyczność rozkładu, efekty grupowania wariancji, efekt dźwigni, efekt autokorelacji stóp zwrotu), których uwzględnienie stanowi wyzwanie dla osób próbujących je modelować. Dobór odpowiedniego modelu jest zadaniem trudnym, gdyż każdy z nich zawiera pewne ograniczenia. Badania empiryczne pokazują, że modele klasy GARCH(p,q) lepiej dopasowują się do danych empirycznych niż ARCH. Ponadto dodatkową zaletą jest m.in. lepszy od ARCH opis rozkładów o grubych ogonach. Niemniej jednak modele GARCH(p,q) nie uwzględniają asymetrycznego wpływu dobrych i złych wiadomości na rynek i nie dają możliwości modelowania efektu dźwigni. Aby zaradzić niektórym wadom modelu GARCH(p,q) skonstruowano nowe klasy – np. modele EGARCH, GARCH-in-mean, TGARCH.

1. Metodologia

1.1. Skład zrównoważonego portfela akcji

Powszechne Towarzystwa Emerytalne (PTE) są instytucjami funkcjonującymi w formie spółek akcyjnych, powoływanymi przez akcjonariuszy w celu zarządzania Otwartymi Funduszami Emerytalnymi (OFE). Zarządzanie otwartym funduszem emerytalnym przez powszechne towarzystwo emerytalne sprowadza się do organizacji systemu ewidencji członków OFE i aktywów funduszu oraz zarządzania procesem inwestowania aktywów funduszu. Ubezpieczony samodzielnie wybiera Powszechne Towarzystwo Emerytalne, które będzie prowadzić jego otwarty fundusz emerytalny (na który trafia od 2014 r. 2,92% składki na ubezpieczenie emerytalne). Aktualnie na polskim rynku funkcjonuje 12 towarzystw emerytalnych, które zarządzają środkami w łącznej wysokości 172,6 mld zł (wg danych KNF na dzień 31.05.2017 r.). PTE inwestują te środki w instrumenty finansowe dopuszczone przez obowiązujące przepisy i w granicach dopuszczonych limitów. PTE osiągają różne wyniki inwestycyjne, w zależności od prowadzonej przez siebie strategii inwestycyjnej i jej skuteczności. Zgodnie z opublikowaną przez KNF Informacją dotyczącą otwartych funduszy emerytalnych z dnia 27.04.2017 r., najwyższą średnią ważoną stopę zwrotu w okresie od 31.03.2016 r. do 31.03.2017 r. osiągnął MetLife OFE (19,5%), najniższą zaś AEGON OFE (14,2%). MetLife OFE jest również liderem spośród wszystkich pozostałych PTE dla okresu trzyletniego (średnia ważona stopa zwrotu na poziomie 17,4%).
Urząd Komisji Nadzoru Finansowego co roku na swojej stronie publikuje szczegółowe dane dotyczące wszystkich PTE. Do wiadomości publicznej podane są nazwy spółek, w które inwestują poszczególne towarzystwa emerytalne, wraz z wielkością tych zaangażowań. W niniejszej pracy wykorzystano dane UKNF za rok 2015. Na ich podstawie oraz biorąc pod uwagę skład portfela benchmarkowego OFE (tj. przyjęto, iż jest to portfel będący sumą wszystkich portfeli akcji PTE), zbudowano portfel zrównoważony. Zdecydowano, że w jego skład wejdzie 20 spółek, w które wszystkie PTE łącznie najwięcej zainwestowały (tj. spółki, które mają największy udział w portfelu benchmarkowym). Tym sposobem wybrano akcje następujących spółek: PKO BP, PKN Orlen, PZU, PEKAO, BZ WBK, ING Bank Śląski, LPP, PGNIG, Mbank, PGE, KGHM, Grupa Kęty, Bank Millenium, Asseco Poland, Grupa Azoty, Cyfrowy Polsat, Amrest Holdings, CCC, Inter Cars oraz Bank Handlowy. Łącznie wartość inwestycji w te spółki przez wszystkie PTE wyniosła 64,7% całkowitej wartości aktywów OFE.

lista <- matrix(c("pko", "pkn", "pzu", "peo", "bzw", "ing", "lpp", "pgn",
                  "mbk", "pge", "kgh", "kty", "mil", "acp", "att", "cps",
                  "eat", "ccc", "car", "bhw" ),ncol=1,byrow=TRUE)
colnames(lista) <- c("s_stock")

1.2. Dane i okres analizy

Dane dotyczące cen zamknięcia ww. spółek zostały pobrane ze strony stooq.com. Następnie zostały obliczone ciągłe stopy zwrotu dla 20 wybranych podmiotów. Stopa zwrotu dla portfela zrównoważonego w danym dniu była równa 1/20 sumy stóp zwrotu wszystkich spółek w danym dniu. W kolejnym kroku zdefiniowano zbiór dla okresu in-sample (od 01.01.2010 r. do 31.12.2014 r.) oraz out-sample (od 01.01.2015 r. do 23.06.2017 r.).

url1 <- "http://stooq.com/q/d/l/?s="
url2<-lista[1,]
url3<-"&i=d"
url<-paste(url1, url2, url3, sep = "")
total <- read.csv(url,
                  header = TRUE,
                  sep = ",",
                  dec = ".",
                  stringsAsFactors = F)
total$Date <- as.Date(total$Date, format = "%Y-%m-%d")
total <- total[, c("Date", "Close")]
colnames(total) <- c("Date", lista[1,])

 for(i in 2:nrow(lista)){
    url1 <- "http://stooq.com/q/d/l/?s="
    url2<-lista[i,]
    url3<-"&i=d"
    url<-paste(url1, url2, url3, sep = "")
    stock <- read.csv(url,
                      header = TRUE,
                      sep = ",",
                      dec = ".",
                      stringsAsFactors = F)
    stock$Date <- as.Date(stock$Date, format = "%Y-%m-%d")
    stock <- stock[, c("Date", "Close")]
    colnames(stock) <- c("Date", lista[i,])
    
    total<-merge(total,stock,by="Date",all=TRUE)

 }
for(i in 1:nrow(lista)){
  zestaw[,i+1]<- diff.xts(log(zestaw[,i+1]))
}
zestaw$bench<-rowSums(zestaw[,2:ncol(zestaw)],na.rm=TRUE)/(ncol(zestaw)-1)


# In sample
zestaw_in<-zestaw[!(zestaw$Date<"2010-01-01"),]
zestaw_in<-zestaw_in[!(zestaw_in$Date>"2014-12-31"),]

#out sample
zestaw_out<-zestaw[!(zestaw$Date<"2015-01-01"),]
zestaw_out<-zestaw_out[!(zestaw_out$Date>"2017-06-23"),]

Poniższy wykres przedstawia stopy zwrotu portfela akcji. Można zauważyć efekty grupowania wariancji, który zostanie bardziej przeanalizowany w dalszej części pracy.

1.3. Charakterystyka badanego szeregu czasowego

1.3.1. Normalnolność rozkładu

Kurtoza wynosząca 4,59 wskazuje, że występuje duża koncentracja wokół średniej, a zatem rozkład najprawdopodobniej nie będzie rozkładem normalnym. Przypuszczenia potwierdza histogram wraz z zaznaczoną krzywą gaussowską, która znacznie odbiega od otrzymanego empirycznego rozkładu stóp zwrotu (jest on smuklejszy od krzywej rozkładu normalnego). Brak normalności rozkładu potwierdza formalnie przeprowadzony test Jarque-Bera. Na poziomie istotności 5% odrzucamy hipotezę zerową zakładającą, że rozkład stóp zwrotu portfela akcji jest normalny.

##             X..zestaw_in.bench
## nobs               1249.000000
## NAs                   0.000000
## Minimum              -0.061634
## Maximum               0.046154
## 1. Quartile          -0.004330
## 3. Quartile           0.006248
## Mean                  0.000527
## Median                0.000803
## Sum                   0.658372
## SE Mean               0.000294
## LCL Mean             -0.000049
## UCL Mean              0.001103
## Variance              0.000108
## Stdev                 0.010374
## Skewness             -0.720391
## Kurtosis              4.574632

## 
##  Jarque Bera Test
## 
## data:  zestaw_in$bench
## X-squared = 1203.2, df = 2, p-value < 2.2e-16

1.3.2. Efekt grubych ogonów

Na poniższym wykresie kwantylowym można zaobserwować, że prawdopodobieństwa wystąpienia obserwacji znacznie oddalonych od średniej wartości (nietypowych) są znacznie wyższe niż dla rozkładu normalnego. Wykres stanowi graficzny test zgodności z przeprowadzonym w poprzednim podpunkcie testem na brak rozkładu normalnego stóp zwrotu portfela.

1.3.3. Leptokurtyczność

Występowanie wyższego szczytu empirycznej funkcji gęstości niż dla rozkładu normalnego (większe skupienie wyników wokół średniej), a także zjawiska grubych ogonów wskazuje na leptokurtyczność.

1.3.4. Efekt grupowania wariancji

Jak zostało wcześniej wspomniane, zaprezentowany wykres stóp zwrotu z portfela akcji wskazuje na występowanie efektów grupowania wariancji (po sobie następują okresy nasilonej zmienności i względnie stabilnych). W sytuacji, gdy w modelu występuje zjawisko grupowania wariancji, to wariancja składnika losowego dla różnych obserwacji jest ze sobą powiązana. Statystyka Durbina-Watsona dla kwadratów stóp zwrotu dla pięciu opóźnień jest mniejsza niż 1,8, zaś p-value=0. Zatem na poziomie istotności 5% należy odrzucić hipotezę zerową, zakładającą brak autokorelacji wśród kwadratów stóp zwrotu. Ponadto przeprowadzony test LM silnie odrzuca hipotezę zerową zakładającą brak występowania efektów ARCH. Wyniki dla stóp zwrotu:

##  lag Autocorrelation D-W Statistic p-value
##    1    0.1024183837      1.793750   0.000
##    2   -0.0452243433      2.088187   0.116
##    3   -0.0696575803      2.136864   0.018
##    4   -0.0391582232      2.075121   0.184
##    5    0.0305017320      1.934071   0.336
##    6    0.0030627307      1.985633   0.880
##    7    0.0753078335      1.839749   0.004
##    8   -0.0337568303      2.055497   0.234
##    9    0.0002542923      1.986361   0.976
##   10    0.0028840194      1.980602   0.964
##  Alternative hypothesis: rho[lag] != 0

Wyniki dla kwadratów stóp zwrotu:

##  lag Autocorrelation D-W Statistic p-value
##    1       0.2449761      1.509840   0.000
##    2       0.1745892      1.650556   0.000
##    3       0.2055989      1.588349   0.000
##    4       0.1721980      1.655076   0.000
##    5       0.2713925      1.456436   0.000
##    6       0.1018981      1.794834   0.018
##    7       0.1750262      1.648483   0.000
##    8       0.1675293      1.663270   0.000
##    9       0.1462243      1.705772   0.006
##   10       0.1222245      1.753620   0.012
##  Alternative hypothesis: rho[lag] != 0
## 
##  ARCH LM-test; Null hypothesis: no ARCH effects
## 
## data:  zestaw_in$bench
## Chi-squared = 170.48, df = 5, p-value < 2.2e-16

Poniższe wykresy przedstawiają wartości funkcji ACF odpowiednnio dla stóp zwrotu i kwadratów stóp zwrotu. Analiza graficzna jest zgodna z powyższymi wynikami na obecność autokorelacji wśród stóp zwrotu i ich kwadratów.

2. Model GARCH(p,q)

2.1. Dobór parametrów modelu

W pierwszym kroku modelowania GARCH(p,q) sprawdzono dopasowanie najbardziej podstawowego modelu GARCH(1,1). Na poziomie istotności 5% wszystkie oszacowane parametry są istotne (dla każdego współczynnika odrzucamy hipotezę zerową zakładającą jego nieistotność). Oszacowana w modelu parametr \(\alpha_1\) wyniósł‚ 0,07, zaś współczynnik \(\beta_1\) 0,90. Wartość kryteriów informacyjnych AIC i BIC wyniosła odpowiednio -6,51 i -6,50. Dla wizualizacji różnicy w dopasowaniu modelu przedstawiono wykres zwrotów oraz warunkowe odchylenie standardowe.

Widzimy, że model obejmuje prawie wszystkie obserwacje. Kolejnym krokiem było dokonanie analizy standaryzowanych reszt modelu. Na podstawie testu Jarque-Bera odrzucamy hipotezę zerową zakładającą, że rozkład standaryzowanych reszt jest normalny. Sprawdzono również efekt autokorelacji wśród reszt.

## 
## Title:
##  GARCH Modelling 
## 
## Call:
##  garchFit(formula = ~garch(1, 1), data = na.omit(zestaw_in$bench), 
##     cond.dist = "norm", include.mean = F, trace = FALSE) 
## 
## Mean and Variance Equation:
##  data ~ garch(1, 1)
## <environment: 0x000000001bb51828>
##  [data = na.omit(zestaw_in$bench)]
## 
## Conditional Distribution:
##  norm 
## 
## Coefficient(s):
##      omega      alpha1       beta1  
## 2.1706e-06  7.8317e-02  9.0090e-01  
## 
## Std. Errors:
##  based on Hessian 
## 
## Error Analysis:
##         Estimate  Std. Error  t value Pr(>|t|)    
## omega  2.171e-06   7.838e-07    2.769  0.00562 ** 
## alpha1 7.832e-02   1.294e-02    6.053 1.42e-09 ***
## beta1  9.009e-01   1.649e-02   54.617  < 2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Log Likelihood:
##  4069.762    normalized:  3.258416 
## 
## Description:
##  Thu Jun 29 23:15:13 2017 by user: Ania 
## 
## 
## Standardised Residuals Tests:
##                                 Statistic p-Value     
##  Jarque-Bera Test   R    Chi^2  235.3052  0           
##  Shapiro-Wilk Test  R    W      0.9805204 6.268967e-12
##  Ljung-Box Test     R    Q(10)  18.80091  0.04286555  
##  Ljung-Box Test     R    Q(15)  24.0383   0.06444597  
##  Ljung-Box Test     R    Q(20)  26.48293  0.1504479   
##  Ljung-Box Test     R^2  Q(10)  13.68066  0.1880624   
##  Ljung-Box Test     R^2  Q(15)  16.13449  0.3731703   
##  Ljung-Box Test     R^2  Q(20)  19.72959  0.4749551   
##  LM Arch Test       R    TR^2   16.57005  0.1664989   
## 
## Information Criterion Statistics:
##       AIC       BIC       SIC      HQIC 
## -6.512028 -6.499706 -6.512040 -6.507396

Na podstawie zaprezentowanych powyżej wyników można przyjąć, że nie występuje istotna autokorelacja reszt, ani ich kwadratów. Poniższa tabela prezentuje wartości kryteriów informacyjnych dla różnych modeli GARCH(p,q).

##               Akaike     Bayes   Shibata Hannan-Quinn
## GARCH(1,1) -6.513714 -6.501392 -6.513725    -6.509081
## GARCH(1,2) -6.522379 -6.505950 -6.522400    -6.516203
## GARCH(1,3) -6.519172 -6.498635 -6.519204    -6.511451
## GARCH(2,1) -6.512019 -6.495589 -6.512039    -6.505842
## GARCH(2,2) -6.520778 -6.500241 -6.520810    -6.513057
## GARCH(2,3) -6.517571 -6.492926 -6.517616    -6.508306
## GARCH(3,1) -6.510374 -6.489837 -6.510406    -6.502653
## GARCH(3,2) -6.519097 -6.494452 -6.519142    -6.509832
## GARCH(3,3) -6.516137 -6.487386 -6.516200    -6.505328

Wszystkie kryteria informacyjne wskazują, że najbardziej odpowiednim modelem jest jednak model GARCH(1,2) (najmniejsze wartości kryteriów informacyjnych). Wobec powyższego postanowiono sprawdzić dopasowanie modelu GARCH(1,2) do danych.

2.2. Model GARCH(1,2)

Na poziomie istotności 5% tylko parametr \(\beta_1\) jest nieistotny. Oszacowany w modelu parametr \(\alpha_1\) wyniósł‚ 0,13, zaś współczynnik \(\beta_2\)=0,73. Równanie warunkowej wariancji modelu GARCH(1,2) jest wobec powyższego następujące:

\[ \sigma^2_t=\omega+ 0.13 u^2_{t-1} + 0.73 \sigma^2_{t-2} \]

## 
## Title:
##  GARCH Modelling 
## 
## Call:
##  garchFit(formula = ~garch(1, 2), data = na.omit(zestaw_in$bench), 
##     cond.dist = "norm", include.mean = F, trace = FALSE) 
## 
## Mean and Variance Equation:
##  data ~ garch(1, 2)
## <environment: 0x000000001d9bb7c0>
##  [data = na.omit(zestaw_in$bench)]
## 
## Conditional Distribution:
##  norm 
## 
## Coefficient(s):
##      omega      alpha1       beta1       beta2  
## 4.0027e-06  1.3109e-01  9.8301e-02  7.3058e-01  
## 
## Std. Errors:
##  based on Hessian 
## 
## Error Analysis:
##         Estimate  Std. Error  t value Pr(>|t|)    
## omega  4.003e-06   1.338e-06    2.992  0.00277 ** 
## alpha1 1.311e-01   2.014e-02    6.508  7.6e-11 ***
## beta1  9.830e-02   6.283e-02    1.565  0.11767    
## beta2  7.306e-01   6.239e-02   11.710  < 2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Log Likelihood:
##  4076.287    normalized:  3.26364 
## 
## Description:
##  Thu Jun 29 23:15:13 2017 by user: Ania 
## 
## 
## Standardised Residuals Tests:
##                                 Statistic p-Value     
##  Jarque-Bera Test   R    Chi^2  142.6487  0           
##  Shapiro-Wilk Test  R    W      0.9850023 4.662234e-10
##  Ljung-Box Test     R    Q(10)  17.60315  0.06203839  
##  Ljung-Box Test     R    Q(15)  23.07116  0.08264061  
##  Ljung-Box Test     R    Q(20)  25.50957  0.1826269   
##  Ljung-Box Test     R^2  Q(10)  6.89571   0.7352525   
##  Ljung-Box Test     R^2  Q(15)  9.637079  0.8419266   
##  Ljung-Box Test     R^2  Q(20)  13.95751  0.8326438   
##  LM Arch Test       R    TR^2   9.114172  0.693149    
## 
## Information Criterion Statistics:
##       AIC       BIC       SIC      HQIC 
## -6.520875 -6.504446 -6.520896 -6.514699

Dla wizualizacji różnicy w dopasowaniu modelu przedstawiono wykres zwrotów oraz warunkowe odchylenie standardowe.

Kolejnym krokiem było dokonanie analizy standaryzowanych reszt modelu. Na podstawie testu Jarque-Bera odrzucamy hipotezę zerową zakładającą, że rozkład standaryzowanych reszt jest normalny. Sprawdzono również efekt autokorelacji wśród reszt.

Na podstawie zaprezentowanych powyżej wyników można przyjąć, że nie występuje istotna autokorelacja reszt, ani ich kwadratów. Wśród reszt nie ma efektów ARCH.

##                  stdres
## nobs        1249.000000
## NAs            0.000000
## Minimum       -4.700939
## Maximum        3.391247
## 1. Quartile   -0.521080
## 3. Quartile    0.670073
## Mean           0.058283
## Median         0.088025
## Sum           72.796025
## SE Mean        0.028248
## LCL Mean       0.002865
## UCL Mean       0.113702
## Variance       0.996629
## Stdev          0.998313
## Skewness      -0.391181
## Kurtosis       1.451453
## 
##  Jarque Bera Test
## 
## data:  stdres
## X-squared = 142.65, df = 2, p-value < 2.2e-16
##  lag Autocorrelation D-W Statistic p-value
##    1     0.078935378      1.840724   0.000
##    2    -0.019061689      2.035669   0.564
##    3    -0.008428382      2.014201   0.794
##    4    -0.004271754      2.005017   0.886
##    5     0.039258114      1.916339   0.188
##  Alternative hypothesis: rho[lag] != 0
## 
##  ARCH LM-test; Null hypothesis: no ARCH effects
## 
## data:  stdres
## Chi-squared = 3.8914, df = 5, p-value = 0.5652

Do dalszej analizy został wybrany model GARCH(1,2).

3. Rozszerzenia modelu GARCH

3.1. Model EGARCH

Autorki postanowiły oszacować model EGARCH, chcąc sprawdzić, czy dla badanych danych wystepuje efekt dźwigni, czyli efekt asymetrycznego wpływu informacji pozytywnych i negatywnych na poziom przyszłej wariancji. Efekt ten objawia się na przykład w ten sposób, że spadek ceny instrumentu powoduje wzrost wariancji stóp zwrotu. Model ten ma pewne zalety względem modelu GARCH, m.in. możliwości występowania ujemnych współczynników i modelowania efektu dźwigni.

Oszacowano 9 modeli EGARCH, po czym porównano dla nich wartości kryteriów informacyjnych, które przedstawione są w poniższej tabeli.

##                Akaike     Bayes   Shibata Hannan-Quinn
## EGARCH(1,1) -6.513714 -6.501392 -6.513725    -6.509081
## EGARCH(1,2) -6.522379 -6.505950 -6.522400    -6.516203
## EGARCH(1,3) -6.519172 -6.498635 -6.519204    -6.511451
## EGARCH(2,1) -6.512019 -6.495589 -6.512039    -6.505842
## EGARCH(2,2) -6.520778 -6.500241 -6.520810    -6.513057
## EGARCH(2,3) -6.517571 -6.492926 -6.517616    -6.508306
## EGARCH(3,1) -6.510374 -6.489837 -6.510406    -6.502653
## EGARCH(3,2) -6.519097 -6.494452 -6.519142    -6.509832
## EGARCH(3,3) -6.516137 -6.487386 -6.516200    -6.505328

Na podstawie powyższych wartości kryteriów informacyjnych wybrano model EGARCH(1,2) i go oszacowano. Dla tej wersji modelu funkcja warunkowej wariancji określona jest następującym równaniem:

\[ ln(\sigma^2_t)=\omega+ \alpha_1 ln(\sigma^2_{t-1} ) + \beta_1 g(Z_{t-1}) +\beta_2 g(Z_{t-2}) \] gdzie \[g(Z_{t-1})= \theta Z_{t-1} + \lambda (|Z_{t-1}|-E(|Z_{t-1}|)), \] \[g(Z_{t-2})= \theta Z_{t-2} + \lambda (|Z_{t-2}|-E(|Z_{t-2}|)). \]

Przeprowadzono estymację.

## 
## *---------------------------------*
## *          GARCH Model Fit        *
## *---------------------------------*
## 
## Conditional Variance Dynamics    
## -----------------------------------
## GARCH Model  : eGARCH(1,2)
## Mean Model   : ARFIMA(0,0,0)
## Distribution : norm 
## 
## Optimal Parameters
## ------------------------------------
##         Estimate  Std. Error  t value Pr(>|t|)
## omega   -0.43046    0.119327  -3.6074 0.000309
## alpha1  -0.12704    0.023624  -5.3776 0.000000
## beta1    0.18304    0.020620   8.8767 0.000000
## beta2    0.76930    0.020227  38.0337 0.000000
## gamma1   0.25532    0.034052   7.4979 0.000000
## 
## Robust Standard Errors:
##         Estimate  Std. Error  t value Pr(>|t|)
## omega   -0.43046    0.138858  -3.1000 0.001935
## alpha1  -0.12704    0.034237  -3.7107 0.000207
## beta1    0.18304    0.014361  12.7455 0.000000
## beta2    0.76930    0.015647  49.1646 0.000000
## gamma1   0.25532    0.049230   5.1863 0.000000
## 
## LogLikelihood : 4089.103 
## 
## Information Criteria
## ------------------------------------
##                     
## Akaike       -6.5398
## Bayes        -6.5193
## Shibata      -6.5398
## Hannan-Quinn -6.5321
## 
## Weighted Ljung-Box Test on Standardized Residuals
## ------------------------------------
##                         statistic  p-value
## Lag[1]                      8.652 0.003267
## Lag[2*(p+q)+(p+q)-1][2]     8.972 0.003575
## Lag[4*(p+q)+(p+q)-1][5]     9.508 0.012297
## d.o.f=0
## H0 : No serial correlation
## 
## Weighted Ljung-Box Test on Standardized Squared Residuals
## ------------------------------------
##                          statistic p-value
## Lag[1]                      0.2318  0.6302
## Lag[2*(p+q)+(p+q)-1][8]     2.0645  0.8489
## Lag[4*(p+q)+(p+q)-1][14]    4.0584  0.8746
## d.o.f=3
## 
## Weighted ARCH LM Tests
## ------------------------------------
##             Statistic Shape Scale P-Value
## ARCH Lag[4]    0.6271 0.500 2.000  0.4284
## ARCH Lag[6]    2.0304 1.461 1.711  0.4841
## ARCH Lag[8]    2.5600 2.368 1.583  0.6286
## 
## Nyblom stability test
## ------------------------------------
## Joint Statistic:  0.9747
## Individual Statistics:              
## omega  0.08210
## alpha1 0.43451
## beta1  0.08369
## beta2  0.08373
## gamma1 0.06216
## 
## Asymptotic Critical Values (10% 5% 1%)
## Joint Statistic:          1.28 1.47 1.88
## Individual Statistic:     0.35 0.47 0.75
## 
## Sign Bias Test
## ------------------------------------
##                    t-value   prob sig
## Sign Bias           0.1121 0.9108    
## Negative Sign Bias  1.0040 0.3156    
## Positive Sign Bias  0.8660 0.3867    
## Joint Effect        3.1277 0.3723    
## 
## 
## Adjusted Pearson Goodness-of-Fit Test:
## ------------------------------------
##   group statistic p-value(g-1)
## 1    20     39.09     0.004305
## 2    30     43.69     0.039259
## 3    40     60.82     0.014219
## 4    50     73.70     0.012790
## 
## 
## Elapsed time : 0.1922469

Na podstawie uzyskanych wartości statysk oraz odpowiadających im p-value można dla każdego z parametru odrzucić hipotezę o jego nieistotności. W związku z tym każdy parametr jest oddzielnie istotny.

Na podstawie uzyskanych oszacowań możemy wnioskować, że występuje efekt asymetrycznego wpływu informacji, gdyż \(\alpha\) jest istotna i przyjmuje wartość uejmną, więc mamy do czynienia z ujemnym szokiem.

Otrzymane wyniki można przedstawić na wykresie ACF wystandaryzowanych reszt:

Na podstawie powyższych wyników możemy mówić o asymetrii reakcji funkcji warunowej wariancji na wiadomości napływające na rynek.

3.2. Model GARCH-in-mean

Model MGARCH został stworzony w celu modelowania liniowej zależności pomiędzy oczekiwaną stopą zwrotu a zmienną w czasie warunkową wariancją. Ważne jest by móc badać czy większa warunkowa zmienność powoduje większy zwrot, czyli czy występuje tzw. “premia za ryzyko”.

Oszacowano 9 modeli GARCH-M, po czym porównano dla nich wartości kryteriów informacyjnych, które są przedstawione w poniższej tabeli.

##                 Akaike     Bayes   Shibata Hannan-Quinn
## GARCH-M(1,1) -6.513714 -6.501392 -6.513725    -6.509081
## GARCH-M(1,2) -6.522379 -6.505950 -6.522400    -6.516203
## GARCH-M(1,3) -6.519172 -6.498635 -6.519204    -6.511451
## GARCH-M(2,1) -6.512019 -6.495589 -6.512039    -6.505842
## GARCH-M(2,2) -6.520778 -6.500241 -6.520810    -6.513057
## GARCH-M(2,3) -6.517571 -6.492926 -6.517616    -6.508306
## GARCH-M(3,1) -6.510374 -6.489837 -6.510406    -6.502653
## GARCH-M(3,2) -6.519097 -6.494452 -6.519142    -6.509832
## GARCH-M(3,3) -6.516137 -6.487386 -6.516200    -6.505328

Na podstawie powyższych wartości kryteriów informacyjnych wybrano model GARCH-M(1,2) i go oszacowano. Dla tej wersji modelu funkcja warunkowej wariancji określona jest następującym równaniem:

\[ \sigma^2_t= \alpha_0 + \alpha_1\epsilon^2_{t-1} + \beta_1\sigma^2_{t-1}+ \beta_2\sigma^2_{t-2} \]

Przeprowadzono estymację.

## 
## *---------------------------------*
## *          GARCH Model Fit        *
## *---------------------------------*
## 
## Conditional Variance Dynamics    
## -----------------------------------
## GARCH Model  : sGARCH(1,2)
## Mean Model   : ARFIMA(0,0,0)
## Distribution : norm 
## 
## Optimal Parameters
## ------------------------------------
##         Estimate  Std. Error  t value Pr(>|t|)
## mu      0.000686    0.000966  0.71083  0.47719
## archm   0.014215    0.111037  0.12802  0.89813
## omega   0.000004    0.000006  0.68783  0.49156
## alpha1  0.136182    0.025858  5.26663  0.00000
## beta1   0.082469    0.074997  1.09963  0.27149
## beta2   0.741254    0.060591 12.23377  0.00000
## 
## Robust Standard Errors:
##         Estimate  Std. Error  t value Pr(>|t|)
## mu      0.000686    0.001243 0.552371 0.580694
## archm   0.014215    0.137789 0.103165 0.917832
## omega   0.000004    0.000041 0.097577 0.922268
## alpha1  0.136182    0.122822 1.108781 0.267525
## beta1   0.082469    0.372192 0.221577 0.824644
## beta2   0.741254    0.159858 4.636957 0.000004
## 
## LogLikelihood : 4081.809 
## 
## Information Criteria
## ------------------------------------
##                     
## Akaike       -6.5265
## Bayes        -6.5019
## Shibata      -6.5266
## Hannan-Quinn -6.5173
## 
## Weighted Ljung-Box Test on Standardized Residuals
## ------------------------------------
##                         statistic  p-value
## Lag[1]                      7.526 0.006082
## Lag[2*(p+q)+(p+q)-1][2]     7.791 0.007351
## Lag[4*(p+q)+(p+q)-1][5]     8.372 0.023809
## d.o.f=0
## H0 : No serial correlation
## 
## Weighted Ljung-Box Test on Standardized Squared Residuals
## ------------------------------------
##                          statistic p-value
## Lag[1]                       1.769  0.1835
## Lag[2*(p+q)+(p+q)-1][8]      3.495  0.5950
## Lag[4*(p+q)+(p+q)-1][14]     5.228  0.7408
## d.o.f=3
## 
## Weighted ARCH LM Tests
## ------------------------------------
##             Statistic Shape Scale P-Value
## ARCH Lag[4]    0.7968 0.500 2.000  0.3720
## ARCH Lag[6]    2.0326 1.461 1.711  0.4836
## ARCH Lag[8]    2.2991 2.368 1.583  0.6810
## 
## Nyblom stability test
## ------------------------------------
## Joint Statistic:  1.1393
## Individual Statistics:              
## mu     0.09267
## archm  0.06992
## omega  0.24270
## alpha1 0.05450
## beta1  0.10222
## beta2  0.09351
## 
## 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.7847 0.432751    
## Negative Sign Bias  1.7570 0.079163   *
## Positive Sign Bias  0.9099 0.363063    
## Joint Effect       12.8022 0.005085 ***
## 
## 
## Adjusted Pearson Goodness-of-Fit Test:
## ------------------------------------
##   group statistic p-value(g-1)
## 1    20     24.84     0.166043
## 2    30     42.87     0.046768
## 3    40     65.36     0.005132
## 4    50     75.30     0.009253
## 
## 
## Elapsed time : 0.4645002

Na podstawie uzyskanych wartości statysk oraz odpowiadających im p-value można zauważyć, że wszystkie parametry są nieistotne poza \(\alpha_1\) i \(\beta_2\). Parametr premii za ryzyko wyszedł nieistotny, w związku z czym odrzucamy hipotezę, jakoby większa warunkowa zmienność powodowała większy zwrot, czyli nie występuje tzw. “premia za ryzyko”.

Warto przedstawić także wykres krzywej napływu wiadomości (News Impact curve) oraz wykres ACF wystandaryzowanych reszt.

4. Oszacowania wartości narażonej na ryzyko (VaR)

4.1. Model GARCH(1,2)

W pierwszym kroku, w celu obliczenia wartości narażonej na ryzyko zostały obliczone standaryzowane stopy zwrotu z portfela oraz 1% kwantyl z próbki. Poniższy wykres prezentuje stopy zwrotu z portfela (kolor czerwony) i VaR (kolor zielony). Straty przekroczyły zakładany poziom wartości narażonej na ryzyko w 0,48% przypadków.

4.2. Model EGARCH(1,2)

Procedurę powtórzono dla rozszerzeń modelu GARCH. W modelu EGARCH(1,2) w 0.4% przypadków straty przekroczyły zakładany poziom VaR. Wykres przedstawia stopy zwrotu z portfela (kolor czerwony) i VaR (kolor zielony).

4.3. Model GARCH-M(1,2)

Wykonano także obliczenia dla drugiego rozszerzenia modelu GARCH. Poniżej przedstawiono wykres stopy zwrotu z portfela (kolor czerwony) i VaR (kolor zielony) dla modelu GARCH-M(1,2). W 0.48% przypadków straty przekroczyły zakładany poziom VaR.

Poniższy wykres prezentuje porównanie trzech oszacowanych wartości narażonych na ryzyko. Można zauważyć, że oszacowania VaR modelem GARCH(1,2) i GARCH-in-mean(1,2) są tożsame.

5. Warunkowa wariancja

Kolejnym etapem analizy było obliczenie bezwarunkowych wariancji oraz oszacowanie warunkowej wariancji dla badanych modeli klasy GARCH i dokonanie jej predykcji.

Poniższy wykres prezentuje porównanie warunkowych wariancji dla trzech analizowanych modeli klasy GARCH - GARCH(1,2) (zielony), EGARCH(1,2) (żółty) i GARCH-M(1,2) (niebieski). Warunkowa wariancja oszacowana z modelu GARCH(1,2) jest bardzo zbliżona (niemalże pokrywa się) do warunkowej wariancji modelu GARCH-in-mean(1,2).

5.1. Model GARCH(1,2)

W pierwszej kolejności dokonano obliczenia bezwarunkowej wariancji na podstawie wzoru \[E\epsilon_t^2 = \frac{\alpha_0}{1-\alpha_1-\beta_1-\beta_2}.\] Bezwarunkowa wariancja wyniosła 0,00009998.

Poniższe wykresy przedstawiają prognozy warunkowej wariancji dla 100 i 500 przyszłych okresów oraz wielkosć wariancji bezwarunkowej.

Na podstawie powyzszych wykresów można wywnioskować, że prognoza warunkowej waiancji w długim okresie zbiega do poziomu wariancji bezwarunkowej.

5.2. Model EGARCH(1,2)

Procedurę powtórzono dla wybranych rozszerzeń modeu GARCH. Dla modelu EGARCH(1,2) bezwarunkowa wariancja wyniosła 0.00012.

## unconditional variance 2 
##             0.0001196596

Poniższe wykresy przedstawiają prognozy warunkowej wariancji odpowiednio dla 100 i 500 przyszłych okresów oraz wielkosć wariancji bezwarunkowej.

Na podstawie powyzszych wykresów można wywnioskować, że prognoza warunkowej waiancji w długim okresie zbiega do poziomu wariancji bezwarunkowej.

5.3. Model GARCH-M(1,2)

Dla modelu GARCH-M(1,2) bezwarunkowa wariancja wyniosła 0.0001.

## unconditional variance 3 
##             0.0001008959

Poniższe wykresy przedstawiają prognozy warunkowej wariancji odpowiednio dla 100 i 500 przyszłych okresów oraz wielkosć wariancji bezwarunkowej.

Na podstawie powyzszych wykresów można wywnioskować, że prognoza warunkowej wariancji w długim okresie zbiega do poziomu wariancji bezwarunkowej.

Podsumowanie

Badany szereg (stopy zwrotu z wybranego portfela akcji) charakteryzował się cechami właściwymi dla większości szeregów, tj. brak normalności rozkładu, leptokurtyczność, grupowanie wariancji, efekt grubych ogonów. W związku z tym, w badaniu zastosowano model GARCH i jego rozszerzenia - Exponential GARCH (aby sprawdzić asymetrię informacji) i GARCH-in-Mean (w celu sprawdzenia, czy występuje premia za ryzyko). Na podstawie kryteriów informacyjnych wybrano najlepsze modele w swojej klasie. Okazało się, że model GARCH(1,2) dobrze opisuje zmienność stóp zwrotu, niemniej jednak model EGARCH(1,2) pozwala również na uchwycenie asymetrii wpływu informacji pozytywnych i negatywnych (tzw. “efekt dźwigni”). Jest to ważne w przypadku modelowania zmienności stóp zwrotu, gdyż wraz ze spadkiem ceny instrumentu występuje tendencja do wzrostu wariancji stóp zwrotu. Analiza modelu GARCH-M(1,2) pozwoliła na odrzucenie hipotezy, że występuje premia za ryzyko (większa warunkowa zmienność nie powoduje większego zwrotu), wobec czego uzyskane wyniki (oszacowania warunkowej wariancji i VaR) były bardzo zbliżone w przypadku tego modelu i modelu GARCH(1,2). Najmniejsza liczba przekroczeń wartości zagrożonej była dla modelu EGARCH(1,2). W długim okresie, dla wszystkich trzech modeli, prognoza warunkowej wariancji zbiegała do wariancji bezwarunkowej. Największa wariancja bezwarunkowa była dla modelu EGARCH(1,2).