# Some Sequential Boundary Crossing Results for Geometric Brownian Motion and Their Applications in Financial Engineering

Authors: Tristan Guillaume
Journal: ISRN Applied Mathematics
 (2011)

---

## Abstract

This paper provides new explicit results for some boundary crossing distributions in a multidimensional geometric Brownian motion framework when the boundary is a piecewise constant function of time. Among their various possible applications, they enable accurate and efficient analytical valuation of a large number of option contracts traded in the financial markets belonging to the classes of barrier and look-back options.


---

## Body

## 1. Introduction

The joint law of the maximum (or the minimum) of a real-valued Brownian motion and its endpoint over a finite time interval is a central result in the study of Brownian motion, particularly with regard to the many applications of the theory in finance, medical imaging, robotics, and biology. It can be obtained as a consequence of the “reflection principle,” which derives from the strong Markov property of Brownian motion (Freedman [1]). Application of Girsanov’s theorem easily generalizes this seminal result to the case of a geometric Brownian motion (GBM), a frequently encountered diffusion process that is the building block of financial engineering. Alternatively, the law can be derived by a partial differential equation approach, using Kolmogorov’s equation for the transition density function of a diffusion process. The distribution of the first passage time by a one-dimensional GBM to a one-sided or a two-sided straight boundary then follows. A few cases where the boundary is curved have been handled (Barba Escribá [2]; Salminen [3]; Kunitomo and Ikeda [4]). For a comprehensive source of formulae, one may refer to Borodin and Salminen [5].Various extensions to these results are often needed to solve practical engineering problems. In particular, one may look for joint distributions of the highest or lowest points hit over several time intervals, and one may deal with a multidimensional GBM. There are few known explicit results in these more general settings, partly because of the laborious analytical calculations involved, but also due to numerical obstacles: the rapidly increasing dimension of the boundary crossing problem leads to analytical solutions that are expressed in terms of functions that are hard to compute with accuracy.  As far as one-sided boundaries are concerned, formulae have been published for the joint law of a sequence of maxima or minima of a one-dimensional GBM over several time intervals (Guillaume [6]), as well as for the joint law of the maxima or minima of a two-dimensional GBM over one time interval (Iyengar [7]; He et al. [8]). A formula for the joint law of the exit times of a one-dimensional GBM from two successive two-sided boundaries is also known (Guillaume [9]).This paper focuses on a sequence of two one-sided straight boundaries conditional on two correlated GBMs, while the value of a third correlated GBM is taken into account at the endpoint of the time interval. The state space is thus three-dimensional. The choice of this particular distribution is motivated both by its usefulness in financial engineering applications and by the fact that it leads to tractable analytical solutions that can be computed with great accuracy and efficiency.Section2 presents the two main formulae of this paper. Section 3 deals with applications of these formulae and discusses their numerical implementation.


## 2. Main Formulae

This section contains the two main formulae of the paper. Let{S1(t),t≥0} and {S2(t),t≥0} be two geometric Brownian motions with constant drift coefficients, α1 and α2, respectively, under a probability measure ℙ, constant positive diffusion coefficients σ1 and σ2, respectively, and constant correlation coefficient ρ1·2. In other words, the processes S1(t) and S2(t) evolve in time according to the following dynamics: (2.1)dS1(t)=α1S1(t)dt+σ1S1(t)dB1(t),dS2(t)=α2S2(t)dt+σ2S2(t)dB2(t),
where B1(t) and B2(t) are standard real-valued Brownian motions and d[B1,B2](t)=ρ1·2dt. Let H1, H2, K1, K2, K3 be positive real numbers and P[·](H1,H2,K1,K2,K3,t1,t2) defined as one of the four following cumulative distribution functions, where t2≥t1≥t0=0:   (1)
(2.2)P[up-and-up](H1,H2,K1,K2,K3,t1,t2)≜P(sup⁡0≤t≤t1S1(t)≤H1,S1(t1)≤K1,S2(t1)≤K2,sup⁡t1≤t≤t2S2(t)≤H2,S2(t2)≤K3),(2)
(2.3)P[down-and-down](H1,H2,K1,K2,K3,t1,t2)≜P(inf⁡0≤t≤t1S1(t)≥H1,S1(t1)≥K1,S2(t1)≥K2,inf⁡t1≤t≤t2S2(t)≥H2,S2(t2)≥K3),(3)
(2.4)P[down-and-up](H1,H2,K1,K2,K3,t1,t2)≜P(inf⁡0≤t≤t1S1(t)≥H1,S1(t1)≥K1,S2(t1)≤K2,sup⁡t1≤t≤t2S2(t)≤H2,S2(t2)≤K3),(4)
(2.5)P[up-and-down](H1,H2,K1,K2,K3,t1,t2)≜P(sup⁡0≤t≤t1S1(t)≤H1,S1(t1)≤K1,S2(t1)≥K2,inf⁡t1≤t≤t2S2(t)≥H2,S2(t2)≥K3).Let us introduce the following notations:(2.6)k1=ln⁡(K1S1(t0)),k2=ln⁡(K2S2(t0)),k3=ln⁡(K3S2(t0)),h1=ln⁡(H1S1(t0)),h2=ln⁡(H2S2(t0)),μ1=α1-σ122,μ2=α2-σ222.
The next proposition provides an exact formula for the four above-mentioned cumulative distribution functions.Proposition 2.1.
LetN3[·,·,·;θ1·2,θ1·3,θ2·3] denote the joint trivariate cumulative distribution function of three standard normal random variables X1, X2, X3, where θa·b is the correlation coefficient between Xa and Xb, (a,b)∈{1,2,3}.
Then, for the up-and-up and the down-and-down distributions, one has(2.7)P[⋅](H1,H2,K1,K2,K3,t1,t2)=N3[λ(k1-μ1t1σ1t1),λ(k2-μ2t1σ2t1),λ(k3-μ2t2σ2t2);ρ1⋅2,t1t2ρ1⋅2,t1t2]-exp⁡(2μ2h2σ22)N3[λ(k1-μ1t1σ1t1)+ρ1⋅22μ2t1σ2,λ(k2+μ2t1σ2t1),λ(k3-2h2-μ2t2σ2t2);ρ1⋅2,-ρ1⋅2t1t2,-t1t2]-exp⁡(2μ1h1σ12)×N3[λ(k1-2h1-μ1t1σ1t1),λ(k2-μ2t1σ2t1-ρ1⋅22h1σ1t1),λ(k3-μ2t2σ2t2-ρ1⋅22h1σ1t2);ρ1⋅2,ρ1⋅2t1t2,t1t2]+exp⁡((2μ1σ12-4μ2ρ1⋅2σ1σ2)h1+2μ2h2σ22)×N3[λ(k1-2h1-μ1t1σ1t1)+ρ1⋅22μ2t1σ2,λ(k2+μ2t1σ2t1-ρ1⋅22h1σ1t1),λ(k3-2h2-μ2t2σ2t2+ρ1⋅22h1σ1t2);ρ1⋅2,-ρ1⋅2t1t2,-t1t2],
where
(2.8)λ=1ifP[⋅](H1,H2,K1,K2,K3,t1,t2)=P[up-and-up](H1,H2,K1,K2,K3,t1,t2),λ=-1ifP[⋅](H1,H2,K1,K2,K3,t1,t2)=P[down-and-down](H1,H2,K1,K2,K3,t1,t2).
And, for the up-and-down and the down-and-up distributions, one has
(2.9)P[⋅](H1,H2,K1,K2,K3,t1,t2)=N3[λ(k1-μ1t1σ1t1),λ(-k2+μ2t1σ2t1),λ(-k3+μ2t2σ2t2);-ρ1⋅2,-t1t2ρ1⋅2,t1t2]-exp⁡(2μ2h2σ22)×N3[λ(k1-μ1t1σ1t1)+ρ1⋅22μ2t1σ2,λ(-k2-μ2t1σ2t1),λ(-k3+2h2+μ2t2σ2t2);-ρ1⋅2,t1t2ρ1⋅2,-t1t2]-exp⁡(2μ1h1σ12)×N3[λ(k1-2h1-μ1t1σ1t1),λ(-k2+μ2t1σ2t1+ρ1⋅22h1σ1t1),λ(-k3+μ2t2σ2t2+ρ1⋅22h1σ1t2);-ρ1⋅2,-t1t2ρ1⋅2,t1t2]+exp⁡((2μ1σ12-4μ2ρ1⋅2σ1σ2)h1+2μ2h2σ22)×N3[λ(k1-2h1-μ1t1σ1t1)+ρ1⋅22μ2t1σ2,λ(-k2-μ2t1σ2t1+ρ1⋅22h1σ1t1),λ(-k3+2h2+μ2t2σ2t2-ρ1⋅22h1σ1t2);-ρ1⋅2,t1t2ρ1⋅2,-t1t2],
where
(2.10)λ=1ifP[⋅](H1,H2,K1,K2,K3,t1,t2)=P[up-and-down](H1,H2,K1,K2,K3,t1,t2),λ=-1ifP[⋅](H1,H2,K1,K2,K3,t1,t2)=P[down-and-up](H1,H2,K1,K2,K3,t1,t2).Corollary of Proposition2.1
The four following joint cumulative distribution functions, that will be useful in Section3, are deduced from Proposition 2.1:(1)
(2.11)P(sup⁡0≤t≤t1S1(t)≤H1,S1(t1)≤K1,S2(t1)≤K2,sup⁡t1≤t≤t2S2(t)≤H2,S2(t2)≥K3)=P[up-and-up](H1,H2,K1,K2,∞)-P[up-and-up](H1,H2,K1,K2,K3),(2)
(2.12)P(inf⁡0≤t≤t1S1(t)≥H1,S1(t1)≥K1,S2(t1)≥K2,inf⁡t1≤t≤t2S2(t)≥H2,S2(t2)≤K3)=P[down-and-down](H1,H2,K1,K2,-∞)-P[down-and-down](H1,H2,K1,K2,K3),(3)
(2.13)P(inf⁡0≤t≤t1S1(t)≥H1,S1(t1)≥K1,S2(t1)≤K2,sup⁡t1≤t≤t2S2(t)≤H2,S2(t2)≥K3)=P[down-and-up](H1,H2,K1,K2,∞)-P[down-and-up](H1,H2,K1,K2,K3),(4)
(2.14)P(sup⁡0≤t≤t1S1(t)≤H1,S1(t1)≤K1,S2(t1)≥K2,inf⁡t1≤t≤t2S2(t)≥H2,S2(t2)≤K3)=P[up-and-down](H1,H2,K1,K2,-∞)-P[up-and-down](H1,H2,K1,K2,K3).Proof of Proposition2.1 is provided in the Appendices A, B, and C.The numerical implementation of Proposition2.1 and its corollary is easy using Genz’s algorithm for the computation of trivariate normal cumulative distribution functions (Genz [10]).In the next Proposition, we introduce a third correlated geometric Brownian motion that will serve as the endpoint of the joint distribution, and we show that this can still be analytically valued. Let{S3(t),t≥0} be a geometric Brownian motion driven by the following dynamics under the initial probability measure ℙ:(2.15)dS3(t)=α3S3(t)dt+σ3S3(t)dB3(t),
where α3 and σ3 are real constants (σ3>0), B3(t) is a real-valued standard Brownian motion and correlations are as follows:(2.16)d[B1,B3](t)=ρ1⋅3dt,d[B2,B3](t)=ρ2⋅3dt.
Let K4 be a positive real number and P[·](H1,H2,K1,K2,K3,K4,t1,t2,t3) defined as one of the four following joint cumulative distribution functions, where t3≥t2:   (1)
(2.17)P[up-and-up](H1,H2,K1,K2,K3,K4,t1,t2,t3)≜P(sup⁡0≤t≤t1S1(t)≤H1,S1(t1)≤K1,S2(t1)≤K2,sup⁡t1≤t≤t2S2(t)≤H2,S2(t2)≤K3,S3(t3)≤K4),(2)
(2.18)P[down-and-down](H1,H2,K1,K2,K3,K4,t1,t2,t3)≜P(inf⁡0≤t≤t1S1(t)≥H1,S1(t1)≥K1,S2(t1)≥K2,inf⁡t1≤t≤t2S2(t)≥H2,S2(t2)≥K3,S3(t3)≥K4),(3)
(2.19)P[down-and-up](H1,H2,K1,K2,K3,K4,t1,t2,t3)≜P(inf⁡0≤t≤t1S1(t)≥H1,S1(t1)≥K1,S2(t1)≤K2,sup⁡t1≤t≤t2S2(t)≤H2,S2(t2)≤K3,S3(t3)≤K4),(4)
(2.20)P[up-and-down](H1,H2,K1,K2,K3,K4,t1,t2,t3)≜P(sup⁡0≤t≤t1S1(t)≤H1,S1(t1)≤K1,S2(t1)≥K2,inf⁡t1≤t≤t2S2(t)≥H2,S2(t2)≥K3,S3(t3)≥K4).Let us introduce the following new notations:(2.21)k4=ln⁡(K4S3(t0)),μ3=α3-σ322.
The next proposition provides an exact formula for the four above-mentioned joint cumulative distribution functions.Proposition 2.2.
Let the real functionΦ[b1,b2,b3,b4;θ1·2,θ1·3,θ1·4,θ2·3,θ3·4], where {b1,b2,b3,b4}∈ℝ4  and each of the real numbers θ1·2,θ1·3,θ1·4,θ2·3,θ3·4 is included in ]-1,1[, be defined as follows:
(2.22)Φ[b1,b2,b3,b4;θ1⋅2,θ1⋅3,θ1⋅4,θ2⋅3,θ3⋅4]=∫x1=-∞b1∫x2=-∞b2∫x3=-∞b31ϕ2∣1ϕ3∣28π3exp⁡(-x222-(x1-θ1⋅2x2)22ϕ2∣12-(x3-θ2⋅3x2)22ϕ3∣22)×N[b4-θ1⋅4x1-θ3⋅4∣1((x3-θ1⋅3x1)/1-θ1⋅32)ϕ4∣1⋅3]dx3dx2dx1,
where N[·] is the standard normal cumulative distribution function and the following definitions apply:
(2.23)ϕ2∣1=1-θ1⋅22,ϕ3∣2=1-θ2⋅32,θ3⋅4∣1=θ3⋅4-θ1⋅3θ1⋅4ϕ3∣1,ϕ4∣1⋅3=1-θ1⋅42-θ3⋅4∣12.
Then, for the up-and-up and the down-and-down distributions, one has
(2.24)P[⋅](H1,H2,K1,K2,K3,K4,t1,t2,t3)=Φ[λ(k1-μ1t1σ1t1),λ(k2-μ2t1σ2t1),λ(k3-μ2t2σ2t2),λ(k4-μ3t3σ3t3);θ1⋅2=ρ1⋅2,θ1⋅3=t1t2ρ1⋅2,θ1⋅4=t1t3ρ1⋅3,θ2⋅3=t1t2,θ3⋅4=t2t3ρ2⋅3]-exp⁡(2μ2h2σ22)×Φ[λ(k1-μ1t1σ1t1)+ρ1⋅22μ2t1σ2t1,λ(k2+μ2t1σ2t1),λ(k3-2h2-μ2t2σ2t2),λ(k4-μ3t3σ3t3-θ3⋅4-θ1⋅3θ1⋅41-θ1⋅32(2h2σ2t2+θ1⋅3ρ1⋅22μ2t1σ2t1)+θ1⋅4ρ1⋅22μ2t1σ2t1);θ1⋅2=ρ1⋅2,θ1⋅3=t1t2ρ1⋅2,θ1⋅4=t1t3ρ1⋅3,θ2⋅3=-t1t2,θ3⋅4=t2t3ρ2⋅3]-exp⁡(2μ1h1σ12)×Φ[λ(k1-2h1-μ1t1σ1t1),λ(k2-μ2t1σ2t1-ρ1⋅22h1σ1t1),λ(k3-μ2t2σ2t2-ρ1⋅22h1σ1t2),λ(k4-μ3t3σ3t3-θ1⋅42h1σ1t1-θ3⋅4-θ1⋅3θ1⋅41-θ1⋅32(ρ1⋅22h1σ1t2-θ1⋅32h1σ1t1));θ1⋅2=ρ1⋅2,θ1⋅3=t1t2ρ1⋅2,θ1⋅4=t1t3ρ1⋅3,θ2⋅3=t1t2,θ3⋅4=t2t3ρ2⋅3]+exp⁡((2μ1σ12-4μ2ρ1⋅2σ1σ2)h1+2μ2h2σ22)×Φ[λ(k1-2h1-μ1t1σ1t1)+ρ1⋅22μ2t1σ2t1,λ(k2+μ2t1σ2t1-ρ1⋅22h1σ1t1),λ(k3-2h2-μ2t2σ2t2+ρ1⋅22h1σ1t2),λ(k4-μ3t3σ3t3+θ1⋅4(ρ1⋅22μ2t1σ2-2h1σ1t1)-θ3⋅4-θ1⋅3θ1⋅41-θ1⋅32(2h2σ2t2-ρ1⋅22h1σ1t2+θ1⋅3(ρ1⋅22μ2t1σ2-2h1σ1t1)));θ1⋅2=ρ1⋅2,θ1⋅3=t1t2ρ1⋅2,θ1⋅4=t1t3ρ1⋅3,θ2⋅3=-t1t2,θ3⋅4=t2t3ρ2⋅3],
where
(2.25)λ=1ifP[⋅](H1,H2,K1,K2,K3,K4,t1,t2,t3)=P[up-and-up](H1,H2,K1,K2,K3,K4,t1,t2,t3),λ=-1ifP[⋅](H1,H2,K1,K2,K3,K4,t1,t2,t3)=P[down-and-down](H1,H2,K1,K2,K3,K4,t1,t2,t3).
And, for the up-and-down and the down-and-up distributions, one has
(2.26)P[.](H1,H2,K1,K2,K3,K4,t1,t2,t3)=Φ[λ(-k1+μ1t1σ1t1),λ(k2-μ2t1σ2t1),λ(k3-μ2t2σ2t2),λ(k4-μ3t3σ3t3);θ1⋅2=-ρ1⋅2,θ1⋅3=-t1t2ρ1⋅2,θ1⋅4=-t1t3ρ1⋅3,θ2⋅3=t1t2,θ3⋅4=t2t3ρ2⋅3]-exp⁡(2μ2h2σ22)×Φ[λ(-k1+μ1t1σ1t1)+ρ1⋅22μ2t1σ2t1,λ(k2+μ2t1σ2t1),λ(k3-2h2-μ2t2σ2t2),λ(k4-μ3t3σ3t3-θ3⋅4-θ1⋅3θ1⋅41-θ1⋅32(2h2σ2t2+θ1⋅3ρ1⋅22μ2t1σ2)+θ1⋅4ρ1⋅22μ2t1σ2);θ1⋅2=-ρ1⋅2,θ1⋅3=t1t2ρ1⋅2,θ1⋅4=-t1t3ρ1⋅3,θ2⋅3=-t1t2,θ3⋅4=t2t3ρ2⋅3]-exp⁡(2μ1h1σ12)×Φ[λ(-k1+2h1+μ1t1σ1t1),λ(k2-μ2t1σ2t1-ρ1⋅22h1σ1t1),λ(k3-μ2t2σ2t2-ρ1⋅22h1σ1t2),λ(k4-μ3t3σ3t3-θ1⋅42h1σ1t1-θ3⋅4-θ1⋅3θ1⋅41-θ1⋅32(ρ1⋅22h1σ1t2-θ1⋅32h1σ1t1));θ1⋅2=-ρ1⋅2,θ1⋅3=-t1t2ρ1⋅2,θ1⋅4=-t1t3ρ1⋅3,θ2⋅3=t1t2,θ3⋅4=t2t3ρ2⋅3]+exp⁡((2μ1σ12-4μ2ρ1⋅2σ1σ2)h1+2μ2h2σ22)×Φ[λ(-k1+2h1+μ1t1σ1t1)+ρ1⋅22μ2t1σ2,λ(k2+μ2t1σ2t1-ρ1⋅22h1σ1t1),λ(k3-2h2-μ2t2σ2t2+ρ1⋅22h1σ1t2),λ(k4-μ3t3σ3t3+θ1⋅4(ρ1⋅22μ2t1σ2-2h1σ1t1)-θ3⋅4-θ1⋅3θ1⋅41-θ1⋅32(2h2σ2t2-ρ1⋅22h1σ1t2+θ1⋅3(ρ1⋅22μ2t1σ2-2h1σ1t1)));θ1⋅2=-ρ1⋅2,θ1⋅3=t1t2ρ1⋅2,θ1⋅4=-t1t3ρ1⋅3,θ2⋅3=-t1t2,θ3⋅4=t2t3ρ2⋅3],
where
(2.27)λ=1ifP[⋅](H1,H2,K1,K2,K3,K4,t1,t2,t3)=P[down-and-up](H1,H2,K1,K2,K3,K4,t1,t2,t3),λ=-1ifP[⋅](H1,H2,K1,K2,K3,K4,t1,t2,t3)=P[up-and-down](H1,H2,K1,K2,K3,K4,t1,t2,t3).Corollary of Proposition2.2
The four following cumulative distribution functions, that will be useful in Section3, are deduced from Proposition 2.2:(1)
(2.28)P(sup⁡0≤t≤t1S1(t)≤H1,S1(t1)≤K1,S2(t1)≤K2,sup⁡t1≤t≤t2S2(t)≤H2,S2(t2)≤K3,S3(t3)≥K4)=P[up-and-up](H1,H2,K1,K2,K3,∞)-P[up-and-up](H1,H2,K1,K2,K3,K4),(2)
(2.29)P(inf⁡0≤t≤t1S1(t)≥H1,S1(t1)≥K1,S2(t1)≥K2,inf⁡t1≤t≤t2S2(t)≥H2,S2(t2)≥K3,S3(t3)≤K4)=P[down-and-down](H1,H2,K1,K2,K3,-∞)-P[down-and-down](H1,H2,K1,K2,K3,K4),(3)
(2.30)P(inf⁡0≤t≤t1S1(t)≥H1,S1(t1)≥K1,S2(t1)≤K2,sup⁡t1≤t≤t2S2(t)≤H2,S2(t2)≤K3,S3(t3)≥K4)=P[down-and-up](H1,H2,K1,K2,K3,∞)-P[down-and-up](H1,H2,K1,K2,K3,K4),(4)
(2.31)P(sup⁡0≤t≤t1S1(t)≤H1,S1(t1)≤K1,S2(t1)≥K2,inf⁡t1≤t≤t2S2(t)≥H2,S2(t2)≥K3,S3(t3)≤K4)=P[up-and-down](H1,H2,K1,K2,K3,-∞)-P[up-and-down](H1,H2,K1,K2,K3,K4).Proof of Proposition2.2 is provided in the Appendices A, B, and C.Practical use of Proposition2.2 and its corollary in engineering applications depend on the accuracy and the efficiency with which the function Φ can be numerically implemented. This question is dealt with in Section 3.


## 3. Applications and Numerical Implementation

This Section deals with applications of the results provided in Section1 in financial engineering. Indeed, Propositions 2.1 and 2.2 can be used as the building blocks for the valuation and the risk management of a large number of option contracts. The main class of instruments under consideration will be barrier options. In their standard form, the latter are contracts whereby the holder is entitled (but not obligated) to buy (call option) or sell (put option) an asset at a prespecified future date (the option expiry) at a prespecified price (the strike price), on condition that the asset price has not (knock-out type) or has (knock-in type) crossed a specific upward level (called the up-barrier) or downward level (called the down-barrier) at any given time from valuation date to expiry. In standard option pricing theory, asset prices are modelled as geometric Brownian motions, so that boundary crossing distributions such as those provided in Section 1 are of immediate use in financial engineering.It must be emphasized that the market for barrier options is huge. They are the most actively traded class of nonstandard options (usually referred to as exotic options). Many variations on the standard barrier option payoff have been designed to match investors’ demand more closely. One of them is the partial-time barrier option, which specifies that the barrier level is monitored during only a fraction of the option lifespan. The basic contracts were studied by Heynen and Kat [11], and more general contracts were valued by Armstrong [12] and Guillaume [13]. Another variation is the step barrier option (Guillaume [6, 9]), whereby the barrier evolves as a step function of time. It is also common to encounter outside barrier options, whose payoff is a function of two asset prices: one of them is compared with the strike price at expiry to determine the moneyness of the option, while the other one is monitored up until expiry to check whether the barrier level has been crossed. Outside barrier options were originally valued by Heynen and Kat [14], and they were further studied by Kwok et al. [15] and Wong and Kwok [16].These different features: partial time barrier, step barrier, and outside barrier, often combine to allow for increased flexibility. But then, performing analytical valuation becomes more and more involved and practitioners have to turn to numerical methods that are slow and relatively inaccurate. This is when the results provided in Section1 become valuable. Indeed, they enable to price a large number of barrier option contracts sharing the three above-mentioned innovative features (partial time, step, and outside barrier), that is, options based on two or even three assets that may knock out depending on whether the underlying assets move within a sequence of prespecified ranges of prices over all or part of the option lifespan. More precisely, Proposition 2.1 and its corollary are the building blocks for the valuation of sequential two-asset knock-out calls and puts whose four pay-off structures are defined as follows:(1)
(K-St2(2))+𝕀{sup0≤t≤t1St(1)<H1,supt1≤t≤t2St(2)<H2} for a two-asset up-and-up knock-out put (substitute (St2(2)-K)+ for (K-St2(2))+ for an up-and-up knock-out call),(2)
±(St2(2)-K)+𝕀{inf0≤t≤t1St(1)>H1,inft1≤t≤t2St(2)>H2} for a two-asset down-and-down knock-out call or put,(3)
±(St2(2)-K)+𝕀{sup0≤t≤t1St(1)<H1,inft1≤t≤t2St(2)>H2} for a two-asset up-and-down knock-out call or put,(4)
±(K-St2(2))+𝕀{inf0≤t≤t1St(1)>H1,supt1≤t≤t2St(2)<H2} for a two-asset down-and-up knock-out put or call,where:(i)
𝕀{} is the indicator function,(ii)
{S1(t),t≥0} and {S2(t),t≥0} are the price processes of two risky assets,(iii)
K is the strike price and t2 is the option expiry,(iv)
H1 is a knock-out barrier monitored with respect to {S1(t),t≥0} over the time interval  [t0,t1],(v)
H2 is a knock-out barrier monitored with respect to {S2(t),t≥0} over the time interval  [t1,t2].The no-arbitrage price of these eight different option contracts is now provided by the following proposition.Proposition 3.1.
Let two geometric Brownian motions{S1(t),t≥0} and {S2(t),t≥0} be the price processes of two risky assets with constant volatilities  σ1 and σ2, respectively, and constant correlation coefficient ρ1·2. Let δ1 and δ2 be two constant pay-out rates on assets S1  and S2, respectively, and r be the constant risk-free interest rate. Let υ1=r-δ1+σ1σ2ρ1·2  and υ2=r-δ2+σ22.
The value, at timet0=0, of a two-asset sequential knock-out option with strike price K, with knock-out barriers H1 and H2 monitored with respect to S1 over the time interval [t0,t1]  and S2 over the time interval [t1,t2], respectively, and with maturity t2, t2≥t1≥t0, is given by
(3.1)V(S1(t0),S2(t0),t1,t2,K,H1,H2)=±[Kexp⁡(-rt2)P1-S2(t0)exp⁡(-δ2t2)P2],
where the following specifications hold.(1)
If the option is an up-and-up put,±=1, and P1 is obtained by substituting α1=r-δ1 and α2=r-δ2 in P[up-and-up](H1,H2,H1,H2,K,t1,t2) in Proposition 2.1, and P2 is obtained by substituting α1=υ1 and α2=υ2 in P[up-and-up](H1,H2,H1,H2,K,t1,t2) in Proposition 2.1.(2)
If the option is a down-and-down call,±=-1, and P1 is obtained by substituting α1=r-δ1 and α2=r-δ2 in P[down-and-down](H1,H2,H1,H2,K,t1,t2) in Proposition 2.1, and P2 is obtained by substituting α1=υ1 and α2=υ2 in P[down-and-down](H1,H2,H1,H2,K,t1,t2) in Proposition 2.1.(3)
If the option is an up-and-down call,±=-1, and P1 is obtained by substituting α1=r-δ1 and α2=r-δ2 in P[up-and-down](H1,H2,H1,H2,K,t1,t2) in Proposition 2.1 and P2 is obtained by substituting α1=υ1 and α2=υ2 in P[up-and-down](H1,H2,H1,H2,K,t1,t2) in Proposition 2.1.(4)
If the option is a down-and-up put,±=1, and P1 is obtained by substituting α1=r-δ1 and α2=r-δ2 in P[down-and-up](H1,H2,H1,H2,K,t1,t2) in Proposition 2.1 and P2 is obtained by substituting α1=υ1 and α2=υ2 in P[down-and-up](H1,H2,H1,H2,K,t1,t2) in Proposition 2.1.By using Corollary of Proposition 2.1, one can value up-and-up and down-and-up call options, as well as down-and-down and up-and-down put options.A sketch of proof of Proposition3.1 is provided in the Appendices A, B, and C.Just as Proposition2.1 enables to value two-asset knock-out calls and puts, Proposition 2.2 can be used to value three-asset knock-out calls and puts. More specifically, if {S3(t),t≥0} is the price process of a third risky asset and t3 is the new option expiry, with t3≥t2≥t1≥t0, then the four following pay-off structures can be analytically tackled:(1)
±(K-St3(3))+𝕀{sup0≤t≤t1St(1)<H1,supt1≤t≤t2St(2)<H2} for a three-asset up-and-up knock-out put or call,(2)
±(St3(3)-K)+𝕀{inf0≤t≤t1St(1)>H1,inft1≤t≤t2St(2)>H2} for a three-asset down-and-down knock-out call or put,(3)
±(St3(3)-K)+𝕀{sup0≤t≤t1St(1)<H1,inft1≤t≤t2St(2)>H2} for a three-asset up-and-down knock-out call or put,(4)
±(K-St3(3))+𝕀{inf0≤t≤t1St(1)>H1,supt1≤t≤t2St(2)<H2} for a three-asset down-and-up knock-out put or call.The no-arbitrage prices of the eight types of option under consideration are now provided by the following proposition.Proposition 3.2.
Let three geometric Brownian motions{S1(t),t≥0}, {S2(t),t≥0}, and {S3(t),t≥0} be the price processes of three risky assets with constant volatilities σ1, σ2, and σ3, respectively, and constant pairwise correlation coefficients ρa·b, (a,b)∈{1,2,3}, at any fixed time t≥0. Let δ1, δ2, and δ3 be three constant pay-out rates on assets S1, S2, and S3, respectively, and r be the constant risk-free interest rate. Let υ1=r-δ1+σ1σ3ρ1·3, υ2=r-δ2+σ2σ3ρ2·3, and υ3=r-δ3+σ32.
The value, at timet0=0, of a three-asset sequential knock-out option with strike price K, with knock-out barriers H1 and H2 monitored with respect to S1 over the time interval [t0,t1]  and S2 over the time interval [t1,t2], respectively, with maturity t3, t3≥t2≥t1≥t0, is given by
(3.2)V(S1(t0),S2(t0),S3(t0),t1,t2,t3,K,H1,H2)=±[Kexp⁡(-rt3)P1-S3(t0)exp⁡(-δ3t3)P2],
where the following specifications hold.(1)
If the option is an up-and-up put,±=1, and P1 is obtained by substituting α1=r-δ1, α2=r-δ2 and α3=r-δ3 in P[up-and-up](H1,H2,H1,H2,H2,K,t1,t2,t3) in Proposition 2.2, and P2 is obtained by substituting α1=υ1, α2=υ2 and α3=υ3 in P[up-and-up](H1,H2,H1,H2,H2,K,t1,t2,t3) in Proposition 2.2.(2)
If the option is a down-and-down call,±=-1, and P1 is obtained by substituting α1=r-δ1, α2=r-δ2 and α3=r-δ3 in P[down-and-down](H1,H2,H1,H2,H2,K,t1,t2,t3) in Proposition 2.2, and P2 is obtained by substituting α1=υ1, α2=υ2 and α3=υ3 in P[down-and-down](H1,H2,H1,H2,H2,K,t1,t2,t3) in Proposition 2.2.(3)
If the option is an up-and-down call,±=-1, and P1 is obtained by substituting α1=r-δ1, α2=r-δ2 and α3=r-δ3 in P[up-and-down](H1,H2,H1,H2,H2,K,t1,t2,t3) in Proposition 2.2 and P2 is obtained by substituting α1=υ1, α2=υ2 and α3=υ3 in P[up-and-down](H1,H2,H1,H2,H2,K,t1,t2,t3) in Proposition 2.2.(4)
If the option is a down-and-up put,±=1, and P1 is obtained by substituting α1=r-δ1, α2=r-δ2 and α3=r-δ3 in P[down-and-up](H1,H2,H1,H2,H2,K,t1,t2,t3) in Proposition 2.2 and P2 is obtained by substituting α1=υ1, α2=υ2 and α3=υ3 in P[down-and-up](H1,H2,H1,H2,H2,K,t1,t2,t3) in Proposition 2.2.By using Corollary of Proposition 2.2, one can value up-and-up and down-and-up call options, as well as down-and-down and up-and-down put options.A sketch of proof of Proposition3.2 is provided in the Appendices A, B, and C.It is easy to notice that Proposition3.2 nests Proposition 3.1 under the following three conditions:(1)
υ1=r-δ1+σ1σ2ρ1·2  and υ2=r-δ2+σ22,(2)
takeP[·](H1,H2,H1,H2,K,U,t1,t2,t3) with U→∞ if you are dealing with a put and with U→-∞ if you are dealing with a call,(3)
compute±[Kexp(-rt2)P1-S2(t0)exp(-δ2t2)P2].It is useful, though, to state Proposition3.1, as it uses more standard functions that can be numerically computed with known high accuracy and thus can serve as a benchmark to test the reliability of the numerical implementation of Proposition 3.2.It can also be noticed that Proposition3.2 enables to value sequential two-asset knock-out options having the property that the moneyness of the option at expiry is tested with respect to asset S1 instead of asset S2 as in Proposition 3.1. This is achieved by implementing Proposition 3.2 with the following specifications:      (1)
υ1=r-δ1+σ12  and υ2=r-δ2+σ1σ2ρ1·2,(2)
t3→t2, σ3→σ1, ρ1·3→1 and ρ2·3→ρ1·2,(3)
compute±[Kexp(-rt2)P1-S1(t0)exp(-δ1t1)P2].Similarly, Proposition3.2 can be used to value sequential two-asset knock-out options with expiry t3, so that, overall, this formula alone suffices to compute 40 different types of sequential multiasset knock-out options.Proposition3.2 implies the numerical implementation of Proposition 2.2. The function Φ in Proposition 2.2 involves three-dimensional numerical integration. Given the fact that the integrand is very smooth, a classical Gauss-Legendre quadrature can be used, along with a simple adaptive rule increasing the number of points in case the required accuracy is not reached after one iteration, starting with an 8-point quadrature. To test the efficiency and the accuracy of this numerical implementation, four series of tests have been implemented. In the first three categories of tests, the option parameters are specified so that the option can be valued by a known analytical formula involving numerical integration of lower dimension. More specifically, tests of category 1 reduce a three-asset up-and-up knock-out put to a two-asset up-and-up knock-out put by letting K4→∞ in P[up-and-up](H1,H2,K1,K2,K3,K4,t1,t2,t3); the benchmark is then Proposition 3.1, that involves the computation of trivariate normal cumulative distribution functions. The latter can be obtained with high precision by Genz’s algorithm (Genz [10]), as mentioned earlier.Tests of category 2 reduce the three-asset up-and-up knock-out put to an outside knock-out put by lettingK1=H1, H2→∞, K3→∞, K4→∞, and t3→t2→t1 in P[up-and-up](H1,H2,K1,K2,K3,K4,t1,t2,t3); the benchmark, then, is the formula for an outside barrier option as given by Heynen and Kat [14], that involves the computation of bivariate normal cumulative distribution functions. Again, the latter can be obtained with high precision using Genz [10].Tests of category 3 reduce the three-asset up-and-up knock-out put to a single asset up-and-out put by lettingH2→∞, K2=K1, K3→∞, K4→∞, t3→t2→t1, σ2=σ1 and ρ1·2→1 in P[up-and-up](H1,H2,K1,K2,K3,K4,t1,t2,t3); the benchmark, then, is the formula for a basic barrier option as given, for example, by Haug [17], that involves the computation of univariate normal cumulative distribution functions.In the last set of tests, the three-asset up-and-up knock-out put cannot be valued as a simpler contract, so that no reduction of the dimension of numerical integration is possible. Our goal, then, is to find a robust numerical valuation method and study its convergence pattern with the results obtained by implementing Proposition3.2. Neither using trees nor numerically solving a partial differential equation stands as the most suitable methods in this case, in view of the double difficulty of having three correlated diffusion processes and path-dependent features in the option payoff. On another hand, crude Monte Carlo simulation would be very poor both in terms of accuracy and efficiency, due to the required time discretization. That is why the most suitable numerical method in this case is probably a Brownian Bridge Monte Carlo approximation (Glasserman [18]), as it involves only 6 random draws of the asset prices for each simulation and avoids a discretization bias. This is indeed the method commonly used by market practitioners when they are faced with a pricing problem involving both multiasset and path-dependent features.For each category of tests, a total of 10,000 different option prices have been computed with randomly drawn volatility, correlation, and expiry parameters. Let us now sum up the main findings of these numerical experiments. First, the average computational time when applying Proposition3.2 on an ordinary personal computer (P7350 2 GHz) is 0.6 seconds, which is fast. The longest recorded computational time was 1.4 second. Speed could obviously be greatly increased by involving more hardware resources. In absolute terms, Proposition 3.2 can thus be considered as an efficient way of computing option prices. In relative terms, Proposition 3.2 is even extremely efficient. Indeed, 10,000,000 simulations are required to obtain a modest 10-4 digits convergence of the Brownian Bridge Monte Carlo approximation with Proposition 3.2, which amounts to a computational time of over 4 minutes! Such slowness becomes a serious issue in real-time trading or when large portfolios of options are valued. As far as accuracy is concerned, the results obtained using Proposition 3.2 always matched to at least 10-8 digits the analytical benchmarks in almost all the numerical experiments conducted. The very rare exceptions were observed when the absolute value of one of the correlation coefficients was above 99.7%. It must be noticed that in this couple of cases, the correlation matrix was found to be almost singular. Besides, risky assets with such levels of correlation are never encountered in the markets, so that it does not seem worthwhile to us to search for more complicated integration rules that would be better suited to these negligible subregions of integration. For all practical purposes, the numerical accuracy of our simple implementation of Proposition 3.2 can thus be deemed as very satisfactory. The results did not vary significantly between tests of categories 1, 2, and 3, which means that an increase in the number of option parameters that are given limit values does not entail instability in the quadrature rule.To close this section, it is worth mentioning that the results provided in Section1 can also apply to other popular derivatives contracts. Among them is the look-back option, which allows investors to sell at the highest and buy at the lowest over a time interval, thus optimizing market exit and market entry on a given asset. In its basic form, a fixed strike look-back call pays off (sup0≤t≤TS(t)-K)+ at expiry T. Many variations on this standard structure have been designed. One of them is the double look-back call, written on two assets S1  and S2, which provides investors with the following payoff at expiry T: (a1sup0≤t≤TS1(t)-a2sup0≤t≤TS2(t)-K)+, where a1 and a2 are two positive weights. The valuation of this contract has been studied by He et al. [8]. A sequential double look-back call is a variation on the double look-back call. Its payoff at expiry t2 writes (a1sup0≤t≤t1S1(t)-a2supt1≤t≤t2S2(t)-K)+. Drawing on Proposition 2.1, the no-arbitrage price of this option is equal to(3.3)exp⁡(-rt2)∫h1=0∞∫h2=-∞∞(a1S1(t0)exp⁡(h1)-a2S2(t0)exp⁡(h2)-K)+∂2P[up-and-up](H1,H2,H1,H2,H2)∂h1∂h2dh2dh1,
where α1=r-δ1 and α2=r-δ2 in Proposition 2.1.The partial derivative in (3.3) can be determined by analytical differentiation of P[up-and-up](H1,H2,H1,H2,H2), but this is relatively cumbersome. Given the simplicity of the computation of Proposition 2.1, it is faster and easier to compute those partial derivatives numerically. One can notice that the valuation of sequential double look-back options is only semianalytical since it still requires numerical double integration. For a given degree of the quadrature rule, it will not be as accurate as the valuation of two-asset and three-asset sequential knock-out options because of the inaccuracies arising from the numerical differentiation of Proposition 2.1.


## 4. Conclusion

This paper provides new explicit solutions for some sequential boundary crossing problems in a multidimensional geometric Brownian motion framework. Among their various possible applications, they form the basis for accurate and fast valuation of a large number of option contracts traded in the markets. The methods presented here could be used to solve higher-dimensional problems, especially the introduction of a third boundary monitored with respect to the third process in the model. However, the resulting formulae would become quite cumbersome and, more importantly, the quality of the numerical integration involved should have to be compared with a Brownian Bridge Monte Carlo approximation. Another possible extension would be to study the analytical tractability of introducing curved boundaries, but it should be pointed out that this form of boundary is rarely encountered in the options markets.



---

*Source: 120253-2011-07-11.xml*