1 Introduction

Long-time large deviations for the multi-asset Wishart stochastic volatility model and option pricing Aurélien Alfonsi111Université Paris-Est, Cermics (ENPC), INRIA, F-77455 Marne-la-Vallée, France.

David Krief 222LPSM, Université Paris Diderot, Paris, France

Peter Tankov 333ENSAE ParisTech, Palaiseau, France

In this paper, we prove a large deviations principle for the class of multidimensional affine stochastic volatility models considered in (Gourieroux, C. and Sufana, R., J. Bus. Econ. Stat., 28(3), 2010), where the volatility matrix is modelled by a Wishart process. This class extends the very popular Heston model to the multivariate setting, thus allowing to model the joint behaviour of a basket of stocks or several interest rates. We then use the large deviation principle to obtain an asymptotic approximation for the implied volatility of basket options and to develop an asymptotically optimal importance sampling algorithm, to reduce the number of simulations when using Monte-Carlo methods to price derivatives.

Key words: Large deviations, Wishart process, Importance sampling, Basket options, Implied volatility

MSC2010: 60F10, 91G20, 91G60

## 1. Introduction

The Heston stochastic volatility model [Heston, 1993] is one of the most popular models in quantitative finance for the evolution of a single asset price. The Wishart stochastic volatility model is its natural extension to a basket of assets, since it coincides with the Heston model in dimension  and preserves the affine structure. This model, proposed in [Gourieroux and Sufana, 2010], assumes that under the risk-neutral probability, the vector of asset prices is modelled as an Itô process

 dSt=Diag(St)(r1dt+~X1/2td~Zt), (1.1)

where the volatility matrix follows the Wishart process with dynamics

 d~Xt=(αa⊤a+~b~Xt+~Xt~b⊤)dt+~X1/2td~Wta+a⊤(d~Wt)⊤~X1/2t, (1.2)

where and are independent standard -dimensional and -dimensional Brownian motions, and is the diagonal matrix whose diagonal elements are given by the vector .

The matrix process (1.2) has been introduced by [Bru, 1991] to model the perturbation of experimental biological data. As shown by [Bru, 1991] and [Cuchiero et al., 2011] in a more general framework, for (resp. ), the SDE (1.2) has a unique strong (resp. weak) solution. Furthermore, since is positive semi-definite [Bru, 1991, Prop. 4], Wishart processes turn out to be very suitable processes to model covariance matrices. This, and the affine property of the Wishart process, led several authors to use them in stochastic volatility models for a single asset, such as [Da Fonseca et al., 2008] and [Benabid et al., 2008] and in the Wishart stochastic volatility model for multiple assets (1.1)–(1.2). Subsequently, this model has been extended by [Da Fonseca et al., 2007] to include a constant correlation between and in a way to preserve the affine structure.

By using the affine property, the Laplace transform of the model (1.1)–(1.2) is computed as follows [Da Fonseca et al., 2007].

 E(eθ⊤log(St))=exp(βθ(t)+\textupTr[γθ(t)~X0]+δ⊤θ(t)log(St)), (1.3)

where and satisfy the matrix Riccati equations

 ∂tβθ(t) =rδ⊤θ(t)1+α\textupTr[γθ(t)] ∂tγθ(t) =~b⊤γθ(t)+γθ(t)~b+2γθ(t)a⊤aγθ(t)−12(% Diag(δθ(t))−δθ(t)δ⊤θ(t)) ∂tδθ(t) =0,

with initial conditions , and . Since the Riccati equations can be solved explicitly, the Laplace transform can be expressed explicitly in terms of matrix exponentials and inverses.

The goal of the present paper is to prove a large deviations principle the Wishart stochastic volatility model (1.1)–(1.2) in the large-time asymptotic regime. Since the Laplace transform of the log-price vector in the Wishart model is known explicitly, a natural path towards a large deviations principle is via Gärtner-Ellis theorem. However, despite the explicit form of the Laplace transform, it is not easy to calculate its long-time asymptotics and to check the assumptions of the theorem because of the multi-dimensional setting. In this paper we therefore focus on a (large enough) subclass of the model (1.1)–(1.2) which enables us to obtain a simpler formula for the limiting Laplace transform and then prove a large deviations principle.

Beyond its theoretical interest, knowing that a given model satisfies a large deviations principle, and knowing the explicit form of the rate function, enables one to develop a number of important applications. One can mention e.g., efficient importance sampling methods for Monte Carlo option pricing; asymptotic formulas for option prices and implied volatilities in various asymptotic regimes, approximate evaluation of risk measures, simulation of rare events and others. We refer the reader to [Pham, 2007] for a review of various applications of large deviations methods in finance. In this paper we develop applications to variance reduction of Monte Carlo methods and to the asymtotic computation of implied volatilities far from maturity.

Our variance reduction method follows previous works of [Guasoni and Robertson, 2008], [Robertson, 2010] and [Genin and Tankov, 2016] and uses Varadhan’s lemma of large deviations theory to approximate the optimal measure change in the importance sampling algorithm. Note that since the Laplace tranform is known explicitly, Fourier inversion methods can be used, as explained in [Da Fonseca et al., 2007]. However, these methods are much less competitive than in dimension  since they require to approximate an integral on . When, for complexity reasons, Fourier methods are not an option, the use of a large number of Monte-Carlo simulations is necessary. [Ahdida and Alfonsi, 2013] present an exact simulation method for Wishart processes and a second order scheme for the Gourieroux and Sufana model (1.1)–(1.2). Thus, it is possible to sample efficiently such processes, and it is relevant to develop variance reduction techniques to reduce computational costs.

The approximation of implied volatility far from maturity extends earlier results on the Heston model and the one-dimensional affine stochastic volatility models [Forde and Jacquier, 2011, Jacquier et al., 2013] to the multidimensional setting of Wishart model. Once again, this approach is more relevant in the multidimensional setting, since in one-dimensional affine models the implied volatility may be quickly computed by Fourier inversion.

In this paper, we denote the set of real squared matrices, the set of symmetric matrices and , (resp. ), the sets of symmetric an non-negative (resp.) positive definite. For a Borel set , we denote by the closure of and by the interior of .

The paper is structured as follows. In Section 2, we describe the model, make certain assumptions on the parameters and give some properties of the model. In Section 3, we prove that the asset log-price vector satisfies large deviations principle when maturity goes to infinity. In Section 4, we calculate the asymptotic put basket implied volatility, following the approach of [Jacquier et al., 2013]. In Section 5, we develop the variance reduction method using Varadhan’s lemma. Finally, in Section 6, we test numerically the results of Sections 4 and 5.

## 2. The Wishart stochastic volatility model

In this section we introduce the subclass of the Wishart stochastic volatility models, in which we are interested in the present paper, and compute the Laplace transform of the log stock price process.

Let be a -dimensional vector stochastic process with dynamics

 dSt=Diag(St)(r1dt+a⊤X1/2tdZt),Si0>0, i=1,…,n, (2.1)

where , , is -dimensional standard Brownian motion and the stochastic volatility matrix is a Wishart process with dynamics

 dXt=(αIn+bXt+Xtb)dt+X1/2tdWt+(dWt)⊤X1/2t,X0=x. (2.2)

with , invertible, and is a matrix standard Brownian motion independent of . Note again that [Bru, 1991, Prop. 4]. Let us also assume that is such that .

###### Remark 2.1.

The model defined in (2.1) and (2.2) is a (quite large) subclass of the one defined in (1.1) and (1.2). Indeed, defining , we have , where is another -dimensional standard Brownian motion and

 d~Xt=(αa⊤a+~b~Xt+~Xt~b⊤)dt+~X1/2td~Wta+a⊤(d~Wt)⊤~X1/2t,~X0=a⊤xa,

where and is another -Brownian motion.

###### Remark 2.2.

In dimension one, the model defined by eqs. (2.1) and (2.2) corresponds to the famous Heston model [Heston, 1993] and being negative definite yields the mean reversion property of the stochastic volatility process.

Defining the log-price , a simple application of Itō’s lemma gives

 dYt=(r1−12((a⊤Xta)11,...,(a⊤Xta)nn)⊤)dt+a⊤X1/2tdZt. (2.3)

We are interested in the Laplace transform of . In order to calculate it, we first cite the following proposition.

###### Proposition 2.3.

[Alfonsi et al., 2016, Prop. 5.1.]. Let , , and with dynamics (2.2). Let be such that

 ∃m∈Sn,v2−mb−bm−2m2∈S+nandw2+m∈S+n.

If , then we have for

 E[exp(−12\textupTr[wXt]−12\textupTr[vRt])] =exp(−α2\textupTr[b]t)det[Vv,w(t)]α/2exp(−12\textupTr[(V′v,w(t)V−1v,w(t)+b)x]),

with

 Vv,w(t)=(∞∑k=0t2k+1~vk(2k+1)!)~w+∞∑k=0t2k~vk(2k)!,~v=v+b2and~w=w−b.

If besides, , then

 Vv,w(t)=~v−1/2sinh(~v1/2t)~w+cosh(~v1/2t)

and

 V′v,w(t)=cosh(~v1/2t)~w+sinh(~v1/2t)~v1/2.

The following proposition provides and explicit formula for the Laplace transform of the log stock price in the model (2.1)–(2.2).

###### Proposition 2.4.

Let be the function defined by

 ϕ(θ):=b2+a(Diag(θ)−θθ⊤)a⊤∈Sn, (2.4)

Let , be the set defined by

 U:={θ∈Rn:ϕ(θ)∈S+n}.

Then, for all , the Laplace transform of is

 E(eθ⊤Yt)=eθ⊤Y0+rθ⊤1t−α2\textupTr[b]t−12\textupTr[(b+ϕ1/2(θ))x−exp(−tϕ1/2(θ))(b+ϕ1/2(θ))V−1(t)x]det[V(t)]α/2,

where

 V(t)=cosh(tϕ1/2(θ))−ϕ−1/2(θ)sinh(tϕ1/2(θ))b.
###### Proof.

By conditioning on the trajectory of , we have

 E(eθ⊤Yt)=E(E(eθ⊤Yt∣∣(Xs)s≤t)),

where

 E(eθ⊤Yt∣∣(Xs)s≤t) =eθ⊤Y0+rθ⊤1t−12∫t0θ⊤((a⊤Xsa)11,...,(a⊤Xsa)nn)⊤−θ⊤a⊤Xsaθds =eθ⊤Y0+rθ⊤1t−12∫t0\textupTr[Diag(θ)a⊤Xsa]−\textupTr[θ⊤a⊤Xsaθ]ds =eθ⊤Y0+rθ⊤1t−12\textupTr[a(Diag(θ)−θθ⊤)a⊤Rt].

Let . Then and

 a(Diag(θ)−θθ⊤)a⊤2−mb−bm−2m2=ϕ(θ)2∈S+n.

Therefore, by Proposition 2.3,

 E(eθ⊤Yt) =eθ⊤Y0+rθ⊤1tE(e−12\textupTr[a(Diag(θ)−θθ⊤)a⊤Rt]) (2.5)

where

 {V(t)=cosh(tϕ1/2(θ))−ϕ−1/2(θ)sinh(tϕ1/2(θ))b,V′(t)=sinh(tϕ1/2(θ))ϕ1/2(θ)−cosh(tϕ1/2(θ))b.

Since , we can write , where is diagonal, is orthonormal and .

 ⎧⎪ ⎪ ⎪ ⎪⎨⎪ ⎪ ⎪ ⎪⎩V(t)=P(cosh(tD1/2)+sinh(tD1/2)D−1/2^b)P⊤,V′(t)=P(sinh(tD1/2)D1/2+cosh(tD1/2)^b)P⊤=ϕ1/2(θ)V(t)−exp(−tϕ1/2(θ))(b+ϕ1/2(θ)).

Replacing by the latter expression finishes the proof.∎

###### Remark 2.5.

Note that, when , is not invertible. The notation is therefore abusive and is to be interpreted as the finite limit

 limS+,∗n∋ϕ→ϕ(θ)ϕ−1/2sinh(tϕ1/2)=∞∑k=0ϕ(θ)kt2k+1(2k+1)!.
###### Remark 2.6.

The set is bounded. Indeed, let , with and . Then, letting , we have

 u⊤ϕ(θ)u=∥b(a⊤)−1¯θ∥2+λ¯θ⊤Diag(¯θ)¯θ−λ2≤∥b(a⊤)−1∥2+λ−λ2

It follows that is contained, e.g., in the set with

 λ∗=max{2,∥b(a⊤)−1¯θ∥√2}.

## 3. Long-time large deviations for the Wishart volatility model

In this section, we prove that the Wishart stochastic volatility model satisfies a large deviation principle when time tends to infinity.

### 3.1. Reminder of large deviations theory

Let us recall some standard definitions and results of large deviations theory. For a wider overview of large deviations theory, we refer the reader to [Dembo and Zeitouni, 1998]. We consider a family of random variables on a measurable space , where is a topological space.

###### Definition 3.1 (Rate function).

A rate function is a lower semi-continuous mapping . A good rate function is a rate function such that, for every , is compact.

###### Definition 3.2 (Large deviation principle).

satisfies a large deviation principle with rate function if, for every , denoting and the interior and the closure of ,

 −infx∈∘AΛ∗(x)≤liminfϵ→0ϵlogP(Xϵ∈A)≤limsupϵ→0ϵlogP(Xϵ∈A)≤−infx∈¯AΛ∗(x).
###### Definition 3.3.

Let be a convex function with domain . is called essentially smooth if is differentiable on and for every , .

The following theorem is the celebrated Gärtner-Ellis theorem of the large deviations theory. [Dembo and Zeitouni, 1998] give a version of this theorem for a family of random variables parameterized by an integer number (see paragraph 2.3 in their book), but the version for families parameterized by a real number is easily deduced from the abstract Gärtner-Ellis theorem given in paragraph 4.5.3.

###### Theorem 3.4 (Gärtner-Ellis).

Let be a family of random vectors in . Assume that for each ,

 Λ(λ):=limϵ→0ϵlogE[e⟨λ,Xϵ⟩ϵ] (3.1)

exists as an extended real number. Assume also that 0 belongs to the interior of . Denoting

 Λ∗(x)=supλ∈Rn⟨λ,x⟩−Λ(λ),

the Fenchel-Legendre transform of , the following hold.

1. For any closed set ,

 limsupϵ→0ϵlogP(Xϵ∈F)≤−infx∈FΛ∗(x).
2. For any open set ,

 liminfϵ→0ϵlogP(Xϵ∈G)≥−infx∈G∩FΛ∗(x),

where is the set of exposed points of , whose exposing hyperplane belongs to the interior of .

3. If is an essentially smooth, lower semi-continuous function, then satisfies a large deviations principle with good rate function .

###### Remark 3.5.

The function of (3.1) is a convex function. Indeed, let and . A direct application of Hölder’s inequality yields

Applying the logarithm then proves that and therefore are convex.

###### Theorem 3.6 (Varadhan’s Lemma, extension of [Guasoni and Robertson, 2008]).

Let be a metric space with its Borel -field. Let be a family of -valued random variables that satisfies a large deviations principle with rate function . If is a continuous function which satisfies

 limsupϵ→0ϵlogE[exp(γφ(Xϵ)ϵ)]<∞

for some , then, for any ,

 supx∈A∘{φ(x)−Λ∗(x)}≤liminfϵ→0ϵlog∫A∘exp(φ(z)ϵ)dμϵ(z)≤limsupϵ→0ϵlog∫¯Aexp(φ(z)ϵ)dμϵ(z)=supx∈¯A{φ(x)−Λ∗(x)},

where denotes the law of

### 3.2. Long-time behaviour of the Laplace transform of the log-price

Let and define the transformation , which corresponds to the long-time behaviour of . We are interested in the function

 θ↦limϵ→0ϵlogE[eϵ−1θ⊤YϵT].

We first give the following lemma.

###### Lemma 3.7.

Let such that est invertible for all . Then, is bounded for all sufficiently large .

###### Proof.

Since is invertible, for all ,

 (A+tB)−1tB={I+(t−t0)~B}−1(t−t0)~Btt−t0,

where . Now, the fact that est invertible for means that the eigenvalues of satisfy or for all . This implies for some , and since the adjugate matrix of has coefficients of order , we get that is bounded for . Therefore, is bounded, and as well, whenever is sufficiently large. ∎

We now characterise the asymptotic behaviour of the Laplace transform of .

###### Proposition 3.8.

Define

 Λ(θ):={T(rθ⊤1−α2\textupTr[b+ϕ1/2(θ)]) if θ∈U∞ if θ∉U. (3.2)

For every ,

 limϵ→0ϵlogE[eϵ−1θ⊤YϵT]=Λ(θ).
###### Proof.

Let . By Proposition 2.4,

 ϵlogE[eϵ−1θ⊤YϵT] =ϵlogE[eθ⊤YT/ϵ] (3.3) =ϵ(θ⊤Y0−12\textupTr[(b+ϕ1/2(θ))x]) +12ϵ\textupTr[exp(−T/ϵϕ1/2(θ))(b+ϕ1/2(θ))V−1(T/ϵ)x] +Trθ⊤1−Tα2\textupTr[b]−α2ϵlogdet[V(T/ϵ)].

Write , where is diagonal, is orthonormal and let . Then

 V(t)=P(cosh(tD1/2)+sinh(tD1/2)D−1/2^b)P⊤,

Let and be square matrices with and . We then have

 cosh(tD1/2)=etD1/22(In+e−2tD1/2)=etD1/22(In+E+O(t−1))

and

 sinh(tD1/2)D−1/2=etD1/22D−1/2(In−e−2tD1/2)=etD1/22(~E+2tE+O(t−1)).

Therefore,

 V(t) =12PetD1/2((In+E)+(2tE+~E)^b+O(t−1)\scalebox0.8$O$(t−1))P⊤ (3.4) =−12P(In+E)etD1/2(^b−1+(tE+~E)+O(t−1)\scalebox0.8$O$(t−1))P⊤b

and

 V−1(t)=−2b−1P(^b−1+(tE+~E)+O(t−1))−1e−tD1/2(In−12E)P⊤

where the invertibility of is guaranteed for every by the existence of the Laplace transform. Since and , and is therefore invertible. Hence

 (^b−1+(tE+~E)+O(t−1))=(^b−1+(tE+~E))(In+O(t−1))

and

 V−1(t)=−2b−1P(In+O(t−1))(^b−1+(tE+~E))−1e−tD1/2(In−12E)P⊤.

But

 =t−1(^b−1+(tE+~E))−1tE+(^b−1+(tE+~E))−1(In−E)e−tD1/2,

where is bounded by Lemma 3.7. Therefore,

 (^b−1+(tE+~E))−1e−tD1/2→0

and as . Using (3.4), we find

 ϵlogdet[V(T/ϵ)] =T\textupTr[D1/2]+ϵlogdet[12(In+E)(In+(ϵ−1TE+~E)^b+O(ϵ))] =T\textupTr[ϕ1/2(θ)]+ϵlogdet[ϵ−1TE^b+12(In+E)(In+~E^b)+O(ϵ)] =T\textupTr[ϕ1/2(θ)]−nϵlog(ϵ)+ϵlogdet[TE^b+ϵ2(In+E+~E^b)+O(ϵ2)].

We have , since the latter determinant is a non-zero polynomial of (for the determinant is clearly positive). Thus, by passing to the limit, . Furthermore, since , is bounded. Therefore,

Finally, passing to the limit in (3.3) finishes the proof. ∎

The next proposition proves the essential smoothness of .

###### Proposition 3.9.

The function defined in (3.2) is essentially smooth.

###### Proof.

The function defined in (3.2) is a lower semi-continuous proper convex function with domain . Furthermore, since for every , , is of class on . Only remains to prove that when goes to the boundary of . Let . By Proposition 3.8

 Λ(θ)=T(rθ⊤1−α2\textupTr[b+ϕ1/2(θ)]).

Then for every ,

 ∂θjΛ(θ)=T(r−α2\textupTr[∂θj[ϕ1/2](θ)]),

where satisfies

 ∂θjϕ(θ)=∂θj[ϕ1/2(θ)ϕ1/2(θ)]=ϕ1/2(θ)∂θj[ϕ1/2](θ)+∂θj[ϕ1/2](θ)ϕ1/2(θ).

Multiplying this equation by and using the cyclic property of the trace, we get

 \textupTr[∂θj[ϕ1/2](θ)]=12\textupTr[ϕ−1/2(θ)∂θjϕ(θ)].

and therefore

 ∂θjΛ(θ)=T(r−α2\textupTr[∂θj[ϕ1/2](θ)])=T(r−α4\textupTr[ϕ−1/2(θ)∂θjϕ(θ)]), (3.5)

where

 ∂θjϕ(θ)=a(eje⊤j−θe⊤j−ejθ⊤)a⊤.

We write with diagonal and denote , which is invertible since is orthonormal and . Then

 \textupTr[ϕ−1/2(θ)∂θjϕ(θ)] =\textupTr[D−1/2P⊤∂θjϕ(θ)P] =\textupTr[D−1/2w⊤(eje⊤j−θe⊤j−ejθ⊤)w]