A
numerical method to price European derivatives
based
on the one factor LIBOR Market Model of interest rates
Luca Vincenzo Ballestra, Graziella Pacelli, Francesco
Zirilli
1.
Abstract
5.
References
7.
Contact us
We consider the problem of pricing European
interest rate derivatives based on the LIBOR Market Model (LMM) with one
driving factor. We derive a closed-form approximation of the transition
probability density functions associated to the stochastic dynamical systems
that describe the behaviour of the forward LIBOR interest rates in the LMM.
These approximate formulae are based on a truncated power series expansion of
the solutions of the Fokker-Planck equations associated to the LMM. The
approximate probability density functions obtained are used to price European
interest rate derivatives using the method of discounted expectations. The
resulting integrals are low-dimensional when the most commonly traded European
interest rate derivatives are considered, and they can be computed efficiently
using elementary numerical quadrature schemes (i.e. Simpson's rule). The
algorithm obtained is very well suited for parallel computing and is tested on
the problem of pricing several derivatives including an European swaption and
an interest rate spread option. In both cases, the numerical method proposed
appears to be accurate (i.e.: relative error of order 10-2, 10-3,
or even 10-4) and approximatively between 278 and 63000 times faster
than previous methods based on the Monte Carlo simulation of the LMM stochastic
dynamical systems. This website distributes to the potential users the software
LIBOR_SWAPTION, a Fortran code that implements the numerical method proposed to
compute the price of European swaptions.
Given the current time value T0 = 0, we define the
so-called “tenor structure” as a finite
set of time values Ti, i = 1, 2, … ,n+1, such that T0 <
T1 < T2 < … Tn+1. Moreover we
define the so-called “years fractions” δi
= Ti+1 - Ti , i = 1, 2,
… , n.
For i =
1, 2, … , n+1, let us consider a bond
which matures at time Ti
,whose price at time t, T0 ≤ t ≤ Ti , is denoted by
Pi(t).
The forward LIBOR rate Li(t) is defined as:
(1)
The LIBOR Market Model (LMM), introduced by
Brace, Gatarek and Musiela [3] , Miltersen, Sandmann, Sondermann [8] , and
Jamshidian [6], is the following set of systems of stochastic differential
equations:
(2)
with initial conditions:
(3)
where Wq+1(t) is a p-dimensional standard Wiener process, dWq+1(t) is its stochastic differential, 
is a p-dimensional
vector valued function of time, that models the volatility of the LIBOR
interest rate Li(t), i
= 1, 2, … , n, · dWq+1
denotes the Euclidean scalar product of and dWq+1,
Li,0 is a
positive real number and is the value assumed by the random variable Li(T0) with probability one, i = 1, 2, … , n, q +1 is an integer such that 0 ≤ q ≤ n that indicates the
bond Pq+1(t)
that has been chosen as numeraire asset (see [4] ,
[6] ), and
(4)
where e=Min(i,q+1).
We consider
the LMM model with one driving factor, that is the stochastic dynamical systems
(2)-(3) written when p=1. In
particular we derive an approximation of the probability density functions
associated to the model (2)-(3) when p=1. When p=1 we write respectively dWq+1(t) and 
instead than dWq+1(t) and , i=1,2,…,n.
Let us
consider the change of variables:
(5)
Using the Ito’s lemma and setting p=1, the set of stochastic differential
equations (2)-(3) can be rewritten as:
(6)
(7)
where
(8)
As announced previously we derive a closed-form
approximation of the probability density functions associated to the model
(6)-(8) [1] .
Let h be a positive integer such
that h≤Min(n,q+1),
Zh(t) = (Zh(t), Zh+1(t),…,
Zn(t)), zh(t) = ( zh(t), zh+1(t),…,
zn(t)),
z´h(t)
= (z´h(t),
z´h+1(t),…,
z´n(t)), τ
= t-Th-1. Note that in the previous sentence the variable
τ should
depend on h. We omit this dependence
to keep the notation simple.
Let ph(τ, zh , z´h),
0< τ <Th-Th-1
, denote the probability density function of having Zh(t) = z´h given Zh(Th-1) = zh ,
Th-1< t ≤Th, h = 1,2,…,Min(n,q+1).
The probability density functions ph(τ, zh , z´h),
h = 1,2,…, Min(n,q+1), are the solutions of the
Fokker-Planck equations associated to the stochastic systems of differential
equations (6)-(8) with the necessary initial and boundary conditions.
Let ηj,h
be a suitable positive constant, we consider the problem (6)-(8) when we
substitute to νj(τ) the function
ηj,h + εj,h (νj(τ) - ηj,h), where εj,h is a real parameter, j=1,2,…,n, h=1, 2,…, Min(n,q+1).
We note that when εj,h =0, j=1,2,…,n, h=1,2,…,Min(n,q+1), we have an approximate problem
with time independent volatilities and when εj,h =1, j=1,2,…,n, h=1,2,…,Min(n,q+1), we have the original problem
(6)-(8).
We approximate the probability density function
ph(τ, zh ,
z´h) using a
perturbative power series in the parameters δi,
εj,h, i = h,h+1,…,n, j = h,h+1,…,n, h=1,2,…, Min(n,q+1)
and we truncate this perturbation series at first-order. That is we look for
the coefficients of the asymptotic expansion:
(9)
In [1] we have obtained simple closed-form
expressions involving elementary functions and the Dirac delta function and its
derivatives for p0h(τ, zh, z´h),
pδ,1h,i(τ, zh,
z´h), pε,1h,j(τ, zh,
z´h), i=1,2,…,n, j=1,2,…,n, h=1,2,…,Min(n,q+1). Moreover, due to the fact that the
LMM has been considered with p=1, the
elementary functions appearing in these closed-form expressions depend only on
one variable, and therefore can be used efficiently to price European interest
rate derivatives using the method of discounted expectations, since the
integrals that solve the pricing problems are low-dimensional.
Let us consider the case of an European payer
swaption, whose payoff, paid at time T1,
is given by the following function:
(10)
where
, (11)
and K>0 is the strike price.
The swaption price at time T0, that we denote with Vsw(T0),
is computed using the method of discounted expectations, that is:
, (12)
where
(13)
and
(14)
Due to the fact that the LMM model has been
considered with one driving factor (p
= 1 in (6)), the evaluation of the integral appearing in (12) leads to the
evaluation of a one-dimensional integral, that can be done using an elementary
quadrature rule.
We consider a (5 x10) European payer swaption based on
the 6-month forward LIBOR rate, that is a swaption maturing at time T1 = 5 year on a swap spanning the tenor structure T1, T2,
…, Tn+1 of
length 10 year, where n = 20, Ti+1-Ti = 0.5 year, i = 1,2,…,n.
The volatility of the forward LIBOR rate Li(t) is chosen as suggested in [5] , that is:
(15)
where a
= 0.2934 year-3/2,
b = 1.251 year-1, c =
0.1315, d = 0 year-1/2. Moreover we consider three different
cases: Li,0=
0.06 year-1, i = 1,2,…,n, and K= 0.04 year-1 (Test Case A), Li,0= 0.05 year-1, i = 1,2,…,n, and K= 0.05 year-1 (Test Case B), Li,0= 0.04 year-1, i =
1,2,…,n, and K= 0.06 year-1
(Test Case C), corresponding to a swaption in the money, at the money, and out
of the money, respectively.
Let us begin with Test Case A. We have computed
the price of the swaption at time T0
= 0 using formula (12). The discount factor P1(T0) in (12) has been chosen
equal to 1. The (one-dimensional) integral contained in (12) is evaluated using
Simpson's rule with 400 subdomains. This choice suffices to keep the relative
error due to the numerical integration in formula (12) to be approximately of
order 10-4. The result obtained using the analytical approximation
(12), that we denote with Vsw,an , is Vsw,an
= 0.1538. The computer time needed to evaluate formula (12) with the choices
described previously using a computer with a Pentium 4 Processor 1700 MHZ 256
MB Ram is 0.17 second.
In order to test the accuracy of this
approximation an accurate estimation of the exact option value, that we denote
by Vsw,true , has been obtained with a
We have performed 100.000.000 simulations of
the payoff function (10) where the Euler-Maruyama scheme has been applied with
500 time discretization steps of equal size to discretize the interval [T0, T1]. The “true” option value obtained is Vsw,true = 0.1539. Note that there is empirical evidence
suggesting that the value computed for Vsw,true has
four significant digits correct. Hence we evaluate the relative error of the
approximate formula (12), that we denote by RelErr,
as follows: RelErr = (0.1539 -
0.1538)/0.1539 = 6.5 x 10-4.
Numerical experiments reveal that the Monte
Carlo simulation converges to Vsw,true, up to an error on the fourth significant decimal digit,
also when the Euler-Maruyama scheme is applied with only 80 time discretization
steps of equal size in the interval [T0,
T1], that is with the
previous choices the time discretization error of the stochastic differential
equations (6) is negligible compared to the other errors present in the
computation.
Figure 1 shows the (relative) numerical error
obtained employing both the approximation formula (12) and the
The computer time needed to perform 700.000
Monte Carlo simulations with 80 time
discretization steps in the simulation of (6), (8) is approximatively 583.2
seconds using the Pentium IV processor.
We conclude that the computer time needed to calculate the swaption
price with the accuracy discussed previously using the numerical method
described here and proposed in [1] is about 583.2/0.17 ≈ 3430 times
smaller than the computer time necessary to the
Figure 1, Test Case A, error of the
analytical formula (12) compared to the
error of the
The numerical simulation relative to Test Case
B and Test Case C is performed similarly to the simulation relative to Test
Case A. The results obtained are summarized in Table 1, where we show, for Test
Cases A, B, and C, the swaption price computed using the analytical
approximation (12), Vsw,an, the relative error RelErr
of the approximation (12) compared to a true value obtained with an
accurate Monte Carlo simulation, the computer time necessary to evaluate
formula (12), CPUTimean,
the computer time necessary to evaluate the swaption value using the Monte
Carlo simulation with approximately the same accuracy obtained using the
analytical method, CPUTimeMC.
We note that the relative error obtained using the approximation (12) is in
Test Case B, equal to 3.5 x 10-3, and, in Test Case C, equal to 1.1 x 10-2. Moreover, in Test
Case B, the analytical approximation
(12) is 62.5/0.17≈ 368 times faster than the
Table 1
|
|
Vsw,an |
RelErr |
CPUTimean |
CPUTimeMC |
|
Test Case A |
0.1538 |
6.5 x 10-4 |
0.17 s |
583.2 s |
|
Test Case B |
0.04816 |
3.5 x 10-3 |
0.17 s |
62.5 s |
|
Test Case C |
0.005668 |
1.1 x 10-2 |
0.17 s |
65.8 s |
We have applied the numerical method described
above to price European swaptions considering tenor structures, volatility
functions, strike prices and initial data different from those used in the
simulation presented here. We have also employed the numerical method to the
problem of pricing interest rate spread options. The results obtained are
contained in [1] , and are very similar to those presented here, in fact in the
numerical experience reported in [1] the numerical method proposed appears to
be accurate (i.e.: relative error of order 10-2, 10-3, or
even 10-4) and approximately between 278 and 63000 times faster than
the Monte Carlo simulation when results of approximately the same accuracy are
produced. The results obtained demonstrate the possibility of pricing interest
rates derivatives under the LMM with one driving factor using numerical methods
based on an appropriate combination of mathematical analysis and elementary numerical
methods rather than on relying exclusively on the
[1] L.
V. Ballestra, G. Pacelli, F.Zirilli, A Numerical Method to Price Interest Rate Derivatives
based on the One Factor LIBOR Market Model of Interest Rates,
Nonlinear Analysis: Hybrid
Systems
2 (2008) 568-589.
[2] F. Black, The Pricing of Commodity
Contracts, Journal of Financial Economics
3 (1976) 167-179.
[3] A.
Brace, M. Gatarek, M. Musiela, The Market Model of Interest Rate Dynamics, Mathematical Finance 7 (1997) 127-155.
[4] D. Brigo, F. Mercurio, Interest Rate Models: Theory and Practice,
[5] D.
Brigo, F. Mercurio, M. Morini, The LIBOR Model Dynamics: Approximations,
Calibration and Diagnostics, European
Journal of Operational Research 163
(2005) 30-51.
[6] F. Jamshidian, LIBOR and Swap Market Models
and Measures, Finance and Stochastics
1 (1997) 293-330.
[7] P.
E. Kloeden, E. Platen, Numerical Solution
of Stochastic Differential Equations,
[8] K.
Miltersen, K. Sandmann, D. Sondermann, Closed Form Solutions for Term Structure
Derivatives with Lognormal Interest Rates, Journal
of Finance 52 (1997) 409-430.
6. LIBOR_SWAPTION:
A Fortran code to compute the price of an European swaption using the one
factor
LIBOR Market Model
The codes
made available are: Fmain.f, eigeis.f.
The archive
LIBORS_WAPTION.zip contains the following files:
Fmain.f: a Fortran main code that computes Vsw(T0).
eigeis.f: a Fortran code that computes the eigenvalues and eigenvectors
of a real symmetric matrix. This code is part of the EISPACK library routines,
and has been taken from
http://www.netlib.org/eispack .
click
here to download LIBOR_SWAPTION.zip
