Analiza porównawcza ryzyka w modelach klasy GARCH.

Wstęp

Celem badania jest analiza porównawcza ryzyka, rozumianego jako oszacowanie funkcji warunkowej wariancji w modelach klasy GARCH. W pracy zostały porównane 3 modele: model GARCH, eGARCH oraz t-GARCH, z różnymi parametrami. Z każdego rodzaju został wybrany jeden model a następnie przeprowadzono analizę wariancji warunkowiej w okresie in sample oraz out of sample. Na koniec policzono Value at Risk dla każdego modelu i na tej podstawie zdecydowano, który model daje najmniejsze ryzyko straty.

Badanie składa się z następujących części:

Wstępna analiza portfela oraz zwrotów z porfela
Korelacje w danych - ACF
Modele GARCH, EARCH, t-GARCH
Wariancja warunkowa w okresie in-sample
Value at Risk w okresie in-sample
Wariancja warunkowa w okresie out-of-sample
Value at Risk w okresie out-of-sample
Analiza Wrażliwości
Podsumowanie

1. Wstępna analiza portfela

Aby przeprowadzić badanie zbudowano portfel równoważny (equally-weighted portfolio), składający się z 5 szeregów czasowych. Na podstawie danych z cen zamknięcia 5 aktywów, stworzono portfel, dla którego policzono dzienne stopy zwrotu. W skład portfela wchodzą następujące akcje i indeksy:

Inter Cars
Comarch
Amica
ceny złota
kurs dolara (USD/PLN)

Dane zostały pobrane ze strony https://stooq.com/, obejmują okres od 26.05.2004 r. do 18.08.2017 r. Po połączeniu portfela. liczba obserwacji wynosi 3236.

Próba in-sample obejmuje 2985 obserwacji z okresu 26.05.004 do 18.07.2016
Próba out-of-sample obejmuje 250 obserwacji z okresu 19.08.2016 do 18.08.2017

```

2.Analiza danych:

Wykorzystane dane wyglądają następująco (wartości podane są w PLN):

head(total)

##            Date    CMRCH      CAR  GOLD    AMICA    USD
## 8924 2004-05-26 56.99995 22.50163 388.3 45.50005 3.8282
## 8925 2004-05-27 58.49992 23.00364 396.1 45.50005 3.8040
## 8928 2004-06-01 57.99993 23.00364 395.2 43.89917 3.7920
## 8929 2004-06-02 57.70015 23.00364 392.5 44.29876 3.7940
## 8930 2004-06-03 57.59951 22.80366 388.9 44.39992 3.7945
## 8931 2004-06-04 57.70015 23.00364 391.7 43.89917 3.7924

tail(total)

##             Date CMRCH    CAR   GOLD  AMICA     USD
## 12344 2017-08-10 200.0 289.00 1290.1 166.30 3.63535
## 12345 2017-08-11 198.8 280.05 1294.0 165.05 3.62520
## 12346 2017-08-14 200.0 278.00 1290.4 168.65 3.63913
## 12348 2017-08-16 200.0 277.00 1282.9 166.00 3.62595
## 12349 2017-08-17 213.0 280.00 1292.4 164.55 3.64065
## 12350 2017-08-18 209.0 278.10 1291.6 163.60 3.63076

Stworzony został również indeks zwrotów oraz średnia stopa zwrotu r ze wszystkich indeksów.

Wykres logarymicznych stóp zwrotów:

plot(total.diff$Date, total.diff$r,
     type = "l", col="red", lwd = 1,
     main = "Logarytmiczne stopy zwrotu portfela")

Na wykresie widać efekt grupowania wariancji. Występują okresy nasilonej zmienności i względnej stabilności.Jeśli spojrzymy na statystyki opisowe, dane zawierają wartości z przedziału [-0.463884;0.108505]. Ujemny współczynnik skośności wskazuje na skośność lewostronną rozkładu, co dodatkowo potwierdza, że mediana jest większa od średniej. Widać to również na wykresie, chociaż różnica jest nie wielka. Z analizy wykresu oraz kurtoza o wartości 7.35 wskazują na leptokurtyczność rozkładu.

Podstawowe statystyki opisowe danego zbioru

basicStats(total.diff$r)

##             X..total.diff.r
## nobs            3235.000000
## NAs                0.000000
## Minimum           -0.077423
## Maximum            0.051075
## 1. Quartile       -0.004326
## 3. Quartile        0.005444
## Mean               0.000386
## Median             0.000447
## Sum                1.248459
## SE Mean            0.000164
## LCL Mean           0.000065
## UCL Mean           0.000707
## Variance           0.000087
## Stdev              0.009315
## Skewness          -0.547366
## Kurtosis           5.248092

Hisogram stóp zwrotu

Również na histogramie stóp zwrotu, z nałożoną linią gęstości rozkładu normalnego, rozkład szeregu wygląda na leptokurtyczny:

hist(total.diff$r, prob = T, breaks = 90, main = "Histogram stóp zwrotu")
curve(dnorm(x,
            mean = mean(total.diff$r, na.rm = T),
            sd = sd(total.diff$r, na.rm = T)),
      col = "darkblue", lwd = 2, add = TRUE)

Wykres Q-Q (Quantile-Quantile)

Wykres Q-Q obrazuje, że prawdopodobieństwo wystąpienia obserwacji oddalonych od średniej (nietypowych zmian) są znacznie wyższe niż w rozkłądzie normalnym.

qqnorm(total.diff$r)
qqline(total.diff$r, col = 2)

Test Jarque Bera
Ostatecznie wynik testu Jarque Bera sprawia, że odrzucamy hipotezę zerową o normalności rozkładu badanej próby.

jarque.bera.test(total.diff$r)

## 
##  Jarque Bera Test
## 
## data:  total.diff$r
## X-squared = 4000, df = 2, p-value <0.0000000000000002

Test Ljunga - Boxa

Kolejnym etapem jest sprawdzenie, czy występują korelacje

Test Ljung - Boxa dla zwrotów:

Box.test(total.diff$r, lag = 5, type = 'Ljung-Box')

## 
##  Box-Ljung test
## 
## data:  total.diff$r
## X-squared = 50, df = 5, p-value = 0.000000001

Test Ljunga - Boxa dla kwadratów zwrotów:

Box.test(total.diff$r^2, lag = 5, type = 'Ljung-Box')

## 
##  Box-Ljung test
## 
## data:  total.diff$r^2
## X-squared = 200, df = 5, p-value <0.0000000000000002

Wyniki testu wskazują na odrzucenie H₀ o braku autokorelacji.

Test Durbina - Watsona

Przeprowadzamy również test Durbina - Watsona dla 5 opóźnień zwrotów. Wartość p-value jest poniżej przyjętego poziomu istotności 5%, a więc odrzucamy hipotezę o braku autokorelacji 1 i 3 opóźnienia. W przypadku kwadratów zwrotów, dla wszystkich 5 opóźnień odrzucamy hipotezę zerową.

durbinWatsonTest(lm(formula = total.diff$r ~ 1),
                 max.lag = 5) # sprawdźmy autokorelację dla 5 pierwszych opóźnień

##  lag Autocorrelation D-W Statistic p-value
##    1          0.1083           1.8   0.000
##    2          0.0301           1.9   0.060
##    3          0.0511           1.9   0.002
##    4          0.0118           2.0   0.502
##    5          0.0013           2.0   0.906
##  Alternative hypothesis: rho[lag] != 0

durbinWatsonTest(lm(formula = total.diff$r ^ 2 ~ 1),
                 max.lag = 5) # sprawdźmy autokorelację dla 5 pierwszych opóźnień kwadratów

##  lag Autocorrelation D-W Statistic p-value
##    1           0.131           1.7   0.000
##    2           0.112           1.8   0.000
##    3           0.096           1.8   0.004
##    4           0.102           1.8   0.000
##    5           0.116           1.8   0.004
##  Alternative hypothesis: rho[lag] != 0

Test na występowanie efektów ARCH

Ostatecznym potewierdzeniem występowania efektów ARCH w badanym szeregu jest test na występowanie efektów ARCH, w którym H₀ o braku efektów jest silnie odrzucana, co mówi, że wśród zwrotów występują efekty ARCH.

## 
##  ARCH LM-test; Null hypothesis: no ARCH effects
## 
## data:  total.diff$r
## Chi-squared = 100, df = 5, p-value <0.0000000000000002

Wykresy ACF zwrotów oraz kwadratów zwrotów

Wykres ACF zwrotów oraz kwadratów zwrotów obrazuje występowanie autokorelacji zwrotów oraz silnej autokorelacji kwadratów zwrotów logarytmicznych.

acf(total.diff$r,
    lag.max = 36,
    na.action = na.pass,
    ylim = c(-0.05, 0.15),
    col = "darkblue",
    lwd = 7,
    main = "Wykres ACF zwrotów")

acf(total.diff$r^2,
    lag.max = 36,
    na.action = na.pass,
    ylim = c(-0.05, 0.15),
    col = "darkblue",
    lwd = 7,
    main = "Wykres ACF kwadratów zwrotów")

Wyniki powyższej analizy - występowanie autkorelacji dla logarytmicznych stóp zwrotu, wskazuje, że można estymować wartości oczekiwanych stóp zwrotu na podstawie ich przeszłych wartości. W tym celu wykorzystamy modele GARCH.

3.Modele GARCH, EARCH, t-GARCH

Modele GARCH

W celu wybrania najlepiej dopasowanego modelu GARCH wykonano porównanie modeli z różnymi wartościami parametrów p i q:

Model GARCH(1,1)
Model GARCH(1,2)
Model GARCH(1,3)
Model GARCH(2,1)

Model Garch(1,1)

garch11 <- garchFit(formula = ~ garch(1, 1),
                    data = na.omit(total.diff$r),
                    include.mean = F,
                    cond.dist = "norm",
                    trace = FALSE)
 
summary(garch11)

## 
## Title:
##  GARCH Modelling 
## 
## Call:
##  garchFit(formula = ~garch(1, 1), data = na.omit(total.diff$r), 
##     cond.dist = "norm", include.mean = F, trace = FALSE) 
## 
## Mean and Variance Equation:
##  data ~ garch(1, 1)
## <environment: 0x0000000029d13a88>
##  [data = na.omit(total.diff$r)]
## 
## Conditional Distribution:
##  norm 
## 
## Coefficient(s):
##         omega         alpha1          beta1  
## 0.00000086661  0.03869889438  0.95115715106  
## 
## Std. Errors:
##  based on Hessian 
## 
## Error Analysis:
##            Estimate   Std. Error  t value             Pr(>|t|)    
## omega  0.0000008666 0.0000002949    2.939              0.00329 ** 
## alpha1 0.0386988944 0.0072435457    5.343         0.0000000917 ***
## beta1  0.9511571511 0.0097745463   97.310 < 0.0000000000000002 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Log Likelihood:
##  10781    normalized:  3.3 
## 
## Description:
##  Sun Sep 10 23:05:00 2017 by user: Michał Bylicki 
## 
## 
## Standardised Residuals Tests:
##                                 Statistic p-Value 
##  Jarque-Bera Test   R    Chi^2  1755      0       
##  Shapiro-Wilk Test  R    W      0.97      0       
##  Ljung-Box Test     R    Q(10)  39        0.000022
##  Ljung-Box Test     R    Q(15)  47        0.000038
##  Ljung-Box Test     R    Q(20)  57        0.000024
##  Ljung-Box Test     R^2  Q(10)  8         0.63    
##  Ljung-Box Test     R^2  Q(15)  10        0.79    
##  Ljung-Box Test     R^2  Q(20)  14        0.81    
##  LM Arch Test       R    TR^2   9         0.7     
## 
## Information Criterion Statistics:
##  AIC  BIC  SIC HQIC 
## -6.7 -6.7 -6.7 -6.7

sum(garch11@fit$coef)

## [1] 0.99

W modelu GARCH(1,1) spełnione zostały wymagane warunki, które świadczą o poprawności wybranego modelu:

Parametry modelu są istotne oraz dodatnie, a ich suma jest mniejsza od 1.
P-value dla testu Ljung-Boxa dla kwadratów wystandaryzowanych reszt wskazuje na brak podstaw do odrzucenia hipotezy, że kwadraty wystandaryzowanych reszt są białym szumem (0.63>0.05).
P-value testu LM ARCH wskazuje na brak występowania efektów ARCH wśród reszt modelu.

Analogicznie wykonano analizę dla pozostałych 3 modeli:

W przypadku modelu GARCH(1,2) nie został spełniony warunek o istotności parametru $\beta_2$.
Dla modelu GARCH(1,3) nie został spełniony warunek o istotności parametrów $\beta_1$ i $\beta_3$.
W modelu GARCH(2,1) nie jest spełniony warunek o istotności parametru $\alpha_2$.

Porównanie modeli - kryteria informacyjne

AIC	BIC	SIC	HQIC	model
-6.66308	-6.65744	-6.66308	-6.66106	garch11
-6.66370	-6.65618	-6.66370	-6.66100	garch12
-6.66471	-6.65531	-6.66471	-6.66134	garch13
-6.66253	-6.65501	-6.66253	-6.65984	garch21

Na podstawie kryteriów informacyjnych wybieramy model GARCH(1,1). Kryteria dla modelu GARCH(1,3) są w większości niższe, co mogłoby wskazywać na lepsze dopasowanie modelu, jednak w GARCH(1,3) mamy do czynienia z nieistotnością dwóch parametrów, dlatego wybrany został GARCH(1,1).

Modele e-GARCH

Analiza modeli e-GARCH i proces wyboru najlepszego dopasowania parametrów p i q przebiegał podobnie, jak
w przypadku modeli GARCH. Wybrano kilka modeli, które zostały na końcu porównane pod względem właściwości statystycznych dopasowania oraz kryteriów informacyjnych.

Model e-GARCH(1,1)
Model e-GARCH(1,2)
Model e-GARCH(2,1)
Model e-GARCH(2,2)

Model e-GARCH(1,1)

spec = ugarchspec(variance.model = list(model ="eGARCH",
                  garchOrder = c(1, 1)),
                  mean.model = list(armaOrder = c(0, 0), include.mean = F),
                  distribution.model = "norm")
 
egarch11 <- ugarchfit( spec = spec, data = na.omit(total.diff$r))
egarch11

## 
## *---------------------------------*
## *          GARCH Model Fit        *
## *---------------------------------*
## 
## Conditional Variance Dynamics    
## -----------------------------------
## GARCH Model  : eGARCH(1,1)
## Mean Model   : ARFIMA(0,0,0)
## Distribution : norm 
## 
## Optimal Parameters
## ------------------------------------
##         Estimate  Std. Error    t value Pr(>|t|)
## omega  -0.076178    0.000847   -89.9495 0.000000
## alpha1 -0.023021    0.005065    -4.5454 0.000005
## beta1   0.991364    0.000092 10821.3498 0.000000
## gamma1  0.086569    0.003755    23.0569 0.000000
## 
## Robust Standard Errors:
##         Estimate  Std. Error  t value Pr(>|t|)
## omega  -0.076178    0.002684  -28.378 0.000000
## alpha1 -0.023021    0.008926   -2.579 0.009907
## beta1   0.991364    0.000160 6204.491 0.000000
## gamma1  0.086569    0.003998   21.654 0.000000
## 
## LogLikelihood : 10783.5 
## 
## Information Criteria
## ------------------------------------
##                     
## Akaike       -6.6643
## Bayes        -6.6568
## Shibata      -6.6643
## Hannan-Quinn -6.6616
## 
## Weighted Ljung-Box Test on Standardized Residuals
## ------------------------------------
##                         statistic      p-value
## Lag[1]                      21.68 0.0000032150
## Lag[2*(p+q)+(p+q)-1][2]     22.94 0.0000007529
## Lag[4*(p+q)+(p+q)-1][5]     26.61 0.0000002834
## d.o.f=0
## H0 : No serial correlation
## 
## Weighted Ljung-Box Test on Standardized Squared Residuals
## ------------------------------------
##                         statistic p-value
## Lag[1]                      6.257 0.01237
## Lag[2*(p+q)+(p+q)-1][5]     7.991 0.02963
## Lag[4*(p+q)+(p+q)-1][9]     9.164 0.07521
## d.o.f=2
## 
## Weighted ARCH LM Tests
## ------------------------------------
##             Statistic Shape Scale P-Value
## ARCH Lag[3]    0.4892 0.500 2.000  0.4843
## ARCH Lag[5]    2.4090 1.440 1.667  0.3877
## ARCH Lag[7]    2.6913 2.315 1.543  0.5736
## 
## Nyblom stability test
## ------------------------------------
## Joint Statistic:  0.7689
## Individual Statistics:              
## omega  0.26209
## alpha1 0.07969
## beta1  0.24855
## gamma1 0.09189
## 
## 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.4108 0.68125    
## Negative Sign Bias  1.9955 0.04607  **
## Positive Sign Bias  0.9473 0.34358    
## Joint Effect        6.4233 0.09274   *
## 
## 
## Adjusted Pearson Goodness-of-Fit Test:
## ------------------------------------
##   group statistic       p-value(g-1)
## 1    20     99.44 0.0000000000006775
## 2    30    118.43 0.0000000000009082
## 3    40    141.98 0.0000000000001266
## 4    50    145.94 0.0000000000143728
## 
## 
## Elapsed time : 0.935375

W modelu e-GARCH(1,1) nie został spełniony warunek o braku autokoreclacji kwadratów standaryzowanych reszt.

Pozostałe wymagania:

istotność parametrów oraz ich suma mniejsza od 1,
brak występowania efetów ARCH zostały spełnione.

Następnie sprawdzono pozostałe modele:

Model e-GARCH(1,2): spełnione założenie o braku autokorelacji kwadratów reszt oraz dla parametrów, ale dla lag=4 występują efekty ARCH (p-value 0.047)
Model e-GARCH(2,1): brak spełnionych warunków o autokorelacji kwadratów reszt oraz występują efekty ARCH.
Model e-GARCH(2,2): nieistoty parametr $\alpha_2$

Tabala wartości kryteriów informacyjnych

Akaike	Bayes	Shibata	Hannan-Quinn
-6.6073	-6.5992	-6.6073	-6.6044	egarch11
-6.6092	-6.5991	-6.6092	-6.6055	egarch12
-6.6115	-6.5993	-6.6115	-6.6071	egarch21
-6.6108	-6.5967	-6.6108	-6.6057	egarch22

Do dalszej analizy wybrany został model e-GARCH(1,2). Na model ten wskazują kryteria informacyjne. Model ten nie spełnia wszystkich wymaganych warunków (widoczne efekty ARCH dla lag=4), jednak p-value dla testu ARCH jest zbliżone do przyjętego poziomu istotności $= 5% $. Te dwa fakty pozwoliły nam zdecydować się na wybór modelu e-GARCH(1,2) do dalszego badania i uznać go za najlepszy spośród porównywanych modeli e-GARCH.

Modele GARCH-t

Podobnie jak dla dwóch poprzednich klas modeli wyboru modelu GARCh-t dokonano za pomocą porównania 4 modeli z różnymi kombinacjami parametrów p i q:

Model GARCH-t(1,1)
Model GARCH-t(1,2)
Model GARCH-t(2,1)
Model GARCH-t(2,2)

Aby model był poprawny muszą zostać spełnione warunki:

spec <- ugarchspec(
 
        # równanie warunkowej wariancji
        variance.model = list(model = "sGARCH",
                              garchOrder = c(1, 1)),
        # równanie warunkowej wartości oczekiwanej
        mean.model = list(armaOrder = c(0, 0),
                          include.mean = F),
        # zakładany rozkład warunkowy reszt
        distribution.model = "std") # rozkład t-Studenta
 
 
garcht11 <- ugarchfit(spec=spec, data=na.omit(total.diff$r))
garcht11

## 
## *---------------------------------*
## *          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   2.3236 0.020148
## alpha1  0.054029    0.007999   6.7547 0.000000
## beta1   0.928304    0.009499  97.7250 0.000000
## shape   5.340667    0.457022  11.6858 0.000000
## 
## Robust Standard Errors:
##         Estimate  Std. Error  t value Pr(>|t|)
## omega   0.000002    0.000002  0.65293  0.51380
## alpha1  0.054029    0.039474  1.36871  0.17109
## beta1   0.928304    0.043136 21.52064  0.00000
## shape   5.340667    0.714628  7.47336  0.00000
## 
## LogLikelihood : 10914.8 
## 
## Information Criteria
## ------------------------------------
##                     
## Akaike       -6.7455
## Bayes        -6.7380
## Shibata      -6.7455
## Hannan-Quinn -6.7428
## 
## Weighted Ljung-Box Test on Standardized Residuals
## ------------------------------------
##                         statistic      p-value
## Lag[1]                      21.07 0.0000044240
## Lag[2*(p+q)+(p+q)-1][2]     22.80 0.0000008217
## Lag[4*(p+q)+(p+q)-1][5]     27.13 0.0000002024
## d.o.f=0
## H0 : No serial correlation
## 
## Weighted Ljung-Box Test on Standardized Squared Residuals
## ------------------------------------
##                         statistic p-value
## Lag[1]                      1.063  0.3025
## Lag[2*(p+q)+(p+q)-1][5]     1.412  0.7616
## Lag[4*(p+q)+(p+q)-1][9]     2.240  0.8746
## d.o.f=2
## 
## Weighted ARCH LM Tests
## ------------------------------------
##             Statistic Shape Scale P-Value
## ARCH Lag[3]    0.1459 0.500 2.000  0.7025
## ARCH Lag[5]    0.6693 1.440 1.667  0.8328
## ARCH Lag[7]    0.9411 2.315 1.543  0.9230
## 
## Nyblom stability test
## ------------------------------------
## Joint Statistic:  273.994
## Individual Statistics:              
## omega  27.1202
## alpha1  0.7753
## beta1   0.6812
## shape   0.2451
## 
## 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.34917 0.7270    
## Negative Sign Bias 1.62939 0.1033    
## Positive Sign Bias 0.05883 0.9531    
## Joint Effect       5.40519 0.1444    
## 
## 
## Adjusted Pearson Goodness-of-Fit Test:
## ------------------------------------
##   group statistic p-value(g-1)
## 1    20     23.90      0.20003
## 2    30     36.54      0.15838
## 3    40     50.92      0.09569
## 4    50     50.73      0.40502
## 
## 
## Elapsed time : 0.379973

Model GARCH-t(1,1): spełnione wszystkie warunki.
Model GARCH-t(1,2): nie jest istotny jeden z parametrów - $\omega$.
Model GARCH-t(2,1): nie jest istotny jeden z parametrów - $\alpha_1$.
Model GARCH-t(2,2): parametr $\alpha_2$ nei jest istotny.

Tabala wartości kryteriów informacyjnych

Akaike	Bayes	Shibata	Hannan-Quinn
-6.6927	-6.6846	-6.6927	-6.6898	garcht11
-6.6928	-6.6827	-6.6928	-6.6891	garcht12
-6.6921	-6.6820	-6.6921	-6.6884	garcht21

Kryteria informacyjne wskazują, że należy wybrać model GARCH-t(1,1) lub GARCH-t(1,2). Dla modelu GARCHt-(1,1) spełnione zostały wszystkie warunki (istotne paramerty, ich suma mniejsza od 1, brak efektów ARCH oraz brak autokorelacji wśród kwadratów reszt), dlatego model GARCH-t(1,1) został wykorzystany w dalszej analizie.

4. Wariancja warunkowa w okresie in-sample

Wrunkowa wariancja w okresie in-sample w omdelu GARCH(1,1)

plot( total.diff$Date[0:2985], garch11@sigma.t[0:2985], type='l', col='blue',xlab='Data', ylab='warunkowa wariancja')

W modelu eGarch(1,2), wariancja warunkowa w okresie in-sampe wygląda następująco:

plot( total.diff$Date[0:2985], egarch12@fit$sigma[0:2985], type='l', col='blue', xlab='Data', ylab='warunkowa wariancja')

W modelu tGarch(1,2), wariancja warunkowa w okresie in-sampe wygląda następująco:

plot( total.diff$Date[0:2985], garcht11@fit$sigma[0:2985], type='l', col='blue', xlab='Data', ylab='warunkowa wariancja')

Widzimy, że wszystkie trzy wykresy zachowują się w podobny sposób. Obserwujemy wysoką wartość waruancji warunkowej w latach 2007-2010, co jest wynikiem kryzysu. Po roku 2012 następuje jednak ograniczenie wahań i spadek warunkowej wariancji. W kolejnym kroku będziemy sprawdzać, jak wygląda wartość narażona na ryzyko.

5. Szacowanie wartości narażonej na ryzyko

Aby oszacować wartość narażoną na ryzyko, należy policzyć 1 % kwantyl empiryczny standaryzowanego rozkładu stopy zwrotu oraz warunkowej wariancji w modelu.

total.in <- total.diff[total.diff$Date <= "2016-08-18",]
total.in$rstd <- (total.in$r - mean(total.in$r, na.rm=T)) / sd(total.in$r ,na.rm = T)
q01 <- quantile(total.in$rstd, 0.01, na.rm = T)
q01

##       1% 
## -2.89709

GARCH(1,1)

total.in$VaR <- q01 * garch11@sigma.t[1:2985]

plot(total.in$Date, total.in$r, col = "red", lwd = 1, type = 'l', main='VaR GARCH(1,1)',
     ylim = c(-0.1, 0.1))
abline(h = 0, lty = 2)
lines(total.in$Date, -total.in$VaR, type = 'l', col = "green")
lines(total.in$Date, total.in$VaR, type = 'l', col = "green")

Trend krzywej VaR zmienia się poprawnie, dostosowując się do zmian wartości zwrotów, zwiększonej niepewności.
Na wykresie widoczne są pewne wykroczenia zwrotów poza przedział wartości narażonej na ryzyko. W tym przypadku dla modelu GARCH(1,1) zakładany poziom wartości narażonej na ryzyko został przekroczony w 0.77% przypadków.

sum(total.in$r< total.in$VaR) / length(total.in$VaR)

## [1] 0.00770519

e-GARCH(1,2)

plot(total.in$Date, total.in$r, col = "red", lwd = 1, type = 'l',
     ylim = c(-0.1, 0.1), main='VaR e-GARCH(1,2)')
abline(h = 0, lty = 2)
lines(total.in$Date, -total.in$VaR1, type = 'l', col = "green")
lines(total.in$Date, total.in$VaR1, type = 'l', col = "green")

Dla modelu e-GARCH(1,2) również obserwujemy przekroczenia warotści narażonej na ryzyko.

sum(total.in$r< total.in$VaR1) / length(total.in$VaR1)

## [1] 0.00670017

Tutaj zakładany poziom wartości narażonej na ryzyko został przekroczony w 0.67% przypadków.

GARCH-t(1,1)

total.in$VaR2 <- q01 * garcht11@fit$sigma[1:2985]

plot(total.in$Date, total.in$r, col = "red", lwd = 1, type = 'l', main='VaR GARCH-t(1,1)',
     ylim = c(-0.1, 0.1))
abline(h = 0, lty = 2)
lines(total.in$Date, -total.in$VaR2, type = 'l', col = "green")
lines(total.in$Date, total.in$VaR2, type = 'l', col = "green")

sum(total.in$r< total.in$VaR2) / length(total.in$VaR1)

## [1] 0.00737018

Dla modelu GARCH-t(1,1) zakładany poziom wartości narażonej na ryzyko został przekroczony w 0.74% przypadków.

Z powyższej analizy wnioskujemy, że dla wybranego portfela i założonej wartości VaR, najlepszym modelem pod kątem ograniczenia ryzyka jest e-GARCH(1,2), gdzie procent przekroczenia przedziału wartości narażonej na ryzyko jest najmniejszy i wynosi 0.67%.

6. Wariancja warunkowa w okresie out-of-sample

W modelu GARCH(1,1), wariancja warunkowa w okresie out-of-sampe wygląda następująco:

var_uncondt <- garch11@fit$matcoef[1] / (1 - garch11@fit$matcoef[2]
                                          - garch11@fit$matcoef[3])
names(var_uncondt) <- "unconditional variance"
var_uncondt

## unconditional variance 
##           0.0000854314

# annualizacja
sqrt(var_uncondt*252)

## unconditional variance 
##               0.146727

tmp1.fore100t <- predict(garch11, n.ahead = 100)
head(tmp1.fore100t)

##   meanForecast  meanError standardDeviation
## 1            0 0.00682156        0.00682156
## 2            0 0.00685042        0.00685042
## 3            0 0.00687887        0.00687887
## 4            0 0.00690691        0.00690691
## 5            0 0.00693456        0.00693456
## 6            0 0.00696182        0.00696182

# wykres oszacowań i prognoz warunkowej wariancji w długim okresie
plot(tmp1.fore100t[,3]^ 2, type = "l",   ylim = c(0.000048, 0.000085))
# dodajmy jeszcze poziom wariancji bezwarunkowej
abline(h = var_uncondt, col = "red", lty = 2, ylim=c(0,0.00008)) # ?? nie pojawia się linia
title(main = "Warunkowa i bezwarunkowa wariancja zwrotów portfela GARCH(1,1)")

W modelu e-GARCH(1,2), wariancja warunkowa w okresie out-of-sampe wygląda następująco:

var_uncondt <- uncvariance(egarch12)

# annualizacja
sqrt(var_uncondt*252)

## [1] 0.196558

tmp1.fore100t <- ugarchforecast(egarch12, n.ahead = 100)
head(tmp1.fore100t@forecast$sigmaFor)

##     1978-11-10 01:00:00
## T+1          0.00722022
## T+2          0.00714393
## T+3          0.00724777
## T+4          0.00722833
## T+5          0.00729305
## T+6          0.00730005

# wykres oszacowań i prognoz warunkowej wariancji w długim okresie
plot(tmp1.fore100t@forecast$sigmaFor ^ 2, type = "l",   ylim = c(0.000048, 0.000085))
# dodajmy jeszcze poziom wariancji bezwarunkowej
abline(h = var_uncondt, col = "red", lty = 2, ylim=c(0,0.00008)) # ?? nie pojawia się linia
title(main = "Warunkowa i bezwarunkowa wariancja zwrotów portfela eGARCH(1,2)")

W modelu GARCH-t(1,1), wariancja warunkowa w okresie out-of-sampe wygląda następująco:

var_uncondt <- garcht11@fit$matcoef[1] / (1 - garcht11@fit$matcoef[2]
                                          - garcht11@fit$matcoef[3])
names(var_uncondt) <- "unconditional variance"
var_uncondt

## unconditional variance 
##           0.0000853966

# annualizacja
sqrt(var_uncondt*252)

## unconditional variance 
##               0.146697

tmp1.fore100t <- ugarchforecast(garcht11, n.ahead = 100)
head(tmp1.fore100t@forecast$sigmaFor)

##     1978-11-10 01:00:00
## T+1          0.00720640
## T+2          0.00724730
## T+3          0.00728726
## T+4          0.00732630
## T+5          0.00736445
## T+6          0.00740173

# wykres oszacowań i prognoz warunkowej wariancji w długim okresie
plot(tmp1.fore100t@forecast$sigmaFor ^ 2, type = "l",   ylim = c(0.000048, 0.000085))
# dodajmy jeszcze poziom wariancji bezwarunkowej
abline(h = var_uncondt, col = "red", lty = 2, ylim=c(0,0.00008)) # ?? nie pojawia się linia
title(main = "Warunkowa i bezwarunkowa wariancja zwrotów portfela GARCH-t(1,1)")

7. Szacowanie wartości narażonej na ryzyko w okresie out-of-sample

GARCH(1,1)

DO obliczenia VaR dla modelu GARCH(1,1) w okresie out-of-sample wykonujemy poniższą pętlę:

for (k in start:finish) {
  tmp.data <- total.diff[total.diff$obs <= (k - 1), ]
  tmp.data$rstd <- (tmp.data$r - mean(tmp.data$r, na.rm = T)) /
    sd(tmp.data$r, na.rm = T)
  q01 <- quantile(tmp.data$rstd, 0.01, na.rm = T)
  
  t.garch11 <- garchFit(formula = ~ garch(1, 1),
                        data = na.omit(tmp.data$r),
                        include.mean = F,
                        cond.dist = "norm",
                        trace = FALSE)
  
  sigma.forecast  <- predict(t.garch11, n.ahead = 1)
  sigma.forecast2 <- sigma.forecast$standardDeviation
  VaR[k - start + 1] <- q01 * sigma.forecast2
}
total.for$VaR1 <- VaR

Wrunkowa wariancja w okresie out-of-sample w modelu GARCH(1,1)

plot(total.for$Date, total.for$r, col = "red", lwd = 1, type = 'l',
     ylim = c(-0.15, 0.15), main="Zwroty portfela i VaR dla GARCH(1,1) w out-of-sample")
abline(h = 0, lty = 2)
lines(total.for$Date, total.for$VaR1, type = 'l', col = "green")

Na wykresie widoczne są pewne wykroczenia zwrotów poza przedział wartości narażonej na ryzyko. W tym przypadku dla modelu GARCH(1,1) zakładany poziom wartości narażonej na ryzyko został przekroczony w 1.2 % przypadków.

sum(total.for$r < total.for$VaR1) / length(total.for$VaR1)

## [1] 0.004

e-GARCH(1,2)

plot(total.for$Date, total.for$r, col = "red", lwd = 1, type = 'l',
     ylim = c(-0.15, 0.15), main="Zwroty portfela i VaR dla e-GARCH(1,2) w out-of-sample")
abline(h = 0, lty = 2)
lines(total.for$Date, total.for$VaR2, type = 'l', col = "green")

Dla modelu e-GARCH(1,2) nie zaobserwowaliśmy przekroczenia wartości narażonej na ryzyko.

sum(total.for$r < total.for$VaR2) / length(total.for$VaR2)

## [1] 0

Tutaj zakładany poziom wartości narażonej na ryzyko został przekroczony w 1.2% przypadków.

GARCH-t(1,1)

VaR <- rep(NA, times = finish - start + 1)
sigma.forecast <- ugarchforecast(garcht11, n.ahead = 1)
sigma.forecast@forecast$sigmaFor

##     1978-11-10 01:00:00
## T+1           0.0072064

sigma.forecast2 <- sigma.forecast@forecast$sigmaFor[1, 1]
q01 * sigma.forecast2

##         1% 
## -0.0210296

time1 <- Sys.time()
for (k in start:finish) {
  tmp.data <- total.diff[total.diff$obs <= (k - 1), ]  
  tmp.data$rstd <- (tmp.data$r - mean(tmp.data$r, na.rm = T)) /
    sd(tmp.data$r, na.rm = T)
  q01 <- quantile(tmp.data$rstd, 0.01, na.rm = T)
  spec <- ugarchspec(variance.model = list(model = "sGARCH",
                                           garchOrder = c(1, 1)),
                     mean.model = list(armaOrder = c(0, 0),
                                       include.mean = T),
                     distribution.model = "std")
  tmp.garcht11 <- ugarchfit(spec = spec, data = na.omit(tmp.data$r))
  sigma.forecast  <- ugarchforecast(tmp.garcht11, n.ahead = 1)
  sigma.forecast2 <- sigma.forecast@forecast$sigmaFor[1, 1]
  VaR[k - start + 1] <- q01 * sigma.forecast2
}

total.for$VaR3 <- VaR

plot(total.for$Date, total.for$r, col = "red", lwd = 1, type = 'l',
     ylim = c(-0.15, 0.15), main="Zwroty portfela i VaR dla GARCH-t(1,1) w out-of-sample")
abline(h = 0, lty = 2)
lines(total.for$Date, total.for$VaR3, type = 'l', col = "green")

sum(total.for$r < total.for$VaR3) / length(total.for$VaR3)

## [1] 0.004

Z powyższej analizy wnioskujemy, że dla wybranego portfela i założonej wartości VaR, najlepszym modelem pod kątem ograniczenia ryzyka jest GARCH-t(1,1). Zwroty z portfela przekraczają Value at Risk w 1,2% przypadków dla modelu GARCH(1,1) oraz 0,4% dla modelu GARCH-t(1,1) procent, jedynie dla modelu e-GARCH(1,2) nie zaobserwowaliśmy przekroczenia przedziału wartości narażonej na ryzyko. Na tej podstawie możemy stwierdzić, iż model e-GARCH(1,2) będzie prawdopodobnie najlepszy do wykonywania prognoz.

8. Analiza Wrażliwości przy zmianie długości okna do estymacji

Dla porównania wykonano analizę VaR dla okresu out-of-sample dla 6 miesięcy w modelu eGARCH(1,2)

e-GARCH(1,2)

plot(total.for1$Date, total.for1$r, col = "red", lwd = 1, type = 'l',
     ylim = c(-0.15, 0.15), main="Zwroty portfela i VaR dla e-GARCH(1,2) w out-of-sample")
abline(h = 0, lty = 2)
lines(total.for1$Date, total.for1$VaR1, type = 'l', col = "green")

Dla modelu e-GARCH(1,2) nie zaobserwowaliśmy przekroczenia wartości narażonej na ryzyko.

sum(total.for1$r < total.for1$VaR1) / length(total.for1$VaR1)

## [1] 0

Podsumowanie

W powyższej pracy dokonano analizy ryzyka pod porfela dla zrównoważonego portfela składającego się z 5 aktywów. W ramach tej analizy porównano trzy modele - GARCH(1,1), e-GARCH(1,2) oraz GARCH-t(1,1). Otrzymane wyniki wskazują na duże podobieństwo między modelami, Wartść narażona na ryzyko w wybranym porfelu waha się między 0,67% (w modelu e-GARCH a 0,77% (GARCH).Niewielkie różnice występują również w warunkowej wariancji. Sprawdzono również, czy znaczenie na uzyskane wyniki ma długość okresu in sample. Wykonana została analiza wrażliwości przy zmianie okna do estymacji dla modelu eGARCH(1,2), który według analizy in-sample okazał się być prawdopodobnie najlpeszym modelem do wykonania prognozy. Jednak, ponieważ w badanym przykładzie, nie zaobserwowano przekroczeń VaR, zgodnie z oczekiwaniami, skrócenie okresu prognozy nie pogorszyło wyników estymacji.

Analiza porównawcza ryzyka w modelach klasy GARCH.

Projekt zaliczeniowy ZASC

Karolina Tylak, Michał Bylicki

Wrzesień 2017

Wstęp

Badanie składa się z następujących części:

1. Wstępna analiza portfela

2.Analiza danych:

3.Modele GARCH, EARCH, t-GARCH

Modele GARCH

Modele e-GARCH

Modele GARCH-t

4. Wariancja warunkowa w okresie in-sample

5. Szacowanie wartości narażonej na ryzyko

6. Wariancja warunkowa w okresie out-of-sample

7. Szacowanie wartości narażonej na ryzyko w okresie out-of-sample

8. Analiza Wrażliwości przy zmianie długości okna do estymacji

Podsumowanie