Hypothesis testing in the CER Model

Eric Zivot

July 30, 2015

Set options and load packages

options(digits=3, width=70)
# install IntroCompFinR package from R-forge
# use install.packages("IntroCompFinR", repos="http://R-Forge.R-project.org")
library(IntroCompFinR)
## Loading required package: xts
## Loading required package: zoo
## 
## Attaching package: 'zoo'
## 
## The following objects are masked from 'package:base':
## 
##     as.Date, as.Date.numeric
library(mvtnorm)
library(PerformanceAnalytics)
## 
## Attaching package: 'PerformanceAnalytics'
## 
## The following object is masked from 'package:graphics':
## 
##     legend
library(tseries)
library(zoo)
Sys.setenv(TZ="UTC")

Get data from package IntroCompfinR

# get data from IntroCompFin package
data(msftDailyPrices, sp500DailyPrices, sbuxDailyPrices)
msftPrices = to.monthly(msftDailyPrices, OHLC=FALSE)
sp500Prices = to.monthly(sp500DailyPrices, OHLC=FALSE)
sbuxPrices = to.monthly(sbuxDailyPrices, OHLC=FALSE)
# set sample to match book chapter
smpl = "1998-01::2012-05"
msftPrices = msftPrices[smpl]
sp500Prices = sp500Prices[smpl]
sbuxPrices = sbuxPrices[smpl]
# calculate returns
msftRetS = Return.calculate(msftPrices, method="simple")
sp500RetS = Return.calculate(sp500Prices, method="simple")
sbuxRetS = Return.calculate(sbuxPrices, method="simple")
msftRetS = msftRetS[-1]
sp500RetS = sp500RetS[-1]
sbuxRetS = sbuxRetS[-1]
msftRetC = log(1 + msftRetS)
sp500RetC = log(1 + sp500RetS)
sbuxRetC = log(1 + sbuxRetS)
# merged data set
cerRetS = merge(msftRetS, sbuxRetS, sp500RetS)
cerRetC = merge(msftRetC, sbuxRetC, sp500RetC)
colnames(cerRetS) = colnames(cerRetC) = c("MSFT", "SBUX", "SP500")
head(cerRetC, n=3)
##             MSFT   SBUX   SP500
## Feb 1998 0.12786 0.0793 0.06808
## Mar 1998 0.05421 0.1343 0.04874
## Apr 1998 0.00689 0.0610 0.00904

Testing \(H_0: \mu = 0\) vs. \(H_1: \mu \ne 0\)

Brute force calculation of t-statistic

  1. Get estimates of \(\mu\)
n.obs = nrow(cerRetC)
muhat.vals = apply(cerRetC, 2, mean)
muhat.vals
##    MSFT    SBUX   SP500 
## 0.00413 0.01466 0.00169
  1. Calculate SEs
# calculate standard errors
sigmahat.vals = apply(cerRetC, 2, sd)
se.muhat = sigmahat.vals/sqrt(n.obs)
se.muhat
##    MSFT    SBUX   SP500 
## 0.00764 0.00851 0.00370
  1. Calculate t-statistics
t.stats = muhat.vals/se.muhat
abs(t.stats)
##  MSFT  SBUX SP500 
## 0.540 1.722 0.457
  1. Apply decision rule: all t-statistics have absolute value less than 2, so do not reject null hypothesis at 5% level of significance.

P-value calculation

P-value = significance level at which null hypothesis is just rejected

# calculate 2-sided p-value
2*(1-pnorm(abs(t.stats)))
##   MSFT   SBUX  SP500 
## 0.5893 0.0851 0.6480

All p-values are less than 5%, so do not reject null hypothesis at 5% level.

Use 95% confidence interval to perform hypothesis test

\(H_0: \mu = 0\) vs. \(H_1: \mu \ne 0\) is rejected at the 5% level if \(\mu = 0\) is not in the 95% confidence interval.

# calculate 95% confidence interval
lower = muhat.vals - 2*se.muhat
upper = muhat.vals + 2*se.muhat
cbind(lower, upper)
##          lower   upper
## MSFT  -0.01116 0.01941
## SBUX  -0.00237 0.03168
## SP500 -0.00570 0.00908

Here, \(\mu = 0\) is in all of the 95% confidence intervals so we do not reject \(H_0: \mu = 0\) vs. \(H_1: \mu \ne 0\) at the 5% level for any asset.

Test for normal distribution

Graphical descriptive statistics suggest that the distribution of MSFT returns have fatter tails than the normal distribution.

fourPanelPlot(cerRetC[, "MSFT"])

However, graphical diagnostics are not a formal statistical test.

Jarque-Bera test for normality

Excess kurtosis value indicates non-normality.

msft.skew = skewness(cerRetC[, "MSFT"])
msft.ekurt = kurtosis(cerRetC[, "MSFT"])
msft.skew
## [1] -0.0901
msft.ekurt
## [1] 2.08

JB statistic confirms non-mormality

JB = n.obs*(msft.skew^2 + 0.25*msft.ekurt^2)/6
JB
## [1] 31.2

JB > 6 so reject \(H_0: r_t \sim N(\mu, \sigma^2)\) at 5% level.

R package tseries function jarque.bera.test()

library(tseries)
jarque.bera.test(cerRetC[, "MSFT"])
## 
##  Jarque Bera Test
## 
## data:  cerRetC[, "MSFT"]
## X-squared = 30, df = 2, p-value = 2e-07

Here, “X-squared” denotes the JB statistic. The small p-value indicates rejection of null hypothesis at any reasonable significance level.

Test for no serial correlation

The R function acf() computes the sample autocorrelations \(\hat{\rho_j}\), plots them and shows the critical values for rejecting the null hypothesis \(H_0: \rho_j = 0\)

acf(cerRetC[, "MSFT"], lwd=2)

Here, \(\hat{rho_1}\) and \(\hat{rho_3}\) look to be statistically differnet from zero at the 5% level.

Compute rolling means

Use the zoo function rollapply() to compute rolling mean estimates.

# compute rolling means over 24 month windows
roll.muhat = rollapply(cerRetC[, "MSFT"], width=24,
                       FUN=mean, align="right")
class(roll.muhat)
## [1] "xts" "zoo"
# first 24 values are NA
roll.muhat[23:25]
##            MSFT
## Dec 1999     NA
## Jan 2000 0.0402
## Feb 2000 0.0311

Plot rolling means with returns

plot.zoo(merge(roll.muhat,cerRetC[, "MSFT"]), plot.type="single",
     main="24 month rolling means for MSFT",ylab="returns",
     lwd=c(2,2), col=c("blue","orange"))
abline(h=0)
grid()
legend(x="bottomleft",legend=c("Rolling mean","Monthly returns"),
       lwd=c(2,2), col=c("blue","orange"))

24-month rolling means for S&P500 index

roll.muhat = rollapply(cerRetC[, "SP500"], width=24,
                       FUN=mean, align="right")
plot.zoo(merge(roll.muhat,cerRetC[, "SP500"]), plot.type="single",
     main="24-month rolling means for SP500",ylab="returns",
     lwd=c(2,2), col=c("blue","orange"))
abline(h=0)
grid()
legend(x="bottomleft",legend=c("Rolling mean","Monthly returns"),
       lwd=c(2,2), col=c("blue","orange"))

Compute rolling volatility (\(\sigma\)) estimates

Use the zoo function rollapply() to compute rolling \(\sigma\) estimates.

# compute rolling sd over 24 month windows
roll.sigmahat = rollapply(cerRetC[, "MSFT"], width=24,
                       FUN=sd, align="right")
class(roll.sigmahat)
## [1] "xts" "zoo"
# first 24 values are NA
roll.sigmahat[23:25]
##           MSFT
## Dec 1999    NA
## Jan 2000 0.124
## Feb 2000 0.125

Plot rolling volatility with returns

plot.zoo(merge(roll.sigmahat,cerRetC[, "MSFT"]), plot.type="single",
     main="24-month rolling volatility for MSFT",ylab="returns",
     lwd=c(2,2), col=c("blue","orange"))
abline(h=0)
grid()
legend(x="bottomleft",legend=c("Rolling sd","Monthly returns"),
       lwd=c(2,2), col=c("blue","orange"))

24-month rolling volatility for S&P500 index

roll.sigmahat = rollapply(cerRetC[, "SP500"], width=24,
                       FUN=sd, align="right")
plot.zoo(merge(roll.sigmahat,cerRetC[, "SP500"]), plot.type="single",
     main="24-month rolling volatility for SP500",ylab="returns",
     lwd=c(2,2), col=c("blue","orange"))
abline(h=0)
grid()
legend(x="bottomleft",legend=c("Rolling sd","Monthly returns"),
       lwd=c(2,2), col=c("blue","orange"))

Compute rolling correlations

To compute rolling correlations, you need to write a function to return the pairwise correlation etween two returns.

rhohat = function(x) {
    cor(x)[1,2]
}

Then use rollapply() with this function.

roll.rhohat = rollapply(cerRetC[,c("MSFT","SP500")],
                       width=24,FUN=rhohat, by.column=FALSE,
                       align="right")
class(roll.rhohat)
## [1] "xts" "zoo"
roll.rhohat[23:25]
##           [,1]
## Dec 1999    NA
## Jan 2000 0.652
## Feb 2000 0.649

Plot rolling correlations (\(\rho_{ij}\))

plot.zoo(roll.rhohat, main="24-month rolling correlations b/w MSFT and SP500",
     ylim=c(0,1), lwd=2, col="blue", ylab="rho.hat")
abline(h=0)   
grid()

Summary of Hypothesis Testing in CER Model