#************************************************************ 
# GARCH Type Models
# install.packages("tseries")
# install.packages("fGarch")
# install.packages("AER")
# install.packages("fBasics")
#------------------------------------------------------------

library(quantmod)
library(fBasics)
library(TSPred)
library(tseries)
library(fGarch)
library(forecast)
# Yahoo Data in R

getSymbols("AAPL",from="2000-01-02", to="2017-06-30")
getSymbols("IBM",from="2000-01-02", to="2017-06-30")
getSymbols("^GSPC",from="2000-01-02", to="2017-06-30")
getSymbols("^VIX",from="2000-01-02", to="2017-06-30")

chartSeries(GSPC, theme="white") 
chartSeries(AAPL, theme="white") 
chartSeries(IBM, theme="white",subset='2007-01::2010-01') 
chartSeries(VIX, theme="white") 

AAPL.rtn = 100*diff(log(AAPL$AAPL.Adjusted))
IBM.rtn = 100*diff(log(IBM$IBM.Adjusted))
GSPC.rtn = 100*diff(log(GSPC$GSPC.Adjusted))

chartSeries(AAPL.rtn, theme="white",subset='2007::2010')
chartSeries(IBM.rtn, theme="white",subset='2007::2010')
chartSeries(GSPC.rtn, theme="white",subset='2007::2010')

par(mfcol=c(2,1))
plot(AAPL$AAPL.Adjusted,col="blue",main="Apple Stock Price")
plot(AAPL.rtn,col="blue",main="Apple Stock Daily Return")


par(mfcol=c(2,1))
plot(GSPC$GSPC.Adjusted,col="blue",main="SP500 Stock Index")
plot(GSPC.rtn,col="blue",main="SP500 Daily Return")


par(mfcol=c(2,1))
     plot(GSPC.rtn,col="blue",main="SP500 Daily Return")
     plot(VIX$VIX.Adjusted,col="blue",main="VIX")
head(GSPC.rtn)
nob <- nrow(IBM.rtn)
IBM_return <- as.ts(IBM.rtn[2:nob],frequency(365)  )    
head(IBM_return)

par(mfcol=c(2,1))
acf(IBM_return )
pacf(IBM_return )

dy_square <- IBM_return^2
acf(dy_square )
pacf(dy_square )


# Here is a fund part to estimate various GARCH models 
y <- IBM_return
auto.arima(y, ic = "aic", trace = TRUE)


g1 <- garchFit(~arma(3,2) + garch(5,0), data = y, trace = FALSE)
summary(g1)
g1_ic <- as.matrix(g1@fit$ics)
colnames(g1_ic) <- c('ARCH')
# plot(g0)

g2 <- garchFit(~arma(3,2) + garch(1,1), data = y, trace = FALSE)
summary(g2)
g2_ic <- as.matrix(g2@fit$ics)
colnames(g2_ic) <- c('GARCH')

# plot(g1)

# Taylor-Schwert ARCH (compare Ding, Granger, Engle, eq. (16))
g3 <- garchFit(~ arma(3,2) + aparch(1,1), delta = 1,
                   data = y, trace = FALSE)
summary(g3)
g3_ic <- as.matrix(g3@fit$ics)
colnames(g3_ic) <- c('APARCH1')

# Threshold GARCH (TGARCH)
g4 <- garchFit(~ arma(3,2) + garch(1,1), delta = 1,
                  leverage = TRUE, data = y, trace = FALSE)
summary(g4)
g4_ic <- as.matrix(g4@fit$ics)
colnames(g4_ic) <- c('TGARCH')

# GJR-GARCH
g5 <- garchFit(~ arma(3,2) + garch(1,1), delta = 2,
                    leverage = TRUE, data = y, trace = FALSE)
summary(g5)
g5_ic <- as.matrix(g5@fit$ics)
colnames(g5_ic) <- c('GRJ_GARCH')

# APARCH
g6 <- garchFit(~arma(3,2) + aparch(1,1), data = y, trace = FALSE)
summary(g6)
g6_ic <- as.matrix(g6@fit$ics)
colnames(g6_ic) <- c('APARCH')
#GARCHM
g7 <- garchFit(~arma(3,2) + garch(1,1), data = y, trace = FALSE,leverage=TRUE)
summary(g7)
g2_ic <- as.matrix(g2@fit$ics)
colnames(g2_ic) <- c('GARCH')


IC_table <- cbind(g1_ic,g1_ic,g3_ic,g4_ic,g5_ic,g6_ic)
t(IC_table)