### Pricing Discrete Barrier Options

*• Barrier options whose barrier is monitored only at*
discrete times are called discrete barrier options.

*• They are more common than the continuously*
monitored versions.

*• The main diﬃculty with pricing discrete barrier options*
lies in matching the monitored times.

*• Here is why.*

### Pricing Discrete Barrier Options (continued)

*• Suppose each period has a duration of Δt and the*

* > 1 monitored times are*

*t*_{0} *= 0, t*_{1}*, t*_{2}*, . . . , t*_{}*= T.*

*• It is unlikely that all monitored times coincide with the*
*end of a period on the tree, or Δt divides t*_{i}*for all i.*

*• The binomial-trinomial tree can handle discrete options*
with ease, however.

*• Simply build a binomial-trinomial tree from time 0 to*
*time t*_{1}*, followed by one from time t*_{1} *to time t*_{2}, and so
*on until time t* .

*****
**j***
**j***
**j***
**j**

*****
**j***
**j***
**j**

**:**
**j**

**R**
**z**

*****

**R**
**z**

*****

**R**
**z**

*****

**R**
**z**

*****

**R**
**z**

*****

*t*

**-**

^{Δt}^{}^{1} ^{Δt}**-**^{1} ^{Δt}**-**^{1} ^{Δt}^{}^{2}**-**

*t*
*2σ**√*

*Δt*1

*2σ**√*
*Δt*2

### Pricing Discrete Barrier Options (concluded)

*• This procedure works even if each t** _{i}* is associated with

*a distinct barrier or if each window [ t*

_{i}*, t*

*) has its own continuously monitored barrier or double barriers.*

_{i+1}### Options on a Stock That Pays Known Dividends

*• Many ad hoc assumptions have been postulated for*
option pricing with known dividends.^{a}

1. The one we saw earlier (p. 305) models the stock price minus the present value of the anticipated dividends as following geometric Brownian motion.

2. One can also model the stock price plus the forward values of the dividends as following geometric

Brownian motion.

aFrishling (2002).

Options on a Stock That Pays Known Dividends (continued)

*• Realistic models assume:*

**– The stock price decreases by the amount of the**
dividend paid at the ex-dividend date.

**– The dividend is part cash and part yield (i.e.,**
*α(t)S*_{0} *+ β(t)S** _{t}*), for practitioners.

^{a}

*• The stock price follows geometric Brownian motion*
between adjacent ex-dividend dates.

*• But they result in binomial trees that grow*
exponentially (recall p. 304).

*• The binomial-trinomial tree can avoid this problem in*
most cases.

Options on a Stock That Pays Known Dividends (continued)

*• Suppose that the known dividend is D dollars and the*
*ex-dividend date is at time t.*

*• So there are m* *= t/Δt periods between time 0 and the*^{Δ}
ex-dividend date.^{a}

*• To avoid negative stock prices, we need to make sure the*
*lowest stock price at time t is at least D, i.e.,*

*Se*^{−(t/Δt)σ}^{√}^{Δt}*≥ D.*

**– Or,**

*Δt* *≥*

*tσ*
*ln(S/D)*

_{2}
*.*

aOr simply assume *m is an integer input and Δt* =^{Δ} *t/m.*

Options on a Stock That Pays Known Dividends (continued)

*• Build a CRR tree from time 0 to time t as before.*

*• Subtract D from all the stock prices on the tree at time*
*t to represent the price drop on the ex-dividend date.*

*• Assume the top node’s price equals S** ^{}*.

**– As usual, its two successor nodes will have prices**
*S*^{}*u and S*^{}*u** ^{−1}*.

*• The remaining nodes’ successor nodes will have prices*
*S*^{}*u*^{−3}*, S*^{}*u*^{−5}*, S*^{}*u*^{−7}*, . . . ,*

same as the CRR tree.

### A Stair Tree

*D*

0

1

2

3

4

*S''*
*S*

*D*
*Su*^{2}−
*D*
*S*−

*D*
*Sd*^{2}−

Options on a Stock That Pays Known Dividends (continued)

*• For each node at time t below the top node, we build*
the trinomial connection.

*• Note that the binomial-trinomial structure remains valid*
*in the special case when Δt*^{}*= Δt on p. 703.*

**– And even with the displacements ±2σ**√

*Δt (as on p.*

731).

Options on a Stock That Pays Known Dividends (concluded)

*• Hence the construction can be completed.*

*• From time t + Δt onward, the standard binomial tree*
will be used until the maturity date or the next

ex-dividend date when the procedure can be repeated.

*• The resulting tree is called the stair tree.*^{a}

aDai (B82506025, R86526008, D8852600) & Lyuu (2004); Dai (B82506025, R86526008, D8852600) (2009).

Other Applications of Binomial-Trinomial Trees

*• Pricing guaranteed minimum withdrawal beneﬁts.*^{a}

*• Option pricing with stochastic volatilities.*^{b}

*• Eﬃcient Parisian option pricing.*^{c}

*• Option pricing with time-varying volatilities and*
time-varying barriers.^{d}

*• Defaultable bond pricing.*^{e}

aH. Wu (R96723058) (2009).

bC. Huang (R97922073) (2010).

cY. Huang (R97922081) (2010).

dC. Chou (R97944012) (2010); C. Chen (R98922127) (2011).

eDai (B82506025, R86526008, D8852600), Lyuu, & C. Wang

### General Properties of Trees

^{a}

*• Consider the Ito process,*

*dX = a(X, t) dt + σ dW,*
*where a(X, t) = O(1) and σ is a constant.*

*• The mean and volatility of the next move’s size are*
*O(Δt) and O(√*

*Δt), respectively.*

*• Note that* *√*

*Δt* * Δt.*

*• The tree spacing must be in the order of σ√*

*Δt if the*
variance is to be matched.^{b}

aChiu (R98723059) (2012); C. H. Wu (R99922149) (2012).

bLyuu & C. Wang (F95922018) (2009, 2011); Lyuu & Wen (D94922003) (2012).

### Merton’s Jump-Diﬀusion Model

*• Empirically, stock returns tend to have fat tails,*

inconsistent with the Black-Scholes model’s assumptions.

*• Stochastic volatility and jump processes have been*
proposed to address this problem.

*• Merton’s (1976) jump-diﬀusion model is our focus.*

### Merton’s Jump-Diﬀusion Model (continued)

*• This model superimposes a jump component on a*
diﬀusion component.

*• The diﬀusion component is the familiar geometric*
Brownian motion.

*• The jump component is composed of lognormal jumps*
driven by a Poisson process.

**– It models the sudden changes in the stock price**

because of the arrival of important new information.

### Merton’s Jump-Diﬀusion Model (continued)

*• Let S*_{t}*be the stock price at time t.*

*• The risk-neutral jump-diﬀusion process for the stock*
price follows

*dS*_{t}

*S*_{t}*= (r* *− λ¯k) dt + σ dW*_{t}*+ k dq*_{t}*.* (96)

*• Above, σ denotes the volatility of the diﬀusion*
component.

### Merton’s Jump-Diﬀusion Model (continued)

*• The jump event is governed by a compound Poisson*
*process q*_{t}*with intensity λ, where k denotes the*
*magnitude of the random jump.*

**– The distribution of k obeys***ln(1 + k)* *∼ N*

*γ, δ*^{2}
with mean ¯*k* *= E (k) = e*^{Δ} ^{γ+δ}^{2}^{/2}*− 1.*

*• The model with λ = 0 reduces to the Black-Scholes*
model.

### Merton’s Jump-Diﬀusion Model (continued)

*• The solution to Eq. (96) on p. 738 is*

*S*_{t}*= S*_{0}*e*^{(r−λ¯}^{k−σ}^{2}^{/2) t+σW}^{t}*U (n(t)),* (97)
where

*U (n(t)) =*

*n(t)*

*i=0*

*(1 + k*_{i}*) .*

**– k**_{i}*is the magnitude of the ith jump with*
*ln(1 + k** _{i}*)

*∼ N(γ, δ*

^{2}).

**– k**_{0} = 0.

**– n(t) is a Poisson process with intensity λ.**

### Merton’s Jump-Diﬀusion Model (concluded)

*• Recall that n(t) denotes the number of jumps that*
*occur up to time t.*

*• As k > −1, stock prices will stay positive.*

*• The geometric Brownian motion, the lognormal jumps,*
and the Poisson process are assumed to be independent.

### Tree for Merton’s Jump-Diﬀusion Model

^{a}

*• Deﬁne the S-logarithmic return of the stock price S** ^{}* as

*ln(S*

^{}*/S).*

*• Deﬁne the logarithmic distance between stock prices S*^{}*and S as*

*| ln(S** ^{}*)

*− ln(S) | = | ln(S*

^{}*/S)|.*

aDai (B82506025, R86526008, D8852600), C. Wang (F95922018), Lyuu,

& Y. Liu (2010).

### Tree for Merton’s Jump-Diﬀusion Model (continued)

*• Take the logarithm of Eq. (97) on p. 740:*

*M** _{t}* = ln

^{Δ}

*S*_{t}*S*_{0}

*= X*_{t}*+ Y*_{t}*,* (98)
where

*X** _{t}* =

^{Δ}

*r* *− λ¯k −* *σ*^{2}
2

*t + σW*_{t}*,* (99)

*Y** _{t}* =

^{Δ}

*n(t)*
*i=0*

*ln (1 + k*_{i}*) .* (100)

*• It decomposes the S*_{0}*-logarithmic return of S** _{t}* into the

*diﬀusion component X*

_{t}*and the jump component Y*

*.*

_{t}### Tree for Merton’s Jump-Diﬀusion Model (continued)

*• Motivated by decomposition (98) on p. 743, the tree*
construction divides each period into a diﬀusion phase
followed by a jump phase.

*• In the diﬀusion phase, X** _{t}* is approximated by the
BOPM.

*• So X*_{t}*makes an up move to X*_{t}*+ σ√*

*Δt with*
*probability p*_{u}*or a down move to X*_{t}*− σ√*

*Δt with*
*probability p** _{d}*.

### Tree for Merton’s Jump-Diﬀusion Model (continued)

*• According to BOPM,*

*p** _{u}* =

*e*

^{μΔt}*− d*

*u*

*− d*

*,*

*p*

*= 1*

_{d}*− p*

_{u}*,*

*except that μ = r*

*− λ¯k here.*

*• The diﬀusion component gives rise to diﬀusion nodes.*

*• They are spaced at 2σ√*

*Δt apart such as the white*
nodes A, B, C, D, E, F, and G on p. 746.

( − Δ*1) t* Δ* ^{t}* ( + Δ

*1) t*

*q*1

*q*_{−}1

*p**u*

*p**d*

*q*0

2 2

*h*= γ +δ

2σ Δ*t*

White nodes are *diﬀusion nodes.*

Gray nodes are *jump nodes.* In
the diﬀusion phase, the solid black
lines denote the binomial structure
of BOPM; the dashed lines denote
the trinomial structure. Here *m is*
set to one for simplicity. Only the
double-circled nodes will remain af-
ter the construction. Note that a
and b are diﬀusion nodes because
no jump occurs in the jump phase.

### Tree for Merton’s Jump-Diﬀusion Model (concluded)

*• In the jump phase, Y** _{t+Δt}* is approximated by moves

*from each diﬀusion node to 2m jump nodes that match*
*the ﬁrst 2m moments of the lognormal jump.*

*• The m jump nodes above the diﬀusion node are spaced*
*at h* =^{Δ}

*γ*^{2} *+ δ*^{2} apart.

*• The same holds for the m jump nodes below the*
diﬀusion node.

*• The gray nodes at time Δt on p. 746 are jump nodes.*

*• The size of the tree is O(n** ^{2.5}*).

### Multivariate Contingent Claims

*• They depend on two or more underlying assets.*

*• The basket call on m assets has the terminal payoﬀ*
max

_{m}

*i=1*

*α*_{i}*S*_{i}*(τ )* *− X, 0*

*,*
*where α*_{i}*is the percentage of asset i.*

*• Basket options are essentially options on a portfolio of*
stocks; they are index options.

*• Option on the best of two risky assets and cash has a*
*terminal payoﬀ of max(S*_{1}*(τ ), S*_{2}*(τ ), X).*

### Multivariate Contingent Claims (concluded)

^{a}

Name Payoﬀ

Exchange option max(*S*1(*τ) − S*2(*τ), 0)*
Better-oﬀ option max(*S*1(*τ), . . . , S** _{k}*(

*τ), 0)*Worst-oﬀ option min(

*S*1(

*τ), . . . , S*

*(*

_{k}*τ), 0)*

Binary maximum option *I{ max(S*1(*τ), . . . , S** _{k}*(

*τ)) > X }*Maximum option max(max(

*S*1(

*τ), . . . , S*

*(*

_{k}*τ)) − X, 0)*Minimum option max(min(

*S*1(

*τ), . . . , S*

*(*

_{k}*τ)) − X, 0)*Spread option max(

*S*1(

*τ) − S*2(

*τ) − X, 0)*

Basket average option max((*S*1(*τ) + · · · + S** _{k}*(

*τ))/k − X, 0)*Multi-strike option max(

*S*1(

*τ) − X*1

*, . . . , S*

*(*

_{k}*τ) − X*

_{k}*, 0)*

Pyramid rainbow option max(*| S*1(*τ) − X*1 *| + · · · + | S** _{k}*(

*τ) − X*

_{k}*| − X, 0)*

Madonna option max(

(*S*1(*τ) − X*1)^{2} + *· · · + (S** _{k}*(

*τ) − X*

*)*

_{k}^{2}

*− X, 0)*

aLyuu & Teng (R91723054) (2011).

### Correlated Trinomial Model

^{a}

*• Two risky assets S*_{1} *and S*_{2} follow
*dS*_{i}

*S*_{i}*= r dt + σ*_{i}*dW*_{i}*in a risk-neutral economy, i = 1, 2.*

*• Let*

*M** _{i}* =

^{Δ}

*e*

^{rΔt}*,*

*V** _{i}* =

^{Δ}

*M*

_{i}^{2}

*(e*

^{σ}

^{i}^{2}

^{Δt}*− 1).*

**– S**_{i}*M*_{i}*is the mean of S*_{i}*at time Δt.*

**– S**_{i}^{2}*V*_{i}*the variance of S*_{i}*at time Δt.*

### Correlated Trinomial Model (continued)

*• The value of S*_{1}*S*_{2} *at time Δt has a joint lognormal*
*distribution with mean S*_{1}*S*_{2}*M*_{1}*M*_{2}*e*^{ρσ}^{1}^{σ}^{2}^{Δt}*, where ρ is*
*the correlation between dW*_{1} *and dW*_{2}.

*• Next match the 1st and 2nd moments of the*

approximating discrete distribution to those of the continuous counterpart.

*• At time Δt from now, there are 5 distinct outcomes.*

### Correlated Trinomial Model (continued)

*• The ﬁve-point probability distribution of the asset prices*
is

Probability Asset 1 Asset 2
*p*_{1} *S*_{1}*u*_{1} *S*_{2}*u*_{2}
*p*_{2} *S*_{1}*u*_{1} *S*_{2}*d*_{2}
*p*_{3} *S*_{1}*d*_{1} *S*_{2}*d*_{2}
*p*_{4} *S*_{1}*d*_{1} *S*_{2}*u*_{2}

*p*_{5} *S*_{1} *S*_{2}

*• As usual, impose u*_{i}*d** _{i}* = 1.

### Correlated Trinomial Model (continued)

*• The probabilities must sum to one, and the means must*
be matched:

1 = *p*_{1} *+ p*_{2} *+ p*_{3} *+ p*_{4} *+ p*_{5}*,*

*S*_{1}*M*_{1} = *(p*_{1} *+ p*_{2}*) S*_{1}*u*_{1} *+ p*_{5}*S*_{1} *+ (p*_{3} *+ p*_{4}*) S*_{1}*d*_{1}*,*
*S*_{2}*M*_{2} = *(p*_{1} *+ p*_{4}*) S*_{2}*u*_{2} *+ p*_{5}*S*_{2} *+ (p*_{2} *+ p*_{3}*) S*_{2}*d*_{2}*.*

### Correlated Trinomial Model (concluded)

*• Let R* *= M*^{Δ} _{1}*M*_{2}*e*^{ρσ}^{1}^{σ}^{2}* ^{Δt}*.

*• Match the variances and covariance:*

*S*1^{2}*V*^{1} = *(p*^{1} *+ p*^{2}*)((S*^{1}*u*^{1})^{2} *− (S*1*M*^{1})^{2}*) + p*^{5}*(S*1^{2} *− (S*1*M*^{1})^{2})
*+(p*3 *+ p*4*)((S*1*d*1)^{2} *− (S*1*M*1)^{2}*),*

*S*_{2}^{2}*V*2 = *(p*1 *+ p*4*)((S*2*u*2)^{2} *− (S*2*M*2)^{2}*) + p*5*(S*_{2}^{2} *− (S*2*M*2)^{2})
*+(p*2 *+ p*3*)((S*2*d*2)^{2} *− (S*2*M*2)^{2}*),*

*S*^{1}*S*^{2}*R* = *(p*^{1}*u*^{1}*u*^{2} *+ p*^{2}*u*^{1}*d*^{2} *+ p*^{3}*d*^{1}*d*^{2} *+ p*^{4}*d*^{1}*u*^{2} *+ p*^{5}*) S*^{1}*S*^{2}*.*

*• The solutions appear on p. 246 of the textbook.*

### Correlated Trinomial Model Simpliﬁed

^{a}

*• Let μ*^{}_{i}*= r*^{Δ} *− σ*_{i}^{2}*/2 and u*_{i}*= e*^{Δ} ^{λσ}^{i}

*√**Δt* *for i = 1, 2.*

*• The following simpler scheme is good enough:*

*p1* = 1

4

1
*λ2* +

*√**Δt*
*λ*

*μ1*

*σ1* + *μ2*
*σ2*

+ *ρ*

*λ2*

*,*

*p2* = 1

4

1
*λ2* +

*√**Δt*
*λ*

*μ1*

*σ1* *−* *μ2*
*σ2*

*−* *ρ*
*λ2*

*,*

*p3* = 1

4

1
*λ2* +

*√**Δt*
*λ*

*−**μ1*

*σ1* *−* *μ2*
*σ2*

+ *ρ*

*λ2*

*,*

*p4* = 1

4

1
*λ2* +

*√**Δt*
*λ*

*−**μ1*

*σ1* + *μ2*
*σ2*

*−* *ρ*
*λ2*

*,*

*p5* = *1 −* 1
*λ2* *.*

aMadan, Milne, & Shefrin (1989).

### Correlated Trinomial Model Simpliﬁed (continued)

*• All of the probabilities lie between 0 and 1 if and only if*

*−1 + λ**√*
Δ*t*

*μ*^{}_{1}
*σ*1

+ *μ*^{}_{2}
*σ*2

* ≤ ρ ≤ 1 − λ**√*
Δ*t*

*μ*^{}_{1}
*σ*1

*−* *μ*^{}_{2}
*σ*2

*,(101)*

1 *≤ λ* (102)

*• We call a multivariate tree (correlation-) optimal if it*
guarantees valid probabilities as long as

*−1 + O(√*

*Δt) < ρ < 1* *− O(√*

*Δt),*
such as the above one.^{a}

aW. Kao (R98922093) (2011); W. Kao (R98922093), Lyuu, & Wen (D94922003) (2014).

### Correlated Trinomial Model Simpliﬁed (concluded)

*• But this model cannot price 2-asset 2-barrier options*
accurately.^{a}

*• Few multivariate trees are both optimal and able to*
handle multiple barriers.^{b}

*• An alternative is to use orthogonalization.*^{c}

aSee Y. Chang (B89704039, R93922034), Hsu (R7526001, D89922012),

& Lyuu (2006); W. Kao (R98922093), Lyuu, & Wen (D94922003) (2014) for solutions.

bSee W. Kao (R98922093), Lyuu, & Wen (D94922003) (2014) for one.

cHull & White (1990); Dai (B82506025, R86526008, D8852600), Lyuu,

& C. Wang (F95922018) (2012).

### Extrapolation

*• It is a method to speed up numerical convergence.*

*• Say f(n) converges to an unknown limit f at rate of*
*1/n:*

*f (n) = f +* *c*

*n* *+ o*

1
*n*

*.* (103)

*• Assume c is an unknown constant independent of n.*

**– Convergence is basically monotonic and smooth.**

### Extrapolation (concluded)

*• From two approximations f(n*_{1}*) and f (n*_{2}) and
ignoring the smaller terms,

*f (n*_{1}) = *f +* *c*
*n*_{1} *,*
*f (n*_{2}) = *f +* *c*

*n*_{2} *.*

*• A better approximation to the desired f is*
*f =* *n*_{1}*f (n*_{1}) *− n*_{2}*f (n*_{2})

*n*_{1} *− n*_{2} *.* (104)

*• This estimate should converge faster than 1/n.*^{a}

*• The Richardson extrapolation uses n*_{2} *= 2n*_{1}.

aIt is identical to the forward rate formula (21) on p. 140!

### Improving BOPM with Extrapolation

*• Consider standard European options.*

*• Denote the option value under BOPM using n time*
*periods by f (n).*

*• It is known that BOPM convergences at the rate of 1/n,*
consistent with Eq. (103) on p. 758.

*• But the plots on p. 288 (redrawn on next page)*

demonstrate that convergence to the true option value
*oscillates with n.*

*• Extrapolation is inapplicable at this stage.*

5 10 15 20 25 30 35 n

11.5 12 12.5 13

Call value

0 10 20 30 40 50 60 n

15.1 15.2 15.3 15.4 15.5

Call value

### Improving BOPM with Extrapolation (concluded)

*• Take the at-the-money option in the left plot on p. 761.*

*• The sequence with odd n turns out to be monotonic*
and smooth (see the left plot on p. 763).^{a}

*• Apply extrapolation (104) on p. 759 with n*_{2} *= n*_{1} + 2,
*where n*_{1} is odd.

*• Result is shown in the right plot on p. 763.*

*• The convergence rate is amazing.*

*• See Exercise 9.3.8 of the text (p. 111) for ideas in the*
general case.

5 10 15 20 25 30 35 n

12.2 12.4 12.6 12.8 13 13.2 13.4

Call value

5 10 15 20 25 30 35 n

12.11 12.12 12.13 12.14 12.15 12.16 12.17

Call value

*Numerical Methods*

All science is dominated by the idea of approximation.

— Bertrand Russell

### Finite-Diﬀerence Methods

*• Place a grid of points on the space over which the*
desired function takes value.

*• Then approximate the function value at each of these*
points (p. 767).

*• Solve the equation numerically by introducing diﬀerence*
equations in place of derivatives.

0 0.05 0.1 0.15 0.2 0.25 80

85 90 95 100 105 110 115

### Example: Poisson’s Equation

*• It is ∂*^{2}*θ/∂x*^{2} *+ ∂*^{2}*θ/∂y*^{2} = *−ρ(x, y), which describes the*
electrostatic ﬁeld.

*• Replace second derivatives with ﬁnite diﬀerences*
through central diﬀerence.

*• Introduce evenly spaced grid points with distance of Δx*
*along the x axis and Δy along the y axis.*

*• The ﬁnite diﬀerence form is*

*−ρ(x*_{i}*, y** _{j}*) =

*θ(x*

_{i+1}*, y*

*)*

_{j}*− 2θ(x*

_{i}*, y*

_{j}*) + θ(x*

_{i−1}*, y*

*)*

_{j}*(Δx)*

^{2}

+*θ(x*_{i}*, y** _{j+1}*)

*− 2θ(x*

_{i}*, y*

_{j}*) + θ(x*

_{i}*, y*

*)*

_{j−1}*(Δy)*^{2} *.*

### Example: Poisson’s Equation (concluded)

*• In the above, Δx* *= x*^{Δ} _{i}*− x*_{i−1}*and Δy* *= y*^{Δ} _{j}*− y** _{j−1}* for

*i, j = 1, 2, . . . .*

*• When the grid points are evenly spaced in both axes so*
*that Δx = Δy = h, the diﬀerence equation becomes*

*−h*^{2}*ρ(x*_{i}*, y*_{j}*) = θ(x*_{i+1}*, y*_{j}*) + θ(x*_{i−1}*, y** _{j}*)

*+θ(x*

_{i}*, y*

_{j+1}*) + θ(x*

_{i}*, y*

*)*

_{j−1}*− 4θ(x*

_{i}*, y*

_{j}*).*

*• Given boundary values, we can solve for the x** _{i}*s and the

*y*

*s within the square [*

_{j}*±L, ±L ].*

*• From now on, θ** _{i,j}* will denote the ﬁnite-diﬀerence

*approximation to the exact θ(x*

_{i}*, y*

*).*

_{j}### Explicit Methods

*• Consider the diﬀusion equation*
*D(∂*^{2}*θ/∂x*^{2}) *− (∂θ/∂t) = 0, D > 0.*

*• Use evenly spaced grid points (x*_{i}*, t** _{j}*) with distances

*Δx and Δt, where Δx* *= x*^{Δ} _{i+1}*− x*_{i}*and Δt* *= t*^{Δ} _{j+1}*− t** _{j}*.

*• Employ central diﬀerence for the second derivative and*
forward diﬀerence for the time derivative to obtain

*∂θ(x, t)*

*∂t*

*t=t*_{j}

= *θ(x, t** _{j+1}*)

*− θ(x, t*

*)*

_{j}Δ*t* + *· · · ,* (105)

*∂*^{2}*θ(x, t)*

*∂x*^{2}

*x=x*_{i}

= *θ(x*_{i+1}*, t) − 2θ(x*_{i}*, t) + θ(x*_{i−1}*, t)*

(Δ*x)*^{2} + *· · · .(106)*

### Explicit Methods (continued)

*• Next, assemble Eqs. (105) and (106) into a single*
*equation at (x*_{i}*, t** _{j}*).

*• But we need to decide how to evaluate x in the ﬁrst*
*equation and t in the second.*

*• Since central diﬀerence around x**i* is used in Eq. (106),
*we might as well use x*_{i}*for x in Eq. (105).*

*• Two choices are possible for t in Eq. (106).*

*• The ﬁrst choice uses t = t** _{j}* to yield the following
ﬁnite-diﬀerence equation,

*θ*_{i,j+1}*− θ*_{i,j}

*Δt* *= D* *θ*_{i+1,j}*− 2θ*_{i,j}*+ θ*_{i−1,j}

*(Δx)*^{2} *.*

(107)

### Explicit Methods (continued)

*• The stencil of grid points involves four values, θ** _{i,j+1}*,

*θ*

_{i,j}*, θ*

_{i+1,j}*, and θ*

*.*

_{i−1,j}*• Rearrange Eq. (107) on p. 771 as*

*θ**i,j+1* *= DΔt*

*(Δx)*^{2} *θ**i+1,j* +

*1 − 2DΔt*
*(Δx)*^{2}

*θ**i,j* *+ DΔt*

*(Δx)*^{2} *θ**i−1,j**.*

*• We can calculate θ*_{i,j+1}*from θ*_{i,j}*, θ*_{i+1,j}*, θ** _{i−1,j}*, at the

*previous time t*

*(see exhibit (a) on next page).*

_{j}### Stencils

t_{j} t_{j }
x_{i }

x_{i }
x_{i}

t_{j} t_{j }
x_{i }

x_{i }
x_{i}

= >

### Explicit Methods (concluded)

*• Starting from the initial conditions at t*_{0}, that is,
*θ*_{i,0}*= θ(x*_{i}*, t*_{0}*), i = 1, 2, . . . , we calculate*

*θ*_{i,1}*,* *i = 1, 2, . . . .*

*• And then*

*θ*_{i,2}*,* *i = 1, 2, . . . .*

*• And so on.*

### Stability

*• The explicit method is numerically unstable unless*
*Δt* *≤ (Δx)*^{2}*/(2D).*

**– A numerical method is unstable if the solution is**
highly sensitive to changes in initial conditions.

*• The stability condition may lead to high running times*
and memory requirements.

*• For instance, halving Δx would imply quadrupling*
*(Δt)** ^{−1}*, resulting in a running time 8 times as much.

### Explicit Method and Trinomial Tree

*• Recall that*

*θ**i,j+1* *= DΔt*

*(Δx)*^{2} *θ**i+1,j* +

*1 − 2DΔt*
*(Δx)*^{2}

*θ**i,j* *+ DΔt*

*(Δx)*^{2} *θ**i−1,j**.*

*• When the stability condition is satisﬁed, the three*
*coeﬃcients for θ*_{i+1,j}*, θ*_{i,j}*, and θ** _{i−1,j}* all lie between
zero and one and sum to one.

*• They can be interpreted as probabilities.*

*• So the ﬁnite-diﬀerence equation becomes identical to*
backward induction on trinomial trees!

### Explicit Method and Trinomial Tree (concluded)

*• The freedom in choosing Δx corresponds to similar*
freedom in the construction of trinomial trees.

*• The explicit ﬁnite-diﬀerence equation is also identical to*
backward induction on a binomial tree.^{a}

**– Let the binomial tree take 2 steps each of length**
*Δt/2.*

**– It is now a trinomial tree.**

aHilliard (2014).

### Implicit Methods

*• Suppose we use t = t** _{j+1}* in Eq. (106) on p. 770 instead.

*• The ﬁnite-diﬀerence equation becomes*
*θ*_{i,j+1}*− θ*_{i,j}

*Δt* *= D* *θ*_{i+1,j+1}*− 2θ*_{i,j+1}*+ θ*_{i−1,j+1}

*(Δx)*^{2} *.*

(108)

*• The stencil involves θ*_{i,j}*, θ*_{i,j+1}*, θ*_{i+1,j+1}*, and θ** _{i−1,j+1}*.

*• This method is implicit:*

**– The value of any one of the three quantities at t*** _{j+1}*
cannot be calculated unless the other two are known.

**– See exhibit (b) on p. 773.**

### Implicit Methods (continued)

*• Equation (108) can be rearranged as*

*θ*_{i−1,j+1}*− (2 + γ) θ*_{i,j+1}*+ θ** _{i+1,j+1}* =

*−γθ*

_{i,j}*,*

*where γ*

*= (Δx)*

^{Δ}

^{2}

*/(DΔt).*

*• This equation is unconditionally stable.*

*• Suppose the boundary conditions are given at x = x*_{0}
*and x = x** _{N +1}*.

*• After θ*_{i,j}*has been calculated for i = 1, 2, . . . , N , the*
*values of θ*_{i,j+1}*at time t** _{j+1}* can be computed as the
solution to the following tridiagonal linear system,

### Implicit Methods (continued)

⎡

⎢⎢

⎢⎢

⎢⎢

⎢⎢

⎢⎢

⎢⎢

⎢⎢

⎢⎢

⎢⎢

⎣

*a* 1 0 *· · ·* *· · ·* *· · ·* 0

1 *a* 1 0 *· · ·* *· · ·* 0

0 1 *a* 1 0 *· · ·* 0

.. .

...

...

...

...

... .. . ..

.

... ... ... ... ... .. .

0 *· · ·* *· · ·* 0 1 *a* 1

0 *· · ·* *· · ·* *· · ·* 0 1 *a*

⎤

⎥⎥

⎥⎥

⎥⎥

⎥⎥

⎥⎥

⎥⎥

⎥⎥

⎥⎥

⎥⎥

⎦

⎡

⎢⎢

⎢⎢

⎢⎢

⎢⎢

⎢⎢

⎢⎢

⎢⎢

⎢⎢

⎢⎢

⎢⎢

⎢⎣

*θ1,j+1*
*θ2,j+1*
*θ3,j+1*

..
.
..
.
..
.
*θN,j+1*

⎤

⎥⎥

⎥⎥

⎥⎥

⎥⎥

⎥⎥

⎥⎥

⎥⎥

⎥⎥

⎥⎥

⎥⎥

⎥⎦

=

⎡

⎢⎢

⎢⎢

⎢⎢

⎢⎢

⎢⎢

⎢⎢

⎢⎢

⎢⎢

⎢⎢

⎢⎣

*−γθ1,j − θ0,j+1*

*−γθ2,j*

*−γθ3,j*
..
.
..
.

*−γθN−1,j*

*−γθN,j − θN+1,j+1*

⎤

⎥⎥

⎥⎥

⎥⎥

⎥⎥

⎥⎥

⎥⎥

⎥⎥

⎥⎥

⎥⎥

⎥⎦
*,*

*where a* =^{Δ} *−2 − γ.*

### Implicit Methods (concluded)

*• Tridiagonal systems can be solved in O(N) time and*
*O(N ) space.*

**– Never invert a matrix to solve a tridiagonal system.**

*• The matrix above is nonsingular when γ ≥ 0.*

**– A square matrix is nonsingular if its inverse exists.**

### Crank-Nicolson Method

*• Take the average of explicit method (107) on p. 771 and*
implicit method (108) on p. 778:

*θi,j+1 − θi,j*
*Δt*

= 1

2

*D* *θi+1,j − 2θi,j + θi−1,j*

*(Δx)2* *+ D* *θi+1,j+1 − 2θi,j+1 + θi−1,j+1*
*(Δx)2*

*.*

*• After rearrangement,*

*γθi,j+1 −* *θi+1,j+1 − 2θi,j+1 + θi−1,j+1*

2 *= γθi,j +* *θi+1,j − 2θi,j + θi−1,j*

2 *.*

*• This is an unconditionally stable implicit method with*
excellent rates of convergence.

### Stencil

*t*

_{j}*t*

_{j+1}*x*

_{i}*x*

_{i+1}*x*

_{i+1}### Numerically Solving the Black-Scholes PDE (82) on p.

### 629

*• See text.*

*• Brennan and Schwartz (1978) analyze the stability of*
the implicit method.

### Monte Carlo Simulation

^{a}

*• Monte Carlo simulation is a sampling scheme.*

*• In many important applications within ﬁnance and*
without, Monte Carlo is one of the few feasible tools.

*• When the time evolution of a stochastic process is not*
easy to describe analytically, Monte Carlo may very well
be the only strategy that succeeds consistently.

aA top 10 algorithm (Dongarra & Sullivan, 2000).

### The Big Idea

*• Assume X*_{1}*, X*_{2}*, . . . , X** _{n}* have a joint distribution.

*• θ* *= E[ g(X*^{Δ} _{1}*, X*_{2}*, . . . , X*_{n}*) ] for some function g is*
desired.

*• We generate*

*x*^{(i)}_{1} *, x*^{(i)}_{2} *, . . . , x*^{(i)}_{n}

*, 1* *≤ i ≤ N*

independently with the same joint distribution as
*(X*_{1}*, X*_{2}*, . . . , X** _{n}*).

*• Set*

*Y*_{i}*= g*^{Δ}

*x*^{(i)}_{1} *, x*^{(i)}_{2} *, . . . , x*^{(i)}_{n}

*.*

### The Big Idea (concluded)

*• Y*_{1}*, Y*_{2}*, . . . , Y** _{N}* are independent and identically
distributed random variables.

*• Each Y** _{i}* has the same distribution as

*Y*

*= g(X*

^{Δ}

_{1}

*, X*

_{2}

*, . . . , X*

_{n}*).*

*• Since the average of these N random variables, Y ,*
*satisﬁes E[ Y ] = θ, it can be used to estimate θ.*

*• The strong law of large numbers says that this*
procedure converges almost surely.

*• The number of replications (or independent trials), N, is*
called the sample size.

### Accuracy

*• The Monte Carlo estimate and true value may diﬀer*
owing to two reasons:

1. Sampling variation.

2. The discreteness of the sample paths.^{a}

*• The ﬁrst can be controlled by the number of replications.*

*• The second can be controlled by the number of*
observations along the sample path.

aThis may not be an issue if the ﬁnancial derivative only requires discrete sampling along the time dimension, such as the discrete barrier option.

### Accuracy and Number of Replications

*• The statistical error of the sample mean Y of the*
*random variable Y grows as 1/√*

*N .*
**– Because Var[ Y ] = Var[ Y ]/N .**

*• In fact, this convergence rate is asymptotically optimal.*^{a}

*• So the variance of the estimator Y can be reduced by a*
*factor of 1/N by doing N times as much work.*

*• This is amazing because the same order of convergence*
*holds independently of the dimension n.*

aThe Berry-Esseen theorem.

### Accuracy and Number of Replications (concluded)

*• In contrast, classic numerical integration schemes have*
*an error bound of O(N*^{−c/n}*) for some constant c > 0.*

**– n is the dimension.**

*• The required number of evaluations thus grows*

*exponentially in n to achieve a given level of accuracy.*

**– The curse of dimensionality.**

*• The Monte Carlo method is more eﬃcient than*

alternative procedures for multivariate derivatives when
*n is large.*

### Monte Carlo Option Pricing

*• For the pricing of European options on a*

dividend-paying stock, we may proceed as follows.

*• Assume*

*dS*

*S* *= μ dt + σ dW.*

*• Stock prices S*_{1}*, S*_{2}*, S*_{3}*, . . . at times Δt, 2Δt, 3Δt, . . .*
can be generated via

*S*_{i+1}*= S*_{i}*e*^{(μ−σ}^{2}^{/2) Δt+σ}

*√**Δt ξ**, ξ* *∼ N(0, 1).*

(109)

### Monte Carlo Option Pricing (continued)

*• If we discretize dS/S = μ dt + σ dW directly, we will*
obtain

*S*_{i+1}*= S*_{i}*+ S*_{i}*μ Δt + S*_{i}*σ√*

*Δt ξ.*

*• But this is locally normally distributed, not lognormally,*
hence biased.^{a}

*• In practice, this is not expected to be a major problem*
*as long as Δt is suﬃciently small.*

aContributed by Mr. Tai, Hui-Chin (R97723028) on April 22, 2009.

### Monte Carlo Option Pricing (continued)

*• Non-dividend-paying stock prices in a risk-neutral*

*economy can be generated by setting μ = r and Δt = T .*

1: *C := 0;* *{Accumulated terminal option value.}*

2: **for i = 1, 2, 3, . . . , N do**

3: *P := S* *× e*^{(r−σ}^{2}^{/2) T +σ}^{√}^{T ξ}*, ξ* *∼ N(0, 1);*

4: *C := C + max(P* *− X, 0);*

5: **end for**

6: *return Ce*^{−rT}*/N ;*

### Monte Carlo Option Pricing (concluded)

*• Pricing Asian options is also easy.*

1: *C := 0;*

2: **for i = 1, 2, 3, . . . , N do**

3: *P := S; M := S;*

4: **for j = 1, 2, 3, . . . , n do**

5: *P := P* *× e*^{(r−σ}^{2}*/2)(T /n)+σ**√*

*T /n ξ*;

6: *M := M + P ;*

7: **end for**

8: *C := C + max(M/(n + 1)* *− X, 0);*

9: **end for**

10: *return Ce*^{−rT}*/N ;*

### How about American Options?

*• Standard Monte Carlo simulation is inappropriate for*
American options because of early exercise (why?).

**– Given a sample path S**_{0}*, S*_{1}*, . . . , S** _{n}*, how to decide

*which S*

*is an early-exercise point?*

_{i}**– What is the option price at each S*** _{i}* if the option is
not exercised?

*• It is diﬃcult to determine the early-exercise point based*
on one single path.

*• But Monte Carlo simulation can be modiﬁed to price*
American options with small biases (pp. 847ﬀ).^{a}

aLongstaﬀ & Schwartz (2001).