ZEngineer ADEC 7460: Discussion Wk 5
MODELING GARCH FOR STOCK EPAM
library(forecast)
library(xts)
library(quantmod)
library(tseries)
library(TTR)
library(magrittr)
library(rugarch)
library(fGarch)
library(PerformanceAnalytics)- Downlaod 5 years of stock prices
G <- new.env()
getSymbols("EPAM", env = G, src = "yahoo",
from = as.Date("2015-07-30"), to = as.Date("2020-07-29"),return.class="ts")## [1] "EPAM"e<- G$EPAM
e.ts<-ts(e,frequency = 252, start = 2015) # no. of days its traded in 1 yr.
autoplot(e.ts[,"EPAM.Adjusted"])
- Use STL to decompose the Adjusted Price
e.Adj<- e.ts[,6]
e.Adj%>%
stl(t.window=13, s.window="periodic", robust=TRUE) %>%
autoplot()
There is a lot of seasonal variability in the data, especially as the trend has increased.
- Create training and validation sets
e.train<-e.ts[1:1228,6]
e.test<-e.ts[1229:1258,6]
vol.train<-e.ts[1:1228,5]
vol.test<-e.ts[1229:1258,5]Model 1: ARIMA
model1<-auto.arima(e.train)
summary(model1)## Series: e.train
## ARIMA(1,1,0) with drift
##
## Coefficients:
## ar1 drift
## -0.2667 0.1315
## s.e. 0.0275 0.0728
##
## sigma^2 estimated as 10.45: log likelihood=-3179.51
## AIC=6365.02 AICc=6365.04 BIC=6380.36
##
## Training set error measures:
## ME RMSE MAE MPE MAPE MASE
## Training set 0.0003433461 3.228165 1.82056 -0.07550476 1.467754 1.001429
## ACF1
## Training set 0.01133663##RMSE: 3.2281ARIMA(1,1,0) with drift
f1<- forecast(model1,h=30)
autoplot(f1)
## Model 2: ARIMA with Volume as Regressor
model2<-auto.arima(e.train,xreg = vol.train)
summary(model2)## Series: e.train
## Regression with ARIMA(1,1,0) errors
##
## Coefficients:
## ar1 drift xreg
## -0.2709 0.1316 0
## s.e. 0.0276 0.0000 0
##
## sigma^2 estimated as 10.43: log likelihood=-3178.19
## AIC=6364.38 AICc=6364.42 BIC=6384.83
##
## Training set error measures:
## ME RMSE MAE MPE MAPE MASE
## Training set 0.0002340523 3.224689 1.827992 -0.07588145 1.475724 1.005517
## ACF1
## Training set 0.01132387##RMSE: 3.2247
f2<- forecast(model2,h=30,xreg = e.test)
autoplot(f2)
Model 3 GARCH
#install.packages("betategarch")
library(betategarch)
model3<- tegarch(e.train)
summary(model3)## $date
## [1] "Thu Jul 30 12:44:21 2020"
##
## $par
## omega phi1 kappa1 kappastar df skew
## -0.01460799 0.91304542 -0.04638260 -0.40142716 9.89901114 0.53629958
##
## $objective
## [1] -7037.25
##
## $convergence
## [1] 1
##
## $iterations
## [1] 12
##
## $evaluations
## function gradient
## 47 73
##
## $message
## [1] "false convergence (8)"f3<- predict.tegarch(model3, n.ahead = 30)
plot(f3)
Accuracy measures
Model 1: AIC=6365.02 AICc=6365.04 BIC=6380.36
Model 2: AIC=6364.38 AICc=6364.42 BIC=6384.83
Model 3: Log-likelihood: -7620.20876 BIC: 12.45566
It is not clear how to interpret the GARCH model results.
- Here the GARCH model is applied to returns
Model GARCH for returns (Close Price)
returns<- CalculateReturns(e[,4])
returns<- returns[-1]
chart.Histogram(returns,
methods = c("add.density","add.normal"),
colorset = c("blue","green","red"))
s<- ugarchspec(mean.model = list(armaOrder = c(1,0)),
variance.model = list(model = "sGARCH"),
distribution.model = "norm")
model4<- ugarchfit(data = returns, spec = s)
model4##
## *---------------------------------*
## * GARCH Model Fit *
## *---------------------------------*
##
## Conditional Variance Dynamics
## -----------------------------------
## GARCH Model : sGARCH(1,1)
## Mean Model : ARFIMA(1,0,0)
## Distribution : norm
##
## Optimal Parameters
## ------------------------------------
## Estimate Std. Error t value Pr(>|t|)
## mu 0.001867 0.000499 3.7382 0.000185
## ar1 -0.013593 0.032327 -0.4205 0.674120
## omega 0.000050 0.000014 3.5839 0.000338
## alpha1 0.142397 0.032264 4.4134 0.000010
## beta1 0.736541 0.058747 12.5376 0.000000
##
## Robust Standard Errors:
## Estimate Std. Error t value Pr(>|t|)
## mu 0.001867 0.000542 3.44188 0.000578
## ar1 -0.013593 0.030045 -0.45243 0.650961
## omega 0.000050 0.000032 1.56457 0.117683
## alpha1 0.142397 0.081203 1.75359 0.079500
## beta1 0.736541 0.136779 5.38489 0.000000
##
## LogLikelihood : 3190.009
##
## Information Criteria
## ------------------------------------
##
## Akaike -5.0676
## Bayes -5.0472
## Shibata -5.0677
## Hannan-Quinn -5.0600
##
## Weighted Ljung-Box Test on Standardized Residuals
## ------------------------------------
## statistic p-value
## Lag[1] 0.04861 0.8255
## Lag[2*(p+q)+(p+q)-1][2] 0.10355 0.9999
## Lag[4*(p+q)+(p+q)-1][5] 0.25462 0.9994
## d.o.f=1
## H0 : No serial correlation
##
## Weighted Ljung-Box Test on Standardized Squared Residuals
## ------------------------------------
## statistic p-value
## Lag[1] 0.002864 0.9573
## Lag[2*(p+q)+(p+q)-1][5] 0.099460 0.9981
## Lag[4*(p+q)+(p+q)-1][9] 1.046498 0.9843
## d.o.f=2
##
## Weighted ARCH LM Tests
## ------------------------------------
## Statistic Shape Scale P-Value
## ARCH Lag[3] 0.1119 0.500 2.000 0.7380
## ARCH Lag[5] 0.1902 1.440 1.667 0.9678
## ARCH Lag[7] 1.2954 2.315 1.543 0.8613
##
## Nyblom stability test
## ------------------------------------
## Joint Statistic: 0.9151
## Individual Statistics:
## mu 0.22141
## ar1 0.07477
## omega 0.20617
## alpha1 0.45329
## beta1 0.34751
##
## 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.5154 0.60637
## Negative Sign Bias 2.2888 0.02226 **
## Positive Sign Bias 0.5780 0.56338
## Joint Effect 11.0750 0.01133 **
##
##
## Adjusted Pearson Goodness-of-Fit Test:
## ------------------------------------
## group statistic p-value(g-1)
## 1 20 106.1 4.145e-14
## 2 30 125.2 6.478e-14
## 3 40 150.1 5.932e-15
## 4 50 145.6 1.624e-11
##
##
## Elapsed time : 0.257777THe information we get from the model is: 1. Returns equation is: R.t = mu + e.t R.t = 0.001867 + error.t
Variability (sigma.t^2): omega + alpha(e.t-1)^2 + betasigma.t-1^2 0.000050 + 0.1424(e.t-1)^2 + 0.7365sigma.t-1^2
Ljung test: p-values are not significant. Therefore we cannot reject the null hypothesis, indicating there is no evidence of serial correlation in teh standardized residuals or squared residuals. Therefore they are behaving like a white noise process.
The Goodness of Fit shows very small p-values, which means this model for teh residual using Normal Distribution is not a good choice and there is scope for improvement in the model.
Forecast
plot(model4, which = "all")##
## please wait...calculating quantiles...
f4<- ugarchforecast(fitORspec = model4,
n.ahead = 30)
#Variability Plot
plot(sigma(f4))
The model expects variability to increase over the test period.