setwd("/Users/Orachai/Documents/BE/4TH YEAR 2ND SEMESTER/EE435")
cat(rep("\n",50))
install.packages("quantmod")  
install.packages("fBasics") 
install.packages("sn")  
install.packages("PerformanceAnalytics") 
install.packages("car") 
install.packages("tseries")  
install.packages("forecast")

library(quantmod) 
library(fBasics)
library(sn)
library(PerformanceAnalytics)
library(car)
library(tseries)
library(forecast)
getSymbols("CAT",from="2000-01-03",to="2021-01-31")
getSymbols("AOT.BK",from="2000-01-03",to="2021-01-31")

da1=CAT
da2=AOT.BK
chartSeries(CAT,theme="white")
chartSeries(AOT.BK,theme="white")

price1=da1[,6] 
logprice1=log(price1) 
plot(logprice1,type='l')
logreturn1=diff(log(price1)) 
simplereturn1 <-exp(logreturn1)-1 

#1 Plot the series of log return and simple return
par("mar")
par(mar=c(1,1,1,1))
par(mfrow=c(2,1))
plot(logreturn1,type = "l", col = "green", lwd = 2,
     main = "Log and Simple Return of CAT")
lines(simplereturn1,col = "red", lwd = 2)

newlogreturn1 <- logreturn1[2:nrow(logreturn1),] 
newsimplereturn1 <- simplereturn1[2:nrow(logreturn1),]
par(mfrow=c(2,1))
plot(newlogreturn1, type = "l", col = "green", lwd = 2,
     main = "Logs and Simple Return of CAT")
lines(newsimplereturn1, col = "red", lwd = 2)

price2=da2[,6] 
logprice2=log(price2) 
plot(logprice2,type='l')
logreturn2=diff(log(price2)) 
simplereturn2 <-exp(logreturn2)-1

newlogreturn2 <- logreturn2[2:nrow(logreturn2),] 
newsimplereturn2 <- simplereturn2[2:nrow(logreturn2),]
par("mar")
par(mar=c(1,1,1,1))
par("mar")
par(mfrow=c(2,1))
plot(newlogreturn2, type = "l", col = "green", lwd = 2,
     main = "Logs and Simple Return of AOT.BK")
lines(newsimplereturn2, col = "red", lwd = 2)

table.Stats(simplereturn1)
table.Stats(simplereturn2)

par(mfrow=c(1,1))
hist(simplereturn1, breaks=100, col="slateblue")
chart.Histogram(simplereturn1,methods = c("add.normal"))
hist(simplereturn2, breaks=100, col="slateblue")
#test for normality
par(mfrow=c(1,1))
qqnorm(newsimplereturn1) 
qqline(newsimplereturn1, col = 2)
jarque.bera.test(newsimplereturn1) 

#4 Test mean = 0
t.test(newlogreturn1)
t.test(newlogreturn2)

#Empirical density plots of daily log returns
par(mfrow=c(1,1))
hist(logreturn1, breaks=100, col="slateblue")
hist(logreturn2, breaks=100, col="grey")



#5 Test Skewness CAT = 0
T=length(newlogreturn1)
s3c=skewness(newlogreturn1)
tst1 = s3c/sqrt(6/T)  
tst1
s3c
pv1 = 2*pnorm(tst1)
pv1

#6 Test excess kurtosis CAT =0
k4c = kurtosis(newlogreturn1)
tst1 = k4c/sqrt(24/T)  
tst1
pv1 = 2*(1-pnorm(tst1)) 
pv1

# Test skewness AOT.BK
T=length(newlogreturn2)
s3a=skewness(newlogreturn2)
tst2 = s3a/sqrt(6/T)  
tst2
s3a
pv2 = 2*(1-pnorm(tst2))
pv2

# Test excess Kurtosis AOT.BK = 0 
k4a = kurtosis(newlogreturn2)
tst2 = k4a/sqrt(24/T)  
tst2
pv2 = 2*(1-pnorm(tst2)) 
pv2