library(quantmod)
## Загрузка требуемого пакета: xts
## Загрузка требуемого пакета: zoo
##
## Присоединяю пакет: 'zoo'
## Следующие объекты скрыты от 'package:base':
##
## as.Date, as.Date.numeric
## Загрузка требуемого пакета: TTR
## Registered S3 method overwritten by 'quantmod':
## method from
## as.zoo.data.frame zoo
library(forecast)
library(tseries)
library(rugarch)
## Загрузка требуемого пакета: parallel
##
## Присоединяю пакет: 'rugarch'
## Следующий объект скрыт от 'package:stats':
##
## sigma
library(ggplot2)
library(stargazer)
##
## Please cite as:
## Hlavac, Marek (2022). stargazer: Well-Formatted Regression and Summary Statistics Tables.
## R package version 5.2.3. https://CRAN.R-project.org/package=stargazer
library(lmtest)
library(tseries)
library(urca)
library(FinTS)
##
## Присоединяю пакет: 'FinTS'
## Следующий объект скрыт от 'package:forecast':
##
## Acf
SIRI <- getSymbols("SIRI", from = "2010-01-01", to = "2022-12-31", auto.assign = FALSE)
Price <- SIRI$SIRI.Close
chartSeries(SIRI)

adf.test(Price)
##
## Augmented Dickey-Fuller Test
##
## data: Price
## Dickey-Fuller = -3.0082, Lag order = 14, p-value = 0.1514
## alternative hypothesis: stationary
#Стационарны ли первые разности?
adf.test(na.omit(diff(Price))) #ряд стационарен
## Warning in adf.test(na.omit(diff(Price))): p-value smaller than printed p-value
##
## Augmented Dickey-Fuller Test
##
## data: na.omit(diff(Price))
## Dickey-Fuller = -15.827, Lag order = 14, p-value = 0.01
## alternative hypothesis: stationary
summary(ur.df(na.omit(diff(Price))))
##
## ###############################################
## # Augmented Dickey-Fuller Test Unit Root Test #
## ###############################################
##
## Test regression none
##
##
## Call:
## lm(formula = z.diff ~ z.lag.1 - 1 + z.diff.lag)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.78135 -0.03273 0.00258 0.03870 0.49704
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## z.lag.1 -1.04272 0.02563 -40.678 <2e-16 ***
## z.diff.lag -0.02906 0.01749 -1.662 0.0967 .
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.0795 on 3267 degrees of freedom
## Multiple R-squared: 0.5374, Adjusted R-squared: 0.5371
## F-statistic: 1897 on 2 and 3267 DF, p-value: < 2.2e-16
##
##
## Value of test-statistic is: -40.6782
##
## Critical values for test statistics:
## 1pct 5pct 10pct
## tau1 -2.58 -1.95 -1.62
dim(SIRI$SIRI.Close %>% diff())
## [1] 3272 1
SIRI$SIRI.Close %>% diff() %>% ggtsdisplay(main="", theme = theme_classic())
## Warning: Removed 1 rows containing missing values (`geom_point()`).

SIRI$SIRI.Close %>% diff() %>% ggtsdisplay(main="", theme = theme_classic(), plot.type="scatter")
## Warning: Removed 1 rows containing missing values (`geom_point()`).
## Warning: Removed 1 rows containing missing values (`geom_point()`).

#значения ACF и PACF
acf(na.omit(SIRI$SIRI.Close %>% diff()), plot = FALSE)
##
## Autocorrelations of series 'na.omit(SIRI$SIRI.Close %>% diff())', by lag
##
## 0 1 2 3 4 5 6 7 8 9 10
## 1.000 -0.074 0.034 -0.014 -0.005 0.029 -0.058 0.031 -0.075 0.004 0.009
## 11 12 13 14 15 16 17 18 19 20 21
## -0.012 0.016 -0.021 0.007 -0.025 -0.018 -0.002 -0.011 0.001 -0.010 0.018
## 22 23 24 25 26 27 28 29 30 31 32
## -0.034 0.017 0.011 -0.002 -0.012 -0.002 -0.015 -0.022 -0.001 -0.048 -0.004
## 33 34 35
## 0.006 -0.006 0.000
pacf(na.omit(SIRI$SIRI.Close %>% diff()), plot = FALSE)
##
## Partial autocorrelations of series 'na.omit(SIRI$SIRI.Close %>% diff())', by lag
##
## 1 2 3 4 5 6 7 8 9 10 11
## -0.074 0.029 -0.009 -0.007 0.029 -0.054 0.021 -0.068 -0.009 0.012 -0.009
## 12 13 14 15 16 17 18 19 20 21 22
## 0.009 -0.013 -0.004 -0.021 -0.026 -0.007 -0.007 -0.005 -0.006 0.013 -0.034
## 23 24 25 26 27 28 29 30 31 32 33
## 0.009 0.011 -0.001 -0.016 0.000 -0.021 -0.022 -0.008 -0.049 -0.010 0.006
## 34 35
## -0.008 -0.005
#построение моделей
mod1 <- auto.arima(Price) #модель - ARIMA(0,1,2)
mod2 <- arima(Price, order = c(0,1,1))
mod3 <- arima(Price, order = c(0,1,2)) #для сохранения через stargazer
mod4 <- arima(Price, order = c(0,1,3))
mod1
## Series: Price
## ARIMA(0,1,2)
##
## Coefficients:
## ma1 ma2
## -0.0710 0.0324
## s.e. 0.0175 0.0174
##
## sigma^2 = 0.006318: log likelihood = 3642.49
## AIC=-7278.98 AICc=-7278.97 BIC=-7260.7
mod2
##
## Call:
## arima(x = Price, order = c(0, 1, 1))
##
## Coefficients:
## ma1
## -0.0697
## s.e. 0.0169
##
## sigma^2 estimated as 0.00632: log likelihood = 3640.75, aic = -7277.51
mod3
##
## Call:
## arima(x = Price, order = c(0, 1, 2))
##
## Coefficients:
## ma1 ma2
## -0.0710 0.0324
## s.e. 0.0175 0.0174
##
## sigma^2 estimated as 0.006314: log likelihood = 3642.49, aic = -7278.98
mod4
##
## Call:
## arima(x = Price, order = c(0, 1, 3))
##
## Coefficients:
## ma1 ma2 ma3
## -0.0714 0.0341 -0.0164
## s.e. 0.0175 0.0175 0.0184
##
## sigma^2 estimated as 0.006312: log likelihood = 3642.89, aic = -7277.77
AIC(mod1)
## [1] -7278.978
AIC(mod2)
## [1] -7277.509
#первые 2000 наблюдений меньше волатильности, чем на последних 1000
Price1 <- na.omit(SIRI$SIRI.Close %>% diff())
var.test(SIRI$SIRI.Close[c(1:2000),], SIRI$SIRI.Close[c(2000:3000),]) #да это так
##
## F test to compare two variances
##
## data: SIRI$SIRI.Close[c(1:2000), ] and SIRI$SIRI.Close[c(2000:3000), ]
## F = 5.9722, num df = 1999, denom df = 1000, p-value < 2.2e-16
## alternative hypothesis: true ratio of variances is not equal to 1
## 95 percent confidence interval:
## 5.358728 6.642889
## sample estimates:
## SIRI.Close
## SIRI.Close 5.972171
## attr(,"names")
## [1] "ratio of variances"
var.test(Price1[c(1:2000),], Price1[c(2000:3000),])
##
## F test to compare two variances
##
## data: Price1[c(1:2000), ] and Price1[c(2000:3000), ]
## F = 0.1937, num df = 1999, denom df = 1000, p-value < 2.2e-16
## alternative hypothesis: true ratio of variances is not equal to 1
## 95 percent confidence interval:
## 0.1738000 0.2154493
## sample estimates:
## SIRI.Close
## SIRI.Close 0.1936959
## attr(,"names")
## [1] "ratio of variances"
#проверка на автокорреляцию
checkresiduals(mod3)

##
## Ljung-Box test
##
## data: Residuals from ARIMA(0,1,2)
## Q* = 30.572, df = 8, p-value = 0.0001674
##
## Model df: 2. Total lags used: 10
#видна автокорреляция, кластеры волатильности
#есть ли ARCH структара?
ArchTest(resid(mod3))
##
## ARCH LM-test; Null hypothesis: no ARCH effects
##
## data: resid(mod3)
## Chi-squared = 564.34, df = 12, p-value < 2.2e-16
#да
#нормальны ли остатки?
jarque.bera.test(resid(mod3))
##
## Jarque Bera Test
##
## data: resid(mod3)
## X-squared = 19991, df = 2, p-value < 2.2e-16
#нет
#попробуем добавить гетероскедтичность
e <- resid(mod3)
acf(e^2)

pacf(e^2)

#похоже на GARCH структуру
e^2 %>% ggtsdisplay(main="", theme = theme_classic())

#GARCH модель
garch_spec1 <- ugarchspec(variance.model=list(model="sGARCH", garchOrder=c(1,1)), mean.model=list(armaOrder=c(0,2)))
fit_garch1 <- ugarchfit(spec = garch_spec1, data = na.omit(diff(Price)))
fit_garch1
##
## *---------------------------------*
## * GARCH Model Fit *
## *---------------------------------*
##
## Conditional Variance Dynamics
## -----------------------------------
## GARCH Model : sGARCH(1,1)
## Mean Model : ARFIMA(0,0,2)
## Distribution : norm
##
## Optimal Parameters
## ------------------------------------
## Estimate Std. Error t value Pr(>|t|)
## mu 0.002258 0.000948 2.3815 0.017244
## ma1 -0.033328 0.018349 -1.8164 0.069311
## ma2 -0.018796 0.018573 -1.0120 0.311515
## omega 0.000015 0.000005 2.9999 0.002700
## alpha1 0.037831 0.002614 14.4715 0.000000
## beta1 0.961169 0.001844 521.2528 0.000000
##
## Robust Standard Errors:
## Estimate Std. Error t value Pr(>|t|)
## mu 0.002258 0.000918 2.46093 0.013858
## ma1 -0.033328 0.019796 -1.68356 0.092267
## ma2 -0.018796 0.020683 -0.90878 0.363464
## omega 0.000015 0.000011 1.35050 0.176857
## alpha1 0.037831 0.005945 6.36334 0.000000
## beta1 0.961169 0.002621 366.74926 0.000000
##
## LogLikelihood : 4213.048
##
## Information Criteria
## ------------------------------------
##
## Akaike -2.5723
## Bayes -2.5612
## Shibata -2.5723
## Hannan-Quinn -2.5683
##
## Weighted Ljung-Box Test on Standardized Residuals
## ------------------------------------
## statistic p-value
## Lag[1] 0.0005998 0.9805
## Lag[2*(p+q)+(p+q)-1][5] 0.2048564 1.0000
## Lag[4*(p+q)+(p+q)-1][9] 1.9190880 0.9877
## d.o.f=2
## H0 : No serial correlation
##
## Weighted Ljung-Box Test on Standardized Squared Residuals
## ------------------------------------
## statistic p-value
## Lag[1] 3.169 0.07506
## Lag[2*(p+q)+(p+q)-1][5] 6.047 0.08770
## Lag[4*(p+q)+(p+q)-1][9] 6.826 0.21419
## d.o.f=2
##
## Weighted ARCH LM Tests
## ------------------------------------
## Statistic Shape Scale P-Value
## ARCH Lag[3] 0.1833 0.500 2.000 0.6686
## ARCH Lag[5] 0.3019 1.440 1.667 0.9398
## ARCH Lag[7] 0.7671 2.315 1.543 0.9482
##
## Nyblom stability test
## ------------------------------------
## Joint Statistic: 0.9908
## Individual Statistics:
## mu 0.11258
## ma1 0.09702
## ma2 0.03604
## omega 0.21126
## alpha1 0.27390
## beta1 0.25369
##
## 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 2.1107 0.0348732 **
## Negative Sign Bias 0.5489 0.5830977
## Positive Sign Bias 3.3564 0.0007986 ***
## Joint Effect 12.2711 0.0065099 ***
##
##
## Adjusted Pearson Goodness-of-Fit Test:
## ------------------------------------
## group statistic p-value(g-1)
## 1 20 198.0 8.409e-32
## 2 30 321.0 5.564e-51
## 3 40 307.8 1.725e-43
## 4 50 315.9 1.090e-40
##
##
## Elapsed time : 0.328692
garch_spec2 <- ugarchspec(variance.model=list(model="sGARCH", garchOrder=c(1,2)), mean.model=list(armaOrder=c(0,2)))
fit_garch2 <- ugarchfit(spec = garch_spec2, data = na.omit(diff(Price)))
fit_garch2
##
## *---------------------------------*
## * GARCH Model Fit *
## *---------------------------------*
##
## Conditional Variance Dynamics
## -----------------------------------
## GARCH Model : sGARCH(1,2)
## Mean Model : ARFIMA(0,0,2)
## Distribution : norm
##
## Optimal Parameters
## ------------------------------------
## Estimate Std. Error t value Pr(>|t|)
## mu 0.002230 0.000944 2.3614 0.018205
## ma1 -0.033175 0.018692 -1.7748 0.075937
## ma2 -0.019223 0.018372 -1.0463 0.295436
## omega 0.000019 0.000008 2.4538 0.014137
## alpha1 0.053922 0.007111 7.5831 0.000000
## beta1 0.533068 0.015769 33.8042 0.000000
## beta2 0.412010 0.015189 27.1263 0.000000
##
## Robust Standard Errors:
## Estimate Std. Error t value Pr(>|t|)
## mu 0.002230 0.000917 2.43159 0.015033
## ma1 -0.033175 0.019659 -1.68749 0.091510
## ma2 -0.019223 0.019915 -0.96522 0.334435
## omega 0.000019 0.000020 0.94823 0.343011
## alpha1 0.053922 0.026360 2.04557 0.040798
## beta1 0.533068 0.014061 37.91068 0.000000
## beta2 0.412010 0.013269 31.05067 0.000000
##
## LogLikelihood : 4216.169
##
## Information Criteria
## ------------------------------------
##
## Akaike -2.5736
## Bayes -2.5606
## Shibata -2.5736
## Hannan-Quinn -2.5690
##
## Weighted Ljung-Box Test on Standardized Residuals
## ------------------------------------
## statistic p-value
## Lag[1] 0.0008563 0.9767
## Lag[2*(p+q)+(p+q)-1][5] 0.2239654 1.0000
## Lag[4*(p+q)+(p+q)-1][9] 1.9154667 0.9878
## d.o.f=2
## H0 : No serial correlation
##
## Weighted Ljung-Box Test on Standardized Squared Residuals
## ------------------------------------
## statistic p-value
## Lag[1] 1.700 0.1923
## Lag[2*(p+q)+(p+q)-1][8] 5.740 0.2622
## Lag[4*(p+q)+(p+q)-1][14] 6.875 0.5221
## d.o.f=3
##
## Weighted ARCH LM Tests
## ------------------------------------
## Statistic Shape Scale P-Value
## ARCH Lag[4] 0.1487 0.500 2.000 0.6997
## ARCH Lag[6] 0.2904 1.461 1.711 0.9476
## ARCH Lag[8] 0.7675 2.368 1.583 0.9553
##
## Nyblom stability test
## ------------------------------------
## Joint Statistic: 1.0935
## Individual Statistics:
## mu 0.11142
## ma1 0.08874
## ma2 0.03336
## omega 0.43602
## alpha1 0.25858
## beta1 0.24545
## beta2 0.24304
##
## 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 2.0879 0.036887 **
## Negative Sign Bias 0.2337 0.815249
## Positive Sign Bias 3.0325 0.002445 ***
## Joint Effect 9.9740 0.018788 **
##
##
## Adjusted Pearson Goodness-of-Fit Test:
## ------------------------------------
## group statistic p-value(g-1)
## 1 20 203.1 8.319e-33
## 2 30 320.2 7.987e-51
## 3 40 304.1 8.667e-43
## 4 50 322.6 6.172e-42
##
##
## Elapsed time : 0.333338
garch_spec3 <- ugarchspec(variance.model=list(model="sGARCH", garchOrder=c(2,2)), mean.model=list(armaOrder=c(0,2)))
fit_garch3 <- ugarchfit(spec = garch_spec3, data = na.omit(diff(Price)))
fit_garch3
##
## *---------------------------------*
## * GARCH Model Fit *
## *---------------------------------*
##
## Conditional Variance Dynamics
## -----------------------------------
## GARCH Model : sGARCH(2,2)
## Mean Model : ARFIMA(0,0,2)
## Distribution : norm
##
## Optimal Parameters
## ------------------------------------
## Estimate Std. Error t value Pr(>|t|)
## mu 0.002230 0.000945 2.361068 0.018222
## ma1 -0.033158 0.018692 -1.773860 0.076086
## ma2 -0.019217 0.018373 -1.045949 0.295585
## omega 0.000019 0.000008 2.380103 0.017308
## alpha1 0.053915 0.013367 4.033564 0.000055
## alpha2 0.000000 0.013625 0.000001 0.999999
## beta1 0.533382 0.015707 33.958187 0.000000
## beta2 0.411703 0.014776 27.863772 0.000000
##
## Robust Standard Errors:
## Estimate Std. Error t value Pr(>|t|)
## mu 0.002230 0.000913 2.443937 0.014528
## ma1 -0.033158 0.019670 -1.685737 0.091846
## ma2 -0.019217 0.019920 -0.964707 0.334692
## omega 0.000019 0.000022 0.891565 0.372626
## alpha1 0.053915 0.023910 2.254918 0.024138
## alpha2 0.000000 0.021766 0.000001 0.999999
## beta1 0.533382 0.013326 40.026957 0.000000
## beta2 0.411703 0.016010 25.715706 0.000000
##
## LogLikelihood : 4216.169
##
## Information Criteria
## ------------------------------------
##
## Akaike -2.5730
## Bayes -2.5581
## Shibata -2.5730
## Hannan-Quinn -2.5677
##
## Weighted Ljung-Box Test on Standardized Residuals
## ------------------------------------
## statistic p-value
## Lag[1] 0.0008092 0.9773
## Lag[2*(p+q)+(p+q)-1][5] 0.2236244 1.0000
## Lag[4*(p+q)+(p+q)-1][9] 1.9151226 0.9878
## d.o.f=2
## H0 : No serial correlation
##
## Weighted Ljung-Box Test on Standardized Squared Residuals
## ------------------------------------
## statistic p-value
## Lag[1] 1.700 0.1923
## Lag[2*(p+q)+(p+q)-1][11] 6.187 0.4110
## Lag[4*(p+q)+(p+q)-1][19] 7.838 0.6938
## d.o.f=4
##
## Weighted ARCH LM Tests
## ------------------------------------
## Statistic Shape Scale P-Value
## ARCH Lag[5] 0.05333 0.500 2.000 0.8174
## ARCH Lag[7] 0.63293 1.473 1.746 0.8607
## ARCH Lag[9] 0.85835 2.402 1.619 0.9495
##
## Nyblom stability test
## ------------------------------------
## Joint Statistic: 3.8964
## Individual Statistics:
## mu 0.11135
## ma1 0.08858
## ma2 0.03335
## omega 0.43500
## alpha1 0.25827
## alpha2 0.30143
## beta1 0.24514
## beta2 0.24273
##
## Asymptotic Critical Values (10% 5% 1%)
## Joint Statistic: 1.89 2.11 2.59
## Individual Statistic: 0.35 0.47 0.75
##
## Sign Bias Test
## ------------------------------------
## t-value prob sig
## Sign Bias 2.0878 0.036890 **
## Negative Sign Bias 0.2336 0.815330
## Positive Sign Bias 3.0323 0.002446 ***
## Joint Effect 9.9729 0.018798 **
##
##
## Adjusted Pearson Goodness-of-Fit Test:
## ------------------------------------
## group statistic p-value(g-1)
## 1 20 202.9 8.848e-33
## 2 30 320.0 8.614e-51
## 3 40 303.8 1.018e-42
## 4 50 322.0 8.016e-42
##
##
## Elapsed time : 0.4039321
infocriteria(fit_garch1)
##
## Akaike -2.572331
## Bayes -2.561155
## Shibata -2.572338
## Hannan-Quinn -2.568329
infocriteria(fit_garch2)
##
## Akaike -2.573628
## Bayes -2.560589
## Shibata -2.573637
## Hannan-Quinn -2.568959
infocriteria(fit_garch3)
##
## Akaike -2.573017
## Bayes -2.558115
## Shibata -2.573029
## Hannan-Quinn -2.567680
infocriteria(fit_garch1) > infocriteria(fit_garch2)
##
## Akaike TRUE
## Bayes FALSE
## Shibata TRUE
## Hannan-Quinn TRUE
infocriteria(fit_garch1) > infocriteria(fit_garch3)
##
## Akaike TRUE
## Bayes FALSE
## Shibata TRUE
## Hannan-Quinn FALSE
infocriteria(fit_garch2) < infocriteria(fit_garch3)
##
## Akaike TRUE
## Bayes TRUE
## Shibata TRUE
## Hannan-Quinn TRUE
#AIC меньший в модели fit_garch2
#визуальная диагностика
#прогнозы: 2SD
plot(fit_garch2, which = 1)

#VAR
plot(fit_garch2, which = 2)
##
## please wait...calculating quantiles...

#условная дисперсия
plot(fit_garch2, which = 3)

#распределение остатков
plot(fit_garch2, which = 8)
