Pricing of options and other financial instruments is one of the most important problems in finance. To price an option, one needs first to choose the model of underlying asset price dynamics, usually in the form of a stochastic diferential equation depending on market parameters such as interest rate, asset volatility, etc., then find a solution of the chosen model in some form, and finally obtain the required option price. The asset volatility can be constant (Black-Scholes model), depend on the value of the underlying asset and time (local volatility model), or satisfy some other stochastic diferential equation (stochastic volatility model). Exact analytical formula for option price was obtained for Black-Scholes model and a few other models with local and stochastic volatility, but in general case one has to use approximations or numerical methods to evaluate options. We study the problem of European option pricing when volatility is a function of underlying asset value and time. To solve the problem, we use Monte Carlo method to construct the empirical distribution density of underlying asset and then determine the option value using Feynman-Katz formula. The issues of parallelization of the option pricing algorithm and its implementation on computers with multicore processors are discussed. Possible applications of local volatility models with parallel computing include modeling and management of large equity portfolios, assessing and managing market and credit risks.
The increasing amount of data being processed and active use of machine learning methods lead to a growing demand for high-performance computing (HPC) solutions in modern financial industry [1]. The traditional applications of HPC in finance include pricing of financial instruments (derivative securities), algorithmic high-frequency trading, market and credit risk
© 2020 Copyright for this paper by its authors.
CEUR Workshop Proceedings htp:/ceur-ws.org IS N1613-073
CEUR Workshop Proceedings (CEUR-WS.org) Use permitted under Creative Commons License Attribution 4.0 International (CC BY 4.0). However, the assumption of constant volatility and based on this assumption Black-Scholes model fail to explain some important efects observed in financial markets, such as volatility smiles [4] and fat tails of financial data distributions [ 5]. Therefore, alternative stochastic models [6] with volatility depending on the value of the underlying asset and time (local volatility models) should be examined to determine the most appropriate model for available ifnancial data. However, for a general local volatility model, exact closed form formula for the option price has not been derived and miscellaneous numerical methods have to be used.
The paper considers the application of high-performance computing in calculating the value of options under assumption that underlying asset price follows stochastic diferential equation (SDE) with volatility, being a function of asset price and time. 2. Local volatility models In risk-neutral framework the price of underlying asset is supposed to follow SDE where > 0 is the risk-free interest rate, (, ) > 0 is a volatility function of underlying asset value and time , is a Wiener random process, with initial condition = +
(, ) , ( 0) = 0 > 0. where
We intend to determine the fair value of a derivative (option or other derivative instrument) ( 0, 0), being a function of the current price of the underlying asset 0 and current time 0, with payment function (, ) = ( ) at the time > 0 of execution (expiration) of the derivative.
Known local volatility models of the form (
To find the value of an option (derivative) in model (
0, 0) ( ) , +∞ −∞ where ( , ; 0, 0) is the probability density of the transition from state ( 0, 0) to state ( , ).
In Black–Scholes model [2, 3] the transition probability density function (pdf) is equal to
(
(arsinh ( √2 ) − arsinh ( √2 √2 √ 2 + 2 2√ − 0 √ 2 + 2 and the transition pdf is
2 0)) }
In normal (Ornstein–Uhlenbeck) model [8] the volatility function is and the transition pdf
In hyperbolic sine model [9] the volatility function is 2
√2 √ − 0 ( − )
also assuming negative underlying asset values.
3. Monte Carlo method for option pricing in local volatility
The Feynman–Katz formula (
0, 0), which is the distribution density of the underlying asset ( ) at
derivative execution (expiration) time for fixed initial values 0 and 0. For a general local
volatility model (
For approximate determination of the distribution density of the random variable ( ), one
can use Monte Carlo method [10], which implies generation of a large number of realizations of
the random variable ( ), namely, a large bundle of trajectories of SDE (
The time segment [ 0, ] is divided into equal parts of length △ = − 0 .
Let = 0 + △ , =
1, , then the simplest numerical approximation of the trajectories of
SDE (
+1 ≈ + △ + ( , ) √△ , where = ( ), ∼
(
Let 1 = min {
parts of length ℎ = ( ) } =1 , 2 = max {
} . The segment [ 1, 2] is divided into equal into the last segment [ 1 + ( − 1 )ℎ, 1 + ℎ ] for = [ 1 + ( − 1 )ℎ, 1 + ℎ ) for = 1, − 1 .
Empirical distribution density of the random variable ( ) can be obtained by the formula
2− 1 and let be the number of values from the set {
and half-intervals of the form
( )
}
=1
that fall
The simulation of trajectories of SDE (
1 ⎧0, ⎩0,
1 + ( − 1 )ℎ ⩽ < 1 + ℎ, = 1, − 1, ⎨⎪ , 1 + ( − 1 )ℎ ⩽ ⩽ 1 + ℎ, ( )
=1 < 1 > 2
. =1 ∑ = [ − ℎ 1 + 1 ] 2 1 2 1 2
(
Then, using the empirical density (
For a European call option with strike price the payment function ( ) is equal to max( −
, 0) , so the price of call option is
( 0, 0) ≈ − ( − 0) 2 − 1
( 1 + ( −
) ℎ − ) ,
where the square brackets in (
The values ( )
, = 1, can be calculated in parallel, and this makes it possible to use parallel 4. Parallel computing in Python There are diferent ways to organize parallel computing in Python [ 14], including multithreading (module t h r e a d i n g ) and multi-processor code execution (method f o r k of module o s or module m u l t i p r o c e s s i n g ). If parallel processes do not interact with each other, then parallel computations are better done using the m u l t i p r o c e s s i n g module API, which is a part of the standard Python library.
To execute multiple calls of the function l v m _ s d e , which returns the value l v m _ s d e is launched for parallel execution using method m a p , which applies the function l v m _ s d e ( ), in parallel, to the list of input parameters i n p u t _ l i s t : i m p o r t m u l t i p r o c e s s i n g a s m p l v m _ p o o l = m p . P o o l ( ) r e s u l t = l v m _ p o o l . m a p ( l v m _ s d e , i n p u t _ l i s t )
The number of calls of the function l v m _ s d e corresponds to the length of the i n p u t _ l i s t list.
The function l v m _ s d e for calculating ( ) according to formula (
SDE simulation parameters (variables s 0 , t _ s t a r t , t _ s t o p , N and lambda expressions d r i f t ,
d i f f u s i o n ) are passed to the function l v m _ s d e as global variables.
5. Empirical distribution construction in local volatility models
Empirical pdf for random variable ( ) can be obtained using Monte Carlo method for various
local volatility models of the form (
Simulation parameters for SDE (
Theoretical probability density functions in local volatility models under consideration are shown in Figure 1.
Lambda expressions for simulation of SDE (
The simulation of SDE (
1.0 1.5 2.0 2.5 Ornstein-Uhlenbeck model
Empirical pdf Theoretical pdf 3.0 0.0 0.5 1.0 1.5 2.0 2.5 Hyperbolic sine model
Empirical pdf Theoretical pdf 3.0 0.0 0.5 1.0 1.5 2.0 2.5 3.0 0.0 0.5 1.0 1.5 2.0 2.5 3.0
Parallel computing using module m u l t i p r o c e s s i n g on MacBook Pro with 8-Core Intel Core i9 processor can significantly accelerate the execution of program code (see Table 1) with almost full load of MacBook Pro i9 processor cores (see Figure 3).
Further speedup of calculation of empirical distribution density can be gained using GPU resources [15]. When using GPU, program code can be accelerated many times depending on the task being solved. For instance, on MacBook Pro with Intel Core i9 processor and discrete graphics card AMD Radeon Pro 5500M (1536 shader processors) the program code can be accelerated hundreds or even thousands of times compared to sequential execution on CPU.
Parallel GPU computing can be implemented in a Python program using the PyCUDA library (for Nvidia graphics cards) or the PyOpenCL library (for AMD graphics cards) [16]. The idea of using Monte Carlo method on a GPU for financial simulation was studied in a number of papers [17], including applications to financial risk measurement [ 18]. Calculating multiple trajectories of a multidimensional system of SDE in parallel using OpenCL is a subject of research in [19].
OpenCL framework lacks a standard random number generator, so in order to take advantages of PyOpenCL library to build empirical distributions and evaluate options with a GPU, one has to implement a pseudorandom number generator. The recognized random number generators implemented on GPU are collected in Random123 library [20].
The publication has been prepared with the support of the “RUDN University Program 5–100”. The research was funded by RFBR, grant No. 19-08-00261.