#
# Written by:
# -- 
# John L. Weatherwax                2009-04-21
# 
# email: wax@alum.mit.edu
# 
# Please send comments and especially bug reports to the
# above email address.
# 
#-----

save_plots = F

set.seed(0) 

#
# Cointegration of Midcap Prices 
#

# P1 EPage 420
# 
library(fEcofin)
data(midcapD.ts)
x = midcapD.ts[,2:11] 
prices = exp(apply(x,2,cumsum))
options(digits=3)

library(urca)
summary(ca.jo(prices))

#
# Cointegration Analysis of Yields:
#

# P2-3 EPage 421
# 
library(fEcofin)
library(urca)
data(mk.zero2)
# the maturities considered: 
mk.maturity[2:11,] * 12 
summary(ca.jo(mk.zero2[,2:11]))

#
# Simulation 
#

# P4-7 EPage 421
# 
mu = 1030000
phi = 0.99

niter = 1.e4 # number of simulations of our price process

profit_target = 1.e6 + 20000.

stop_loss = 1.e6 - 50000.

above = rep(0,niter) # indices as to whether we have to liquidate our portfolio during our trading period 
below = rep(0,niter) # indices as to whether we have to liquidate our portfolio during our trading period 
profit = rep(0,niter) # dollar profit made / lost on this run 
strat_holding_period = rep(0,niter)
set.seed(2009)
for( i in 1:niter ){
  t = 0
  P_t = 1.e6 # initial value of 1M

  # Walk forward in time until we liquidate 
  while( TRUE ){
    # Simulate the portfolios new value: 
    P_tp1 = mu + phi * ( P_t - mu ) + rnorm( 1, mean=0, sd=5000 )
    t = t+1 

    if( P_tp1 >= profit_target ){
      above[i] = 1
      profit[i] = 20000.
      strat_holding_period[i] = t
      break
    }
    if( P_tp1 <= stop_loss ){
      below[i] = 1
      profit[i] = -50000.
      strat_holding_period[i] = t
      break
    }
    
    P_t = P_tp1 
  }
}
strat_return = ( profit / 50000. ) / strat_holding_period # compute daily return
stopifnot( abs( mean(above)+mean(below) - 1) < 1.e-6 )

print(mean(profit))
print(mean(below))
print(mean(strat_holding_period))
print(mean(strat_return) * 253) # expected yearly return