Transcription
Applied Mathematics and Computation 293 (2017) 461–479Contents lists available at ScienceDirectApplied Mathematics and Computationjournal homepage: www.elsevier.com/locate/amcOn a one time-step Monte Carlo simulation approach of theSABR model: Application to European optionsÁlvaro Leitao a,b, , Lech A. Grzelak c,a, Cornelis W. Oosterlee b,aaDelft Institute of Applied Mathematics, TU Delft, Delft, The NetherlandsCWI-Centrum Wiskunde & Informatica, Amsterdam, The NetherlandscQuantitative Analytics, ING, Amsterdam, The Netherlandsba r t i c l ei n f oKeywords:Computational financeStochastic-local volatility modelsSABR modelCopulasa b s t r a c tIn this work, we propose a one time-step Monte Carlo method for the SABR model. Webase our approach on an accurate approximation of the cumulative distribution function ofthe time-integrated variance (conditional on the SABR volatility), using Fourier techniquesand a copula. Resulting is a fast simulation algorithm which can be employed to priceEuropean options under the SABR dynamics. Our approach can thus be seen as an alternative to Hagan’s analytic formula for short maturities that may be employed for modelcalibration purposes. 2016 Elsevier Inc. All rights reserved.1. IntroductionThe Stochastic Alpha Beta Rho (SABR) model [1] is an established SDE system which is often used for interest rates and FXmodeling in practice. The model belongs to the so-called stochastic local volatility (SLV) models. The idea behind SLV modelsis that the modeling of volatility is partly done by a local volatility and partly by a stochastic volatility contribution, aimingto preserve the advantages and minimize the disadvantages of the individual models.In the original paper [1], the authors have provided a closed-form approximation formula for the implied volatility under SABR dynamics. This is important for practitioners, as it can be used highly efficiently within the calibration exercise.However, the closed-form expression is derived by perturbation theory and its applicability is thus not general. The formulais for example not always accurate for small strike values or for long time to maturity options. In [2], Oblój provided acorrection to the Hagan formula but the same drawbacks remained. In order to improve the approximation for the SABRimplied volatility, different techniques have been proposed in the literature, based on either perturbation theory [3,4], heatkernel expansions [5,6] or model map techniques [7].The objective of a bank is typically to employ the SABR model for pricing exotic derivatives. Before this pricing can takeplace, the model needs to be calibrated to European options, and, if possible, to some exotic options as well. It is marketpractice to use Hagan’s formula to calibrate the SABR model, however, since it is only an approximation it will not alwaysresult in a sufficient fit to the market. In other words, a SABR Monte Carlo simulation may give different implied volatilitiesthan the Hagan formula.In this work, we aim to develop a one time-step Monte Carlo method, which is accurate for European derivative contractsup to two years. This kind of contract is often traded in foreign-exchange (FX) markets. The present approach may thus Corresponding author at: Delft Institute of Applied Mathematics, TU Delft, Delft, The Netherlands.E-mail addresses: A.LeitaoRodriguez@tudelft.nl (Á. Leitao), L.A.Grzelak@tudelft.nl (L.A. Grzelak), c.w.oosterlee@cwi.nl (C.W. .0300 096-30 03/ 2016 Elsevier Inc. All rights reserved.
462Á. Leitao et al. / Applied Mathematics and Computation 293 (2017) 461–479be used to improve the calibration of the SABR model. In particular, we focus on the efficient simulation of SABR’s timeintegrated variance, conditional on the volatility dynamics, as this conditional distribution is not available in closed-form.We approximate the distribution by the two marginal distributions involved, i.e. the volatility and time-integrated variancedistributions, by using a copula. The copula methodology has been used intensively in the financial world, see, for example,[8], in risk management and option pricing. In the field of risk management, some relevant works are Li [9], Embrechtset al. [10], Cherubini and Luciano [11] or Rosenberg and Schuermann [12]. In the case of pricing multi-asset options, someexamples of the use of copulas are Cherubini and Luciano [13] or Rosenberg [14].The approximation of the cumulative distribution function (CDF) of the time-integrated variance is obtained by applying arecursive procedure, as described in [15], which was originally developed to price arithmetic Asian options. The derivationis based on the definition of the characteristic function of the time-integrated variance and the use of Fourier inversiontechniques. We adapt the procedure to our problem and improve it in terms of computational costs. Once both marginaldistributions are available, we employ the copula technique to get a multivariate distribution which approximates the conditional distribution of the time-integrated variance given the volatility. Several copulas are considered and an analysis of theirperformance is carried out. The conditional time-integrated variance approximation will be employed in the asset simulationfor the complete SABR model.The paper is organized as follows. In Section 2, the SABR model and the simulation steps are briefly introduced.Section 3 describes the procedure to derive the CDF of SABR’s time-integrated variance. The copula approach to approximatethe conditional time-integrated variance distribution is explained in Section 4. In Section 5, some numerical experiments arepresented. We conclude in Section 6.The one-step methodology restricts the use of the proposed method to option maturities up to two years and Europeantype options. In a forthcoming paper, we will generalize the methodology to the multiple time-step case.2. The SABR modelThe SABR model is based on a parametric local volatility component in terms of a model parameter, β . The formaldefinition of the SABR model readsdS(t ) σ (t )Sβ (t )dWS (t ),S(0 ) S0 exp (rT ),σ ( 0 ) σ0 .dσ (t ) ασ (t )dWσ (t ),(1)where S(t ) S̄ (t ) exp (r (T t ) ) denotes the forward value of the underlying asset S̄ (t ), with r the interest rate, S0 the spotprice and T the contract’s final time. Quantity σ (t) denotes the stochastic volatility, with σ (0 ) σ0 , WS (t) and Wσ (t) are twocorrelated Brownian motions with constant correlation coefficient ρ (i.e. WSWσ ρ t). The open parameters of the model areα 0 (the volatility of the volatility), 0 β 1 (the elasticity) and ρ (the correlation coefficient).A useful aspect for practitioners is that each model parameter can be assigned to a specific feature of the market observed implied volatility when it is represented against the strike prices, commonly known as a volatility smile or skew. Theα parameter mainly affects the smile curvature. Exponent β is connected to the forward process distribution, being normalfor β 0 and log-normal for β 1. In terms of the volatility smile, β also affects the curvature. Correlation parameter ρrotates the implied volatility curve around the at-the-money point, obtaining “more pronounced smile” or “more pronouncedskew” shapes.2.1. SABR Monte Carlo simulationIn the first equation of the SABR model (1), the forward asset dynamics are defined as a constant elasticity of variance(CEV) process [16]. Based on the works of Schroder [17] and Islah [18], an analytic expression for the CDF of the SABRconditional process has been obtained. For some S(0) 0, the conditional CDF of S(t) with an absorbing boundary at S(t ) t0, and given the volatility, σ (t), and the conditional time-integrated variance, 0 σ 2 (s )ds σ (t ), reads Pr S(t ) X S(0 ) 0, σ (t ),wherea 1ν (t )b 2 χ 2 (x; S ( 0 ) 1 β (1 β ) 0t σ 2 (s )ds 1 χ 2 (a; b, c ), 2ρσ(t) σ(0)() ,α1 2β ρ 2 ( 1 β ),(1 β )(1 ρ 2 )(2)X 2 ( 1 β ),(1 β )2 ν (t ) tν (t ) (1 ρ 2 ) σ 2 (s )ds,c 0andδ , λ) is the non-central chi-square CDF.It should be noted that this formula is exact for the case ρ 0 and results in an approximation otherwise. A shiftedprocess with an approximated initial condition is employed in the derivation for ρ 0. Based on Eq. (2), an “exact” MonteCarlo simulation scheme for the SABR model can be defined based on inverting the conditional SABR CDF, when the conditional time-integrated variance and the volatility is already simulated. This forms the basis of the one time-step SABR
Á. Leitao et al. / Applied Mathematics and Computation 293 (2017) 461–479463approach, where exact simulation is required for accuracy reasons when using only one large time-step (meaning no time discretization). Furthermore, it is well known that the convergence ratio of the Monte Carlo method is 1/ n, with n thenumber of paths or scenarios.According to Eq. (2), in order to apply an “exact” one time-step Monte Carlo simulation for the SABR dynamics, we needto perform the following steps (with the terminal time T): Simulation of SABR’s volatility. From Eq. (1), we observe that the volatility process of the SABR model is governed by alog-normal distribution. As is well-known, the solution follows a geometric Brownian motion, i.e. the exact simulation ofthe volatility at terminal time, σ (T), reads12σ (T ) σ (0 ) exp(αWˆ σ (T ) α 2 T ),(3)ˆ σ (T ) is an independent Brownian motion.where W T Simulation of SABR’s time-integrated variance, conditional on the terminal value of the volatility, i.e., 0 σ 2 (s )ds σ (T ). Exactsimulation of the time-integrated variance in the SABR model is not directly possible since its conditional distributiongiven the volatility is not known in closed-form. One possibility is to simulate it by a Monte Carlo method, but thisimplies the use of a nested simulation which can be very expensive or even unaffordable. Another possibility is to find Tan approximation for the conditional distribution of SABR’s time-integrated variance given σ (T). The integral 0 σ 2 (s )dscan also be approximated by some quadrature rule once σ (0) and σ (T) (or intermediate values of the volatility process)are available. Simulation of SABR’s forward asset process. We can simulate the forward dynamics by inverting the CDF in Eq. (2). Forthis, the conditional SABR’s time-integrated variance (besides the volatility) is required. Concerning the forward assetsimulation, there is no analytic expression for the inverse SABR distribution so the inversion has to be calculated bymeans of some numerical approximation. Another possibility is to employ SDE discretization schemes, like Taylor-basedschemes (Euler, log-Euler, Milstein or full-truncation) or more involved ones (low-bias or QE schemes). However, the useof such schemes gives rise to several issues that have to be managed, like the discretization error, negativity in the assetvalues and computational cost.In Algorithm 1, a schematic representation of the complete one time-step Monte Carlo simulation procedure for the SABRmodel is presented. We include references to the sections of this paper where the respective parts of the SABR simulationare discussed.Algorithm 1: SABR’s one time-step Monte Carlo simulation.Data: S0 , σ0 , α , β , ρ , T , n.for Path p 1. . .n do ˆ σ (t ) N (0, T ).Simulation of WSimulation of σ (T ) σ (0 ) (see Eq. (3)). TSimulation of 0 σ 2 (s )ds σ (T ) (see Section 4.3). TSimulation of S(T ) S(0 ), σ (T ), 0 σ 2 (s )ds (see Section 4.4).In the above steps of the SABR model exact simulation, the challenging part is the simulation of the time-integrated Tvariance, 0 σ 2 (s )ds, conditional on σ (T) (and σ (0)). Sometimes moment-matching techniques can be used. Here, however,we propose a computationally efficient approximation, based on a copula multi-variate distribution, to simulate the condi Ttional distribution of 0 σ 2 (s )ds given the volatility σ (T). Simulation by means of a copula technique requires the CDF of Tthe involved marginal distributions. Although the distribution of σ (T) is known, the distribution of 0 σ 2 (s )ds needs to beapproximated. We propose a method based on a Fourier technique, which was already employed for Asian options in [15].In Section 3, this derivation is presented in detail. THereafter, for notational convenience, we will use Y (T ) : 0 σ 2 (s )ds.3. CDF of SABR’s time-integrated varianceIn this section, we present a procedure to approximate the CDF of the time-integrated variance, Y(T). Since we will workin the log-space, an approximation of the CDF of log Y(T), Flog Y(T) , will be derived. First of all, we approximate Y(T) by itsdiscrete analog, i.e. T0σ 2 ( s )d s M j 1σ 2 (t j ) t : Yˆ (T ),(4)
464Á. Leitao et al. / Applied Mathematics and Computation 293 (2017) 461–479where M is the number of discrete time points1 , i.e., t j j t, j 1, . . . , M and t to the logarithmic domain, where we aim to find an approximation of Flog Yˆ (T ) , i.e.Flog Yˆ (T ) (x ) x TM.Yˆ (T ) is subsequently transformedflog Yˆ (T ) (y )dy,(5)with flog Yˆ (T ) the probability density function (PDF) of log Yˆ (T ). flog Yˆ (T ) is, in turn, found by approximating the associatedcharacteristic function, φlog Yˆ (T ) , and applying a Fourier inversion procedure. The characteristic function and the desired PDFof log Yˆ (T ) form a so-called Fourier pair. Based on the works in [15,19], we develop a recursive procedure to recover thecharacteristic function of flog Yˆ (T ) . We start by defining the sequence, σ 2 (t j )R j logσ 2 (t j 1 ) log σ 2 (t j ) log σ 2 (t j 1 ) ,(6)where Rj is the logarithmic increment of σ 2 (t) between tj and t j 1 , j 1, . . . , M. As mentioned, in the SABR model contextthe volatility process follows log-normal dynamics. Since increments of Brownian motion are independent and identicallyddistributed, the Rj are also independent and identically distributed, i.e. R j R. In add
number of paths or scenarios. According to Eq. (2), in order to apply an “exact”one time-step Monte Carlo simulation for the SABR dynamics, we need to perform the following steps (with the terminal time T): Simulation of SABR’s volatility. From Eq. (1), we observe that the volatility process of the SABR model is governed by a
