A Finite Element Framework for Option Pricing with the Bates Model Submitted for publication to Computing and Visualization in Science

A Finite Element Framework for Option Pricing with the Bates Model thanks: Submitted for publication to Computing and Visualization in Science

E. Miglio,C. Sgarra

In the present paper we present a finite element approach for option pricing in the framework of a well-known stochastic volatility model with jumps, the Bates model. In this model the asset log-returns are assumed to follow a jump-diffusion model where the jump component consists of a Lévy process of compound Poisson type, while the volatility behavior is described by a stochastic differential equation of CIR type, with a mean-reverting drift term and a diffusion component correlated with that of the log-returns. Like in all the Lévy models, the option pricing problem can be formulated in terms of an integro-differential equation: for the Bates model the unknown (the option price) of the pricing equation depends on three independent variables and the differential operator part turns out to be of parabolic kind, while the nonlocal integral operator is calculated with respect to the Lévy measure of the jumps. In this paper we will present a variational formulation of the problem suitable for a finite element approach. The numerical results obtained for european options will be compared with those obtained with different methods.

MOX– Modellistica e Calcolo Scientifico

[0pt] Dipartimento di Matematica “F. Brioschi”

[0pt] Politecnico di Milano

[0pt] via Bonardi 9, 20133 Milano, Italy

[0pt] edie.miglio@polimi.it

Dipartimento di Matematica “F. Brioschi”

[0pt] Politecnico di Milano

[0pt] via Bonardi 9, 20133 Milano, Italy

[0pt] carlo.sgarra@polimi.it

Keywords: Option Pricing, Stochastic Volatility Models, Lévy Processes, Partial-Integro-Differential-Equations, Finite Element Methods

AMS Subject Classification: 62P05, 60G35, 45K05, 65L60

1 Introduction

A huge effort has been made in the last few years in order to overcome the intrinsic limitations of the Black-Scholes model. Although it has been a great success as a first attempt to provide an evaluation for financial derivatives, it was soon clear that it’s description of financial market behavior was not satisfactory. Very well known observed empirical features of the log prices were not correctly described by this model: heavy tails, volatility clustering, aggregational gaussianity are some peculiarities that cannot be explained on the basis of the lognormal assumption on which the Black-Scholes model stands. The volatility smile is another relevant phenomenon that cannot be explained on the basis of a Black-Scholes description. Several different approaches have been exploited in order to give a more satisfactory description of financial markets, but the main contributions in this direction can be grouped in two different classes of models, the so called stochastic volatility models and the models with jumps. An extended literature is available on both kind of approaches and they both give a more realistic description of the prices evolution in financial markets, but separately considered they perform very well only in some situations. While jump models can succesfully reproduce the volatility smiles on short term maturity ranges, stochastic volatility models give a better description of the same phenomenon on long maturity terms. This has naturally led to the introduction of more complicated, but more realistic models in which both features of stochastic volatility and jumps can be present. The three more popular models in which this integration of jumps and stochastic volatility has been performed are the BNS model introduced by O. Barndorff-Nielsen and N. Shepard [3], [4], the model introduced by Bates in [5], and the time-changed Lévy models introduced by Carr, D. Madan, H. Geman and M. Yor in [18]. While in the former the volatility dynamics is driven by a positive Lévy process correlated with the jump process in the log-price of the asset, in the latter the volatility dynamics is governed by a time-changed Lévy process. We shall concentrate in the present work on the second model we have just mentioned, the Bates model in which a Merton jump-diffusion model is combined with a stochastic volatility model of the Heston type. As R. Cont and P. Tankov have pointed out, the performance of the time changed Lévy models in calibrating to market option prices are usually much better than those obtained in a BNS framework [14]: in the last model, in fact, the possible implied volatility patterns are restricted by the requirement that the same parameter characterize both the jumps in the returns and the volatility; on the other hand the calibration performances of the Bates model are comparable to those of the time changed Lévy processes, ”Thus the Bates model appears to be at the same time the simplest and the most flexible of the models” [17], pag. 495. In the framework of option pricing via PDE (PIDE for models with jumps) several different approaches have been exploited both for stochastic volatility models and models with jumps. As far as finite element methods are concerned, we just quote the papers by Y. Achdou and N. Tchou [1] and by N. Hilber, A.-M. Matache and C. Schwab [14] for the first class of models and the papers by A.-M. Matache, T. von Petersdorff and C. Schwab [11] and by A.-M. Matache, P.-A. Nitsche and C. Schab [12] for the second. For models including both features the numerical methods proposed until now are much less. Some authors have considered finite-difference schemes for these models. D. Hilber, A.-M. Matache and Schwab have provided a finite-element approach to a large class of stochastic volatility models without jumps in [15], including the Heston model. In the present paper we shall present a finite-element approach to option pricing in a Bates model framework. In the next section we’ll recall the basic features of the Bates model, while in section 3 we’ll provide the PIDE formulation of the option pricing in this model; in section 4 we shall present the variational formulation of the problem. In section 5 we’ll describe the numerical method proposed and in section 6 we’ll expose some comments on the results obtained. In section 7 we’ll outline some conclusions and some futures perspectives of the present work.

2 The Bates model

In the Bates model the asset price evolution is given by:


where the log-returns and its volatility  satisfy the following stochastic differential equations:


where and are two standard Wiener processes with correlation and is a Lévy process of compound Poisson type . Let’s assume for the parameters the following restrictions:


The last requirement is in order to assure that the volatility process will never hit zero. We’ll assume moreover , this implying that the Lévy-Khinchine representation formula for the process will be of the following type:


where and is the Lévy measure of . We’ll denote by the jump measure of  and by its predictable compensator. We’ll have moreover . We can write then:


The stochastic differential equation for the price will be then:

Remark 1

In the original model of [5], the process  is a compound Poisson process,


where is a standard Poisson process with intensity  and are iid , with . The corresponding cumulant function is in that case

Remark 2

If then we obtain Heston’s stochastic volatility model from [8]. If and we obtain Merton’s jump diffusion model from [13] with normal jumps. Consequently we might consider the Bates model as an extension of a Merton jump-diffusion model with stochastic volatility, or as an extension of a Heston volatility model with jumps in the returns.

Lemma 1

The dynamics of the asset price process is given by


In particular if


the process is a local martingale.

Proof: This follows immediately from Itô’s formula for general semimartingales.

As in other affine stochastic volatility models with and without jumps, it is possible to obtain the characteristic function of the log-price in closed form. This characteristic function has been calculated by D. Bates in [5]; a detailed computation is provided also in [17]; it is given by the following expression:




Once the characteristic function of the log-price process is known in a closed form, the valuation problem for vanilla options can be easily solved by an FFT-related technique like that provided in the paper by P. Carr and D. Madan [6].

A quadratic approach to option hedhing in the Bates model has been suggested in [10].

3 Option pricing via PIDE approach

By applying Ito lemma, and introducing the market price of risk , we obtain the following partial integro-differential equation for the price of a European call option  in the framework of the Bates model :


with the following final condition at :


and the following boundary conditions both in and :


where are the Black-Scholes values of the call options at time for underlying price and strike with parameters


By assuming lognormal jumps we have


where are the expected value and variance respectively of the normal distribution describing the jumps’ size.

The variables and can assume values in the following domains: , and .

Remark 3

The market price of risk can be obtained in different ways in the frame of general equilibrium models; consumption-based models give a risk premium proportional tu . In the following, we’ll assume without any lost of generality, that the market price of risk associated to the volatility is zero. The method can be generalized to a different choice in a straightforward way.

Remark 4

By the commonly used substitution ,  the previous equation becomes:


is the Lévy density of the jumps. In the case of gaussian jumps for (i.e. lognormal for ) it will be of the following form:


It will be then a gaussian density with expected value and variance multiplied by the intensity of the Poisson process.

Remark 5

The boundary conditions for the new unknown are now the following. For :


The final condition is now:


4 Variational Formulation

The integro-differential equation given before can be written in the following ”divergence form”:


where the symbol denotes the total derivative, the vector is given by:


and the matrix by:


A variational formulation for the PDE arising in the Heston model has been given in [15], while the variational formulation for the PIDE in an Exponential Lévy framework have been provided in [11], [12], [21], where existence and uniqueness of the solution of the variational problem associated with the differential and integro-differential equations respectively have been proved in suitable weighted Sobolev spaces and detailed analyses of both localization and discretization errors have been provided. Without performing the same analysis here we are going to proceed with the variational formulation for the present problem following the same line.

Introducing the following bilinear form:


it is possible to define a suitable discretization for our integro-differential problem 25.

Remark 6

The equation 17 has degenerate coefficients both in the and in the variables; while the substitution removes the degeneracy in the variable, this is still present for the variable. In particular we want to discuss the boundary . If no restrictions on the model’s parameters would be present, a condition on this boundary should be imposed in order to have a well-posed problem and the correct condition is given by 28 . On the other hand, if the restriction on the parameters of the model given in section 2 5 holds, the variable never hits that boundary. A closer inspection to the bilinear form associated to the problem suggests that when Green’s lemma is applied to the l.h.s. of equation 25 an ”inflow” condition for our backward parabolic equation appears if the restrictions on the parameters hold and this in turn implies that the condition on need not to be imposed in order to have a well-posed problem.

5 Numerical Approach

In this section we will introduce a Finite Element Discretization of the above described PIDE. A similar approach for the Merton’s and Kou’s model has been introduced in [2].

5.1 Temporal Discretization

The integro-differential equation can be written in the following form


where the symbol denotes the total derivative.

The time interval will be discretized using a time step such that ; moreover we will use the following notation . Starting from this formulation we can obtain the following temporal discretization


where is the value of the price evaluated at the foot of the characteristic line at time and is the solution of the following initial value problem


This last ODE can be solved using either the implicit Euler method or a more accurate Runge-Kutta method.

The characteristic Galerkin method, described above, is an Eulerian-Lagrangian approach and it is stable with a mild stability criterion allowing therefore the use of a large time step when appropriate.

5.2 Localization

The infinite domain has to be reduced to a finite one: in particular we will consider a rectangular domain . Using this reduced domain the boundary conditions have to be modifed accordingly in the following way:


Moreover the extrema of the integral term have to be reduced to finite values; to this aim we can use




5.3 Spatial Discretization

The domain will be discretized using an unstructured triangular mesh.

The discrete weak formulation reads as follows: find such that


where is suitable functional space. In the present paper the unknown will be approximated using (linear) finite element i.e.


where is the number of nodes of the triangulation and .

The problem has been solved using FreeFEM++.

6 Numerical Results

In this section we shall provide some comments on the numerical results obtained in order to assess the effectiveness of the proposed numerical method. In particular we’ll present the implicit volatility surfaces for the following sets of parameters taken from [16], where a suitable calibration methodology has been developed for a large class of stochastic volatility models with jumps. Moreover we’ll provide a graphical comparison between the solution obtained with our finite element approach and that obtained by the usual method proposed by P. Carr and D. Madan based on the Fast Fourier Transform [6].

The range of the strike prices for the european call option considered is for Fig.1, 3, 4, while for Fig. 2 is . The maturities are between 0 and 3 years. The initial value of the underlying has been set .

Set 0.21568 0.04937 0.23828 -0.44793 0.33502 0.033582 0.26969 -0.42404 0.13279 0.18193 0.37518 -0.59722 0.48443 0.022097 0.21903 -0.40066 Set -0.11889 0.17189 0.13674 -0.077973 0.11048 0.33785 0.080396 0.057373 0.05218 -0.12938 0.16878 0.15977

Figure 1: Implicit Volatility Surface for parameter set S1.
Figure 2: Implicit Volatility Surface for parameter set S2.
Figure 3: Implicit Volatility Surface for parameter set S3.
Figure 4: Implicit Volatility Surface for parameter set S4.

All the volatility surfaces obtained exhibit both smiles and skews for short maturities. While the skew is quite pronounced for the set of parameters denoted by , it seems less evident for the other sets , . The smiles appear less and less pronounced for longer maturities for all sets. The explanation of the different behavior exhibited by the volatility surfaces seems to be related more to the different values of the diffusion coefficient of the volatility dynamics than to the other parameters, while the ”leverage” coefficient seems to be responsible of both the skews and the smiles appearing in the volatility surfaces. The behavior of the smiles in correspondence to the ”at the money values” of the call option looks moreover quite realistic.

The next two figures represent the solution obtained with the present finite element method versus the solution obtained via a Fast Fourier transform method. While the latter is represented by the continuous line, the former is indicated by the dots. The set of parameters characterizing the model are those denoted before by , respectively. The range of prices is . The calculation has been performed for , .

A very good agreement between the solutions can be immediately recognized, although it looks quite evident that the solution obtained via the finite element method slightly overprices the call option for higher values of the underlying at time 0.

In order to check the robustness of the present method with respect to parameters changed, several trials have been performed corresponding to other sets of parameters belonging to an enlarged range of parameters and the results obtained look very close to those just presented.

The CPU time required for a calculation of the call option price for a single value of the underlying is about 3 minutes, while that required for the volatility surfaces shown here is about 3 hours on a 2Ghz Centrino Duo with 1MB RAM. The 2D mesh used in the computation is composed by 2634 triangles and 1398 nodes.

Figure 5: Comparison between the FFT solution and the FEM solution for set : and .
Figure 6: Comparison between the FFT solution and the FEM solution for set : and .

7 Concluding Remarks

We have presented a finite element approach to the european option valuation problem formulated via a partial integro-differential equation. Several choices of the functional space suitable for the spatial discretization are possible and we have made the most simple choice in order to obtain a fast and efficient algorithm. Other choices have been made by some authors in similar contexts, like in [9], [12], where Wavelets have been used in an extensive way to produce an accurate algorithm, but the convergence of the method with this choice seems to be not very fast. In the present context, with a single underlying, the choice of piece-wise polynomial function spaces seems to perform slightly better.

A natural development of the present analysis will be the variational formulation and the implementation of a finite element method for American option pricing with the Bates model. This problem can be formulated as a free boundary problem for the same Partial Integro-Differential equation given before. The finite element approach seems to be the most promising way to obtain fast and accurate algorithms for this problem. This will be the subject of future investigations.


  • [1] Achdou, Y., Tchou, N.: Variational analysis for the Black-Scholes equation with stochastic volatility, ESAIM: Mat. Mod. and Num. Analysis, 36/3, 373-395 (2002).
  • [2] Almendral, A., OOsterlee, C.W.: Numerical valuation of options with jumps in the underlying, Appl. Num. Math., 53, 1-18 (2005).
  • [3] Barndorff-Nielsen, O.E., Shepard, N.: Ornstein-Uhlenbeck-based models and some of their uses in financial economics (with discussion), J. R. Stat. Soc. B, 63, 167-241 (2001).
  • [4] Barndorff-Nielsen, O.E.: Modelling by Lévy processes for financial econometrics, in Lévy Processes-Theory and Applications, 283-318. Eds. O. Barndorff-Nielsen, T. Mikosch, S. Resnick, Birkauser (2001).
  • [5] Bates, D.: Jumps and stochastic volatility: the exchange rate processes implicit in Deutschemark opions, Rev. Fin. Studies, 9, 69-107 (1996).
  • [6] Carr, P., Madan, D.: Option valuation using the fast Fourier transform, J. of computational finance, 2, 61-73 (1998).
  • [7] Cont, R., Voltchkova, E.: A Finite Difference Scheme for Option Pricing in Jump-Diffusion and Exponential Levy Models, SIAM J. Numer. Anal. 43/4, 1596-1626 (2005).
  • [8] Heston, S.L.: A Closed-Form Solution for Options with Stochastic Volatility with Applications to Bond and Currency Options, Rev. Fin. Studies, 6/2, 327-343 (1993).
  • [9] Hilber, N., Reich, N., Schwab, C., Winter, C.: Finite Element based derivative pricing on general price processes, Presented at Risk Day 2006, ETH Zuerich.
  • [10] Hubalek, F., Sgarra, C.: Quadratic Hedging for the Bates Model, Int. Journ. of Theor. Appl. Finance, 10/5, 873-885 (2007).
  • [11] Matache, A.-M., von Petersdorff, T., Schwab, C.: Fast Deterministic Pricing of Options on Lévy Driven Assets, ESAIM: Mathematical Modelling and Numerical Analysis, 38/1, 37-71 (2003).
  • [12] Matache, A.-M., Nitsche, P.-A., Schwab, C.: Wavelet Galerkin Pricing of American Options on Lévy Driven Assets, Quantitative Finance, 5/4, 403-424 (2005).
  • [13] Merton, R.: Option pricing when underlying stock returns are discontinuous, J. Fin. Economics, 3, 125-144 (1976).
  • [14] Schoutens, W.: LEVY PROCESSES IN FINANCE, Wiley (2003).
  • [15] Hilber, Matache, A.-M., Scwab, C.: Sparse wavelet methods for option pricing under stochastic volatility, J. Comp. Finance, 8/4, 1-42 (2005).
  • [16] Balmisse-Thomsen, J.C.: Model Risk in pricing of Exotics, Master Thesis, Aarhus University (2006)
  • [17] Cont, R., Tankov, P.: FINANCIAL MODELLING WITH JUMPS, Chapman/Hall-CRC (2004).
  • [18] Carr, P., Geman, H., Madan, D., Yor, M.: Stochastic Volatility for Lévy Processes, Math. Finance, 13, 3455-382 (2003).
  • [19] Venardos, E.: Derivatives Pricing and Ornstein-Uhlenbeck Type Stochastic Volatility, Ph.D. Thesis Nuffield College of Oxford (2001).
  • [20] Wu, X.: Stochastic Volatility with Lévy Processes: Calibration and Pricing, Ph.D. Thesis, Dept. of Finance, University of Maryland (2005).
  • [21] Zhang, X.: Analyse Numeriques des Options Américaines dans un Modèle de Diffusion avec des Sauts, Ph.D. Thesis, CERMA (1994).
Comments 0
Request Comment
You are adding the first comment!
How to quickly get a good reply:
  • Give credit where it’s due by listing out the positive aspects of a paper before getting into which changes should be made.
  • Be specific in your critique, and provide supporting evidence with appropriate references to substantiate general statements.
  • Your comment should inspire ideas to flow and help the author improves the paper.

The better we are at sharing our knowledge with each other, the faster we move forward.
The feedback must be of minimum 40 characters and the title a minimum of 5 characters
Add comment
Loading ...
This is a comment super asjknd jkasnjk adsnkj
The feedback must be of minumum 40 characters
The feedback must be of minumum 40 characters

You are asking your first question!
How to quickly get a good answer:
  • Keep your question short and to the point
  • Check for grammar or spelling errors.
  • Phrase it like a question
Test description