FINC621 Winter 2018-19 Lab Worksheet 04

Stochastic Modeling & Numerical Differentiation (finc621-lab04)

About

This worksheet looks at simulating, both arithmetics and geometric, Brownian motion for stock prices, pricing options, and calculating Greeks using numerical differentiation.

Setup

Remember to always set your working directory to the source file location. Go to ‘Session’, scroll down to ‘Set Working Directory’, and click ‘To Source File Location’. Read carefully the below and follow the instructions to complete the tasks and answer any questions. Submit your work to RPubs as detailed in previous notes.

Note

Always read carefully the instructions on Sakai. For clarity, tasks/questions to be completed/answered are highlighted in red color (visible in preview) and numbered according to their particular placement in the task section. Quite often you will need to add your own code chunk.

Execute all code chunks, preview, publish, and submit link on Sakai follwoing the naming convention. Make sure to add comments to your code where appropriate. Use own language!

---

Task 1: Simulating Brownian Motion

This task follows the two examples in the book R Example 5.1/p 148 for simulating an arithmetic Brownian motion and R Example 5.2/p152 for simulating a geometric motion. Two new packages will be required for this worksheet. They are included in the code chunk below.

#Install package quantmod 
if(!require("quantmod",quietly = TRUE))
  install.packages("quantmod",dependencies = TRUE, repos = "https://cloud.r-project.org")

#Install package sde for Stochastic Differential Equation
if(!require("sde",quietly = TRUE))
  install.packages("sde",dependencies = TRUE, repos = "https://cloud.r-project.org")

#Install package fOptions for pricing and evaluating basic options
if(!require("fOptions",quietly = TRUE))
  install.packages("fOptions",dependencies = TRUE, repos = "https://cloud.r-project.org")

##### 1A) Follow example in book to simulate and plot example of an arithmetic Brownian motion. Consider three values for n (number of points): \(2^2, 2^5, 2^{12}\). Explain behavior as n increases.

##inputs:
alpha=0; sigma=1; T=1; n=2^(2); X0=0.1;
##Generate 1 trajectory
dt=T/n
t=seq(0,T,by=dt)
x=c(X0,alpha*dt+sigma*sqrt(dt)*rnorm(n,mean=0,sd=1))
Xt=cumsum(x)
plot(t,Xt,type='l',xlab="time")

##inputs:
alpha=0; sigma=1; T=1; n=2^(5); X0=0.1;
##Generate 1 trajectory
dt=T/n
t=seq(0,T,by=dt)
x=c(X0,alpha*dt+sigma*sqrt(dt)*rnorm(n,mean=0,sd=1))
Xt=cumsum(x)
plot(t,Xt,type='l',xlab="time")

##inputs:
alpha=0; sigma=1; T=1; n=2^(12); X0=0.1;
##Generate 1 trajectory
dt=T/n
t=seq(0,T,by=dt)
x=c(X0,alpha*dt+sigma*sqrt(dt)*rnorm(n,mean=0,sd=1))
Xt=cumsum(x)
plot(t,Xt,type='l',xlab="time")

As the value of n increases, the time steps T/n becomes smaller and thus will improve the quality of the approximate sample path to a true Brownian path. Therefore, the higher the n, the closer it is to Brownian motion.

##### 1B) Write the mathematical equation representing the values along the y-axis in the above plot. Separately express what the value of each variable in the equation is.

The value along the Y axis can be mathematically written as:

\(X(t) = X(t-1) +\alpha\frac{T}{n} + \sigma\varepsilon_{t}\sqrt{\frac{T}{n}}\)

The value \(\sigma\) measures the expected dri where, X(t) is an arithmetic Browninan, n is the number of steps/periods in the time interval, t is the steps in the time series. \(\frac{T}{n}\) represents the length of t, \(\alpha\) is the drift of the brownian motion.\(\varepsilon_{t}\) s a randomly generated number from a normal distribution of mean 0 and variance 1, \(\sigma\) is the volatility of the Brownian motion

##### 1C) Follow example in book to simulate a geometric Brownian motion.

library(sde)

Loading required package: MASS
Loading required package: stats4
Loading required package: fda
Loading required package: splines
Loading required package: Matrix

Attaching package: ‘fda’

The following object is masked from ‘package:graphics’:

    matplot

sde 2.0.15
Companion package to the book
‘Simulation and Inference for Stochastic Differential Equations With R Examples’
Iacus, Springer NY, (2008)
To check the errata corrige of the book, type vignette("sde.errata")

mu=0.16; sigma=0.2; P0=40;
T = 1/12 ##1 month
nt=50; n=2^(8)
##Generate nt trajectories
dt=T/n; t=seq(0,T,by=dt)
X=matrix(rep(0,length(t)*nt), nrow=nt)
for (i in 1:nt) {X[i,]= GBM(x=P0,r=mu,sigma=sigma,T=T,N=n)}
##Plot
ymax=max(X); ymin=min(X) 
#bounds for simulated prices
plot(t,X[1,],t='l',ylim=c(ymin, ymax), col=1,ylab="Price P(t)",xlab="time t")
for(i in 2:nt){lines(t,X[i,], t='l',ylim=c(ymin, ymax),col=i)}

##### 1D) Write the mathematical equation representing the values along the y-axis in the above plot. Separately express what the value of each variable in the equation is.

This can be mathematically written as:

\(P_T =((u- \frac{\sigma^2}{2})T + \sigma\varepsilon_{t}\sqrt{\Delta t})\)

where,\(P_T\), follows a geometric Brownian motion, the return \(r_t\) is normally distributed with mean \((u- \frac{\sigma^2}{2})T\) is the expected return and variance \(\sigma^2\) is the volatility of the price, \(\varepsilon_{t}\) is a randomly generated number from a normal distribution, t is the time series

##### 1E) How is the geometric Brownian simulation different from the arithmetic simulation? Elaborate.

A limitation of the arithmetic Brownian motion when used for modelling the price of a stock is the possibility of getting negative prices. However, this can be fixed by considering prices to be log normally distributed instead. Using log normally distributed data we generate the geometric Brownian motion and this keeps the price of stock postive.

Task 2: Option Pricing, Greeks, and Numerical Differentiation

Follows the example in the book R Example 5.3/p 156 and R Example 5.5/p 157

##### 2A) Calculate the Call and Put price of the given European Option. Explain what the variables r and b in the function call GBSOption represent.

library(fOptions)
GBSOption(TypeFlag = "c", S = 60, X = 65, Time = 1/4, r = 0.08,b = 0.08, sigma = 0.30)

library(fOptions)
GBSOption(TypeFlag = "p", S = 60, X = 65, Time = 1/4, r = 0.08,b = 0.08, sigma = 0.30)

Title:
 Black Scholes Option Valuation 

Call:
 GBSOption(TypeFlag = "p", S = 60, X = 65, Time = 1/4, r = 0.08, 
     b = 0.08, sigma = 0.3)

Parameters:
          Value:
 TypeFlag p     
 S        60    
 X        65    
 Time     0.25  
 r        0.08  
 b        0.08  
 sigma    0.3   

Option Price:
 5.846286 

Description:
 Tue Dec 11 21:11:56 2018 

The price for the European call option is 2.133372 and for the European put option is 5.846286. The variable “r” in GSB option represents the annualized rate of interest, a numeric value and “b” is the annualized cost of carry rate, a numeric value.

##### 2B) Calculate the particular Greeks Delta, Gamma, and Vega for the above corresponding Call option.

library(fOptions)
GBSGreeks(Selection="delta",TypeFlag="c",S=60,X=65,Time=1/4,r=0.08,b=0.08,sigma=0.30)

[1] 0.3724829

library(fOptions)
GBSGreeks(Selection="gamma",TypeFlag="c",S=60,X=65,Time=1/4,r=0.08,b=0.08,sigma=0.30)

[1] 0.04204276

library(fOptions)
GBSGreeks(Selection="vega",TypeFlag="c",S=60,X=65,Time=1/4,r=0.08,b=0.08,sigma=0.30)

[1] 11.35154

##### 2C) Calculate the Delta, Gamma, and Vega for same option using instead numerical differentiation. Write the mathematical equations corresponding to each numerical differentiation.

Using numerical differentiation, Delta \(\Delta\) = \(\Phi\)(\(d_{1}\)), \(d_{1}\) = \(\frac{ln (P_t/K) +(r+\sigma^2/2)(T-t)}{\sigma \sqrt{T-t}}\)

To get the values for \(d_{1}\), plug in the values in the above equation with the values of the call option

\(d_{1}\) = \(\frac{ln (60/65) +(0.08+0.3^2/2)(0.25)}{0.3 \sqrt{0.25}}\)

s= 60; k =65; r=0.08; sigma = 0.3; t= 0.25
d1 = (log(s/k)+ (r+sigma^2/2)*(t))/(sigma*sqrt(t))
d1

[1] -0.3252847

###value of delta
dnorm(d1)

[1] 0.3783848

Gamma \(\Gamma\) = \(\frac{\phi (d_{1})}{P \sigma \sqrt{T-t}}\) \(\Gamma\) = \(\frac{\phi (-0.3253)}{60 *0.3 \sqrt{0.25}}\)

dnorm (d1)/(s*sigma*sqrt(t))

[1] 0.04204276

Vega \(\nu\) = \(P\sqrt{T-t}\phi(d-{1})\)

\(\nu\) = \(60\sqrt{0.25}\phi(-0.32)\)

s*sqrt(t)*dnorm(d1)

[1] 11.35154

2D) Compare results from 2B and 2C. Share observations.

The values for Delta, Gamma and Vega calculated using the GSB Option and numerical differentiation are the same.

*http://computationalfinance.lsi.upc.edu