#### Chapter 5 ####
da=read.table("d-c2c-0110.txt",header=T) # Load data
csco=log(da$CSCO+1)
cat=log(da$CAT+1)
library(fGarch)
m1=garchFit(~1+garch(1,1),data=csco,trace=F)
summary(m1)
m2=garchFit(~1+garch(1,1),data=cat,trace=F)
summary(m2)
vcsco=volatility(m1)
vcat=volatility(m2)
xp=csco+cat
xm=csco-cat
m3=garchFit(~1+garch(1,1),data=xp,trace=F)
summary(m3)
vxp=volatility(m3)
m4=garchFit(~1+garch(1,1),data=xm,trace=F)
summary(m4)
vxm=volatility(m4)
CoV=(vxp^2-vxm^2)/4
COR=CoV/(vcat*vcsco)
source("EWMAvol.R")
M1=EWMAvol(rtn)
tdx=c(1:2515)/252+2001
par(mfcol=c(2,1))
cr2=M1[,3]/sqrt(M1[,1]*M1[,2])
range(cr2,COR)
plot(tdx,COR,xlab='year',ylab='cor',ylim=c(-0.35,1),type='l')
title(main='(a) GARCH based')
plot(tdx,cr2,xlab='year',ylab='cor',ylim=c(-0.35,1),type='l')
title(main='(b) EWMA with theta = 0.94') 
#### time-varying betas
da=read.table("d-sp500-0110.txt")
sp5=da[,1]
da=read.table("d-c2c-0110.txt",header=T)
cat=log(da$CAT+1)
xp=cat+sp5
xm=cat-sp5
m1=garchFit(~1+garch(1,1),data=xp,trace=F)
summary(m1)
m2=garchFit(~1+garch(1,1),data=xm,trace=F)
summary(m2)
m3=garchFit(~1+garch(1,1),data=sp5,trace=F)
summary(m3)
vxp=volatility(m1)
vxm=volatility(m2)
vsp5=volatility(m3)
beta=(vxp^2-vxm^2)/(4*vsp5^2)
tdx=c(1:2515)/252+2001
m4=lm(cat~sp5)
summary(m4)
plot(tdx,beta,xlab='year',ylab='beta',type='l')
abline(h=c(1.146))
idx=c(1:2515)[beta==max(beta)]  #Locate the outlier 
idx
### Global minimum variance portfolio
library(fGarch)
rtn=cbind(ba,cat,ibm,msft,pg)
V1=cov(rtn[1:756,])
V1inv=solve(V1)
One=matrix(1,5,1)
Wgt=V1inv#*#One
D=sum(Wgt*One)
Wgt=Wgt/D
print(Wgt)
print(1/sqrt(D))
print(sqrt(diag(V1)))
rtn=cbind(abt,ibm,wmt)
source("GMVP.R")
M2=GMVP(rtn,start=2011)
names(M2)
wgt=M2$weights
range(wgt)
prtn=M2$returns
mean(prtn)
sqrt(var(prtn))
Mean=apply(rtn[2012:2515,],2,mean)
Mean
v1=sqrt(apply(rtn[2012:2515,],2,var))
print(v1)
minV=sqrt(M2$minVariance)
Vol=sqrt(M2$variances)
range(minV,Vol)
tdx=c(1:505)/2515+2009
plot(tdx,wgt[1,],xlab='year',ylab='weights',type='l',ylim=c(-.75,1.5))
lines(tdx,wgt[2,],lty=2)
lines(tdx,wgt[3,],lty=3) 
plot(tdx,Vol[,1],xlab='year',ylab='vol',type='l',ylim=c(0,0.04))
lines(tdx,Vol[,2],lty=2)
lines(tdx,Vol[,3],lty=3)
lines(tdx,minV,lty=4)
###
library(fGarch)
da=read.table("w-petroprice.txt",header=T)
price=ts(pet,frequency=52,start=c(1997,1))
dp=ts(diff(price),frequency=52,start=c(1997,2))
par(mfcol=c(2,1))
plot(price,xlab='year',ylab='price')
plot(dp,xlab='year',ylab='changes')
cprice=diff(price)
m2=arima(cprice,order=c(3,0,0),seasonal=list(order=c(2,0,0),period=5))
m2=arima(cprice,order=c(3,0,0),seasonal=list(order=c(2,0,0),period=5),include.mean=F)
m2
length(cprice)
m2=arima(cprice,seasonal=list(order=c(2,0,0),period=5),include.mean=F)
m2
adjcp=cprice[11:716]-0.0983*cprice[6:711]-0.1152*cprice[1:706]
acf(adjcp)
pacf(adjcp)
m3=garchFit(~arma(3,0)+garch(1,1),data=adjcp,trace=F,include.mean=F)
summary(m3)
plot(m3)
m4=garchFit(~arma(3,0)+garch(1,1),data=adjcp,trace=F,include.mean=F,cond.dist="std")
summary(m4)
plot(m4)
m5=garchFit(~arma(1,0)+garch(1,1),data=adjcp,trace=F,include.mean=F,cond.dist="sstd")
summary(m5)
plot(m5)
# Backtesting 
M3=arima(adjcp,order=c(3,0,0),include.mean=F)
source("backtest.R")
M3F=backtest(M3,adjcp,650,2,inc.mean=F)
source("backtestGarch.R")
M4F=backtestGarch(adjcp,650,2,inc.mean=F,cdist="sstd")