# 1) CAPM Functions #Import dataset# —> students should change the import path according to their desktop(or save CAPM text file to downloads) library(readr) CAPM <- read_delim("Downloads/CAPM.txt", "\t", escape_double = FALSE, trim_ws = TRUE) View(CAPM) #Optional# da <- CAPM head(da) #Create t-value intercept Function# tintercept <- function (y,x){ lm_result <- summary(lm(y~x)) print(lm_result$coefficients[1,3]) } #Example of Function Output# tintercept(da$RtnInf, da$RtnSp500) #Create t-value slope Function# tslope <- function (y,x){ lm_result <- summary(lm(y~x)) print(lm_result$coefficients[2,3]) } #Example of Function Output# tslope(da$RtnInf, da$RtnSp500) #Verify Function Outputs# summary(lm(da$RtnInf~da$RtnSp500)) # 2) CAPM Data Analysis, Visualisation and Processing (Tsay, 2014 - Chapter 1 R Scripts) #Import dataset# —> students should change the import path according to their desktop(or save CAPM text file to downloads) library(readr) CAPM <- read_delim("Downloads/CAPM.txt", "\t", escape_double = FALSE, trim_ws = TRUE) View(CAPM) #Install and load required packages# install.packages("fBasics") install.packages("rmutil") library(fBasics) library(rmutil) #View data# da = CAPM head(da) #Arrange plot space layout# par(mfrow = c(1,2)) #S&P Global Infrastructure Index# SPGTIND = da$RtnInf hist(SPGTIND, nclass = 30, col = "blue") # Histogram# d1 = density(SPGTIND) # Obtain density estimate of SPGTIND# d1 range(SPGTIND) # Range of SPGTIND# x1 = seq(-.15, .15,.001) # Create a sequence of x with increment 0.001# y1 = dnorm(x1, mean(SPGTIND), stdev(SPGTIND)) plot(d1$x, d1$y, xlab = 'SPGTIND', ylab = 'density', type = 'l', col = "deepskyblue") lines(x1, y1, lty = 2) #S&P 500 Index# SP500 = da$RtnSp500 hist(SP500, nclass = 30, col = "blue") # Histogram# d2 = density(SP500) # Obtain density estimate of SP500# range(SP500) # Range of SP500# x2 = seq(-.15,.1,.001) # Create a sequence of x with increment 0.001# y2 = dnorm(x2, mean(SP500), stdev(SP500)) plot(d2$x, d2$y, xlab = 'SP500', ylab = 'density', type = 'l', col = "deepskyblue") lines(x2, y2, lty = 2) #S&P U.S. Treasury Bond Current 10-Year Index# US_10 = da$Tb10y hist(US_10, nclass = 30, col = "blue") # Histogram# d3 = density(US_10) # Obtain density estimate of US_10# range(US_10) # Range of US_10# x3 = seq(-.03,.03,.001) # Create a sequence of x with increment 0.001# y3 = dnorm(x3, mean(US_10), stdev(US_10)) plot(d3$x, d3$y, xlab = 'US_10', ylab = 'density', type = 'l', col = "deepskyblue") lines(x3, y3, lty = 2) #Difference Between Log returns of S&P Global Infrastructure Index and U.S. Treasury Bond Index# ER1 = da$InfrMTb10 hist(ER1, nclass = 30, col = "blue") # Histogram# d4 = density(ER1) # Obtain density estimate of ER1# range(ER1) # Range of ER1# x4 = seq(-.15,.15,.001) # Create a sequence of x with increment 0.001# y4 = dnorm(x4, mean(ER1), stdev(ER1)) plot(d4$x, d4$y, xlab = 'ER1', ylab = 'density', type = 'l', col = "deepskyblue") lines(x4, y4, lty = 2) #Difference Between Log returns of S&P 500 Index and U.S. Treasury Bond Index# ERP = da$Sp500MTb10 hist(ERP, nclass = 30, col = "blue") # Histogram# d5 = density(ERP) # Obtain density estimate of ERP# range(ERP) # Range of ERP# x5 = seq(-.2,.1,.001) # Create a sequence of x with increment 0.001# y5 = dnorm(x5, mean(ERP), stdev(ERP)) plot(d5$x, d5$y, xlab = 'ERP', ylab = 'density', type = 'l', col = "deepskyblue") lines(x5, y5, lty = 2) #Re-arrange plot space layout# par(mfrow = c(2,3)) #Data Visualisation - Histograms + Laplace Distribution + Normal Distribution# hist(SPGTIND, nclass = 30, freq = F, col = "blue", xlim = c(-.15,.15)) y6 = dlaplace(x1, median(SPGTIND), 0.008446) lines(x1, y6, col = "red") # Add the Laplace distribution line# lines(x1,y1, col = "green") # Add the normal distribution line# legend("topright", legend = c("Laplace","Normal"), lty = 1, col = c("red","green"), cex = 0.5, bty = "n") hist(SP500, nclass = 30, freq = F, col = "blue", xlim = c(-.15,.1)) y7 = dlaplace(x2, median(SP500), 0.009430) lines(x2, y7, col = "red") #Add the Laplace distribution line# lines(x2,y2, col = "green") # Add the normal distribution line# legend("topright", legend = c("Laplace","Normal"), lty = 1, col = c("red","green"), cex = 0.5, bty = "n") hist(US_10, nclass = 30, freq = F, col = "blue", xlim = c(-.03,.03)) y8 = dlaplace(x3, median(US_10), 0.004596) lines(x3, y8, col = "red") #Add the Laplace distribution line# lines(x3,y3, col = "green") # Add the normal distribution line# legend("topright", legend = c("Laplace","Normal"), lty = 1, col = c("red","green"), cex = 0.5, bty = "n") hist(ER1, nclass = 30, freq = F, col = "blue", xlim = c(-.15,.15)) y9 = dlaplace(x4, median(ER1), 0.010309) lines(x4, y9, col = "red") #Add the Laplace distribution line# lines(x4,y4, col = "green") # Add the normal distribution line# legend("topright", legend = c("Laplace","Normal"), lty = 1, col = c("red","green"), cex = 0.5, bty = "n") hist(ERP, nclass = 30, freq = F, col = "blue", xlim = c(-.2,.1)) y10 = dlaplace(x5, median(ERP), 0.012310) lines(x5, y10, col = "red") #Add the Laplace distribution line# lines(x5,y5, col = "green") # Add the normal distribution line# legend("topright", legend = c("Laplace","Normal"), lty = 1, col = c("red","green"), cex = 0.5, bty = "n") #Re-arrange plot space layout# par(mfrow = c(1,3)) #Log Returns# rt <- seq(1,nrow(da),1) plot(rt, da$RtnInf, type = "l", xlab = "No. data") title(main = "Log Return of Infrastructure Index") plot(rt, da$RtnSp500, type = "l", xlab = "No. data") title(main = "Log Return of S&P 500 Index") plot(rt, da$Tb10y, type = "l", xlab = "No. data") title(main = "Log Return of U.S. Treasury Bond Index") #One Sample t-tests# t.test(da$RtnInf) t.test(da$RtnSp500) t.test(da$Tb10y) #Skewness tests# s1 = skewness(da$RtnInf) s2 = skewness(da$RtnSp500) s3 = skewness(da$Tb10y) N = nrow(da) st1 = s1/sqrt(6/N) st2 = s2/sqrt(6/N) st3 = s3/sqrt(6/N) print(st1) print(st2) print(st3) #Computation of p-value# pp1 = 2*(1-pnorm(st1)) pp2 = 2*(1-pnorm(st2)) pp3 = 2*(1-pnorm(st3)) #Kurtosis tests# k1 = kurtosis(da$RtnInf) k2 = kurtosis(da$RtnSp500) k3 = kurtosis(da$Tb10y) kt1 = k1/sqrt(24/N) kt2 = k2/sqrt(24/N) kt3 = k3/sqrt(24/N) print(kt1) print(kt2) print(kt3) #Normality Tests - Jarque–Bera tests# normalTest(da$RtnInf, method = "jb") normalTest(da$RtnSp500, method = "jb") normalTest(da$Tb10y, method = "jb") #Data Processing - As covered in Excel# View(CAPM) #Import unprocessed dataset# —> students should change the import path according to their desktop(or save CAPM(raw) text file to downloads) library(readr) CAPM_raw_ <- read_delim("Downloads/CAPM(raw).txt", "\t", escape_double = FALSE, trim_ws = TRUE) View(CAPM_raw_) daraw = CAPM_raw_ RtnInf_raw = diff(log(daraw$Pinfr), lag = 1) RtnSp500_raw = diff(log(daraw$SP500), lag = 1) Tb10y_raw = diff(log(daraw$PTB10), lag = 1) RtnInf = c(0,RtnInf_raw) RtnSp500 = c(0, RtnSp500_raw) Tb10y = c(0, Tb10y_raw) InfrMTb10 = RtnInf - Tb10y Sp500MTb10 = RtnSp500 - Tb10y dafinal = cbind(daraw, RtnInf, RtnSp500, Tb10y, InfrMTb10, Sp500MTb10) View(dafinal)