Swing Options Valuation: a BSDE with Constrained Jumps Approach Marie Bernhart∗ Huyˆen Pham † Peter Tankov‡ Xavier Warin§ 1 1 January 6, 2011 0 2 n a J Abstract 5 We introduce a new probabilistic method for solving a class of impulse control ] P problemsbasedontheirrepresentationsasBackwardStochasticDifferentialEquations C (BSDEsforshort)withconstrainedjumps. Asanexample,ourmethodisusedforpric- ing Swing options. We dealwith the jump constraintby a penalization procedure and . n apply a discrete-time backward scheme to the resulting penalized BSDE with jumps. i f We study the convergenceof this numericalmethod, with respectto the main approx- - q imation parameters: the jump intensity λ, the penalization parameter p > 0 and the [ time step. In particular, we obtain a convergencerate of the error due to penalization 1 oforder(λp)α−21,∀α∈ 0,12 . Combiningthis approachwithMonte Carlotechniques, v we then work out the valuation problem of (normalized) Swing options in the Black (cid:0) (cid:1) 5 and Scholes framework. We present numerical tests and compare our results with a 7 classical iteration method. 9 0 KeywordsBackwardstochasticdifferentialequationswithconstrainedjumps,Impulse . 1 control problems, Swing options, Monte Carlo methods 0 1 1 : v i X r a ∗Laboratoire de Probabilit´es et Mod`eles Al´eatoires, Universit´es Paris 6-Paris 7, CNRS UMR 7599 and EDF R&D,92141 Clamart, France. Email: [email protected] †Laboratoire de Probabilit´es et Mod`eles Al´eatoires, Universit´es Paris 6-Paris 7, CNRS UMR 7599, France. Email: [email protected] ‡Centre de Math´ematiques Appliqu´ees, Ecole Polytechnique, 91128 Palaiseau, France. Email: [email protected] §EDF R&D, 92141 Clamart, France and Laboratoire de Finance des March´es de l’Energie, Universit´e Paris Dauphine. Email: [email protected] 1 1 Introduction In this report, we introduce a new probabilistic method for solving impulse control problemsbasedontheirrepresentationsasBackwardStochasticDifferentialEquations (BSDEs for short) with constrained jumps. As an example, our method is used for pricing Swing options in the Black and Scholes framework. BSDEs provide alternative characterizations of the solution to multiple-obstacle, optimal switching (see among others [19, 8, 20, 26, 13]) and more generally impulse control problems: Kharroubi et al. [21] recently introduced a family of BSDEs with constrained jumps providing a representation of the solution to such problems. A challenging question is that of the numerical approximation of this kind of BSDEs with constrained jumps. A discrete-time backward scheme for solving BSDEs with jumps (without con- straint)has been introducedby Bouchardand Elie [3]. Inour case, the maindifficulty comesfromthe constraint,whichconcernsthe jumpcomponentofthe solution. These BSDEsdonotaprioriinvolveanySkorohodtypeminimalitycondition. Inconsequence, classicalapproachesbyprojectedschemes(discretelyreflectedbackwardschemes)used for example by [2] and [10] cannot be used. We consider a penalization procedure to deal with the constraint on jumps and provide a convergencerate of the penalized solutionto the exact solution. This allows us to establish a convergence rate of the error between the solution of the considered impulse control problem and the numerical approximation given by the discrete-time solution to the penalized BSDE with jumps, as the penalization coefficient and the number of time steps go to infinity. The rest of the report is structured as follows: in Section 2, we set the considered impulse control problem in the mathematical framework of BSDEs with constrained jumps. We present in Section 3 our penalization approach and provide a global con- vergence rate of our approximation. In Section 4, our method is used for pricing multi-exercise options, so-called (normalized) Swing options. This multiple optimal stopping time problem leads to a particularly degenerate three-dimensional impulse controlproblem. We combineourBSDE-basedapproachwithMonteCarlotechniques and deal with Swing options with a small maximal number of exercises rights, due to large computational times. We compare our pricing results with those obtained by a classical iteration-based approach proposed for example by [9]. 2 BSDE Representation for Impulse Control Prob- lems Let T be a given time horizon. We work in a complete probability space (Ω,F,P), on which is defined a d-dimensional Brownian motion W and a Poisson process N with intensity λ>0. We denote by F=(F ) , the augmentation of the natural filtration t t≥0 generated by W and N, by FW = (FW) the one generated by W, and by P, the t t≥0 σ-algebra of predictable sub-sets of Ω×[0,T]. 2 Notation Throughout this report, the euclidean norm defined on Rd or on R will be indis- criminately denoted by |·|. The matrix transposition is denoted by ⊥. In addition, unlessspecifiedotherwise,C willdenoteastrictlypositiveconstantdependingonlyon Lipschitz constants of the coefficients of the problem, see assumptions (H) and (H′) below, and constants T, |b(0)|, |σ(0)|, |γ(0)|, |f(0)|, |κ(0)| and |g(0)|. Besides, we shall use the standard notations: • S2, the set of real-valued ca`dla`g adapted processes Y =(Y ) such that t 0≤t≤T 1 2 kYk := E sup |Y |2 <∞ . S2 t (cid:18) (cid:20)0≤t≤T (cid:21)(cid:19) • A2, the sub-set of S2 such that A2 := K ∈S2 :(K ) nondecreasing , K =0 . t 0≤t≤T 0 (cid:8) (cid:9) • L2([0,T]), the set of real-valued adapted processes (φ ) such that F t 0≤t≤T T E |φ |2dt <∞ . t "Z0 # • L2(W), the set of real-valued P-measurable processes Z =(Z ) such that t s≤t≤r 1 T 2 kZk := E |Z |2dt <∞ . L2(W) t "Z0 #! • L2(N), the set of real-valued P-measurable processes V =(V ) such that t s≤t≤r 1 T 2 kVk := E |V |2λdt <∞ . L2(N) t "Z0 #! • V denotesthesetofP-measurableessentiallyboundedprocesses,valuedin(0,∞) and Vp ={νp ∈V :νp ≤p a.s.}. t 2.1 A Class of Impulse Control Problems We consider the class of impulse control problems whose value function is defined by: T v(t,x)= sup E g(Xt,x,u)+ f(Xt,x,u)ds+ κ(Xt,x,u) . (1) T s τ− u=(τk)k≥1∈U(t,T] Zt t<Xkτ≥k≤1T k 3 An impulse strategy u=(τ ) is said to be admissible for problem (1) (and belongs k k≥1 to U ) if it is a non-decreasing sequence of FW-stopping times valued in (t,T] (we (t,T] set by convention τ =t) such that, if 0 nu :=♯{k ≥1:t<τ ≤T} (t,T] k denotes the (random) number of interventions of the strategy u before time T, then 2 E nu <C, (2) (t,T] (cid:12) (cid:12) for some universal constant C > 0. (cid:12)The c(cid:12)ontrolled state variable Xt,x,u is a ca`dla`g (cid:12) (cid:12) process with dynamics s s Xt,x,u =x+ b(Xt,x,u)dr+ σ(Xt,x,u)dW + γ(Xt,x,u), ∀s≥t. (3) s r r r τ− Zt Zt t<Xτk≤s k Between two successive intervention times τ and τ , the state variable evolves as a k k+1 diffusion process and the controller makes an integral profit f. At each decided inter- ventiontime τ ,hegivesanimpulseto thesystem: the stateprocessjumps withasize k Xu −Xu =γ(Xu ) and he obtains the intervention gain κ. τk τ− τ− k k We consider standard assumptions on the coefficients of the problem: (H) b:Rd 7→Rd, σ :Rd 7→Rd×d and γ :Rd 7→Rd are Lipschitz continuous and γ is uniformly bounded. f :Rd 7→R, κ:Rd 7→R and g :Rd 7→R are Lipschitz continous. (H′) The maps b, σ, γ, and g belong to C1(Rd) and have Lipschitz continuous b derivatives. A straightforwardcomputation using (H), (2) and Gronwall’s lemma shows that ∀(t,x)∈[0,T]×Rd, sup E sup Xt,x,u 2 <∞. (4) s u∈U(t,T] (cid:20)t≤s≤T (cid:21) (cid:12) (cid:12) Finally, we will assume the existence of an optimal(cid:12)strateg(cid:12)y u∗ = (τ∗) ∈ U k k≥1 (t,T] to problem (1). We refer for example to [4] and [25] in the infinite horizon case, for specific conditions on the coefficients of the problem which ensures such an existence. 2.2 Link to BSDEs with Constrained Jumps Let us consider the BSDE with constrained jumps T T T T Y =g(X )+ f(X )ds− Z dW − V dN + dK , ∀0≤t≤T t T t s t s s t s s t s (Vt+κ(Xt−)≤R0, ∀0≤t≤RT R R (5) where X is the (uncontrolled) jump diffusion process with dynamics dX =b(X )dt+σ(X )dW +γ(X )dN . (6) t t t t t− t 4 Under(H),thisSDEadmitsanuniquesolutioninS2 anditisshownin[21]thatunder the additional assumption (H ) given below, (5) admits a unique minimal solution 1 (Y,Z,V,K)∈S2×L2(W)×L2(N)×A2 with K predictable. The solution (Y,Z,V,K) is said to be minimal if and only if it has the smallest component Y in the (infinite) class of solutions to (5). (Y ) is called the value t t≥0 process and jumps with a size V =Y −Y . t t t− (H ) There exists a solution (Y˜,Z˜,K˜)∈S2×L2(W)×A2 to 1 T T T T Y =g(X )+ f(X )ds− Z dW + κ(X )dN + dK . t T s s s s− s s Zt Zt Zt Zt (H′) (H ) holds and Y˜ =v˜(t,X ),∀0≤t≤T for some v˜ with linear growth. 1 1 t t (H ) There exists a non negative function ϕ∈C2(Rd) and a constant ρ>0 s.t. 2 Lϕ+f ≤ρϕ, ϕ−Hϕ>0, ϕ≥g, lim ϕ(x) =∞, |x|→∞ 1+|x| in which L is the local component of the generator of the process X and H, the intervention operator: 1 Lv(t,x)=b(x)·D v(t,x)+ Tr σσ⊥(x)D2v(t,x) , x 2 x Hv(t,x)=v(t,x+γ(x))+κ(x) .(cid:0) (cid:1) Let (Yt,x,Zt,x,Vt,x,Kt,x) be the solution to (5) when X ≡ (Xt,x) is s s s s t≤s≤T s t≤s≤T the solutionstartingatx int toSDE (6). Under assumptions(H), (H′)and(H ), [21] 1 2 show that the solution to impulse control problem (1) coincides with initial value of component Yt,x: Yt,x =v(t,x) (7) t andisequaltothe(unique)solutionwithlineargrowthtoquasi-variationnalinequality min −∂v(t,x)−Lv(t,x)−f(t,x); ∂t v(t,x)−Hv(t,x)} = 0, ∀(t,x)∈[0,T)×Rd, (8) min{v(T−,x)(cid:8)−g(x);v(T−,x)−Hv(T−,x)} = 0, ∀x∈Rd. Let us mention that in general, the terminal condition v(T−,·) = g is irrelevant, becauseofthepossiblediscontinuityofY inT−duetoconstraints: therelaxedterminal condition in (8) expresses the possibility of a jump at time T−. Remark 1. Forabetterintuition,thefollowinginterpretationtosolution(Y,Z,V,K) holds when assuming v ∈C1,2([0,T],Rd): ∀0≤t≤T, Y = v(t,X ) t t Z = σ(t,X )D v(t,X ) t t− x t− V = v(t,X +γ(X ))−v(t,X ) t t− t− t− = Hv(t,X )−v(t,X )−κ(X ) t− t− t− K = t −∂v −Lv−f (s,X )ds. t 0 ∂t s R (cid:0) (cid:1) 5 The constraint in (5) means thus that the obstacle condition is satisfied, namely v(t,X )−Hv(t,X )≥0. t− t− 3 Convergence of the Numerical Approximation by Penalization ItdoesnotseempossibletousetheminimalityconditionofthesolutiontoBSDEwith constrained jumps (5) directly in a numerical scheme. We thus propose an approach by penalization of the jump constraint. The penalized constraint is introduced in the BSDE driver: when the constraint is fulfilled, this penalization term disappears, and otherwise penalizes the driver with an exploding factor p. In Theorem 1, we provide an explicit rate of convergence of our approximation withrespecttotheparametersintroduced: namely,thejumpintensityλ,thepenaliza- tion coefficient p and the time step. Such an error estimate is essential for numerical purposes (understanding of the numerical impact of those parameters) and allows to adjust in practice the fineness of the time grid in relation to (λ,p). 3.1 Approximation by Penalization Given a parameter value p>0, the penalized BSDE is: T Yp =g(X )+ f(X )+p(Vp+κ(X ))+λ ds (9) t T s s s− Zt h i T T − ZpdW − VpdN , ∀0≤t≤T s s s s Zt Zt which admits an unique solution (Yp,Zp,Vp) ∈ S2 ×L2(W)×L2(N) from the clas- sical theory of BSDEs with jumps. In addition, the sequence of penalized solutions (Yp,Zp,Vp) tendsinL2([0,T])×L2(W)×L2(N)totheminimalsolution(Y,Z,V)to p F (5) as p goes to infinity, see [21]. Besides, the convergence of (Yp) to Y is monotone p and increasing. Let (Yp,t,x,Zp,t,x,Vp,t,x) be the solution to (9) when X ≡ (Xt,x) . We con- s t≤s≤T sider the following error introduced by this penalization procedure: Ep := sup v(t,x)−Yp,t,x . (10) t 0≤t≤T (cid:12) (cid:12) (cid:12) (cid:12) For any t<η ≤T, let us introduce T vη(t,x)= sup E g(Xt,x,u)+ f(Xt,x,u)ds+ κ(Xt,x,u) (11) T T s τ− u=(τk)k≥1∈U(t,T−η] Zt t<Xkτ≥k≤1T k which corresponds to initial problem (1) restricted to the sub-set of strategies taking values in (t,T −η]. We shall denote by uη∗ = (τkη∗)k≥1 an η21-optimal strategy to 6 problem (11) (the existence of an optimal strategy is not ensured) and by nη∗ be the number of impulses in strategy uη∗ that is: nη∗ :=♯ k ≥1:t<τη∗ ≤T −η . k (cid:8) (cid:9) We will use the following additional assumptions: (Hn) There exists some n¯ ∈N∗ such that ∀j ≥n¯, P(nη∗ ≥j)≤l(j) for some map l such that l(j)≤e−Cj for some some constant C >0. (H∗) There exists a map h such that h(ε)=Oε→0(ε21) and ∀ε>0, P min τη∗ −τη∗ ≤ε ≤h(ε). k+1 k k≥1 (cid:18) (cid:19) (cid:12) (cid:12) (cid:12) (cid:12) Remark 2 (Assumptions (Hn) and (H∗)). Both assumptions (Hn) and (H∗) are di- rectlysatisfiedfortheproblemofSwingoptionsvaluationsincethenumberofexercises right is almost surely bounded by some n and there is some fixed time delay δ >0 max between two consecutive interventions. More generally, (Hn) is intuitively satisfied as soon as the controlled state variable is constrained almost surely and admits jumps of constant sign, see Example 1. (H∗) is verified for sufficiently smooth problems, see for example the case of optimal forest management studied in [1]. Example 1. Let us assume that the state variable defined in (3) is such that • b is uniformly bounded and σ >0 constant, • for some constant c>0, sup γ(x)≤−c, x∈Rd and that the optimal strategy u∗ implies Xu∗ ≥ 0 a.s. Then a straightforward com- T putation shows that the (random) number n∗ of optimal impulses before time T (0,T] satisfies, for any a>0, P n∗ >n =O e−an as n→+∞. (0,T] (cid:16) (cid:17) (cid:0) (cid:1) Proposition 1. Let assumptions (H), (Hn), (H∗), (H′) and (H ) be satisfied. Then 1 2 the penalization error in (10) admits the following bound as p goes to infinity: n¯C¯n¯ 1 Ep ≤C , ∀α∈ 0, . (cid:18)(λp)12−α(cid:19) (cid:18) 2(cid:19) for some constants C >0 and C¯ >1, which do not depend either on λ, p, n¯ or α. 7 Proof. Weprovidethemainargumentsoftheproofandreferthereaderto[1]formore details. The main idea comes from the following explicit functional representation available for Yp,t,x, see [21]: T T Yp,t,x =esssupEνp g(Xt,x)+ f(Xt,x)ds+ κ(Xt,x)dN , (12) t T s s− s νp∈Vp " Zt Zt # where Eνp denotes the expectation under the probability measure Pνp equivalent to P on (Ω,F ) with Radon-Nikodym density T dPνp dP =e−R0T(νsp−1)λdseR0Tln(νsp)dNs. (cid:12)FT (cid:12) (cid:12) Thespecificityofsuchachang(cid:12)eofmeasureisthatitimpactsonlythejumppartsofthe processes: under Pνp, the Brownian motion W remains unchanged whereas N has a (stochastic) intensity (λνp) . We shall denote by Np, the doubly stochastic Poisson s s≥0 process (Cox process) with intensity (λνp) under P, by (τp) the sequence of its s s≥0 k k≥1 jump dates and by Xp the solution to (6) driven by Np. In view of (11) and (12), we introduce a convenient measure change, which intu- itively forces the penalized solution to jump as soon as possible after that an optimal impulse occurs: p if 1 6=Np, ∀s≥0, νp = k≥1 {τkη∗<s} s (13) s (0 elsPe. Notice that νp is a P-measurable process bounded by p a.s. By definition of the counting process Np, ∀s≥0, P τp−τη∗ >s| 1 =k−1 =e−λps. k k {τjp≤τkη∗} j≥1 X Inotherwords,the incrementτp−τη∗ hasanexponentialdistributionwithparameter k k (λp), conditionally to the fact that Xp has jumped one time less than Xuη∗. This allows us to compute an estimate of the distance sup vη(t,x)−Yp,t,x , T t 0≤t≤T (cid:12) (cid:12) which holds for Ep by sending η to 0(cid:12)together with a c(cid:12)ontinuity argumentof the value function in its maturity variable. 3.2 Convergence Rate of the Numerical Scheme Givena regulartime gridπ ={t =0,t ,...,t =T}, we assume that the solutionX 0 1 N to(6)canbesimulatedonπeitherperfectlyorbyusinganEulerschemeanddenoteits discrete-timeversionbyXπ. Alongthe linesof[3],weconsiderthefollowingbackward 8 discrete-time scheme for numerically solving the penalized BSDE with jumps (9), Yp,π =g(Xπ ) tN tN ∀t ∈π,t <T : n n VZttppnn,,ππ == ∆λ∆t1nt1+n1+E1EYhYtpnt,+pnπ,1+π∆1∆WN˜tntn++11|F|Ftntni (14) Yp,π =E Yp,πh|F i where∆t =ttn −+t h,hf∆(tXWn+tπn1)+tni(cid:16)sipth(cid:0)eVtBpn,rπo+wnκi(aXntπinn)c(cid:1)r+em−enVttpno,πn(cid:17)[λti,∆ttn+]1and∆N˜ n+1 n+1 n tn+1 n n+1 tn+1 the compensated version of the Poisson increment ∆N on [t ,t ). tn+1 n n+1 We consider the classical discretizationerror between the continuous-time solution (Yp,Zp,Vp)in (9) andits discrete-time approximation(Yp,π,Zp,π,Vp,π) in(14), that is 1 2 Eπ(Yp):= max sup E Yp−Yp,π 2 n<N−1tn≤t≤tn+1 t tn ! (cid:12) (cid:12) (cid:12) (cid:12)1 Eπ(Zp):= N−1 tn+1E Zp−Zp,π 2dt 2 n=0Ztn t tn ! X (cid:12) (cid:12) (cid:12) (cid:12) 1 Eπ(Vp):= N−1 tn+1E Vp−Vp,π 2λdt 2 . n=0Ztn t tn ! X (cid:12) (cid:12) (cid:12) (cid:12) Because of the lack of first order regularity of the driver of the penalized BSDE fp(x,v):=f(x)+ p(v+κ(x))+−v λ, ∀(x,v)∈Rd×R, (15) classicalregularizationargumen(cid:0)ts for the FBSDE(cid:1)coefficients and Malliavin differenti- ation representations allow us to provide an explicit convergencerate of order |π|12 for errors Eπ(Yp) and Eπ(Vp), but only of order |π|14 for error Eπ(Zp), see Proposition2. Theimpactofthepenalizationcoefficientpontheconvergenceofbackwarddiscrete- time schemesis wellknowninpractice,evenifthere does notexist, to ourbest knowl- edge, anyexplicit computationinthe literature(see forexample the numericalexperi- ments of [23]for the resolutionby penalizationofa BSDE with one reflectingbarrier). Basically, as p increases at fixed discrete-time step, the quantity fp(·)∆t explodes, n+1 leading to a numerical explosion of the approximate values Yp,π, see (14). We show rigorously in Proposition 2 that the discretization error grows exponen- tially with (λp2). This is due to the linear dependence in λ and p of the BSDE driver fp, see (15), and estimate computations based on the use of Gronwall’s lemma. Proposition 2. Assume (H). Then, as soon as 1 |π|=O , (16) λp2 (cid:18) (cid:19) 9 we get the following bounds as p goes to infinity: Eπ(Yp)≤C (1+λ)2λpC¯λp2|π|12 Eπ(Vp)≤C(cid:16)(1+λ)2λ32p2C¯λp2|π(cid:17)|12 . (cid:16) (cid:17) Under (H′), there exists a version of Zp such that, Eπ(Zp)≤C (1+λ)2λ52p3C¯λp2|π|14 , (cid:16) (cid:17) for some constants C >0 and C¯ >1, which do not depend either on λ, p or |π|. Proof. This follows from the same arguments as [15]: computations using Itˆo and Gronwall’s lemma and regularization and Malliavin differentiation arguments applied to the penalized BSDE with jumps (9). We refer the reader to [1] for a detailed proof. Propositions 1 and 2 enable us to establish a global convergence rate of the error introduced by our approximationby penalization. Theorem 1. Let the assumptions of Proposition 1 be satisfied. Then 1 1 Ep+Eπ(Yp)≤C +(1+λ)2λpC¯λp2|π|12 , ∀α∈ 0, (17) (cid:18)(λp)12−α (cid:19) (cid:18) 2(cid:19) for some constants C >0 and C¯ >1, which do not depend either on λ, p, |π|, n¯ or α. Thus, for a sufficiently small time step |π| with respect to λ and p, the global error is such that [Ep+Eπ(Yp)]∗ =O 1 , ∀α∈ 0,1 . (λp)12−α 2 (cid:18) (cid:19) (cid:0) (cid:1) Remark 3 (Global convergence rate). At fixed time step |π|, the convergence rate strongly deteriorates as λ or p increases, see (17). The numerical method is more sensitive to p than to λ according to (16) and (17). (16) constitutes a necessary condition for the convergence of the backward discrete-time scheme. In practice, the penalization parameter will need to be chosen relatively small and the time step |π| very small to avoid multiple jump times on each time step (otherwise, this introduces a bias). 4 Application to Swing Options Valuation Inthissection,ourmethodisappliedforthevaluationofSwingoptionsintheBlackand Scholes framework. This constitutes a multiple optimal stopping time problem which can be reformulated as a particularly degenerate three-dimensional impulse control problem. We have been able to achieve convergence for a small maximal number of exercises rights (n ≤2) due to the slow computational speed of our method. max 10