setwd("/Users/CHATKEAWPAISAL/Desktop/4 term 2/EE435/6")
cat(rep("\n", 50)) 

library(quantmod) 
library(fBasics)
library(sn)
library(PerformanceAnalytics)
library(car)
library(tseries)
library(forecast)
library(fGarch)
#install.packages("fUnitRoots")
library(fUnitRoots)
#install.packages('MTS')
require(MTS)

da=read.table("m-unrate-MIILIN.txt",header=T)
head(da)
x=cbind(da[,2],da[,3],da[,4])
x=log(x)
zt=diffM(x)*100
colnames(zt)=c("MI","IL","IN")
tdx=da[,1]/12
dim(da)
par(mar=c(1,1,1,1))
MTSplot(zt,tdx[2:492])
VARorder(zt)
m1=VAR(zt,2)
m2=refVAR(m1,thres=1.645)
MTSdiag(m1)
#Answer 1.1 The final fitted model

require(vars)
da=read.table("m-unrate-MIILIN.txt",header=T)
head(da)
x=cbind(da[,2],da[,3],da[,4])
x=log(x)
rt=diff(x)*100
colnames(rt)=c("MI","IL","IN")
varfit=VAR(rt,p=2)
summary(varfit)
impresp=irf(varfit)
plot(impresp)
#Answer 1.2 

require(MTS)
da1=read.table("m-m1cnwti.txt")
zt1=cbind(da1$V1,da1$V2)
colnames(zt1)=c("M1g","oilg")
dim(da1)
MTSplot(zt1)
#Answer 2.1

VARorder(zt1)
m3=VAR(zt1,3)
#Answer 2.2

require(vars)
varfit1=VAR(zt1,p=3)
summary(varfit1)
impresp1=irf(varfit1)
plot(impresp1)
#Answer 2.3

#VARpred(m3,6)
#Answer 2.4
fevd(varfit1,n.ahead=6)
#Answer 2.5