#Problem 2
#Introduction
#This problems asks us to estimate ARMA model, GARCH model, TGARCH(m,s). Similarly, the problem tries to make us familiarize with stargazer and rugarch package.
# METHODS
#The solution of this problem needs many packages. Let me start by attaching different packages needed for solving the problem.
library(Quandl)
## Warning: package 'Quandl' was built under R version 3.2.4
## Loading required package: xts
## Warning: package 'xts' was built under R version 3.2.4
## Loading required package: zoo
## Warning: package 'zoo' was built under R version 3.2.4
##
## Attaching package: 'zoo'
## The following objects are masked from 'package:base':
##
## as.Date, as.Date.numeric
library(tseries)
## Warning: package 'tseries' was built under R version 3.2.4
library(forecast)
## Warning: package 'forecast' was built under R version 3.2.4
## Loading required package: timeDate
## This is forecast 6.2
library(urca)
## Warning: package 'urca' was built under R version 3.2.4
library(fGarch)
## Warning: package 'fGarch' was built under R version 3.2.4
## Loading required package: timeSeries
## Warning: package 'timeSeries' was built under R version 3.2.4
##
## Attaching package: 'timeSeries'
## The following object is masked from 'package:zoo':
##
## time<-
## Loading required package: fBasics
## Warning: package 'fBasics' was built under R version 3.2.4
##
## Rmetrics Package fBasics
## Analysing Markets and calculating Basic Statistics
## Copyright (C) 2005-2014 Rmetrics Association Zurich
## Educational Software for Financial Engineering and Computational Science
## Rmetrics is free software and comes with ABSOLUTELY NO WARRANTY.
## https://www.rmetrics.org --- Mail to: info@rmetrics.org
library(rugarch)
## Warning: package 'rugarch' was built under R version 3.2.4
## Loading required package: parallel
library(stargazer)
##
## Please cite as:
## Hlavac, Marek (2015). stargazer: Well-Formatted Regression and Summary Statistics Tables.
## R package version 5.2. http://CRAN.R-project.org/package=stargazer
#After this let me import the data from Quandl. I am going to consider the close stock price of Microsoft. I am using following codes for that.
mct <- Quandl("GOOG/NASDAQ_MSFT", api_key="T8gDpM8iFjXR9pfbXBmx", type="zoo")
mc<-mct$Close
logmct<- log(mc)
diffmct<- diff(logmct,1)
par(mfrow=c(1,2))
plot(logmct,xlab="Time", ylab="Log of MCT")
plot(diffmct,xlab="Time", ylab="Log Differences of MCT")
# From the figures above I saw that original data is not stationary and log differenced data is stationary.
par(mfrow=c(1,2))
acf(coredata(diffmct),type='correlation',na.action=na.pass,lag=24, ,main="ACF of diffMCT")
acf(coredata(diffmct),type='partial',na.action=na.pass,lag=24,main="PACF of diffMCT")
#Let me estimate the best ARMA model
arma <- arima(coredata(diffmct), order=c(5,0,1))
arma
##
## Call:
## arima(x = coredata(diffmct), order = c(5, 0, 1))
##
## Coefficients:
## ar1 ar2 ar3 ar4 ar5 ma1 intercept
## -0.0196 -0.0530 -0.0311 -0.0253 -0.0100 -0.0192 8e-04
## s.e. 0.7175 0.0302 0.0405 0.0257 0.0306 0.7172 2e-04
##
## sigma^2 estimated as 0.0005352: log likelihood = 17801.19, aic = -35586.38
tsdiag(arma,gof.lag=24)
arma1 <- arima(coredata(diffmct), order=c(5,0,2))
arma1
##
## Call:
## arima(x = coredata(diffmct), order = c(5, 0, 2))
##
## Coefficients:
## Warning in sqrt(diag(x$var.coef)): NaNs produced
## ar1 ar2 ar3 ar4 ar5 ma1 ma2
## -0.1053 -0.3851 -0.0482 -0.046 -0.0249 0.0664 0.3290
## s.e. NaN 0.2365 NaN NaN NaN NaN 0.2278
## intercept
## 8e-04
## s.e. 2e-04
##
## sigma^2 estimated as 0.0005351: log likelihood = 17802.04, aic = -35586.07
tsdiag(arma1,gof.lag=24)
arma2 <- arima(coredata(diffmct), order=c(4,0,4))
arma2
##
## Call:
## arima(x = coredata(diffmct), order = c(4, 0, 4))
##
## Coefficients:
## Warning in sqrt(diag(x$var.coef)): NaNs produced
## ar1 ar2 ar3 ar4 ma1 ma2 ma3 ma4
## 0.1526 0.0432 -0.0249 -0.188 -0.1899 -0.0911 0.0056 0.176
## s.e. NaN NaN NaN NaN NaN NaN NaN NaN
## intercept
## 8e-04
## s.e. 2e-04
##
## sigma^2 estimated as 0.0005351: log likelihood = 17801.74, aic = -35583.48
tsdiag(arma2,gof.lag=24)
arma3 <- arima(coredata(diffmct), order=c(2,0,2))
arma3
##
## Call:
## arima(x = coredata(diffmct), order = c(2, 0, 2))
##
## Coefficients:
## ar1 ar2 ma1 ma2 intercept
## 0.1336 0.201 -0.1744 -0.2493 8e-04
## s.e. 1.6474 1.231 1.6299 1.2749 2e-04
##
## sigma^2 estimated as 0.0005355: log likelihood = 17799.2, aic = -35586.4
tsdiag(arma3,gof.lag=24)
#As AIC of ARMA(2,0,2) is the lowest. So this will be considered for further analysis. ARMA(2,0,2) seems to be good model among randomly tested models.
#ARCH model
par(mfrow=c(1,2))
acf(arma3$residuals^2,type='correlation',na.action=na.pass,lag=24,main="ACF of Squared Residuals")
acf(arma3$residuals^2,type='partial',na.action=na.pass,lag=24,main="PACF of Squared Residuals")
sarma3 <- arima(arma3$residuals^2, order=c(2,0,2))
tsdiag(sarma3,gof.lag= 36)
arch <-garchFit(~arma(2,2) + garch(7,0), data=diffmct, trace=FALSE)
## Warning in arima(.series$x, order = c(u, 0, v), include.mean =
## include.mean): possible convergence problem: optim gave code = 1
## Warning in sqrt(diag(fit$cvar)): NaNs produced
par(mfrow=c(2,2))
plot(arch@residuals,type="l", xlab="",ylab="",main="Residual")
plot(arch@sigma.t,type="l", xlab="",ylab="",main="Conditional Deviation")
plot(arch@residuals/arch@sigma.t, type="l", xlab="",ylab="",main="Standard Residuals")
plot((arch@residuals/arch@sigma.t)^2, type="l", xlab="",ylab="",main="Squared Standard Residuals")
check1 <- arch@residuals/arch@sigma.t
check2 <- (arch@residuals/arch@sigma.t)^2
par(mfcol=c(2,2), cex=0.6, mar=c(2,2,3,1))
acf(check1, type="correlation", lag.max= 36, ylab="", main="ACF for standardized residuals")
acf(check1, type="partial", lag.max= 36, ylab="", main="PACF for standardized residuals")
acf(check2, type="correlation", lag.max= 36, ylab="", main="ACF for squared standardized residuals")
acf(check2, type="partial", lag.max=36, ylab="", main="PACF for squared standardized residuals")
#For Q(10)
scheck1 <- arima(check1, order=c(10,0,0))
tsdiag(scheck1,gof.lag=36)
#For Q(20)
scheck2 <- arima(check2, order=c(20,0,0))
tsdiag(scheck2,gof.lag=36)
normalTest(diffmct, method="jb")
##
## Title:
## Jarque - Bera Normalality Test
##
## Test Results:
## STATISTIC:
## X-squared: 66354.9472
## P VALUE:
## Asymptotic p Value: < 2.2e-16
##
## Description:
## Wed Mar 23 10:37:55 2016 by user: t420
#GARCH Model with t-student innovations.
garch1<-garchFit(~garch(2,2), data=coredata(diffmct), trace=FALSE, cond.dist="sstd")
## Warning in sqrt(diag(fit$cvar)): NaNs produced
par(mfrow=c(2,2))
plot(garch1@residuals,type="l", xlab="",ylab="",main="Residuals")
plot(garch1@sigma.t,type="l", xlab="",ylab="",main="Conditional Standard Deviation")
plot(garch1@residuals/garch1@sigma.t, type="l", xlab="",ylab="",main="Standardized Residuals")
plot((garch1@residuals/garch1@sigma.t)^2, type="l", xlab="",ylab="",main="Squared Standardized Residuals")
check1 <- garch1@residuals/garch1@sigma.t
check2 <- (garch1@residuals/garch1@sigma.t)^2
par(mfcol=c(2,2), cex=0.6, mar=c(2,2,3,1))
acf(check1, type="correlation", lag.max=36, ylab="", main="ACF for residuals")
acf(check1, type="partial", lag.max=36, ylab="", main="PACF for residuals")
acf(check2, type="correlation", lag.max=36, ylab="", main="ACF fo standardized residuals")
acf(check2, type="partial", lag.max=36, ylab="", main="PACF for squared standardized residuals")
scheck2 <- arima(check2, order=c(3,0,0))
tsdiag(scheck2,gof.lag=36)
# GARCH(m,s) model
tgarch<- garchFit(~garch(2,2), data=coredata(diffmct), trace=FALSE, leverage=TRUE)
## Warning in sqrt(diag(fit$cvar)): NaNs produced
par(mfrow=c(2,2))
plot(tgarch@residuals,type="l", xlab="",ylab="",main="Residuals")
plot(tgarch@sigma.t,type="l", xlab="",ylab="",main="Conditional Standard Deviation")
plot(tgarch@residuals/tgarch@sigma.t, type="l", xlab="",ylab="",main="Standardized Residuals")
plot((tgarch@residuals/tgarch@sigma.t)^2, type="l", xlab="",ylab="",main="Squared Standardized Residuals")
scheck2 <- arima(check2, order=c(10,0,0))
tsdiag(scheck2,gof.lag=36)
# GARCH model with t-student innovations
ah <- garchFit(~ garch(7,0), data = diffmct, trace = F, cond.dist = "sstd")
ah
##
## Title:
## GARCH Modelling
##
## Call:
## garchFit(formula = ~garch(7, 0), data = diffmct, cond.dist = "sstd",
## trace = F)
##
## Mean and Variance Equation:
## data ~ garch(7, 0)
## <environment: 0x00000000230509a0>
## [data = diffmct]
##
## Conditional Distribution:
## sstd
##
## Coefficient(s):
## mu omega alpha1 alpha2 alpha3 alpha4
## 0.00092692 0.00012196 0.17098959 0.17446704 0.09603789 0.15396924
## alpha5 alpha6 alpha7 skew shape
## 0.10359261 0.09792705 0.14606113 1.04751371 4.09515436
##
## Std. Errors:
## based on Hessian
##
## Error Analysis:
## Estimate Std. Error t value Pr(>|t|)
## mu 9.269e-04 1.941e-04 4.776 1.79e-06 ***
## omega 1.220e-04 9.821e-06 12.418 < 2e-16 ***
## alpha1 1.710e-01 2.366e-02 7.226 4.98e-13 ***
## alpha2 1.745e-01 2.537e-02 6.878 6.07e-12 ***
## alpha3 9.604e-02 2.110e-02 4.551 5.34e-06 ***
## alpha4 1.540e-01 2.373e-02 6.488 8.71e-11 ***
## alpha5 1.036e-01 2.115e-02 4.897 9.72e-07 ***
## alpha6 9.793e-02 2.144e-02 4.568 4.93e-06 ***
## alpha7 1.461e-01 2.288e-02 6.385 1.71e-10 ***
## skew 1.048e+00 1.462e-02 71.626 < 2e-16 ***
## shape 4.095e+00 2.031e-01 20.167 < 2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Log Likelihood:
## 19182.47 normalized: 2.529667
##
## Description:
## Wed Mar 23 10:38:23 2016 by user: t420