---
title: "CER Model"
author: "Eric Zivot"
date: "Thursday, April 30, 2015"
output: slidy_presentation
---

## Set options and load packages

```{r}
options(digits=4, width=70)
# install IntroCompFinR package from R-forge
# use install.packages("IntroCompFinR", repos="http://R-Forge.R-project.org")
library(IntroCompFinR)
library(mvtnorm)
library(PerformanceAnalytics)
library(zoo)
Sys.setenv(TZ="UTC")
```

## Get data from package IntroCompfinR

```{r}
# get daily prices on Microsoft from IntroCompFin package
data(msftDailyPrices)
msftPrices = to.monthly(msftDailyPrices, OHLC=FALSE)
# set sample to match book chapter
smpl = "1998-01::2012-05"
msftPrices = msftPrices[smpl]
# calculate returns
msftRetS = Return.calculate(msftPrices, method="simple")
msftRetS = msftRetS[-1]
msftRetC = log(1 + msftRetS)
head(msftRetC, n=3)

```

## Microsoft cc returns to calibrate Monte Carlo simulation

```{r}
plot.zoo(msftRetC, main="", 
         lwd=2, col="blue", ylab="Monthly CC Return")
abline(h=0, lwd=2)
abline(h=(mean(msftRetC)+sd(msftRetC)), lty="dashed", lwd=2)
abline(h=(mean(msftRetC)-sd(msftRetC)), lty="dashed", lwd=2)

```

Mean is `r mean(msftRetC)` and SD is `r sd(msftRetC)`

## Four panel plot of MSFT cc returns

```{r}
fourPanelPlot(msftRetC)

```

Returns have fatter tails than normal and do not show much time dependence.

## Simulating observations from the CER model calibrated to MSFT cc returns

```{r}
mu = 0.004
sd.e = 0.10
n.obs = length(msftRetC)
set.seed(111)
sim.e = rnorm(n.obs, mean=0, sd=sd.e)
sim.ret = mu + sim.e
sim.ret = xts(sim.ret, index(msftRetC))
colnames(sim.ret) = "MSFT.sim"
head(sim.ret, n=3)
```

## Plot simulated returns

```{r}
plot.zoo(sim.ret, main="", 
         lwd=2, col="blue", ylab="Monthly CC Return")
abline(h=0, lwd=2)
abline(h=(mu+sd.e), lty="dashed", lwd=2)
abline(h=(mu-sd.e), lty="dashed", lwd=2)
```

Simulated cc returns look like MSFT cc returns except they have constant vol

## Four panel plot of simulated returns

```{r}
fourPanelPlot(sim.ret)
```

Simulated cc returns share similar stylized facts as MSFT cc returns

## Simulated log prices for Msft

```{r}
sim.p = 2.592 + mu*seq(n.obs) + cumsum(sim.e)
sim.P = exp(sim.p)
# decompose RW into components
sim.sume = zoo(cumsum(sim.e), index(sim.ret))
sim.p = zoo(sim.p, index(sim.ret))
sim.mt = sim.p - sim.sume
sim.P = zoo(sim.P, index(sim.ret))
dataToPlot = merge(sim.sume, sim.mt, sim.p, sim.P)
```

## Plot random walk log prices

```{r}
plot(dataToPlot, screens=c(1,1,1,2), main="", ylab="",
     lwd=2, col=c("red", "green", "blue", "black"),
     lty=c("dotted", "dashed", "solid", "solid"))
legend(x="center", legend=c("p(t)-E[p(t)]", "E[p(t)]", "log p", "P"),
       lwd=2, col=c("red", "green", "blue", "black"),
       lty=c("dotted", "dashed", "solid", "solid"))

```


## Plot MSFT prices with random walk decomposition

```{r}
msftLogPrice = log(msftPrices)
p0.msft = as.numeric(log(first(msftPrices)))
mt = p0.msft+ mu*seq(n.obs+1) 
msft.sume = msftLogPrice - mt
dataToPlot = merge(msft.sume, mt, msftLogPrice, msftPrices)
plot.zoo(dataToPlot, screens=c(1,1,1,2), main="", ylab="",
         lwd=2, col=c("red", "green", "blue", "black"),
         lty=c("dotted", "dashed", "solid", "solid"))
legend(x="center", legend=c("p(t)-E[p(t)]", "E[p(t)]", "log p", "P"),
       lwd=2, col=c("red", "green", "blue", "black"),
       lty=c("dotted", "dashed", "solid", "solid"))

```

Actual log prices look like random walk simulation!

## Data for multivariate simulation

```{r}
# get data from IntroCompFin package
data(sp500DailyPrices, sbuxDailyPrices)
sp500Prices = to.monthly(sp500DailyPrices, OHLC=FALSE)
sbuxPrices = to.monthly(sbuxDailyPrices, OHLC=FALSE)
# set sample to match book chapter
sp500Prices = sp500Prices[smpl]
sbuxPrices = sbuxPrices[smpl]
# calculate returns
sp500RetS = Return.calculate(sp500Prices, method="simple")
sbuxRetS = Return.calculate(sbuxPrices, method="simple")
sp500RetS = sp500RetS[-1]
sbuxRetS = sbuxRetS[-1]
sp500RetC = log(1 + sp500RetS)
sbuxRetC = log(1 + sbuxRetS)
# merged data set
cerRetC = merge(msftRetC, sbuxRetC, sp500RetC)
colnames(cerRetC) = c("MSFT", "SBUX", "SP500")
head(cerRetC, n=3)
```

## Plot 3 assets

```{r}
my.panel <- function(...) {
  lines(...)
  abline(h=0)
}
plot.zoo(cerRetC, col="blue", lwd=2, main="", panel=my.panel)

```

## Plot 3 assets

```{r}
pairs(coredata(cerRetC), col="blue", cex=1.25, pch=16)
```

## Sample descriptive statistics

```{r}
apply(cerRetC, 2, mean)
apply(cerRetC, 2, sd)
cov(cerRetC)
cor(cerRetC)

```

## Multivariate simulation from CER model

```{r}
# set parameters for simulation
muVec = apply(cerRetC, 2, mean)
covMat = cov(cerRetC)
# multivariate simulation
set.seed(123)
returns.sim = rmvnorm(n.obs, mean=muVec, sigma=covMat)
colnames(returns.sim) = colnames(cerRetC)
returns.sim = zoo(returns.sim, index(cerRetC))
head(returns.sim, n=3)
```

## Plot simulated returns

```{r}
plot.zoo(returns.sim, col="blue", lwd=2, main="", panel=my.panel)

```

Simulated returns look like actual returns except they have constant vol

## Plot simulated returns

```{r}
pairs(coredata(returns.sim), col="blue", cex=1.25, pch=16)
```

Pairwise correlations in simulated returns match actual returns