#Set Direction of Work
#Data
#Date Monthly Returns S&P/TSX Composite Index
SPTSX = read.table("S&PTSX Composite Index.txt", header = FALSE, sep = "")%>%
tibble()
names(SPTSX) <- c("Date", "Monthly Returns", "S&P/TSX Composite Index")
#Date Bank of Canada 1 Month Commercial Paper Rates
BankofCanada1MonthCommercialPaperRates = read.table("Bank of Canada 1 Month Commercial Paper Rates.txt", sep = "")%>%
tibble()
names(BankofCanada1MonthCommercialPaperRates) <- c("Date", "Bank of Canada 1 Month Commercial Paper Rates")
#Date Monthly Return S.P.Canada.Aggregate.Bond.Index.Total.Return
SPCanadaAggregateBondIndexTotalReturn = read.table("S&P Canada Aggregate Bond Index Total Return.txt", sep = "")%>%
tibble()
names(SPCanadaAggregateBondIndexTotalReturn) <- c("Date", "Monthly Return", "S.P.Canada.Aggregate.Bond.Index.Total.Return")
#Date Monthly Return Bloomberg Commodity total return index
BloombergCommodityTotalReturnIndex = read.table("Bloomberg Commodity Total Return Index.txt", sep = "")
names(BloombergCommodityTotalReturnIndex) <- c("Date", "Monthly Return", "Bloomberg Commodity total return index")
x_referenceAssets = cbind(SPTSX[,2],SPCanadaAggregateBondIndexTotalReturn[,2],BankofCanada1MonthCommercialPaperRates[,2])
names(x_referenceAssets) <- c("SPTSX Monthly Return", "SPCanda Monthly Return", "Bank of Canada Monthly Return")
i3 = matrix(c(1),nrow=length(x_referenceAssets[,1]), ncol=2)
i2 = matrix(c(1),nrow=length(x_referenceAssets[,1]), ncol=1)
i1 = matrix(c(1),nrow=length(x_referenceAssets[,1]), ncol=3)
lead_x_ref <- dplyr::lead(x_referenceAssets,1)%>%na.omit()
#Power utility function
## MODIFIED BY USING a linear combination TO COMPUTE THE RAs FUNCTIONS
RA2 = function(w, x_referenceAssets){
y = 0.95
Condition1 = (y*(w%*%t((1+lead_x_ref))%>%t())^(-2)*(1+lead_x_ref)-i2)*(1+x_referenceAssets[,1])
Condition2 = (y*(w%*%t((1+lead_x_ref))%>%t())^(-2)*(1+lead_x_ref)-i2)*(1+x_referenceAssets[,2])
Condition3 = (y*(w%*%t((1+lead_x_ref))%>%t())^(-2)*(1+lead_x_ref)-i2)*(1+x_referenceAssets[,3])
E = cbind(Condition1, Condition2, Condition3)
return(E)
}
U2 = gmm(RA2, x_referenceAssets, t0 = c(0.1, 0.1,0.1), type=c("twoStep"),wmatrix = c("ident"))
RA4 = function(w, x_referenceAssets){
y = 0.95
Condition1 = (y*(w%*%t((1+lead_x_ref))%>%t())^(-4)*(1+lead_x_ref)-i2)*(1+x_referenceAssets[,1])
Condition2 = (y*(w%*%t((1+lead_x_ref))%>%t())^(-4)*(1+lead_x_ref)-i2)*(1+x_referenceAssets[,2])
Condition3 = (y*(w%*%t((1+lead_x_ref))%>%t())^(-4)*(1+lead_x_ref)-i2)*(1+x_referenceAssets[,3])
E = cbind(Condition1, Condition2, Condition3)
return(E)
}
U4 = gmm(RA4, x_referenceAssets, t0 = c(0.1,0.1,0.1), type=c("iterative"),wmatrix = c("ident"))
RA6 = function(w, x_referenceAssets){
y = 0.95
Condition1 = (y*(w%*%t((1+lead_x_ref))%>%t())^(-6)*(1+lead_x_ref)-i2)*(1+x_referenceAssets[,1])
Condition2 = (y*(w%*%t((1+lead_x_ref))%>%t())^(-6)*(1+lead_x_ref)-i2)*(1+x_referenceAssets[,2])
Condition3 = (y*(w%*%t((1+lead_x_ref))%>%t())^(-6)*(1+lead_x_ref)-i2)*(1+x_referenceAssets[,3])
E = cbind(Condition1, Condition2, Condition3)
return(E)
}
U6 = gmm(RA6, x_referenceAssets, t0 = c(0.1,0.1,0.1), type=c("iterative"),wmatrix = c("ident"))
RA8 = function(w, x_referenceAssets){
y = 0.95
Condition1 = (y*(w%*%t((1+lead_x_ref))%>%t())^(-8)*(1+lead_x_ref)-i2)*(1+x_referenceAssets[,1])
Condition2 = (y*(w%*%t((1+lead_x_ref))%>%t())^(-8)*(1+lead_x_ref)-i2)*(1+x_referenceAssets[,2])
Condition3 = (y*(w%*%t((1+lead_x_ref))%>%t())^(-8)*(1+lead_x_ref)-i2)*(1+x_referenceAssets[,3])
E = cbind(Condition1, Condition2, Condition3)
return(E)
}
U8 = gmm(RA8, x_referenceAssets, t0 = c(0.1,0.1,0.1), type=c("iterative"),wmatrix = c("ident"))
RA10 = function(w, x_referenceAssets){
y = 0.95
Condition1 = (y*(w%*%t((1+lead_x_ref))%>%t())^(-10)*(1+lead_x_ref)-i2)*(1+x_referenceAssets[,1])
Condition2 = (y*(w%*%t((1+lead_x_ref))%>%t())^(-10)*(1+lead_x_ref)-i2)*(1+x_referenceAssets[,2])
Condition3 = (y*(w%*%t((1+lead_x_ref))%>%t())^(-10)*(1+lead_x_ref)-i2)*(1+x_referenceAssets[,3])
E = cbind(Condition1, Condition2, Condition3)
return(E)
}
U10 = gmm(RA10, x_referenceAssets, t0 = c(0.1,0.1,0.1), type=c("iterative"),wmatrix = c("ident"))
wpower = rbind(coef(U2),coef(U4),coef(U6),coef(U8),coef(U10))
#wpower = coef(U2)
# MODIFIED USING MATRIX ALGEBRA TO COMPUTE THE MUs VARIABLES
MU2 = (wpower[1,] %*% t(1+lead_x_ref)^(-2))
MU4 = (wpower[2,] %*% t(1+lead_x_ref)^(-4))
MU6 = (wpower[3,] %*% t(1+lead_x_ref)^(-6))
MU8 =(wpower[4,] %*% t(1+lead_x_ref)^(-8))
MU10 = (wpower[5,] %*% t(1+lead_x_ref)^(-10))
#####
bloombergF <- dplyr::lead(1+BloombergCommodityTotalReturnIndex[,2],1)-(1+BankofCanada1MonthCommercialPaperRates[,2])
bloombergF <- bloombergF%>%na.omit()
x1 <- dplyr::lead(1+x_referenceAssets[,1],1)-(1+BankofCanada1MonthCommercialPaperRates[,2])
x1 <- na.omit(x1)
x2 <- dplyr::lead(1+x_referenceAssets[,2],1)-(1+BankofCanada1MonthCommercialPaperRates[,2])
x2 <- na.omit(x2)
datareg <- cbind(bloombergF, x1, x2, t(MU2), t(MU4), t(MU6), t(MU8), t(MU10))%>%
na.omit()
names(datareg)<- c("BloombergF", "sptsxF", "spCanadaF", "MU2", "MU4", "MU6", "MU8", "MU10")
# VERSION lm estimation WITH INTERCEPT
BCOMPower_Regression1 = lm(BloombergF ~ . , data = datareg)
BCOMPower_Regression2 = lm(BloombergF ~ sptsxF + spCanadaF , data = datareg)
summary(BCOMPower_Regression1)
##
## Call:
## lm(formula = BloombergF ~ ., data = datareg)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.096442 -0.021455 -0.000708 0.020312 0.130649
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 11.784 8.963 1.315 0.190
## sptsxF -4.931 3.927 -1.256 0.210
## spCanadaF 5.978 5.484 1.090 0.277
## MU2 -28.939 22.262 -1.300 0.195
## MU4 4.725 4.016 1.176 0.240
## MU6 37.064 27.699 1.338 0.182
## MU8 -31.583 23.150 -1.364 0.173
## MU10 7.462 5.342 1.397 0.163
##
## Residual standard error: 0.03715 on 306 degrees of freedom
## Multiple R-squared: 0.2596, Adjusted R-squared: 0.2427
## F-statistic: 15.33 on 7 and 306 DF, p-value: < 2.2e-16
summary(BCOMPower_Regression2)
##
## Call:
## lm(formula = BloombergF ~ sptsxF + spCanadaF, data = datareg)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.127115 -0.022001 -0.000772 0.021824 0.128371
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -0.001803 0.002155 -0.836 0.4036
## sptsxF 0.521050 0.051920 10.036 <2e-16 ***
## spCanadaF -0.338917 0.176972 -1.915 0.0564 .
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.03721 on 311 degrees of freedom
## Multiple R-squared: 0.2453, Adjusted R-squared: 0.2405
## F-statistic: 50.55 on 2 and 311 DF, p-value: < 2.2e-16
L <-matrix(c(1, 0, 0, 0, 0, 0,
0, 0, 0, 0, 0, 0,
0, 0, 0, 0, 0, 0,
0, 1, 0, 0, 0, 0,
0, 0, 1, 0, 0, 0,
0, 0, 0, 1, 0, 0,
0, 0, 0, 0, 1, 0,
0, 0, 0, 0, 0, 1),6,8)
int <- rep(1,314)
X <- cbind(int,as.matrix(datareg[,-1]))
tX <- t(X)
XX <- t(X)%*%X
XXinv <- solve(XX, diag(1,8,8))
s2 <- BCOMPower_Regression1[["residuals"]]^2%>%sum()/307
Rb <- L%*%coef(BCOMPower_Regression1)
RXXinvR <- L%*%XXinv%*%t(L)
sRxR <- s2*RXXinvR
sRxRinv <- solve(sRxR, diag(1,6,6))
W <- t(Rb) %*%sRxRinv%*% Rb
chitst <- c(W, pchisq(W,6, lower.tail = FALSE))
names(chitst) <- c("X2", "P(> X2)")
print(chitst,3)
## X2 P(> X2)
## 6.622 0.357
#Testing MU2, MU6, and MU10 and Newy-West
L1 <-matrix(c(1, 0, 0, 0,
0, 0, 0, 0,
0, 0, 0, 0,
0, 1, 0, 0,
0, 0, 0, 0,
0, 0, 1, 0,
0, 0, 0, 0,
0, 0, 0, 1),4,8)
Rb <- L1%*%coef(BCOMPower_Regression1)
NW <- NeweyWest(BCOMPower_Regression1, lag = 5, prewhite = FALSE)
RNWR <- L1%*%NW%*%t(L1)
RNWRinv <- solve(RNWR, diag(1,4,4))
W <- t(Rb) %*%RNWRinv%*% Rb
chitst <- c(W, pchisq(W,6, lower.tail = FALSE))
names(chitst) <- c("X2", "P(> X2)")
print(chitst,3)
## X2 P(> X2)
## 4.023 0.674
#All Mu variables and Newey-West
L <-matrix(c(1, 0, 0, 0, 0, 0,
0, 0, 0, 0, 0, 0,
0, 0, 0, 0, 0, 0,
0, 1, 0, 0, 0, 0,
0, 0, 1, 0, 0, 0,
0, 0, 0, 1, 0, 0,
0, 0, 0, 0, 1, 0,
0, 0, 0, 0, 0, 1),6,8)
Rb <- L%*%coef(BCOMPower_Regression1)
NW <- NeweyWest(BCOMPower_Regression1, lag = 5, prewhite = FALSE)
RNWR <- L%*%NW%*%t(L)
RNWRinv <- solve(RNWR, diag(1,6,6))
W <- t(Rb) %*%RNWRinv%*% Rb
chitst <- c(W, pchisq(W,6, lower.tail = FALSE))
names(chitst) <- c("X2", "P(> X2)")
print(chitst,3)
## X2 P(> X2)
## 4.13 0.66