Affine POD Galerkin Schemes for Option Pricing in Jump-Diffusion Models

(1)

AFFINE POD GALERKIN SCHEMES FOR OPTION PRICING IN JUMP-DIFFUSION MODELS

by

Jianjie Lu

Submitted to the Department of Mathematics and Statistics in partial fulfilment of the requirements for the degree of

Master of Science at the

Universität Konstanz

January 2014

1st Supervisor: Prof. Stefan Volkwein 2nd Supervisor: Prof. Michael Kupper

Konstanzer Online-Publikations-System (KOPS) URL: http://nbn-resolving.de/urn:nbn:de:bsz:352-0-253591

(2)

Contents I

1 Introduction

The present thesis considers a pricing problem of financial derivatives in jump- diffusion models and focuses on the analysis and implementation of numerical methods for pricing European options. To be general, the numerical methods for option pricing can be grouped into two main categories: the stochastic methods based on Monte Carlo simulations and the deterministic methods based on the numerical solution of the Fokker-Planck-type partial differential equations implied by the pricing problem. We present here the latter class of methods and particularly address the Galerkin finite element method to solve the pricing equations numerically.

Jump-diffusion models, possibly equipped with local volatility terms, have gained their importance in modelling financial derivative pricing problem. It is now recognized that such jump-type models allow for a more dynamic and realistic representation of empirical risky asset price movements, comparing with pure continuous models such as the successful Black-Scholes model.

Given its surprising flexibility and moderate complexity in modelling, models with jumps have become increasingly popular in the option pricing industry.

Unfortunately, from the mathematical point of view, pricing options under the jump-diffusion model faces some challenges. First, on the theoretical side, such models imply an incomplete financial market and standard stochastic analysis has to be adjusted. In particular, option prices have to be obtained by solving partial integro-differential equations (PIDEs), which contain, in addition to the differential operators, non-local integral terms, for which special treatments need to be considered. Second, on the numerical side, by certain proper spatial discretization schema, the non-local integral term leads to dense linear systems of equations and efficient numerical methods are often needed, especially for pricing complex financial derivatives and calibration of model parameters.

Through this thesis, we wish to offer some insightful evidence not only at the theoretical level but also at the numerical level. To carry out, we first focus on the derivation of the PIDE, justification and numerical analysis of

(5)

solution of the PIDE. In particular, we show that the corresponding variational formulation admits a unique solution (or weak solution of the PIDE). To em- phasize the modelling flexibility, an essential derivation of the pricing PIDE is also presented using some basic facts of the exponential Lévy model, a more convenient and general way to characterize the jump-diffusion models.

In the second main part, with the model calibration problem in mind, we address the complexity of solving dense linear systems of equations resulted from the non-local integral term by applying the model order reduction method based on proper orthogonal decomposition (POD). It leads to reduced differential equations. By the finite element method, option prices (as solution of PIDEs) are projected on a sequence of basis functions. The POD provides with a possibility that the solution may be sufficiently represented by a signif- icantly smaller number of basis functions, resulting in small linear systems.

Put it in another way, only the essential information in the (linear) system is retained by POD. We also present a priori estimate of the corresponding approximation error.

While our analysis in the first part could be well understood for readers with a background of computational finance, the application of POD method in pricing financial derivatives is relatively new. For this reason, we first give an illustration of the POD method using image compression, where the application of the POD method is straightforward.

1.1 Proper Orthogonal Decomposition

Our original image is a white and black image of resolution 600×800, shown in figure 1.1-(f ). In the language of mathematical computation, the image can be represented using a matrix of dimension 600×800 with entries between 0 (black) and 1 (white). In total, we have 0.48 million numbers.

While each column of the matrix can be interpreted as a linear combina- tion of unit vectors (basis functions), it is by no means the most efficient way to express the columns. According to the POD method, we know then that there exit more complex basis functions containing the most ‘essential’ information of the image, leading to a more efficient representation of the columns.

More convincingly, the image can be reconstructed in such a way that the information is compressed by representing the columns as a linear combina- tion of these new basis functions. For example, the reconstruction shown in

(6)

1 Introduction 3

(a) 1 (b) 2 (c) 3

(d) 10 (e) 50 (f) original image

Figure 1.1: Illustration of POD method using image compression. The original image is of resolution 600×800 shown in figure (f ). Figures (a)-(e) show the different reconstructions of the original image using 1, 2, 3, 10 and 50 basis functions, respectively.

figure 1.1-(a) uses only the first basis function (600×1 numbers). Together with the weights for each column (800×1 numbers), the corresponding compressed image uses only 0.29% of the original memory storage. Figures 1.1-(b) to 1.1-(e) show different reconstructions of the original image with increasing number of basis functions. The last one with 50 basis functions seems to have already captured the most essential energy of the original image, which still corresponds to only 14.58% of the original memory storage.

1.2 Thesis Organization

In chapter 2, we motivate the choice of jump-diffusion models in option pricing by generalizing the Black-Scholes model. The Black-Scholes model suffers from a few shortcomings so that option prices computed by the model are not consistent with the option prices quoted in the major financial markets. In particular, we point out two well-documented empirical bias of the Black-Scholes model: stochastic jumps underlying the random risky price process and the implied volatility smile. By overlaying the exponential Brownian motion with stochastic jumps, jump-diffusion models have bee proven to be more realistic

(7)

and flexible to represent price dynamics. On the way to pricing options, we apply the classical risk-neutral pricing approach and conclude that option prices in jump-diffusion models can also be computed in a deterministic way through partial integro-differential equations (PIDEs), whose solutions typically require efficient numerical techniques. Without going into details, we also state the calibration problem in jump-diffusion models, which can be seen as a PIDE constrained optimization problem, emphasizing in turn that the design of fast numerical solution of PIDEs is of great importance.

In chapter 3, we first derive the pricing PIDE in jump-diffusion models with constant volatility using the standard arguments of stochastic analysis and show that the result should also hold for general Lévy processes. Trans- lating the PIDE into a local volatility variant , we derive a numerical schema to solve the PIDE in detail using the affine Galerkin finite element method while the schema is illustrated using the Merton’s jump-diffusion model. In particular, the variational formulation involves the weighted Sobolev space and further we show that the variation formulation has a unique solution within an abstract parabolic problem setting, for which fundamental results are already available in the literature. The well-known non-local integral term due to jumps is explicitly addressed not only on the localization stage but also on the numerical computation stage, leaving us with an inconvenient dense matrix after the spatial discretization. The complete schema is then implemented using MATLAB and numerical results are also presented.

In chapter 4, we speed up the Galerkin schema by applying one of the most popular model order reduction techniques, namely the proper orthogonal decomposition (POD) method we illustrated earlier. Still with the Galerkin approach, we replace the affine finite element basis functions with more advanced basis functions. If a set of information about the characteristics of the PIDE is already at hand, these advanced basis functions can be computed using the POD method and can be used as a reduced basis in the Galerkin schema.

We outline the details how to compute the POD basis functions, which are strongly connected to the singular value decomposition (SVD), and therefore efficient numerical computation is possible. With anticipation that the POD method will be more powerful when the PIDE has to be solved repeatedly, for instance, in a calibration model, we investigate the approximation error introduced by the reduced basis and give a priori theoretical error estimates. We test the theoretical results by implement the reduced Galerkin schema, called POD Galerkin schema, using MATLAB. Numerical results are presented.

(8)

1 Introduction 5

In chapter 5, we conclude the thesis by summarizing the observations to the option pricing problem in jump-diffusion models, and provide some recipes for future research.

(9)

2 Basic Theory in Option Pricing

The history of options can be traced back to as early as around 500 B.C. in the ancient Greece. However, for a long time, options were priced based on personal judgements because the market was highly illiquid. Option prices tended to be higher than nowadays quoted market prices. Since options are traded systematically in large exchanges, researchers started to price options using mathematical models. In 1900, Bachelier introduced the first mathematical model to price options, which was based on the arithmetic Brownian motion. Still, it is not until the influential works of Black and Scholes and of Merton in 1973 that investors could get a through understanding how options can be priced based on a mathematical formula. Their approach was a break- through not only in theory but also in practice and has put a huge influence on the way that the modern derivatives market works.

In this chapter, we will follow the modelling direction of Black and Scholes and discuss briefly the general pricing principle behind the Black-Scholes model, namely the “risk-neutral pricing” rule. Our focus will then lay on the Merton’s jump-diffusion model–an extension of the Black-Scholes model which includes stochastic jumps. Option prices can then be computed by solving partial differential equations or, in the case of jump-diffusion models, partial integro-differential equations. We start with a short introduction of the option markets.

2.1 Options and Markets

Options belong to the category of financial derivatives, whose value are derived from some prescribed and often risky assets, known as the underlying assets, or simplyunderlying. In this thesis, we consider the simplest financial options, European call options.

(10)

2 Basic Theory in Option Pricing 7

10 15 20 25 30 35 40 45 50 55 60

0 5 10 15 20 25 30

ST

Payoff

Figure 2.1: Payoff of a European call option with strike priceK =30.

Definition 2.1.1 A European call option is a contract under which its owner may buy a prescribed asset (underlying) at a prescribed time in the future (ma- turity) for a prescribed price (strike price).

Such a contract gives theholderof the option a right, not an obligation, to buy the underlying. However, the other party to the contract, known as thewriter of the option, does have a potential obligation. She must sell the underlying if the holder exercises the option by buying the underlying at the strike price. Since it is the writer that bears the risk associated to the option, she requires compensation for selling such an option, known asoption price.

The price has to be paid by the holder at the time of opening the contract.

Option pricing models are supposed to finding such a price.

At the maturity date, since the price of the underlying can be already read on the market, one knows the option value for sure. The option value at maturity date is called thepayoff of an option. The payoff function of a European call option with strike priceK is given by

H(ST) =max(ST −K, 0),

which is illustrated in Figure (2.1). S_T denotes the underlying price at the maturity dateT. Therefore, the holder of a European call option will exercise the option only ifS_T >K, obtaining the payoff of amountS_T −K.

Besides European call options, there are many other types of options allowing various conditions under which options can be exercised. In opposite to call options, options to sell assets are known asputoptions. While European options give the holders the right to make the transaction only at maturity,

(11)

there are options relaxing the restriction, known as exotic options. For example, American options allow the holders to exercise the options at any time at or before maturity. Moreover, the so-called path dependent options such as Asian, lookback and barrier options are also traded on the market. More details on options can be found, for example, in[17].

The resulted derivative market plays a vital role in the modern finance industry. It has two primary uses, namely speculation and hedging. Options provide with a cheap way of building substantial risk in portfolio strategies, which might produce so-called leveraged profit by taking some extreme po- sition on certain asset, known as speculation with options. In contrast, the subject of hedging goes to the other direction: risk reduction. To hedge potential risks, one takes advantage of the particular correlation between the underlying and the option price movements.

In theory, one can lock in an arbitrage return as follows: sell an option for more than it is worth and hedge away all the potential risk for the rest of the option’s life. This simple arbitrary trading strategy is in fact a foundation in the option pricing theory and has led to elegant financial models.

2.2 Option Pricing Modelling

In this section we deal with option pricing strategies in general and present some central results for pricing European options. The main idea is borrowed from the celebrated Black-Scholes option pricing model. More specifically, we first model the underlying prices by stochastic processes and then option prices are calculated based on the risk-neutral pricing principle. We illus- trate this approach by starting with the classical Black-Scholes model and jump-diffusion models can be viewed as generalization of the Black-Scholes model.

To my knowledge so far, almost all the existing financial models invoke certain assumptions on the financial market. For the classical option pricing theory, the following assumption should come as no surprise.

Assumption 2.2.1 On the financial market, there exists a deterministic risk-free rate r ≥0which delivers a guaranteed return such that money can be deposited and borrowed from a bank account freely on a continuously compounded basis.

That is, 1 currency unit deposited in the account at time t =0will be worth e^{r t} currency units at time t >0. Further, the market is frictionless, that is, there are

(12)

no transaction costs. There are no default risk and no dividend payments. All the market participants are rational such that the market is efficient, that is, there are no instantaneous arbitrage opportunities.

Throughout the entire thesis, we assume that the Assumption 2.2.1 holds.

For stochastic pricing models, it is conventional to characterize the uncer- tainty of the economy by a filtered probability space (Ω,F,(Ft)t≥0,P) with continuous timet ∈[0,T]. Here,Ωis the whole set of possible scenarios which is equipped with aσ-algebraF.(Ft)denotes the filtration generated by the σ-algebra up to time t and for 0≤s ≤t ≤T holds thatFs ⊆ Ft ⊆ FT ⊂ F. Pdenotes thephysical probability measure. Without loss of generality, we considert =0 to be our initial time (today).

2.2.1 A Basic Model for Asset Prices

In the theory of asset pricing, prices of (risky) assets have been often modelled by Markov processes since the influential hypothesis of efficient market was introduced. Basically speaking, asset prices in an efficient market should behave in a random pattern. First, the price today should havefullyreflected the public information in the past history. Second, asset prices change immediately to any new public information.

Consider a underlying whose price is modelled by a Markov process (S_t)t∈[0,T]. In this context, the filtration(Ft)we defined before can be then interpreted as the history of the asset price which contains all the information known up to timet. The most widely used model of asset price movements is given by the followingstochastic differential equation:

dS_t

S_t =µdt +σdX_t. (2.1)

whereµandσare assumed to be constants.¹ The notation d·is used for the small change in any quantity over the time interval of length∆t. By continuous models we mean models in which∆t →0.

Instead of modelling the absolute asset price movements, the equation (2.1) gives a model of the corresponding return, ^dS^t

St , which is a common financial indicator. Essentially, the model decomposes the return into two parts.

On the one hand, there is a predictable and deterministic return of amount

1 Bothµandσcan be functions ofS_t andt, for example, to account for empirical observations of asset price movements.

(13)

µdt, possibly generated by a risk-free bank account. The parameterµis the expected return of the asset. On the other hand, the termσdX_t models the random change in the asset price, possibly as a response to some unexpected new information on the market. Hereσis a variable called the volatility of the asset return withσ²being its variance. The quantity of dX_t is the sample from some specific probability distribution, for example, a normal distribution.

Since dX_t is the only term that contains the randomness, it is certainly of great importance to capture the typical feature of asset prices.

By the fundamental theorem of asset pricing, the condition that our market is arbitrage-free guarantees that there exists at least one risk-neutral probability measureQequivalent to the physical probability measureP.

Definition 2.2.2 Let(Y_t)t∈[0,T]be a stochastic process defined on the probability space(Ω,F,(Ft)t≥0,P). A probability measureQsatisfying E^Q[Yt]<∞for all t ∈[0,T]and

E^Q[Ys|Ft] =Y_t Q−a.s. for s ≥t ≥0 (2.2) is called a risk-neutral measure (or martingale measure). In addition, if

∀A∈ F Q(A) =0⇔P(A) =0

holds, thenQis called an equivalent risk-neutral measure (or equivalent mar- tingale measure). E^Q[·|Ft]denotes the conditional expectation operator with respect to the risk-neutral probability measureQ.

In particular, the no arbitrage condition indicates that the discounted price process must satisfy

E^Q

e^{−r T}S_T|F0

=S₀. (2.3)

The consideration leads us eventually to the so-calledrisk-neutral pricing rule, which is one of the most important results in the option pricing theory.

For more theoretical results on risk-neutral pricing we refer to[12]. The risk- neutral pricing is first applied in the Black-Scholes model.

2.2.2 The Black-Scholes Model

Consider a European option pricing problem in which the underlying price follows a continuous stochastic process of the form (2.1).

(14)

In 1973, BLACK and SCHOLES [4] and MERTON [25] proposed a option pricing model, which suggests that the underlying price movements follow a stochastic differential equation given by

dS_t

S_t =rdt +σdW_t. (2.4)

wherer is the risk-free return andW_t denotes a standard Brownian motion on the physical probability space.

Definition 2.2.3 A stochastic process W = (W_t)t≥0adapted to a filtration(Ft)t≥0

is called a standard Brownian motion on probability space(Ω,F,P) if the follow conditions are satisfied:

1. W₀=0P-almost surely;

2. W has independent increments: W_t −W_s is independent for all0≤s<t ; 3. W has stationary and normally distributed increments: W_t−W_s is nor-

mally distributed with mean zero and variance t −s for all0≤s <t . It is well-known that the solution of (2.4) is given by

S_t =S₀exp

r −σ² 2

t +σWt

. (2.5)

Stochastic processes given by (2.5) are often called exponential Brownian motions, or geometric Brownian motions.

Based on the exponential Brownian motion, the model ultimately yields an analytical formula for pricing the European option with the underlying priceS_t, known as the Black-Scholes formula. The option pricing model built a cornerstone in the pricing of European options. Since then option prices can be computed based on fundamental arguments of hedging and absence of arbitrage, which has become almost a standard paradigm in finance termed as “risk-neutral pricing”.

The derivation of the pricing formula in the Black-Scholes model can be achieved in many different ways. Here we briefly show the way through partial differential equations. Let a European call option with strike priceK and maturityT be given and its price at timet be denoted byC(t,S). By the

(15)

risk-neutral pricing rule, the option price can be calculated as the conditional expectation of the payoff under the risk-neutral probability measureQ:

C(t,S) =e^−r⁽^T^−t⁾E^Q[max(ST −K, 0)]. (2.6) Applying the Itô formula and standard dynamical hedging arguments, one can arrive at the following partial differential equation (PDE), often known as the Black-Scholes differential equation in the literature:

∂C

∂t (t,S) +r S∂C

∂S (t,S) +σ²S² 2

∂²C

∂S²(t,S)−r C_t =0. (2.7) The above partial differential equation has been well studied and is solvable under various terminal conditions, representing different types of options. For example, in the case of a European call option, the terminal condition is given byC(T,S) =max(S−K, 0)(cf. section 2.1). Solving (2.7) for European options yields the Black-Scholes formula. The formula for a European call option takes the form

C^{B S}(t,S) =S N(d1)−K e^−r⁽^T^−t⁾N(d2) (2.8) with

d₁=ln(S/K) + (r +σ²/2)(T −t) σp

T −t ,

d₂=d₁−σp T −t.

Note thatN(·)denotes the cumulative distribution function of the standardized normal distribution.

Being perhaps the most successful model in the history of option pricing, the Black-Scholes model remains its popularity over the last decades. However, it is built on critical assumptions, leading to serious pricing problems.

One of the critical points in the Black-Scholes derivation is that the underlying asset price dynamics must have a continuous sample path, or Brownian motion. Presence of any possible discontinuity, such as stochastic jumps, makes the Black-Scholes model invalid. By examining time series plots of asset prices, large price changes do happen in financial assets such as stocks and currencies, possibly due to some announcement or financial crisis. More- over, for a more subtle jump structure which might be invisible for “untrained”

eyes, CARR and WU [7]develop more advanced methodology and find evidence that there are jump components in the U.S. stock market, especially for short-maturity options.

(16)

Another well-documented problem is the implied volatility smile. Ever since large amount of option price data is available on the market, traders start to compute the implied volatility by equalling the Black-Scholes formula to the observed market prices. According to the Black-Scholes formula (2.8), such implied volatility should be a flat line if plotted against different strike prices.

However, the observed implied volatility curves against the strike prices tend to be downward sloping, which look like “smile” curves. This phenomenon is documented in finance as the volatility smile. The presence of volatility smiles is more obvious for options with short maturities. When time to maturities increases, the smiles tend to flatten out.

The volatility smile suggests that the risk-neutral distribution of asset returns is left-skewed instead of normal and consequently out-of-the-money options are underpriced by the Black-Scholes formula.² In extreme market conditions, investors are willing to pay extra premia for out-of-the-money options which in effect provide insurance against market crashes. Further evidence shows that put options are underpriced by the Black-Scholes formula more severely than call options, indicating that investors do expect jumps on the market, particularly, big down jumps. We refer to[18]for more discussions on the volatility smile.

With these phenomena in mind, studies point to various generalizations of the Black-Scholes model. One possible alternative is to overlay continuous Brownian motion with jumps. Such models are becoming powerful and fashionable in the option pricing theory because: (1) the richness of jump structure can lead to a variety of volatility smile patterns; (2) the models re- tain the mathematical tractability of Brownian motions and the computing complexity can be controlled in a similar level of evaluating the Black-Scholes formula, at least, for European options.

2.2.3 The Merton’s Jump-Diffusion Model

The first jump-diffusion model was introduced by MERTON[24]in 1976. The risk-neutral dynamics of the underlying asset price is modelled by a so-called exponential jump-diffusionprocess given by

S_t =S₀exp

r t +Y_t , (2.9)

2 While an out-of-the-money call option is referred to the option whose strike price is above the current price of the underlying, an out-of-the-money put option is the option whose strike price is below the current price of the underlying.

(17)

(a) Compounded Poisson Process

(b) Exponential Brownian Motion vs. Exponential Jump−Diffusion Process

Brownian Motion Jump−Diffusion

Figure 2.2: Simulation examples of compound Poisson process, exponential Brownian motion and exponential jump-diffusion process, respectively.

whereY_t is a time-homogeneous jump-diffusion process

Y_t =−σ²

2 t +σW_t +

Nt

X

i=1

Z_i. (2.10)

N_t is a Poisson process with intensityλgenerating the jumps andZ_i is independent, identically distributed random variable which controls the jump size.

In the Merton’s model,Z_i is assumed to be normally distributed with density function

f(z) = 1 p2πσJ

exp

−(z −µJ)² 2σ²_J

. (2.11)

Note that the compound Poisson processPN_t

i=1Z_i with intensityλallows only finite number of jumps. Figure (2.2) shows typical trajectories of a compound Poisson process, an exponential Brownian motion and an exponential jump- diffusion process. Between jumps, the price process evolves like an exponential Brownian motion, and after each jumps, the price is multiplied bye^Zⁱ.

(18)

Remark 2.2.4 Given the jump-diffusion model specified by (2.9) and (2.11), price of a European option can be computed as a series of Black-Scholes prices.

For a European call option with strike price K and maturity T , its price at time t =0with current underlying price S₀is given by

C^M(0,S₀) = X∞

n=0

e⁻^λ^˜^T(λT˜ )ⁿ

n! C^{B S}(0,S₀;r_n,σn), (2.12) whereλ˜ =λe^q,σ²_n =σ²+nσ²_J/T , rn =r −λ(e^q−1) +^{n q}_T , C^{B S}(·,·;r,σ)is the Black-Scholes option price with volatilityσand risk-free rate r and q=µJ+¹₂σ²_J.

In contrast with the Black-Scholes model, pricing formula like (2.12) is hard to be established in jump-type models based on hedging arguments.³ Due to the jumps, the standard Itô formula has to be adjusted. As jumps are not tradable on the market, jump-diffusion models imply an incomplete financial market. Therefore, developing hedging strategies needs more rigorous considerations. We refer to [9]for various hedging strategies designed for jump-type processes.

In case of jump-type models, we will arrive at a so-called partial integro- differential equation (PIDE), which describes the option price over time. In contrast to a regular partial differential equation, a PIDE involves an additional integral term in the equation, which leads to difficulties in pricing options.

Often, solving such pricing problems with PIDEs can not be done analytically and requires numerical methods.

Under some regularity conditions, Fourier transform methods are proven to be efficient, for applications, see[6, 9, 21]. More generally, numerical solutions of such PIDEs can be obtained using the finite difference method (FDM) [2, 10]and the finite element method (FEM)[23, 29]. An introduction of the FDM and FEM for option pricing can be found in[1].

2.3 Model Calibration in Option Pricing

As one might have known, practical computation of option prices relies heavily on the estimation of model parameters. For jump-diffusion models introduced earlier, parameters such as the volatilityσ, the jump intensityλand the density

3 Derivation of the formula (2.12) can be found in[24], however, it is done by verifying the formula is indeed a solution of the PIDE type pricing equation.

(19)

functionf(z), have to estimated based on the current market information.⁴ Identifying such model parameters to reproduce the traded option prices on the market is often known asmodel calibration.

Specifically, suppose there are call option prices quoted on the market at timet =0 (which is often the case in practice), the model calibration problem can be formally formulated as follows: given market pricesC_ifor a set of bench- mark call options with different strikesK_i and maturitiesT_i (i ∈ {1, 2,· · ·,M}), find a risk-neutral modelQ, i.e., the model parameters under the risk-neutral probability measure, which prices these options correctly:

C_i=e^{−r T}E^Q

max(ST_i−K_i, 0)

∀i∈ {1, 2,· · ·,M}. (2.13) Once a risk-neutral models is calibrated, the fitted model might then be used to price some exotic, illiquid or over-the-count (OTC) options.⁵

If option prices are computed based on jump-diffusion models, we obtain a constrained calibration problem as a result. Often, such a problem does not have a solution which exactly reproduces quoted option prices on the market.

Moreover, since market option prices are typically given up to bid-ask spreads, exact calibration makes little sense. Thus, the essence of a calibration problem is to find the best approximation of market option prices based on the option pricing model at hand, which can be achieved, for instance, by minimizing the quadratic error between model prices and market prices. Practically, one has to solve the following (constrained) optimization problem: given a risk-neutral jump-diffusion model and the quoted pricesC_i of call options for maturities T_i and strikesK_i (i∈ {1, 2,· · ·,M}), find

θ^∗=arg min XM

i=1

ωi|C^θ(T_i,K_i)−C_i|², (2.14)

whereC^θ(·,·)is the model price computed using the jump-diffusion model with parameterθ ={σ,λ,f(·)}andωi (i =1, 2,· · ·,M)denotes some appropriate weights.

We will show that the model priceC^θ can be formulated as solution of PIDEs in case of jump-diffusion models (cf. chapter 3). As a result, the calibration problem is nothing but a PIDE constrained optimization problem. Solving

4 The risk free rater can be deduced from certain market information independent of the underlying, for instance, the bond prices.

5 There is no guarantee that all such options can be priced through a deterministic P(I)DE approach. In some cases, other techniques such as Monte-Carlo-Simulation have to be applied.

(20)

such an optimization problem needs some advanced technical consideration and is beyond the scope of this thesis. We refer to[2, 9, 29]for more details.

Nevertheless, it is important to note that the proposed POD Galerkin schema should be viewed as effort to solve the calibration problem efficiently.

(21)

3 Affine Galerkin Pricing of Options in Jump-diffusion Models

This chapter deals with option pricing problem in the Merton’s jump-diffusion model based on the PIDE approach. We start with a derivation of the pricing PIDE in jump-diffusion models using classical stochastic analysis, while the jump-diffusion processes are recharacterized as a special case under the general Lévy class, i.e., we consider jump-diffusion processes as Lévy processes which take the form of (2.10).

Equipped with a local volatility term, the PIDE plays also an important role for the calibration problem. Focusing on the European vanilla options, we obtain a log-moneyness version of the PIDE which can be then posed on the whole real line. As the main part of the chapter, numerical analysis based on the PIDE is carried out in detail. In particular, we develop an appropriate variational formulation and show that it admits a unique solution in an abstract parabolic functional setting. Due to the non-local integral term, the artificial boundary conditions have to be taken account beyond the truncated domain, leading to a dense matrix by applying a spatial discretization. We propose a pricing schema using the Galerkin finite element method while the time domain is discretized by the implicit Crank-Nicolson method.

Note that the efficiency of the proposed pricing schema will be addressed later in chapter 4 and here we are more interested in the numerical perfor- mance of the pricing schema, which is investigated through our numerical experiment at the end of this chapter. In particular, option prices can be obtained in general jump models using Fourier techniques if the characteristic function of the process is available and numerically can be computed very efficiently. In comparison with the fast Fourier methods, the Galerkin FEM is more general and can be used for exotic options such as American and barrier

(22)

3 Affine Galerkin Pricing of Options in Jump-diffusion Models 19

options[22]. Moreover, combined with efficient optimization algorithms, it is also well-suited in the context of calibration models[29].

Let us comment briefly on the characterization of jump-diffusion models in the Lévy class, which is a common way to study jump-diffusion models in the literature. First, the characterization of jump-diffusion processes as a special case of Lévy process is mathematically advantageous and, more importantly, demonstrates the flexibility of Lévy models in financial modelling, which has lead to some elegant models, particularly, for option pricing. Second, the Galerkin pricing recipe turns out to be not restrictive to the jump-diffusion models. The PIDE representation for option pricing can be easily extended to more general Lévy models and the complexity of the corresponding numerical analysis based on Lévy models does not seem to be explosive and has been already addressed by some recent works[23].

3.1 Exponential Lévy Processes

In this section, we recall some elementary results of exponential Lévy models.

The results are stated without proof. For readers interested in more details we refer to[3].

Under exponential Lévy models, the risk-neutral dynamics of the asset priceS_t is modelled by the exponential of a Lévy process:

S_t =S₀exp(r t +Y_t),

wherer is the constant risk-free rate andY_t follows a Lévy process.

Definition 3.1.1 A cadlag (right-continuous with left limits) stochastic process Y = (Y_t)t≥0with starting value equal to zero is said to be a Lévy process if it has the following properties:

1. Independent increments: for every 0 ≤ t₀ ≤ t₁ ≤ t₂ ≤ · · · ≤ t_n ≤ t , the random variables Y_t₀,Y_t₁−Y_t₀, . . . ,Y_t_n−Y_t_n−1are independent;

2. Stationary increments: the distribution of Y_t_+h−Y_t does not depend on t ; 3. Stochastic continuity:lim_h→0P(|Y_t₊_h−Y_t| ≥") =0 ∀".

The most special example of Lévy processes is the Brownian motion, which has continuous path. All the other Lévy processes have discontinuous

(23)

paths with presence of jumps, i.e., the continuity of the paths is weakened to stochastic continuity, allowing for jumps at random time. To describe the jump structure, it is convenient to associate a Poisson random measure J_Y to the Lévy process:

J_Y(ω;t,B) =X

s≤t

1_B(∆Ys(ω)) ∀B ⊂R, (3.1) which counts the number of jumps occurring between time 0 tot with jump size belongs to a measurable setB. The intensity ofJ_Y is then given by

ν(B) =E[J_Y(1,B)], (3.2) which is defined onR\{0}and known as theLévy measure. The Lévy measure tells the expected number of jumps during each unit time whose size belongs toB and verifies

Z

|y|≤1

|y|²ν(dy)<∞ and Z

|y|>1

ν(dy)<∞. (3.3)

In general,νis not a finite measure and, as a result, the integralR

ν(dx)is not well-defined.

Further, the distribution of the increments ofY is not fully specified by the definition, in particular, it need not be the normal distribution. In most cases, the distribution of increments for Lévy processes does not have any explicit probability density function and is rather expressed by its characteristic function, which can be viewed as the inverse Fourier transform of the corresponding density function.

By the Lévy-Khintchine representation, the characteristic function ofY reads :

φ(z) =E e^{i z Y}^t

=e^t^ψ(^z⁾ (z ∈R), ψ(z) =iγz −1

2σ²z²+ Z

R

(e^{i z x}−1−i z x1_|x|≤1)ν(dx),

where γ, σ ≥ 0 are constants andν the Lévy measure given by (3.2). The triplet(ν,γ,σ²)is known as thecharacteristic tripletof the Lévy processY and determines uniquely its distribution, which can be viewed as an implication of the following important result, called the Lévy-Itô decomposition theorem.

(24)

Theorem 3.1.2 Let(Yt)t≥0be a Lévy process andνits Lévy measure. Then there exists aγand a Brownian motion(Bt)t≥0with varianceσ² such that

Y_t =γt +B_t +Y_t^l+lim

"→0Y˜_t^", (3.4)

where

Y_t^l = Z

|y|≥1,s∈[0,t]

y J_Y(ds×dy), Y˜_t^"=

Z

"≤|y|<1,s∈[0,t]

y{J_Y(ds×dy)−ν(dy)ds}

and J_Y is a Poisson random measure with intensityν(dy)dt .

Notice that, without the term ˜Y_t^"in (3.4), Lévy processes correspond to jump-diffusion processes. Since the measureν(dy)might approach infinity for dy being small, ˜Y_t^"could be viewed as an infinite superposition of independent small jumps, generated by the centred Poisson processes. Therefore, general Lévy processes may have infinitely many jumps in each time interval, while jump-diffusion processes allow only a finite number of jumps.

Generally speaking, the Lévy-Itô decomposition theorem indicates that each Lévy process can be decomposed into a Brownian motion with drift and a possibly infinite sum of independent compound Poisson processes. In other words, each Lévy process can be approximated with arbitrary precision by jump-diffusion processes, a property which can be useful in practice as well as in theory.

For financial modelling, we have learnt that the discounted price process must form a martingale under the risk-neutral probability. In case of exponential Lévy models, it requirese^Y^t to be a martingale, which can be shown to be the case if and only if

Z

|y|>1

e^yν(dy)<∞, γ+σ² 2 +

Z

R

(e^y−1−y1_|y_|≤1)ν(dy) =0. (3.5)

As mentioned earlier, the standard Itô formula has to adjusted to account for the jump activities. We now state the Itô formula for Lévy process, whose application in the option pricing theory is fundamental.

(25)

Theorem 3.1.3 Let(Yt)t≥0be a Lévy process with the characteristic triplet(ν,γ,σ²).

Then for any C^1,2function f :[0,T]×R→R, the process f(t,Y_t)can be written as

df(t,Y_t) = ∂ f

∂t (t,Y_t₋)dt +∂ f

∂y (t,Y_t₋)dY_t+σ² 2

∂²f

∂ y²(t,Y_t₋)dt +

f(t,Y_t₋+∆Yt)−f(t,Y_t₋)−∆Yt

∂ f

∂y (t,Y_t₋)

,

where∆Yt =Y_t −Y_t₋denotes the ‘jump’ of Y at time t .

A simple application of the Itô formula to the asset priceS_t =S₀exp(r t +Y_t) (cf. section 2,2,3) yields an explicit representation of the risk-neutral price dynamics under the exponential Lévy model, i.e.,

dS_t

S_t₋ = (r +σ²

2 )dt +dY_t + (e^∆^Y^t−1−∆Yt), (3.6) a result which we will need later to derive the pricing PIDE.

3.2 Partial Integro-differential Equations for Option Pricing

If the underlying price process follows the exponential jump-diffusion process stated by (2.9), we show that there exits a PIDE representation for the option pricing problem, which describes the price of options over time.

Let(ν,γ,σ)be the Lévy characteristic triplet of the jump-diffusion process Y whereσis now strictly positive. The Lévy measure for the jump-diffusion process (2.9) takes the form

ν(dy) =λf(dy). (3.7)

That is, jumps occur with intensityλand the jump sizes are normally distributed with the density function f(·)given by (2.11).

(26)

Remark 3.2.1 Since the normal density function f(·)decays exponentially with the term e^−y², the Lévy measureνin (3.7) satisfies the following two conditions:

Z

|y|>1

e^2yν(dy)<∞, (3.8)

Z

R

e^y^+µ|y^|y f(y)dy <∞. (3.9) The condition (3.8) indicates that the exponential jump-diffusion process is always L²(R)-integrable, in particular, the process e^Y will be a square integrable martingale, which turns out to be useful in the proof of the result below. The condition (3.9) will be of importance in the numerical part discussed later on.

Consider a European call option with maturityT and strike priceK. Its payoff function is give byH(ST) =max(ST −K, 0), which is obviously Lipschitz continuous:

|H(x)−H(y)| ≤c|x −y|, ∀x,y ≥0 (3.10) with the Lipschitz constantc =1. Apply the risk-neutral pricing rule to compute the option valueC =C(t,S_t)at timet as

C(t,S) =e^−r⁽^T^−t⁾E^Q[H(S_T)|S_t =S] =e^−r⁽^T^−t⁾E^Q

H(S e^r⁽^T^−t⁾⁺^Y^T^−t) , with the risk-neutral probability measureQ. For simplicity,Qwill be omitted in the sequel. Applying the Itô formula for Lévy process and using the property of martingales, we obtain the following result.

Theorem 3.2.2 Let an exponential Lévy model S_t =S₀exp(r t +Y_t)be given, for which (3.8) is satisfied. Ifσ >0or

∃δ∈(0, 2) such that lim

ε→0infε^−δ Z ε

−ε

|x|²dν(x)>0, (3.11)

then the option value C(t,S)can be computed as C :[0,T]×[0,∞)→R

(t,S)7→C(t,S) =e^−r⁽^T^−t⁾E[H(ST)|S_t =S],

(27)

where C(·,·)∈C⁰([0,T]×[0,∞))∩C^1,2((0,T)×(0,∞))and satisfies the following backward partial integro-differential equation

∂C

∂t (t,S) +r S∂C

∂S(t,S) +σ²S² 2

∂²C

∂S²(t,S)−r C_t +

Z

R

ν(dx)

C(t,S e^x)−C(t,S)−S(e^x−1)∂C

∂S(t,S)

=0

(3.12)

with(t,S)∈[0,T)×(0,∞)and the terminal condition

C(T,S) =H(S) ∀S>0. (3.13) Proof. We restrict ourselves to the jump-diffusion case. First note that the smoothness ofC(t,S)inS is implied by the condition (3.11) and int can be shown by applying Fourier methods.

Notice that the discounted option price ˆC_t =e^{−r t}C(t,S_t)is by construc- tion a martingale. Given the smoothness conditions, we are allowed to apply the Itô formula (cf. Theorem 3.1.3) to ˆC_t to obtain

d ˆC_t =e^{−r t}

−r C_t +∂C

∂t (t,S_t₋) +σ²S_t² 2

∂²C

∂S²(t,S_t₋)

dt +∂C

∂S(t,S_t₋)dS_t

+e^{−r t}

C(t,S_t₋e^∆Y^t)−C(t,S_t₋)−S_t₋(e^∆Y^t −1)∂C

∂S (t,S_t₋)

.

ForS_t following an exponential jump-diffusion process, we use (3.6) with processY_t given by (2.10). Then, it follows that

d ˆC_t =α(t)dt +dM_t, (3.14) where

α(t) =e^{−r t}

−r C_t +∂C

∂t (t,S_t₋) +σ²S_t² 2

∂²C

∂S²(t,S_t₋) +r S_t₋∂C

∂S (t,S_t₋)

+e^{−r t} Z

R

ν(dy)

C(t,S_t₋e^y)−C(t,S_t₋)−S_t₋(e^y −1)∂C

∂S (t,S_t₋)

, dM_t =e^{−r t}∂C

∂t (t,S_t₋)σS_t₋dW_t +

Z

R

e^{−r t}

C(t,S_t₋e^y)−C(t,S_t₋)J˜_Y(dt ×dy)

and the compensated Poisson measure ˜J_Y(dy) = J_Y(dy)−ν(dy).

(28)

Now with conditions (3.8) and (3.10), it is able to show thatM_t is actually a martingale. See, for example,[9, Chapt. 12]. The fact that ˆC_t is a martingale thus implies that the drift termα(t)has to disappear, since ˆC_t−M_t =Rt

0 α(t)dt is continuous with finite variation. That is,α(t) =0 holds almost surely, from which the PIDE (3.12) follows.

Given (3.12), model parameters have to be calibrated in order to compute the option prices. Technically, it can be done by solving the pricing equation (3.12) repeatedly for each possible pair(T,K)and then adjusting the parameters in such a way that the error between the model prices and the quoted market prices is minimized. Formulated as an optimization problem, the evaluation of such error itself would involve solving a whole family of equations like (3.12), which could be a challenge in terms of computational cost (cf.

section 2.3).

In practice, one has the possibility to reparametrize the pricing equation (3.12) and may arrive at the local volatility model proposed by DUPIRE[13] in 1994, which treats volatility as a function of timet and asset priceS, i.e., σ(t,S)instead of a fixedσ. Indeed, there exists an equivalent variant of (3.12) based on the local volatility model by PIRONNEAU[27].

Proposition 3.2.3 Let C(t,S;T,K) be the solution of (3.12) with maturity T and strike price K . Then C(0,S₀;T,K) =w(T,K), where w solves the forward PIDE

∂w

∂ τ(τ,K) +r K ∂w

∂K (τ,K)−σ(τ,K)²K² 2

∂²w

∂K²(τ,K)

− Z

R

ν(dy)e^y

w(τ,K e⁻^y)−w(τ,K) +K(e^y−1)∂w

∂K (τ,K)

=0 (3.15) with(τ,K)∈(0,T)×(0,∞)and initial condition

w(0,K) =H(S0). (3.16) The above formulation makes it possible to compute the values in timet =0 of all the options (with the same underlying) spanning different maturities and strike prices by solving (3.15) only once. For this reason, model calibration problems are usually formulated using the local volatility models.

Further, from a numerical point of view, we introduce the so-called log- moneyness version of equation (3.15), which is posed on the whole real line

(29)

with noK⁰s in the coefficients. The termlog-moneynessof a given option with current underlying priceS₀ is frequently used in option theory and defined as ln(_S^K

0)while the quotient ^K_S

0 is known as themoneynessof the option.

Proposition 3.2.4 Letη:=R

R(e^y −1)f(y)dy . The log-moneyness version of PIDE is given by

∂w

∂ τ(τ,x)−σ²(τ,x) 2

∂²w

∂x²(τ,x) + (r +1

2σ²(τ,x)−λη)∂w

∂x (τ,x) +λ(1+η)w(τ,x)−λ

Z

R

w(τ,x −y)e^yf(y)dy =0

(3.17)

with(τ,x)∈(0,T)×Rand initial condition

w(0,x) =w₀(x) =max(1−e^x, 0). (3.18) Then

wˆ(τ,K) =S₀w(τ, ln(K S₀)) solves (3.15).

Proof. First note that the Lévy measure for jump-diffusion processesν(dy) = λf(dy). Performing the variable transformationx =ln(_S^K

0), we get

∂w

∂K (τ,K) =S₀K⁻¹∂w

∂x (τ,x),

∂²w

∂K²(τ,K) =−S₀K⁻²∂w

∂x (τ,x) +S₀K⁻²∂²w

∂x²(τ,x).

In addition, we have^∂_{∂ τ}^w^ˆ(τ,x) =S₀^∂_{∂ τ}^w(τ,x)and for call optionsH(S0) =S₀max(1− e^x, 0). Using the property of the density functionf(·), we have the desired result after scaling both the PIDE and the initial condition byS₀⁻¹.

3.3 Variational Formulation

Our option pricing strategy is based on the numerical solution of the PIDE (3.17). The proposed and implemented numerical method relies on a variational form of the PIDE (3.17). Although numerical analysis in the sequel is mainly based on European call options, we aim at more complex options such as American and Barrier options, which are treatable using similar methodology[9].