left=1.3in, right=1.3in, top=1in,bottom=1in, includefoot, headheight=13.6pt

# Efficient Pricing and Calibration of High-Dimensional Basket Options

Lech A. Grzelak<EMAIL_ADDRESS>Juliusz Jablecki<EMAIL_ADDRESS>Dariusz Gatarek<EMAIL_ADDRESS>Financial Engineering,
Rabobank, Utrecht, the Netherlands  Mathematical Institute, Utrecht
University, Utrecht, the Netherlands  Financial Risk Management Department,
National Bank of Poland, Warsaw, Poland  Chair of Quantitative Finance,
Faculty of Economic Sciences, University of Warsaw, Warsaw, Poland  Systems
Research Institute, Polish Academy of Sciences, Warsaw, Poland

###### Abstract

This paper studies equity basket options – i.e. multi-dimensional derivatives
whose payoffs depend on the value of a weighted sum of the underlying stocks –
and develops a new and innovative approach to ensure consistency between
options on individual stocks and on the index comprising them. Specifically,
we show how to resolve a well-known problem that when individual constituent
distributions of an equity index are inferred from the single-stock option
markets and combined in a multi-dimensional local/stochastic volatility model,
the resulting basket option prices will not generate a skew matching that of
the options on the equity index corresponding to the basket. To address this
“insufficient skewness”, we proceed in two steps. First, we propose an
“effective” local volatility model by mapping the general multi-dimensional
basket onto a collection of marginal distributions. Second, we build a
multivariate dependence structure between all the marginal distributions
assuming a jump-diffusion model for the effective projection parameters, and
show how to calibrate the basket to the index smile. Numerical tests and
calibration exercises demonstrate an excellent fit for a basket of as many as
30 stocks with fast calculation time.

###### keywords:

Basket Options, Index Skew, Monte Carlo, Local Volatility, Stochastic
Volatility, Collocation Methods

##  1 Introduction

Basket options are contingent claims in which the underlying is a group of
assets – typically equities (single stocks as well as equity indices/ETFs),
commodities or currencies. As with standard options, the holder of a basket
option has the right, but no obligation, to buy (call) or sell (put) the group
at a specified strike price which itself is based on the weighted value of the
component assets. Basket options can be traded on their own, but commonly
feature implicitly in structured products like equity basket linked notes
(ELNs) or equity basket certificates of deposit, where return to investors
depends on the percentage change in the value of stocks or stock indices in
the basket, with partial or total capital protection. For example, on April
16, 2021, JP Morgan priced $\$$9.6 million of notes linked to a basket of 30
unequally weighted U.S. stocks with exposure to infrastructure111Cf. the
prospectus availible in the online records of the U.S. Securities and Exchange
Commission at:
https://www.sec.gov/Archives/edgar/data/19617/000089109221003578/e13291-424b2.htm.
The structure is effectively a combination of a European-style basket option
and a hypothetical zero coupon bond of the issuing entity.

The key feature of basket options, and any structured notes based on them, is
their inherently multidimensional risk profile which presents two sets of
related challenges. First, there is the obvious computational problem of
handling a potentially large number of correlated risk factors driving the
basket, which in turn requires involved Monte Carlo simulation and calls for
the deployment of some dimension-reduction technique, especially in light of
the fact that many such structures exhibit early exercise or other exotic
features (e.g. Asian or lookback options).

Second, and perhaps more importantly, there is the problem of ensuring
consistency between model dynamics chosen for the basket and its constituents,
on the one hand, and market pricing on the other. Typically, stocks or stock
indices/ETFs included in the basket will have their own vanilla options
markets and the basket collection itself can also be traded on its own
account. Thus, ensuring consistency entails not only fitting the implied
volatility smiles/skews for individual assets in the basket – which can be
done with excellent precision using e.g. the local volatility model or the
Heston model – but also their covariance structure, such as to match closely
the implied volatility smile for the market index corresponding to the basket.

Market practice is to use the local volatility model for individual basket
members and historical correlations estimates between underliers to relate the
respective risk factors in a Monte Carlo simulation.222For example, in the
Bloomberg basket options pricing template correlations are, by default,
estimated over a 5 year period, whereby to eliminate noise, a given percentile
of rolling 6-month cross-correlation estimates is chosen in the
parameterization of the full correlation matrix. In addition, to avoid
situations in which the correlation matrix fails to be positive semi-definite,
a flat correlation level is often imposed. However, this procedure does not
automatically ensure that the model price of a basket option is consistent
with the market price of the option on the index corresponding to the basket.
And even if, as per market practice, the correlation coefficient is chosen
specifically to match the price of at-the-money index options, the procedure
generates an implied volatility smile/skew that for the basket that is
significantly less pronounced than that observed for the index options.

Consider for example options on the EURO STOXX Food and Beverage Index
(Bloomberg ticker SX3E). The index currently comprises 11 large cap European
equities from the food and beverage sector, such as French Danone, Dutch
Heineken and Irish Kerry Group, each of which has its own listed options
market. Matching 1Y ATM implied volatility for the index (15.51$\%$) requires
setting $\rho=0.22$ (data as of 1 April 2022). Repricing the basket with this
flat correlation parameter for different moneyness levels generates an implied
volatility skew of about 3.6 vol points,333We define the skew here loosely as
the difference in implied volatilities between the 85-120$\%$ ATM levels.
which is less than half of the skew observed in the market (7.5 vol points;
cf. Figure 1, LHS). Similarly, fitting correlation to 1Y ATM vol does not
allow to reproduce the term structure of index implied volatilities observed
in the market (Figure 1, RHS).

Figure 1: Local volatility model calibrated to Euro STOXX Food and Beverages
Index (SX3E) and a corresponding basket of its constituents; LHS: skew for 1Y
index options; RHS: implied volatility term structure of index options; data
as of 1 April 2022; correlation parameter $\rho=0.22$ calibrated to 1Y ATM
implied SX3E volatility.

The above heuristic argument shows that a covariance structure implied by the
modeling choice, like the multi-dimensional local volatility (or perhaps
Heston) model, will not in general match the market well. This pattern becomes
even more pronounced when handling really large baskets of tens or hundreds of
dimensions, as the generated implied volatilities degenerates to almost a flat
line, resembling implied volatilities produced by a simple Black-Scholes
model. These shortcomings have been known in the literature for quite some
time, as a number of papers have documented that when a constant correlation
is picked to match the price of the ATM implied volatility of the index, it
generates a skew that is smaller by roughly a factor of two than the market
skew – as indeed shown in Figure 1 (see [3], [5], and [4]).

One could argue at this point that since the implied volatilities for the
index are generally available, there is no need to model individual underlying
assets, and the problem of pricing a basket corresponding exactly to the index
is a theoretical exercise with little practical relevance. Note, however, that
by building a model for the basket that matches the covariance structure
implied by index options (while simultaneously calibrating to individual
stocks’ options) we can easily price products whose value depends on any
linear combination of any subset of stocks included in the broader index.
Hence, in essence, we propose a tool for transferring market-implied
covariance patterns between liquid and illiquid baskets. This allows us, for
example, to price a basket of utilities stocks included in FTSE 100 or health
care stocks in Dow Jones, while making sure that the prices of such illiquid
baskets will be consistent with the covariance structures implied by prices of
the respective liquid index options.

Several approaches have recently been proposed in the literature to tackle the
problem of pricing basket options and their skews. Perhaps most notably, so
called local correlation models – which come in various shapes and forms –
allow correlation between stocks to depend on their prices, just like in
standard local volatility models ([14], [15], [21]). Local correlation offered
the first consistent remedy to the problem of insufficient skewness referred
to above (cf. Figure 1; [3]), however not without a cost. Local correlation
models tend to be slow, difficult to implement and lead to the usual problems
of ensuring a positive semi-definite correlation matrix, all of which limits
their practical usefulness in a truly multidimensional setting.

A related family of models proposes stochastic instantaneous correlation,
usually introduced via Jacobi processes either in a pure form or with jumps
(see e.g. [1] and [24]). Such models aim not only to fit the market prices for
correlation products but also, more specifically, to address the so called
correlation skew, i.e. an empirically documented tendency for cross-
correlations to rise when the market falls. Again, however, the added
complexity leads to a further increase in dimensionality and the associated
computational burden.

Still others propose to enforce the desired level of the skew through non-
Gaussian copulas (e.g. [17]), however the resulting model is “black boxy” in
nature and lacks dynamics, which means it can only be applied to products with
a single maturity, as observed in [24]. Moreover, finding a multi-dimensional
copula that isn’t Gaussian can be challenging in practice, and anyway the need
to calibrate the model through brute-force Monte Carlo makes the whole
approach slow and costly.

Finally, [6] proposes an analytical framework to handle stochastic correlation
via Wishart processes while [16] generalizes Merton’s jump diffusion option
pricing approach to the multi-asset case, however without offering a
comprehensive solution that would match both the basket skew as well as
marginals of the constituent stocks.

Against such a background, we try to contribute to the existing literature by
building a highly efficient numerically and robust model capable of handling
even very large baskets while allowing quick and precise calibration. In terms
of structure, the model is a local volatility type and arises by mapping the
basket on a set of proxy variables ensuring that both marginal distributions
for individual stocks are matched and the covariance structure of the basket
corresponds to that of the traded index. Specifically, to ensure excellent
fitting of individual stocks’ marginal distributions, we leverage the so-
called Stochastic Collocation Monte Carlo Sampler (cf. [11]), which is a
computationally cheap method for approximating expensive distributions (we
make no assumptions as to the models driving the respective stocks in the
basket). We then control the overall covariance structure of thus
reconstituted basket via a judicious choice of “kernel” processes, whereby –
similarly as in [16] – we utilize Merton’s jump-diffusion model [19] and show
that the mean and volatility of Poisson jumps, together with cross-
correlations of individual kernel processes, produce a rich enough dependence
structure to fit both ATM index vols and skew with excellent precision. The
proposed approach can be extended into a fully stochastic volatility setting;
however, as we shall see, there appears to be no need for the added complexity
as the classic Merton diffusion already delivers acceptable results.

The rest of the paper proceeds as follows. Section 2 presents the general idea
of the model and discusses the stochastic collocation approach to basket
reconstitution. Section 3 discusses calibration strategy followed by a
discussion of numerical exercises in Section 4. Finally, Section 5 draws
conclusions.

##  2 The Modeling Framework

We consider a collection of $N$ assets444Without loss of generality, we shall
henceforth think of the underlying assets as stocks, however the method
developed below is obviously general and, mutatis mutandis, applies to other
instruments as well. $S_{1}(t),S_{2}(t),\dots,S_{N}(t)$ and define a basket,
$B(\cdot)$, as consisting of $w_{i}$ shares of each individual stock
$S_{n}(t)$, such that at any time $t$ its price is given by:

$\displaystyle
B(t)=\sum_{i=1}^{N}\omega_{i}S_{i}(t),\;\;\;\omega_{j}\in\mathbb{R}.$ (2.1)

In what follows, for modeling purposes, we shall distinguish between a basket
$B(\cdot)$ as a portfolio of stocks, and a stock index, with price process
$I(t)$, representing the underlying instrument for equity index options (e.g.
an ETF), with realistic values of $N$ ranging from 30-50 (DAX Index) up to
hundreds (S&P 500 Index) or even thousands constituents (Russell 3000 Index).

Our approach to valuing structured payoffs on a basket $B(t)$ will be endowed
with three critical features: (i) low dimensionality to facilitate quick and
efficient pricing; (ii) consistency with the values of options on individual
stocks making up the basket; and (iii) consistency with market prices of
options on the stock index corresponding to the basket.

We pursue features (i)-(iii) above by building a one-dimensional local
volatility process 555The proposed framework can also be extended with a
stochastic volatility process. Such an extension is trivial and will, for
simplicity, be omitted.:

$\displaystyle{\rm d}\bar{B}(t)=r\bar{B}(t){\rm
d}t+\sigma_{LV}(t,\bar{B}(t))\bar{B}(t){\rm d}W(t),$ (2.2)

with $r$ being an interest rate and $\sigma_{LV}$ the local volatility
function, defined as follows:

$\displaystyle\small\sigma_{LV}^{2}(t,k)=\frac{\frac{\partial}{\partial
t}\left[{\rm e}^{-rt}\int_{k}^{\infty}(y-k)f_{B(t)}(y){\rm
d}y\right]+rk\left(F_{B(t)}(k)-1\right)}{\frac{1}{2}k^{2}_{i}f_{B(t)}(k)},$
(2.3)

where $F_{B(t)}(\cdot)$ is the CDF of the basket $B(t)$ and $f_{B(t)}(\cdot)$
the corresponding PDF. We shall require that the approximate model for the
basket $\bar{B}(t)$ generates an implied volatility smile/skew matching the
implied volatility surface observed in the market for options on the
corresponding stock index. This will ensure that prices of index options
derived from our basket model will by definition match perfectly (cf. [7],
[8], [12]). Furthermore, it will also guarantee consistency between model-
derived prices of baskets on any subset
$\\{S_{i_{1}},...,S_{i_{K}}\\}\subseteq\\{S_{1},...,S_{N}\\}$ of the $N$
stocks included in the index and the prices of liquid index options.

In practical terms, our approach entails two critical steps. First, for each
$T_{j}$ we reconstitute the basket $B(T_{j})$ by projecting each process
$S_{i}$ – assumed to be known only through its marginal distributions/options
prices – on a polynomial of synthetic variables $g_{i,j}(X_{i}(T_{j}))$, such
that

$\displaystyle S_{i}(T_{j})\approx{\rm
e}^{g_{i,j}(X_{i}(T_{j}))},\;\;\;\text{and}\;\;\;B(T_{j})\approx\sum_{i=1}^{N}{\rm
e}^{g_{i,j}(X_{i}(T_{j}))}.$ (2.4)

The construction of $g_{i,j}(\cdot)$ will be described in detail below, but
its main role will be to ensure the fit of marginal distributions for any
choice of $X_{i}(T_{j})$. The reconstitution of the basket paves the way for
the second stage of the modeling procedure, namely imposing a multivariate
dependence structure on the surrogate kernel processes such as to match the
index implied volatility skew, or equivalently, matching marginal
distributions between the basket and the index for each maturity point
$T_{j}$, i.e. ensuring that: $B(T_{j})\stackrel{{\scriptstyle
d}}{{=}}I(T_{j}).$ This step will be handled through parameters of the
“kernel” processes $X_{j}(t)$, which we assume to be driven by Merton-type
jump diffusions. This procedure will allow us to both calibrate the model to
individual stocks, and control/modify the covariance structure of all the
underlying assets.

###  2.1 Covariance Structure of the Basket and the “Leaking Correlation”
Problem

At this stage it may be worthwhile to pause for a moment and provide some
further motivation for why controlling basket dependence structure is so
tricky and a judicious choice of kernel processes is required. After all, it
could be argued that the model price of a basket option can be fitted to a
desired market level by altering the correlations $\rho_{i,j}$ between the
Brownian motions driving individual stocks $S_{i},S_{j}$. Such correlation
coefficients would then essentially become model inputs describing the
stochastic nature of the underlying assets. However, as discussed in [2],
correlations can be expected to “leak”, i.e. dissolve, in a Monte Carlo
setting, so that effective correlations between the respective basket
constituents ends up weaker than assumed. Clearly, this will limit the extent
to which we may be able to control the basket smile generated by our model,
leading to erroneous pricing of derivatives away from at-the-money levels.

To illustrate the correlation leakage problem, consider two stock exchange
processes $S_{1}(t)$ and $S_{2}(t)$ driven by the Heston model with correlated
Brownian motions666The respective dynamics are given by ($j=1,2$): ${\rm
d}S_{j}(t)=rS_{j}(t){\rm d}t+v_{j}^{1/2}(t)S_{j}(t){\rm d}W_{j,1}(t)$, ${\rm
d}v_{j}(t)=\kappa_{j}(\bar{v}_{j}-v_{j}(t)){\rm
d}t+\gamma_{j}v_{j}^{1/2}(t){\rm d}W_{j,2}(t)$ with correlations ${\rm
d}W_{j,1}(t){\rm d}W_{j,2}(t)=\rho_{j}{\rm d}t$, ${\rm d}W_{1,1}(t){\rm
d}W_{2,1}(t)=\rho_{1,2}{\rm d}t$ and ${\rm d}W_{j,2}(t){\rm
d}W_{k,2}(t)=0\cdot{\rm d}t$. For reference, we set $S_{1}(t_{0})=1,$
$S_{2}(t_{0})=2.5,$ $r=0,$ $\kappa_{1}=1,$ $\kappa_{2}=0.5,$ $\gamma_{1}=1,$
$\gamma_{2}=0.6,$ $\rho_{S_{1},v_{1}}=-0.5$, $\rho_{S_{2},v_{2}}=-0.7$,
$v_{1,0}=0.1$, $v_{2,0}=0.05$, $\bar{v}_{1}=0.1$ and $\bar{v}_{2}=0.05.$. In
the simulation of the underlying processes, we consider a simple Euler
discretization as in [20]. In Figure 2 the cases for the Feller’s condition
are presented for a range of correlation values. We see that the imposed
correlations are not preserved in time, even when Feller’s condition is
satisfied.

Figure 2: Correlation, $\text{corr}(\log S_{1}(t),\log S_{2}(t))$, as a
function of time. Left: Feller condition satisfied; Right: Feller condition
not satisfied.

As a potential remedy to “leaking correlation”, we could attempt to enforce
the desired dependence patterns via copula. This, however, does not guarantee
success either. To see why, consider the previous example of the Heston model
but now – to make the exposition clearer – assume that the two stocks are
independent, i.e. ${\rm d}W_{1,1}(t){\rm d}W_{2,1}(t)=0\cdot{\rm d}t$. In such
a setup the marginal, cumulative distribution function of each stock can be
recovered using Fourier transform (cf. [20]), and the covariance structure can
then be imposed using a copula. The basket, $B(T)$, at given maturity $T$ is
thus computed via:

$\displaystyle
B(T)=F_{S_{1}(T)}^{-1}(U_{1})+F_{S_{2}(T)}^{-1}(U_{2})=S_{1}(T)+S_{2}(T),$

where the uniformly distributed random variables $U_{1}$ and $U_{2}$ are
joined using a Gaussian copula $\mathcal{C}_{\mathcal{R}}.$ Using such a model
we can now price a call option and compute the corresponding implied
volatility for a range of copula parameters.

Figure 3 shows that varying the copula correlation has a significant impact on
the level of corresponding implied volatilities. However, it is also clear
that while varying the covariance structure alone might help in matching, say
ATM index vols, Gaussian copula is unlikely to generate a desired level of
skewness – a problem we already alluded to in the introduction. One potential
remedy would be to consider more elaborate copula functions (as e.g. proposed
in [17]), however we opt for a different approach. As has already been hinted
above, and as we shall explain in detail below, choosing Merton’s jump
diffusion dynamics [18] for the synthetic variables mapping individual stocks
offers great flexibility in fitting both ATM index vols and skew with
excellent precision.

Figure 3: Implied volatilities for a basket $B(T)$ for $d=2$ driven by the
Heston model.

###  2.2 Basket Reconstruction via Projection with Collocation

As discussed before, our ELV model relies on the availability of the marginal
distribution of each asset, $S_{i}(t)$, at every time $T_{j}$. Such a
distribution may be given explicitly via a specific form of the SDE governing
the evolution of $S_{i}$. Alternatively, the relevant distributions can be
inferred from options prices on $S_{i}$ – either obtained from the market or
generated by some unknown model – leaving us completely agnostic as to the
shape or form of the asset-pricing model. In this section we follow the latter
more general approach and derive the distribution through the following well
known relation:

$\displaystyle F_{S_{i}(T_{j})}(y)=1+\frac{\partial
V_{c}(t_{0},T_{j},S_{i}(t_{0}),K))}{\partial
K}\Big{|}_{K=y},\;\;\;f_{S_{i}(T_{j})}(y)=\frac{\partial^{2}V_{c}(t_{0},T_{j},S_{i}(t_{0}),K))}{\partial
K^{2}}\Big{|}_{K=y},$ (2.5)

where $F_{S_{i}(T_{j})}(y)$ is the CDF of $S_{i}(\cdot)$ observed at time
$T_{j}$ and $f_{S_{i}(T_{j})}(y)$ is the corresponding PDF. We shall take such
distributions as given and ultimately try to reconstruct the marginal
distributions of the basket $B(T_{j})$ defined as the weighted sum of all
individual assets. This step will be ultimately performed by “coupling” the
individual marginals via correleted kernel processes.

However, our first goal is to project each basket constituent, $S_{i}(T_{j})$,
on a polynomial, $g_{i,j}^{m}(\cdot)$, of synthetic variables $X_{i}(T_{j})$,
such as to ensure a perfect fit to the respective marginal distributions for
every set of model parameters $X_{i}(T_{j})$. To this end we resort to the so-
called Stochastic Collocation Monte Carlo Sampler (SCMC) discussed in [11]
(cf. also an overview in A). Using the SCMC method we can express $\log
S_{i}(T_{j})$ as:

$\log S_{i}(T_{j})\approx
g^{m}_{i,j}(X_{i}(T_{j}))=\sum_{k=0}^{m-1}s_{i,j,k}\ell_{k}(x_{i,j,k},X_{i}(T_{j})),$
(2.6)

where $\ell_{k}(\cdot)$ are the Lagrange basis function evaluated at the
collocation points $x_{i,j,k}$, $k=0,\dots,m-1$ based on the kernel variable
$X_{i}(T_{j})$. Alternatively, the representation above can be re-expressed
using a polynomial representation,

$g^{m}_{i,j}(X_{i}(T_{j}))=\sum_{k=0}^{m-1}\hat{\alpha}_{i,j,k}X_{i}^{k}(T_{j})=\hat{\alpha}_{i,j,0}+\hat{\alpha}_{i,j,1}X_{i}(T_{j})+...+\hat{\alpha}_{i,j,m-1}X_{i}^{m-1}(T_{j}),$
(2.7)

where coefficients $\hat{\alpha}_{i,j,k}$ are known explicitly [10]. The
expansion in (2.6) ensures that the CDFs of $S_{i}(T_{j})$ and
$g^{m}_{i,j}(X_{i}(T_{j}))$ agree on the so-called collocation points, i.e.,

$F_{\log
S_{i}(T_{j})}(s_{i,j,k})=F_{g_{m}(X_{i}(T_{j}))}(s_{i,j,k}),\;\;\;k=0,\dots,m-1,$

where $s_{i,j,k}=F_{\log S_{i}(T_{j})}^{-1}(F_{X}(x_{i,j,k}))$ and where
$X_{i}$ is the collocating variable.

The choice regarding variables $X_{i}(T_{j})$, although in principle free, is
crucial for fitting to the marginal distribution of $S_{i}(T_{j})$ and
controlling the basket’s covariance structure. The most appealing choice
regarding variables $X_{i}$ is a normal distribution [11], which would imply
sampling from multivariate Gaussian copula, thus a computationally cheap,
numerical procedure. However, although Gaussian random variables are
sufficient for finding optimal $\hat{\alpha}_{i,j,k}$ in (2.6) guaranteeing
excellent fit to marginal distribution of $S_{i}(T_{j})$ such a model will
generate low skew of the basket and therefore the final calibration to the
index will not be satisfactory (see the remarks in Section 2.1). As presented
later in the article, the non-trivial covariance structure of the underlying
assets is necessary for controlling the implied volatility skew. This is
confirmed by empirical studies [3].

In the ELV model presented in this article, the kernel that is used is
Merton’s jump diffusion model [18], which under the risk-neutral measure
consists of a Brownian motion and a compound Poisson process, which is defined
by:

$\displaystyle
X_{i}(t)=X_{i}(t_{0})+\left(r-\frac{1}{2}\sigma_{i}^{2}-\xi_{p,i}\mathbb{E}[{\rm
e}^{J_{i}}-1]\right)t+\sigma_{i}W_{i}(t)+\sum_{k=1}^{X_{\mathcal{P},i}(t)}J_{k,i},$
(2.8)

with $X_{i}(t_{0})=\log S_{i}(t_{0})$, for $\sigma_{i}>0$, Brownian motion,
$W_{i}(t)$, Poisson process $X_{\mathcal{P},i}(t)$, $t\geq 0$ with parameter
$\xi_{p,i}$ such $\mathbb{E}[X_{\mathcal{P},i}(t)]=\xi_{p,i}t$ and where the
jump sizes $J_{k,i}$ are i.i.d. given by
$J_{k,i}\sim\mathcal{N}\left(\mu_{J,i},\sigma_{J,i}^{2}\right)$. For a given
index $i$ the jumps, $J_{k,i}$, Brownian motion, $W_{i}(t)$, and Poisson
process, $X_{\mathcal{P},i}(t)$ are assumed to be independent.

######  Remark (Alternative dynamics for $X_{i}$).

Alternatively, one may also consider the more advanced process structure
$X_{i}(\cdot)$ by adding stochastic volatility. Here, however, we will focus
on a simplified model, and we will show that such a structure is sufficient
for excellent calibration of the market data, even for large baskets
consisting of $30$ assets.

Once the collocating variable has been chosen, we need to establish the
mapping procedure from the asset $\log S_{i}(T_{j})$ onto the surrogate
variable $X_{i}(\cdot).$ In order to achieve accurate approximations one needs
to firstly determine the so-called collocation points, $x_{i,j}$,
$j=1,\dots,N$, of “surrogate” variable $X_{t}(\cdot)$\- these collocation
points are based on moments of $X_{i}(t)$. This ensures that the collocation
points are the zeros of the orthogonal polynomial corresponding to the
distribution of $X_{i}(\cdot)$, and we can establish the connection with the
computation of integrals by Gauss quadrature. Typically, for $N$-collocation
points $2N$ moments of $X_{i}(\cdot)$ are needed. However, for $X_{i}(\cdot)$
defined in (2.8) one may, due to the presence of jumps, expect an effect often
called “moment explosion” which translates to huge moments causing significant
numerical instability of the method. As an alternative approach, we use the
collocation points $z_{k}$, $k=0,\dots,m-1$ corresponding to the standard
normal random variable777The reason why we choose a standard normal
distribution in the alternative approach is twofold. First, even for a
fundamental distribution as the standard normal results are highly accurate –
this is also the case in e.g. [9]. By choosing a different distribution,
results may be further enhanced. Secondly, as mentioned in [11], choosing the
normal distribution is also motivated by the Cameron-Martin Theorem [22],
which states that polynomial chaos approximations based on the normal
distribution converge to any distribution.  $Z\stackrel{{\scriptstyle{\rm
d}}}{{=}}\mathcal{N}(0,1)$. Such a strategy implies, for a given time $T_{j}$,
the following relation between $\log S_{i}(T_{j})$, $X_{i}(\cdot)$ and $Z$:

$\displaystyle x_{i,j,k}$ $\displaystyle=$ $\displaystyle
F^{-1}_{X_{i}(T_{j})}(F_{Z}(z_{k})),$ (2.9) $\displaystyle\log s_{i,j,k}$
$\displaystyle=$ $\displaystyle F^{-1}_{\log
S_{i}(T_{j})}(F_{X_{i}(T_{j})}(x_{i,j,k})),$ (2.10)

which illustrates how to avoid computation of the collocation points of
$X_{i}(T_{j})$ using moments of $X_{i}(T_{j}).$ It is important to note that
although variable $X_{i}(T_{j})$ is not directly used to compute the
collocation points of $S_{i}(T_{j})$ the approximating function
$g^{m}_{i,j}(X_{i}(T_{j}))$ is explicitly defined in terms of $X_{i}$ and
points $x_{i,j,k}$ and $\log s_{i,j,k}$. The details on the optimality of such
“substitution” can be found in [23].

######  Remark (Grid stretching).

When dealing with heavy-tailed distributions (of leptokurtic type or the
distribution is highly skewed), as here proposed for $X_{i}(T_{j})$, to avoid
numerically unstable inversions $F^{-1}_{X_{i}(T_{j})}(a)$ for either
$a\rightarrow 0$ or $a\rightarrow 1$ it is recommended to use the so-called
“grid stretching” technique [11] which facilitates stable inversions,
especially for a high number of the collocation points, $m$.

Availability of the ChF in closed form allows finding, via Fourier inversion
[20], the corresponding CDF, $F_{X_{i}(T_{j})}(x)$, and thus finding of the
corresponding collocation points, $x_{i,j,k}$.

$\phi_{X_{i}(t)}(u)=\mathbb{E}\left[{\rm
e}^{iuX_{i}(t)}\right]=\exp\left(iu(X_{i}(t_{0})+\bar{\mu}_{i}t)-\frac{1}{2}\sigma_{i}^{2}u^{2}t+\xi_{p,i}t{\rm
e}^{iu\mu_{J,i}-\frac{1}{2}u^{2}\sigma_{J,i}^{2}}\right),$ (2.11)

with $\bar{\mu}_{i}=r-\frac{1}{2}\sigma_{i}^{2}-\xi_{p,i}({\rm
e}^{\mu_{J,i}+\frac{1}{2}\sigma_{J,i}^{2}}-1).$

Then, by utilizing the COS method [20], one is able to determine the
corresponding CDF:

$\displaystyle
F_{X_{i}(T_{j})}(y)\approx{\sum_{k=0}^{N_{c}-1}}\bar{F}_{k}\psi(a,b,y),\;\;\;\bar{F}_{k}=\frac{2}{b-a}\Re\left\\{\phi_{X_{i}(T_{j})}\left(\frac{k\pi}{b-a}\right)\cdot\exp\left(-i\frac{ka\pi}{b-a}\right)\right\\},$
(2.12)

with

$\psi(a,b,y)=\left\\{\begin{array}[]{ccc}\frac{b-a}{k\pi}\sin\left(\frac{k\pi(y-a)}{b-a}\right),&\text{for}&k=1,2,\dots,N_{c}-1,\\\
(x-a),&\text{for}&k=0,\end{array}\right.$ (2.13)

where $N_{c}$ indicates the number of the expansion terms, $a$ and $b$ are the
domain parameters, typically determines in terms of cummulants.

In the next step, we build a multivariate dependence structure between all the
marginal distributions using a system of correlated SDEs. The basket,
$B(T_{j})$, yields, in terms of the latent variables, the following form:

$\displaystyle B(T_{j})$ $\displaystyle=$
$\displaystyle\sum_{i=1}^{N}S_{i}(T_{j})\approx\sum_{i=1}^{N}\exp\left(g^{m}_{i,j}(X_{i}(T_{j}))\right)$
(2.14) $\displaystyle=$
$\displaystyle\sum_{i=1}^{N}\exp\left(\sum_{k=0}^{m-1}\hat{\alpha}_{i,j,k}X_{i}^{k}(T_{j})\right),$

with $X_{i}(T_{j})$ defined in (2.8). Although we have initially assumed that
each $X_{i}$ has its own independent jumps $J_{k,i}$ controlled by
$X_{\mathcal{P}_{i}(t)}$, our numerical studies have shown that without
sacrificing much in terms of skew fitting possibilities, one can choose a
sparse for of the model with $X_{\mathcal{P}_{i}(t)}$ being the same for all
the assets, i.e.,

$X_{\mathcal{P}_{i}(t)}\equiv X_{\mathcal{P}(t)}.$

In such form, the dependence between the individual stocks in the basket is
controlled via two main elements: correlated Brownian motions $W_{i}(T_{j})$
and common Poisson process $X_{\mathcal{P}}(T_{j})$ for all Merton’s
processes, $X_{i}.$ Thus, the basket is reconstructed using a combination of
marginal distributions for each asset. The basket PDF, $f_{B(t)}$, can then be
obtained by Monte Carlo. It should be stressed that the sampling from
multivariate normal is computationally inexpensive. The expensive part
involves the inversion of marginal distributions, but as explained above, we
handle this using collocation, mitigating the computational cost quite
considerably.

Importantly, the proposed methodology does not rely on a calibration of the
process for marginal distributions, but solely on the computation of
$\hat{\alpha}_{i,j,k}$, for $i=1,\dots,N$, $j=1,\dots,N_{T},$ $k=1,\dots,m$,
and inexpensive samples from the multidimensional distribution of
$(X_{i}(\cdot),\dots,X_{N}(\cdot))$.

######  Remark (Marginal distribution).

We would like to stress that due to the collocation method and the CDF mapping
procedure, determined coefficients $\hat{\alpha}_{i,j,k}$ will ensure
agreement between the target CDF, $F_{\log S_{i}(T_{j})}(\cdot)$ and the
surrogate, $F_{X_{i}(T_{j})}(\cdot).$ Once the marginal distributions are
calibrated for any configurations of model parameters of $X_{i}(\cdot)$, the
smile/skew of the index needs to be matched. Later, in this article, we will
show that the model proposed allows for intuitive control of different
volatility shapes of the basket, therefore facilitating calibration to the
index.

##  3 Basket Calibration to the Index

Given that the ELV model, by its construction, guarantees fit to marginal
distributions one still needs to calibrate the basket to index volatilities.
The brute force strategy of such optimization can be performed using Monte
Carlo simulation, i.e., for different model parameters, the samples of
$X_{i}(T_{j})$ are fed to basket equation in (2.14) and the option prices are
computed. Such routine is typically sub-optimal as it would require multiple
iterations and samplings for $X_{i}(\cdot).$ To speed up the process and save
computational burden, we opt for a different approach, based on moment
matching.

Recall that in its most general form the model might appear to have an
overwhelming number of parameters to be calibrated. Fortunately, as we have
already hinted above, many simplifications can be resorted to without
sacrificing calibration precision, while improving tractability and sparsity.

Firstly, we note that the ELV model’s heart lies in the idea of construction
of the local volatility function; therefore, only the marginal distribution
for the basket $B(T_{j})$ is relevant. This implies that the calibration can
be performed sequentially, or parallel, for every time $T_{j}$,
$j=1,\dots,N_{T}.$ Secondly, it is in practice safe to assume that the
individual stocks share not only the type of the kernel process but also the
driving parameters. Specifically, we assume henceforth that all the assets
have a common counting process $X_{\mathcal{P}}(T_{j})$, independent i.i.d.
jumps $J\equiv J_{k,i}\sim\mathcal{N}(\mu_{J},\sigma^{2}_{J})$ and correlated
Brownian motions ${\rm d}W_{i}(t){\rm d}W_{j}(t)=\rho{\rm d}t,$ which leads to
the following simplified representation:

$\displaystyle
X_{i}(T_{j})=X_{i}(t_{0})+\big{(}r-\frac{1}{2}\sigma^{2}-\xi_{p}\mathbb{E}[{\rm
e}^{J}-1]\big{)}T_{j}+\sigma
W_{i}(T_{j})+\sum_{k=1}^{X_{\mathcal{P}}(T_{j})}J_{k,i}.$

Thus, in order to calibrate the basket $B(T_{j})$ to an index we need, for
every $T_{j}$, to determine the set of optimal model parameters $\sigma$,
$\xi_{p}$, $\mu_{J}$ and $\sigma_{J}$, driving the model’s covariance
structure.

Each of the jump parameters has a different effect on the shape of the implied
volatility curve, i.e., $\sigma_{J}$ has a significant impact on the
curvature, $\xi_{\tiny p}$ controls the overall level of the implied
volatility, whereas $\mu_{J}$ influences the implied volatility slope (the
skew).

Figure 4: Impact of different jump parameters on the shape of the implied
volatility in Merton’s jump diffusion model. First figure: impact of
$\sigma_{J}$; Second figure: impact of $\xi_{p}$; Third figure: impact of
$\mu_{J}.$

Having already gained much in terms of tractability, we can further save the
computational burden by performing calibration via the moment matching
procedure. Given that the density of the index, $f_{I}(y)$, exists and can be
calculated from options prices, as per Equation (2.5), the objective is to
match moments from the index,$I(t)$, and the basket $B(t)$. It is important to
note that the polynomial coefficients $\alpha_{i,j,k}$ are implicitly a
function of model parameters, i.e.,

$\displaystyle\alpha_{i,j,k}:=\alpha_{i,j,k}(\sigma,\xi_{p},\mu_{J},\sigma_{J}),$
(3.1)

therefore any change in the model parameters requires re-computation of these
coefficients. Since this procedure only relies on the computation of a few
inversions of the CDF in Equation (2.12), it is computationally cheap.

Since the first moment is trivial let us start with the second moment. By
definition we have:

$\displaystyle\mathbb{E}[B^{2}(T_{j})]$ $\displaystyle=$
$\displaystyle\sum_{i_{1}=1}^{N}\sum_{i_{2}=1}^{N}\omega_{i_{1}}\omega_{i_{2}}\mathbb{E}[S_{i_{1}}(T_{j})S_{i_{2}}(T_{j})].$

One the other hand, assuming that we are able to determine the correlation
$\hat{\rho}(S_{i_{1}},S_{i_{2}})$ between $S_{i_{1}}(T_{j})$ and
$S_{i_{2}}(T_{j})$ (this correlation will differ from $\rho_{i_{1},i_{2}}{\rm
d}t={\rm d}W_{i_{1}}(t){\rm d}W_{i_{2}}(t)$ between the corresponding Brownian
motions), we find:

$\hat{\rho}(S_{i_{1}},S_{i_{2}})\stackrel{{\scriptstyle\text{def}}}{{=}}\frac{\mathbb{E}[S_{i_{1}}(T_{j})S_{i_{2}}(T_{j})]-\mathbb{E}[S_{i_{1}}(T_{j})]\mathbb{E}[S_{i_{2}}(T_{j})]}{\sigma_{{i_{1}}}\sigma_{{i_{2}}}},$

where $\sigma_{i,1}$ and $\sigma_{i,2}$ are standard deviations of
$S_{i_{1}}(T_{j})$ and $S_{i_{2}}(T_{j})$ respectively (these quantities can
be calculated directly using (2.5)). Given that the marginal distributions of
$S_{i_{1}}$ and $S_{i_{2}}$ are available, both: first moments and the
standard deviations can be easily obtained, by for example proper integration
of the market-implied densities in (2.5). Then, we have:

$\displaystyle\mathbb{E}[S_{i_{1}}(T_{j})S_{i_{2}}(T_{j})]=\hat{\rho}(S_{i_{1}},S_{i_{2}})\sigma_{{i_{1}}}\sigma_{{i_{2}}}+\mathbb{E}[S_{i_{1}}(T_{j})]\mathbb{E}[S_{i_{2}}(T_{j})],$
(3.2)

therefore the second moment of the basket becomes:

$\displaystyle\mathbb{E}[B^{2}(T_{j})]$ $\displaystyle=$
$\displaystyle\sum_{i_{1}=1}^{N}\sum_{i_{2}=1}^{N}\omega_{i_{1}}\omega_{i_{2}}\Bigg{(}\hat{\rho}(S_{i_{1}},S_{i_{2}})\sigma_{{i_{1}}}\sigma_{{i_{2}}}+\mathbb{E}[S_{i_{1}}(T_{j})]\mathbb{E}[S_{i_{2}}(T_{j})]\Bigg{)}.$
(3.3)

Since the correlation, $\hat{\rho}(S_{i_{1}},S_{i_{2}})$, measures the linear
relationship between assets $S_{i_{1}}$ and $S_{i_{2}}$, we will compute this
coefficient based on a fewer number of the collocation points 888The strategy
proposed in this part does not require “re-calibration” of $\alpha_{i,j,k}$
coefficients, but only neglects the coefficients of higher order., $m=2$ which
yields the following representation:

$\displaystyle\hat{\rho}(S_{i_{1}},S_{i_{2}})$ $\displaystyle\approx$
$\displaystyle\rho\left({\rm
e}^{\hat{\alpha}_{i_{1},j,0}+\hat{\alpha}_{i_{1},j,1}X_{i_{1}}(T_{j})},{\rm
e}^{\hat{\alpha}_{i_{2},j,0}+\hat{\alpha}_{i_{1},j,1}X_{i_{2}}(T_{j})}\right)$
$\displaystyle=$ $\displaystyle\rho\left({\rm
e}^{\hat{\alpha}_{i_{1},j,1}X_{i_{1}}(T_{j})},{\rm
e}^{\hat{\alpha}_{i_{2},j,1}X_{i_{2}}(T_{j})}\right)=:\bar{\rho}(S_{i_{1}},S_{i_{2}}),$

with $X_{i}(\cdot)$ defined in (2.8), with correlated Brownian motions,

$\displaystyle\bar{\rho}(S_{i_{1}},S_{i_{2}})$
$\displaystyle\stackrel{{\scriptstyle\text{def}}}{{=}}$
$\displaystyle\frac{\mathbb{E}[{\rm
e}^{\hat{\alpha}_{i_{1},j,1}X_{i_{1}}(T_{j})+\hat{\alpha}_{i_{2},j,1}X_{i_{2}}(T_{j})}]-\mathbb{E}[{\rm
e}^{\hat{\alpha}_{i_{1},j,1}X_{i_{1}}(T_{j})}]\mathbb{E}[{\rm
e}^{\hat{\alpha}_{i_{2},j,1}X_{i_{2}}(T_{j})}]}{\sqrt{\mathbb{V}\text{ar}[{\rm
e}^{\hat{\alpha}_{i_{1},j,1}X_{i_{1}}(T_{j})}]\mathbb{V}\text{ar}[{\rm
e}^{\hat{\alpha}_{i_{2},j,1}X_{i_{2}}(T_{j})}]}}.$

We immediately note the connection between the moments and the corresponding
ChF of $X_{i}(t)$ defined in (2.11) resulting in,

$\displaystyle\phi_{X_{i}(t)}(u)=\mathbb{E}[{\rm
e}^{iuX_{i}(t)}],\;\;\;\text{thus}\;\;\;\mathbb{E}[{\rm
e}^{aX_{i}(t)}]=\phi_{X_{i}(t)}(-ia),\;\;a\in\mathbb{R},$

resulting in the following representation

$\displaystyle\bar{\rho}(S_{i_{1}},S_{i_{2}})=\frac{\mathbb{E}[{\rm
e}^{c_{1}X_{i_{1}}(T_{j})+c_{2}X_{i_{2}}(T_{j})}]-\phi_{X_{i_{1}}}(-ic_{1})\phi_{X_{i_{2}}}(-ic_{2})}{\sqrt{\left(\phi_{X_{i_{1}}}(-2ic_{1})-\phi_{X_{i_{1}}}^{2}(-ic_{1})\right)\left(\phi_{X_{i_{2}}}(-2ic_{2})-\phi_{X_{i_{2}}}^{2}(-ic_{2})\right)}},$
(3.5)

with $c_{1}=\hat{\alpha}_{i_{1},j,1}$, $c_{2}=\hat{\alpha}_{i_{2},j,1}$ and
$\phi_{X_{i}(T_{j})}$ defined in Equation 2.11. The only term that caries the
correlation between both processes is the expectation involving correlated
Merton’s processes. This expectation is given in closed form and it is
presented in Lemma 3.1 below.

######  Lemma 3.1 (Expectation for a 2D Merton’s model).

For a given 2D Merton’s model defined as:

$\displaystyle
X_{1}(t)=\mu_{1}+\sigma_{1}W_{1}(t)+\sum_{k=1}^{X_{\mathcal{P}}(t)}J_{k,1},\;\;\;X_{2}(t)=\mu_{2}+\sigma_{2}W_{2}(t)+\sum_{k=1}^{X_{\mathcal{P}}(t)}J_{k,2},$
(3.6)

with correlated Brownian motions ${\rm d}W_{1}(t){\rm
d}W_{2}(t)=\rho_{1,2}{\rm d}t,$ common counting Poisson process
$X_{\mathcal{P}}(t)$ with the intensity parameter $\xi_{p}$ and independent
identically distributed jumps, $J\equiv
J_{k,\cdot}\sim\mathcal{N}(\mu_{J},\sigma_{J}^{2})$, the following expression
holds:

$\displaystyle\mathbb{E}[{\rm e}^{aX_{1}(T)+bX_{2}(T)}|\mathcal{F}(t_{0})]$
$\displaystyle=$
$\displaystyle\exp\left[a\mu_{1}+b\mu_{2}+\frac{1}{2}\hat{\sigma}^{2}T+\xi_{p}T\left({\rm
e}^{(a+b)\mu_{J}+\frac{1}{2}(a^{2}+b^{2})\sigma_{J}^{2}}-1\right)\right]$
(3.7) $\displaystyle=:$ $\displaystyle\omega_{X}(a,b),$

where
$\hat{\sigma}^{2}=a^{2}\sigma_{1}^{2}+b^{2}\sigma_{2}^{2}+2\rho_{1,2}ab\sigma_{1}\sigma_{2}$.

###### Proof.

Proof can be found in B. ∎

Utilizing the results above we find the closed form solution for the
expectation and therefore the correlation coefficient,
$\bar{\rho}(S_{i_{1}},S_{i_{2}})$, by simply setting where
$c_{1}=\hat{\alpha}_{i_{1},j,1}$, $c_{2}=\hat{\alpha}_{i_{2},j,1}$,
$\mu_{i}=X_{i}(t_{0})+\big{(}r-\frac{1}{2}\sigma^{2}-\xi_{p}\mathbb{E}[{\rm
e}^{J}-1]\big{)}T$, and $\sigma_{1}=\sigma_{2}=\sigma.$

Once the correlation coefficient, $\bar{\rho}(S_{i_{1}},S_{i_{2}})$, is
determined one is able to perform the ATM calibration of the basket to the
index.

$\displaystyle\min_{\sigma_{J},\mu_{J},\xi_{p},\rho}\left(\mathbb{E}[B^{2}(T_{j})]-\mathbb{E}[I^{2}(T_{j})]\right)^{2},\;\;\;\mathbb{E}[I^{2}(T_{j})]=\int_{\mathbb{R}}x^{2}f_{I(T_{j})}(x){\rm
d}x,$ (3.8)

where the index density $f_{I(T_{j})}(x)$ is implied from the option quotes
for the index (2.5) and $\mathbb{E}[B^{2}(T_{j})]$ is defined in (3.3). It is
important to note that, each change of the model parameters of the kernel
process does impact parameters $\alpha_{i_{1},j,1}$ and $\alpha_{i_{2},j,1}$,
implying that the mapping coefficients, $\alpha_{\cdot}$, need to be
recomputed at each iteration. This however is extremely cheap operation as it
only requires the CDF mappings presented in (2.9).

Following a similar strategy a third moment for the basket can be derived. The
detailed derivations for the approximation results are presented in C. By
definition of the third moment we have:

$\displaystyle\mathbb{E}[B^{3}]$ $\displaystyle=$
$\displaystyle\sum_{i_{1}=1}^{N}\sum_{i_{2}=1}^{N}\sum_{i_{3}=1}^{N}\omega_{i_{1},i_{2},i_{3}}\Big{[}\hat{\rho}({S_{i_{1},i_{2}},S_{i_{3}}})\sigma_{i_{1},i_{2}}\sigma_{{i_{3}}}+\mathbb{E}[S_{i_{1},i_{2}}]\mathbb{E}[S_{i_{3}}]\Big{]}$
$\displaystyle=$
$\displaystyle\sum_{i_{1}=1}^{N}\sum_{i_{2}=1}^{N}\sum_{i_{3}=1}^{N}\omega_{i_{1},i_{2},i_{3}}\left[\hat{\rho}({S_{i_{1},i_{2}},S_{i_{3}}})\sigma_{{i_{1},i_{2}}}\sigma_{{i_{3}}}+\Big{(}\hat{\rho}(S_{i_{1}},S_{i_{2}})\sigma_{{i_{1}}}\sigma_{{i_{2}}}+\mathbb{E}[S_{i_{1}}]\mathbb{E}[S_{i_{2}}]\Big{)}\mathbb{E}[S_{i_{3}}]\right],$

with $\omega_{i_{1},i_{2},i_{3}}:=\omega_{i_{1}}\omega_{i_{2}}\omega_{i_{3}}$,
and $S_{i_{1},i_{2}}:=S_{i_{1}}S_{i_{2}}$. The first moments
$\mathbb{E}[S_{i_{1}}]$, $\mathbb{E}[S_{i_{2}}]$ and $\mathbb{E}[S_{i_{3}}]$
are known explicitly using the market market data for individual assets.
Standard deviations, $\sigma_{i_{1}}$ and $\sigma_{i_{2}}$, of assets
$S_{i_{1}}(T_{j})$ and $S_{i_{2}}(T_{j})$ can be computed utilizing the market
implied density in (2.5). Correlation,
$\hat{\rho}({S_{i_{1},i_{2}},S_{i_{3}}})$, is approximated, as for the second
moment, and yields:

$\displaystyle\hat{\rho}({S_{i_{1},i_{2}},S_{i_{3}}})\approx\frac{\omega_{X}(c_{1},c_{2},c_{3})-\omega_{X}(c_{1},c_{2})\phi_{X_{i_{3}}}(-ic_{3})}{\sqrt{\left(\omega_{X}(2c_{1},2c_{2})-\omega_{X}^{2}(c_{1},c_{2})\right)\left(\phi_{X_{i_{3}}}(-2ic_{3})-\phi_{X_{i_{3}}}^{2}(-ic_{3})\right)}},$

where: $c_{1}:=\hat{\alpha}_{i_{1},j,1}$, $c_{2}:=\hat{\alpha}_{i_{2},j,1}$,
$c_{3}:=\hat{\alpha}_{i_{3},j,1}$ with $\omega_{X}(c_{1},c_{2},c_{3})$ defined
in (C) and $\omega_{X}(c_{1},c_{2})$ in (3.7). Finally, the variance
$\sigma_{{i_{1},i_{2}}}^{2}$ is given explicitly by:

$\displaystyle\sigma_{{i_{1},i_{2}}}^{2}$ $\displaystyle:=$
$\displaystyle\mathbb{V}\text{ar}[S_{i_{1}}S_{i_{2}}]=\mathbb{E}[S_{i_{1}}S_{i_{2}}]-\mathbb{E}[S_{i_{1}}]\mathbb{E}[S_{i_{2}}]$
(3.9) $\displaystyle=$
$\displaystyle\hat{\rho}(S_{i_{1}},S_{i_{2}})\sigma_{{i_{1}}(T_{j})}\sigma_{{i_{2}}(T_{j})},$

where $\hat{\rho}(S_{i_{1}},S_{i_{2}})$ is given in (3) and it is approximated
by $\bar{\rho}(S_{i_{1}},S_{i_{2}})$ and derived in (3.5). Availability of the
third moment allows for an extended optimization problem proposed in (3.8),
i.e.,

$\displaystyle\min_{\sigma_{J},\mu_{J},\xi_{p},\rho}\sum_{k=2}^{3}\left(\mathbb{E}[B^{k}(T_{j})]-\mathbb{E}[I^{k}(T_{j})]\right)^{2},\;\;\;\mathbb{E}[I^{k}(T_{j})]=\int_{\mathbb{R}}x^{k}f_{I(T_{j})}(x){\rm
d}x,$ (3.10)

which will not only allow for the ATM fit but also the level of the implied
volatility skew can be matched.

The quality of the approximating moment formulae will be discussed in the
follow-up section (see Table 5), where a basket consisting of 5 assets will be
considered.

##  4 Numerical Results

This numerical section is dedicated to the numerical aspects and
implementation details of the proposed model. We start with a step-by-step
example outlining the procedure of constructing a 5-dimensional basket
including the calibration of the collocation method to market-implied
volatilities. We follow up with a discussion of the impact of model parameters
on the resulting basket implied volatility and a finally conclude with a
fully-fledged calibration of a high-dimensional basket to market data.

###  4.1 Illustrative 5D Example

This section presents a step-by-step procedure for building up a basket using
the ELV method. We illustrate in detail the complete process, starting from
individual asset calibration to basket construction and the impact of the
model parameters on the basket. In this example, we consider a basket
consisting of 5 assets 9991) UnitedHealth; 2) Home Depot; 3) Goldman Sachs; 4)
Microsoft Corp; 5) salesforce.com Inc that are a part of the DJIA 30. In the
experiment, we take $T_{j}=1y$ as of 8/12/2021. The implied volatilities for
this set of assets are presented in Table 1.

Table 1: Implied volatilities for 5 assets, for $T=1$, observed on 8/12/2021. Spot values for the underlying assets are: $S_{1}(t_{0})=468.86$, $S_{2}(t_{0})=411.25$, $S_{3}(t_{0})=397.32$, $S_{3}(t_{0})=334.97$, $S_{5}(t_{0})=266.31$ asset | 0.8 | 0.85 | 0.9 | 0.95 | 0.975 | 1 | 1.025 | 1.05 | 1.1 | 1.15 | 1.2
---|---|---|---|---|---|---|---|---|---|---|---
stock 1 | 32.20 | 31.00 | 29.98 | 29.07 | 28.65 | 28.25 | 27.87 | 27.51 | 26.85 | 26.28 | 25.83
stock 2 | 29.96 | 28.87 | 27.99 | 27.31 | 27.03 | 26.79 | 26.58 | 26.41 | 26.15 | 25.98 | 25.90
stock 3 | 32.23 | 31.07 | 30.11 | 29.34 | 29.02 | 28.73 | 28.48 | 28.26 | 27.90 | 27.65 | 27.47
stock 4 | 31.86 | 30.56 | 29.39 | 28.35 | 27.88 | 27.43 | 27.02 | 26.64 | 25.99 | 25.47 | 25.09
stock 5 | 34.58 | 33.84 | 33.25 | 32.76 | 32.55 | 32.35 | 32.17 | 32.00 | 31.71 | 31.46 | 31.26

As the first step of our modelling procedure, we parametrize the market-
implied volatilities using the SABR-based parametrization formula [13]. The
parameters obtained in this step are presented in Table 2.

Table 2: Implied volatility parametrization using parameterization in [13]. The meaning of the parameters is as follows: $\beta_{H}$ -“exponent beta”, $\alpha_{H}$\- “initial vol”, $\rho_{H}$\- “correlation” and $\gamma_{H}$\- “vol-vol”. In the calibration procedure $\beta_{H}$ was fixed at 0.9. Interest rates were fixed at $r=0.0.$ asset | $\beta_{H}$ | $\alpha_{H}$ | $\rho_{H}$ | $\gamma_{H}$
---|---|---|---|---
stock 1 | 0.90 | 0.52 | -0.52 | 0.55
stock 2 | 0.90 | 0.48 | -0.25 | 0.65
stock 3 | 0.90 | 0.51 | -0.30 | 0.66
stock 4 | 0.90 | 0.49 | -0.48 | 0.66
stock 5 | 0.90 | 0.56 | -0.27 | 0.46

Once each implied volatility skew/smile is parameterized the next step is use
(2.5) and using the collocation method project the so-called “market
distribution” on $X_{i}$. We set the number of the collocation points to
$m=7$, and as discussed in Remark in Section 2.2 we use the grid stretching
technique to control the tail behaviour. Firstly, we take
$Z~{}\sim\mathcal{N}(0,1)$ and we set the corresponding min/max quantiles at
$q_{1}=0.01$ and $q_{m=7}=0.99$. Once the corresponding collocation points
$z_{k}$, $k=1,\dots,m$ are determined [11] the inverse for $x_{i,k}$ in (2.9)
need to take place. Note that in principle each asset, $S_{i}(T_{j})$, may be
projected on the same variable $X_{i}.$ Even if $X_{i}$’s are correlated (via
Brownian motion) it does not affect marginal distributions for $S_{i}$, but it
will be crucial in the basket construction. Table 3 illustrates the details on
computation of triples $(z_{i,k},x_{i,k},s_{i,k})\equiv(z_{k},x_{k},s_{i,k})$,
$k=1,\dots,7,$ $i=1,\dots,5.$

Table 3: Projection details for the collocation method. Parameters for $X_{i}$ were the following: $\sigma=0.1$, $\sigma_{J}=0.02$, $\mu_{J}=0.0$, $\xi_{p}=0.25$ and $T=0.5$. Note that $T$ here is considered as a parameter therefore it can be fixed for all the assets. In the experiment we have also taken the stretched grid with $p=0.01$, as proposed in [11], implying the the lower and upper bound for the CDF, $F_{Z}(z_{1})$ and $F_{Z}(z_{7})$.  | $k=1$ | $k=2$ | $k=3$ | $k=4$ | $k=5$ | $k=6$ | $k=7$
---|---|---|---|---|---|---|---
$z_{k}$ | -3.7504 | -2.3668 | -1.1544 | 0 | 1.154 | 2.3668 | 3.7504
$F_{Z}(z_{k})$ | 0.0100 | 0.0710 | 0.237 | 0.5000 | 0.7630 | 0.929 | 0.9900
$x_{k}$ | -0.1684 | -0.1074 | -0.0539 | -0.0030 | 0.0478 | 0.1013 | 0.1623
$s_{1,k}$ | 4.9851 | 5.6027 | 5.9341 | 6.1545 | 6.3270 | 6.4891 | 6.6839
$s_{2,k}$ | 4.9591 | 5.5291 | 5.8174 | 6.0101 | 6.1734 | 6.3512 | 6.6100
$s_{3,k}$ | 4.7692 | 5.4420 | 5.7672 | 5.9771 | 6.1500 | 6.3350 | 6.6025
$s_{4,k}$ | 4.5586 | 5.2737 | 5.6126 | 5.8221 | 5.9828 | 6.1404 | 6.3508
$s_{5,k}$ | 4.4752 | 4.9965 | 5.3177 | 5.5592 | 5.7702 | 5.9855 | 6.2538

Given the the projection points presented in Table 3 and the underlying
process $X_{i}(\cdot)$ we are able to replicate the marginal distribution
functions of the assets observable in the market. By simulating a process
$X_{i}(T_{j})$ for every asset $i=1,\dots,N$ with the parameters specified in
Table 3 we evaluate the approximating polynomial for $m=7$ with $T=1y,$

$\log S_{i}(T)\approx
g^{7}_{i}(X_{i}(T))=\sum_{k=0}^{m-1}s_{i,k}\ell_{k}(x_{k},X_{i}(T)),\;\;\;\ell_{k}(x_{k},X_{i}(T))=\prod_{l=1,\\\
l\neq k}^{m}\frac{X_{i}(T)-x_{l}}{x_{k}-x_{l}},$ (4.1)

with points $x_{k}$ and $s_{i,k}$ defined in Table 3 and $X_{i}(\cdot)$
defined in (2.8). Ultimately, all the kernel processes $X_{i}$, $X_{j}$ need
to be correlated with a given correlation coefficient. This, however, is not
necessary if we are only interested in marginal distributions, thus in a fit
of each asset to market implied volatility.

The calibration quality for the considered five assets is presented in Table
4. Again, we report excellent results; the error does not exceed $0.2\%$.

Table 4: Table illustrates the calibration of the collocation method to market-implied volatilities. “mkt” stands for the market implied volatility obtained from the parameterization of implied volatility, “coll” represents the results obtained from the collocation method and “err” is the difference of implied volatilities, in $\%$.  | relative strike |  |
---|---|---|---
asset | IV % | 0.8 | 0.85 | 0.9 | 0.95 | 0.975 | 1 | 1.025 | 1.05 | 1.1 | 1.15 | 1.2
$S_{1}$ | mkt | 32.2 | 31.0 | 30.0 | 29.1 | 28.7 | 28.3 | 27.9 | 27.5 | 26.9 | 26.3 | 25.8
| coll | 32.1 | 31.0 | 30.0 | 29.0 | 28.6 | 28.2 | 27.8 | 27.4 | 26.8 | 26.2 | 25.7
| err | 0.1 | 0.0 | 0.0 | 0.1 | 0.1 | 0.1 | 0.1 | 0.1 | 0.1 | 0.1 | 0.1
$S_{2}$ | mkt | 29.9 | 28.9 | 28.0 | 27.3 | 27.0 | 26.8 | 26.6 | 26.4 | 26.1 | 26.0 | 25.9
| coll | 29.9 | 28.9 | 28.0 | 27.3 | 27.0 | 26.8 | 26.5 | 26.4 | 26.1 | 25.9 | 25.9
| err | 0.0 | 0.0 | 0.0 | 0.0 | 0.0 | 0.0 | 0.0 | 0.0 | 0.0 | 0.0 | 0.1
$S_{3}$ | mkt | 32.2 | 31.1 | 30.2 | 29.4 | 29.0 | 28.7 | 28.5 | 28.2 | 27.9 | 27.6 | 27.5
| coll | 32.0 | 31.0 | 30.0 | 29.3 | 28.9 | 28.7 | 28.4 | 28.2 | 27.8 | 27.6 | 27.4
| err | 0.2 | 0.1 | 0.1 | 0.1 | 0.1 | 0.1 | 0.1 | 0.1 | 0.1 | 0.1 | 0.1
$S_{4}$ | mkt | 31.9 | 30.6 | 29.4 | 28.4 | 27.9 | 27.4 | 27.0 | 26.6 | 26.0 | 25.5 | 25.1
| coll | 31.7 | 30.4 | 29.2 | 28.2 | 27.7 | 27.3 | 26.9 | 26.5 | 25.9 | 25.4 | 25.0
| err | 0.2 | 0.2 | 0.2 | 0.1 | 0.1 | 0.1 | 0.1 | 0.1 | 0.1 | 0.1 | 0.1
$S_{5}$ | mkt | 34.5 | 33.9 | 33.3 | 32.8 | 32.6 | 32.4 | 32.2 | 32.0 | 31.7 | 31.5 | 31.3
| coll | 34.6 | 33.9 | 33.3 | 32.8 | 32.6 | 32.4 | 32.2 | 32.0 | 31.7 | 31.5 | 31.3
| err | -0.1 | 0.0 | 0.0 | 0.0 | 0.0 | 0.0 | 0.0 | 0.0 | 0.0 | 0.0 | 0.0

As the final step of this illustrative example we construct a basket,
$B(T)=\sum_{i=1}^{5}S_{i}(T)$ with $S_{i}(t)$ given in (4.1). For this
purpose, we utilize Equation (2.14), where the processes $X_{i}$ need to be
simulated with Monte Carlo. As discussed before, we consider a standard
Poisson process, $X_{\mathcal{P}}(T)$, for all the underlying assets and
correlated Brownian Motions, $W(t).$ These processes can be pre-simulated and
then fed into the asset Equation (4.1). This procedure will produce Monte
Carlo paths for the basket, $B(t)$.

The next step is to determine optimal kernel model parameters. The objective
is to find the parameters that will generate basket implied volatility as
close as possible to the implied volatility of the index. In this illustrative
experiment, we do not consider an index; however, we can analyze the quality
of the approximating formulae and the impact of the model parameters on the
basket volatilities.

In the calibration procedure, we will follow Equation (3.10), which relies on
fitting based on the variance and the skew. In Table 5 the numerical results
for estimated variance and the skew for varying model parameters are shown.
The results are excellent, especially for the standard deviation,
$\sigma_{B}$, where the error is about $0.1$.

Table 5: Table illustrates the quality of the approximation of the moments derived for a considered basket $B(T)$, $T=1$, and varying set of the model parameters. The first 5 columns of the table illustrate particular parameter configuration (dots indicate unchanged parameters). $\sigma_{B}:=\sqrt{\mathbb{V}\text{ar}[B(T)]}$.  parameters | Monte Carlo | Approximation
---|---|---
$\sigma$ | $\sigma_{J}$ | $\mu_{J}$ | $\xi_{p}$ | $\rho$ | $\sigma_{B}$ | $\mathbb{E}[B^{3}]/\sigma_{B}^{3}$ | $\sigma_{B}$ | $\mathbb{E}[B^{3}]/\sigma_{B}^{3}$
0.2 | 0.02 | 0.0 | 0.25 | 0.5 | 28.17 | 98.42 | 28.1 | 90.84
0.2 | $\cdot$ | $\cdot$ | $\cdot$ | $\cdot$ | 28.19 | 98.17 | 28.17 | 90.17
0.5 | $\cdot$ | $\cdot$ | $\cdot$ | $\cdot$ | 28.18 | 98.25 | 28.19 | 89.95
$\cdot$ | 0.01 | $\cdot$ | $\cdot$ | $\cdot$ | 28.19 | 98.21 | 28.18 | 90.02
$\cdot$ | 0.04 | $\cdot$ | $\cdot$ | $\cdot$ | 27.88 | 101.25 | 27.68 | 94.85
$\cdot$ | $\cdot$ | -0.1 | $\cdot$ | $\cdot$ | 29.77 | 84.48 | 29.61 | 78.13
$\cdot$ | $\cdot$ | -0.25 | $\cdot$ | $\cdot$ | 31.18 | 74.23 | 31.86 | 63.47
$\cdot$ | $\cdot$ | $\cdot$ | 0.1 | $\cdot$ | 28.2 | 98.11 | 28.17 | 90.13
$\cdot$ | $\cdot$ | $\cdot$ | 0.15 | $\cdot$ | 28.19 | 98.22 | 28.14 | 90.47
$\cdot$ | $\cdot$ | $\cdot$ | $\cdot$ | 0.1 | 19.44 | 277.62 | 19.33 | 269.78
$\cdot$ | $\cdot$ | $\cdot$ | $\cdot$ | 0.8 | 33.31 | 62.87 | 33.31 | 56.03

In the final part of this section, we illustrate, in Table 6, the timing
results. As presented in the table, the most time-consuming part is to
parameterize market discrete implied volatilities into a parametric form.
Typically, such a process is performed only once, and only the relevant,
optimized parameters are stored. The collocation method’s calibration and the
collocation points’ generation is a considerably cheap operation. Finally, the
Monte Carlo simulation required sampling of all the underlying processes and
construction of the basket.

Table 6: Timing results in seconds for the basket consisting of 5 assets. The results were performed on a standard home PC where no parallelization was used. Timings are aggregated for all the assets and are reported in seconds. “IV param” indicates the time needed to parameterize the implied volatilities- this process may be considered irrelevant for this framework as it is a part of data processing. “Gen.$x_{i}$” indicates the time needed to generate the collocation points $z_{i}$ and $x_{i}$ based on the process $X_{i}$. “Gen.$y_{i}$” represents the time needed to calibrate the collocation method via the inverse of the corresponding CDFs. “Monte Carlo” corresponds to the complete basket simulation with five assets and 50k paths for each stock.  IV param. (s) |  Gen.$x_{i}$ (s) |  Gen.$y_{i}$ (S) | Monte Carlo (s)
---|---|---|---
0.5765 | 0.040 | 0.0876 | 0.240

In the next section the impact of the model parameters on basket implied
volatility will be studied.

###  4.2 Impact of the Model Parameters on Basket Implied Volatility

This section analyses the impact of the kernel model parameters on the basket
implied volatilities induced by the ELV model. As described earlier, due to
the collocation method, a particular choice of model parameters does not have
a material impact on the fit of individual assets- they, due to the mapping
procedure, remain intact. However, it has an impact on the covariance
structure of the basket. This effect can be measured in terms of implied
volatilities. In the experiment we set the parameters of process $X(t)$ (all
the parameter values are the same for all the assets) $\sigma=0.10$,
$\sigma_{J}=0.02$, $\mu_{J}=0$, $\xi_{p}=0.25$, $T=0.5$ with the correlation
between all the Brownian motions, ${\rm d}W_{i}(t){\rm d}W_{j}(t)=\rho{\rm
d}t$, $\rho=0.5$ and measure the impact on the implied volatility of the
basket, $B(T)$\- the basket consists of 5 assets and it’s construction is
defined in Section 4.1.

Figure 5 and Figure 6 show that101010Results for $\sigma$ are not presented
here as they resembled the impacts of $\xi_{p}$ and $\sigma_{J}$: $\sigma$,
$\sigma_{J}$ and $\xi_{p}$ have a moderate impact on the level of basket
implied volatilities. The two remaining parameters are the most relevant:
$\mu_{J}$ and $\rho$. To a large extent, correlation, $\rho$ determines the
level of basket implied volatilities (higher correlation lower the implied
volatility) while $\mu_{J}$ in a significant way controls the implied
volatility skew, as demonstrated in Figure 7. We considered negative values
for $\mu_{J}$ in the experiment, implying negative shocks to the underlying
process. We observe that a more negative value of $\mu_{J}$ produces more
skew.

Figure 5: Impact of different jump parameters on the shape of the basket
implied volatility in the ELV model. First figure: impact of $\sigma_{J}$;
Second figure: impact of $\mu_{J}$.

Figure 6: Impact of different jump parameters on the shape of the basket
implied volatility in the ELV model. First figure: impact of $\xi_{p}$; Second
figure: impact of $\rho$\- correlation between Brownian motions.

Since both model parameters $\sigma_{J}$ and $\xi_{p}$ impact the basket
implied volatility in a similar fashion, this means that in the calibration
routine, one of them can be effectively fixed, reducing the complexity of the
calibration procedure. First, the calibration process can be performed in an
iterative way, where the $\mu_{J}$ (skew) is chosen; secondly, the ATM implied
volatility is found such that the basked and index ATM volatility matched
perfectly.

Figure 7: Basket implied volatilities with re-calibrated model parameters for
ATM level.

The presented results suggest that the model can only generate implied
volatility skew. The basked implied volatility shapes may also be affected by
the basket’s composition. In Figure 4 we have shown that the kernel process
allows for the generation of a smile, though limited. On the other hand, each
asset’s volatility contributes to the overall basket vol. This will be
confirmed in the numerical experiment where we will consider a DJIA index with
30 underlying assets.

###  4.3 Calibration of a High Dimensional Basket

In this section, we perform a practical experiment where we build a bridge
between a basket and an index. This section considers a DJIA (Dow Jones
Industrial Average) with 30 underlying assets, as of 12/08/2021, with the spot
price of $357.55$. This experiment poses an extension of a case considered in
Section 4.1 where a 5D case was examined. Now, however, we will also perform
the model calibration, allowing us for consistent pricing of both an index and
a combination of basket constitutes.

The essence of the ELV model lies in the construction of a local volatility
model that is constructed of basket marginal distributions that are consistent
with the index. As such, the calibration procedure can be performed separately
for each expiry date $T_{i}$, where we utilize analytic expression for the
basket’s moments (see Equation (3.10) that are matched with the moments
computed from the index. The calibration is straightforward, and as indicated
earlier, there are two dominant model parameters, $\mu_{J}$ and $\rho$, that
play a central role in the skew calibration. Given the closed-form formulas
for the moments, the calibration is swift and straightforward and can be
performed with a simple search algorithm.

In Table 7 the estimated model parameters are presented, while the moments of
the basket compared to the index and “exact” moments obtained from Monte Carlo
are tabulated in Table 8.

Table 7: Calibrated model parameters. $T_{i}$ | $\sigma_{M}$ | $\sigma_{J}$ | $\mu_{J}$ | $\xi_{P}$ | $T$ | $\rho$
---|---|---|---|---|---|---
$1m$ | 0.25 | 0.020 | -0.20 | 0.25 | 0.30 | 0.35
$2m$ | 0.15 | 0.020 | -0.20 | 0.25 | 0.50 | 0.25
$3m$ | 0.15 | 0.020 | -0.20 | 0.25 | 0.50 | 0.25
$6m$ | 0.15 | 0.020 | -0.20 | 0.25 | 0.50 | 0.37
$1y$ | 0.1 | 0.020 | -0.20 | 0.25 | 0.50 | 0.45
$1.5y$ | 0.1 | 0.020 | -0.20 | 0.25 | 0.50 | 0.45
$2y$ | 0.20 | 0.020 | -0.25 | 0.50 | 0.45 | 0.38
Table 8: Moment approximation quality for DJIA with 30 underlying assets, as defined in Table 5.  | Index (Market) | Basket (Monte Carlo) | Basket (Proxy)
---|---|---|---
$T_{i}$ | $\sigma_{B}$ | $\mathbb{E}[B^{3}]/\sigma_{B}^{3}$ | $\sigma_{B}$ | $\mathbb{E}[B^{3}]/\sigma_{B}^{3}$ | $\sigma_{B}$ | $\mathbb{E}[B^{3}]/\sigma_{B}^{3}$
$1m$ | 19.08 | 6,634.4 | 19.0 | 6,714.1 | 18.7 | 6,982.6
$2m$ | 29.07 | 1,895.9 | 29.3 | 1,858.0 | 28.6 | 1,958.6
$3m$ | 36.47 | 969.8 | 36.0 | 1,009.5 | 35.2 | 1,059.5
$6m$ | 53.34 | 319.9 | 54.3 | 304.0 | 53.7 | 302.8
$1y$ | 78.56 | 107.3 | 81.9 | 95.9 | 82.8 | 85.8
$1.5y$ | 94.52 | 65.3 | 99.0 | 58.0 | 99.7 | 50.6
$2y$ | 109.07 | 45.4 | 113.5 | 41.3 | 113.1 | 35.5

Figures 9 and 10 show, for a range of maturities, the implied volatilities
obtained from the calibrated basket and the index. In addition, Figure 8 zooms
in on the fit to marginal distributions of a selected four stocks included in
the index.

Figure 8: Calibration quality for four selected stocks. Red squared indicate
market quote, blue line corresponds to the collocation method.

The results confirm excellent performance of the ELV model in fitting both
individual stocks’ distributions and their covariance structure. For all the
available index quotes, the model is within the bid-ask spreads. We have also
included the ELV model, where jumps are not included for comparison reasons.
Such a case may be considered standard local volatility. The additional degree
of freedom given by the jumps is essential in calibrating the basket to the
index. This is also confirmed by the calibration results shown in the Table 7
where negative $\mu_{J}$ was obtained for all the option expiries.

Figure 9: Implied volatilities for DJIA 30 compared to ELV model with and
without jumps; LHS: $T=2M$, RHS: $T=3M.$

Figure 10: Implied volatilities for DJIA 30 compared to ELV model with and
without Jumps; LHS: $T=1y$, RHS: $T=2y.$

###  4.4 Cost of Basket Miscalibration

In this section, we assess the impact of basket miscalibration. In the
experiment, we will continue with the 30-dimensional basket calibrated in
Section 4.3. We will consider a subset of assets $S_{i_{1}},\dots S_{i_{K}}$
of the basket and price an exotic derivative using two variants of the ELV
framework. In the first case, we take the model with a covariance structure
only driven by Brownian motion, while the alternative is the model where jumps
are included. Our previous experiment has shown that including jumps
facilitates a better fit to the index. Now, we quantify the potential cost of
this mispricing.

In the first step of pricing, the local volatility surface for the ELV
dynamics, $\bar{B}$ in (2.2), needs to be constructed. This is performed based
on the calibrated basket in Equation (2.4). The number of expiries, $N_{T}$,
available for building the local volatility surface depends on the required
accuracy and discretization time-step of (2.2). Although, in practice, a
finite and somewhat limited set of expiries, $T_{j}$, is available in the
market, for the intermediate expiry, $T$, $T_{j}<T<T_{j+1}$ we can impose a
smooth transition between the corresponding model parameters, i.e.,

${\bf p}(T)=\frac{T_{j+1}-T}{T_{j+1}-T_{j}}{\bf
p}(T_{j})+\frac{T-T_{j}}{T_{j+1}-T_{j}}{\bf p}(T_{j+1}),$

where ${\bf p}(T_{j})$, ${\bf p}(T_{j+1})$ represent the parameters of the
kernel process $X(T_{j})$, $X(T_{j+1})$ observed at times $T_{j}$ and
$T_{j+1}$, respectively.

In this experiment we consider a basket consisting of $5$ stocks, already
introduced earlier in Section 4.1 where the subset of DJIA stocks has been
considered, and two exotic derivatives, namely, an arithmetic Asian option,
defined as:

$V_{A}(t_{0})=N{\rm
e}^{-rT}\mathbb{E}\left[\left(\sum_{l=1}^{L}S(T_{l})-K_{A}\right)^{+}|\mathcal{F}(t_{0})\right],$

with $N$ being the notional and where $K_{A}$ is the strike price, and the up-
and-out barrier option:

$V_{Bar}(t_{0})=N{\rm
e}^{-rT}\mathbb{E}\left[(B(T)-K_{B})^{+}1_{\tau(H,B)>T}\Big{|}\mathcal{F}(t_{0})\right],$

where $K_{B}$ is the strike price, $H>B(t_{0})$ is the barrier, and
$\tau(x,B):=\inf\\{t\geq 0:B(t)>x\\}$, and the notional, $N$. In the
experiment we consider the ATM option with $K_{B}=123.82$ and the notional of
$1$. Although relatively standard in the financial world, both derivatives are
sensitive to different characteristics of the underlying model.

The pricing for varying strikes and barrier levels is presented in Figure 11.
The results show that although in the case of the Asian option, the skew-
related miscalibration does not have a significant impact in the case of the
up-and-out barrier option, the skew plays a vital role. As illustrated in
Table 9 the price differences may be even ten-fold.

Table 9: Comparison of Asian and Barrier option pricing for ELV models with and without jumps. Asian option price differences |  |  |  |  |
---|---|---|---|---|---
$K$ | 100.05 | 108.46 | 116.87 | 125.28 | 133.69 | 142.10 | 150.51
diff | 0.0392 | -0.0555 | -0.1439 | -0.0464 | 0.0870 | 0.0602 | 0.0289
Barrier option price differences |  |  |  |  |
H | 123.8225 | 178.1855 | 232.5484 | 286.9113 | 341.2742 | 395.6371 | 450
diff | 0 | -1.881 | -1.2491 | -0.4639 | -0.1986 | -0.0891 | -0.0735

Figure 11: Figure illustrates prices of Asian and Barrier option prices for:
LHS varying strike, $K$, RHS: varying the barrier level, $H$.

##  5 Conclusions

We have proposed above the Effective Local Volatility model capable of
accurate, consistent and swift pricing of basket options – both plain vanilla
and of the more exotic type. The model is accurate, because – as demonstrated
in Section 4.3 – it calibrates to liquid index options with only minimal
error, well within the bid-ask spread. The model is consistent, because it
simultaneously calibrates to individual basket constituents (via the functions
$g_{i,j}(\cdot)$) and their index-implied covariance structure (via the kernel
processes $X_{i}(\cdot)$ which we assume to be of Merton’s jump diffusion
type). Finally, the model works swiftly, because the typically numerically
“expensive” task, like sampling from stocks’ marginal distributions, is
handled via the collocation method, and basket calibration is performed by
means of moment matching technique, leveraging the analytical expressions for
basket moments which we derive in Section 3. Importantly, unlike copula-based
or local-correlation models which tend to be black-boxy and intractable, the
proposed approach allows users to control the basket implied volatility skew
using just a few parameters in a clear and transparent way. Thus, it is
straightforward to control the dependence structure of the basket and trace
how implied volatility smiles transfer from individual assets to the whole
basket (or possibly the index). While interesting in its own right, this
feature makes it possible to price structured products including baskets on
any combination of index constituents stocks: the assets simply “inherit” the
covariance structure implied by the index which we then feed into a one-
dimensional local volatility model, allowing us to price a desired structured
product in a quick and easy way.

While the model as currently proposed strikes a good balance between accuracy
and tractability and works very well for basket dimensions relevant in
practice, the approach can nonetheless be extended in two main directions.
First, one could postulate more sophisticated dynamics for the kernel
processes than the simple jump diffusion model. This might allow for even
better control of the implied volatility shapes, however possibly at some cost
to tractability, slowing down the calibration. Second, one could add a
stochastic component to the Effective Local Volatility model which might allow
for more accurate pricing of callable structures.

## References

  * [1] A. Ahdida and A. Alfonsi.  A mean-reverting sde on correlation matrices.  Stochastic Processes and their Applications, 123(4):1472–2630, 2013\.
  * [2] L. B. G. Andersen.  Simple and efficient simulation of the Heston stochastic volatility model.  Journal of Computational Finance, 11:3:1–22, 2008.
  * [3] K. N. Bakshi, G. and D. Madan.  Stock return characteristics, skew laws, and the differential pricing of individual equity options.  Review of Financial Studies, (16):101–143, 2003.
  * [4] N. P. Bollen and R. E. Whaley.  Does net buying pressure affect the shape of implied volatility functions?  Journal of Finance, 59(2):711–753, 2004.
  * [5] N. Branger and C. Schlag.  Why is the index smile so steep?  Review of Finance, 8(1):109–127, 2004.
  * [6] J. De Fonseca, G. M., and C. Tebaldi.  Option pricing when correlations are stochastic: an analytical framework.  Review of Derivative Research, 10:151–180, 2007.
  * [7] B. Dupire.  Pricing with a smile.  Risk, 7:18–20, 1994.
  * [8] B. Dupire.  A unified theory of volatility.  Paribas working paper, 1996.
  * [9] L. A. Grzelak.  The collocating local volatility framework – a fresh look at efficient pricing with smile.  International Journal of Computer Mathematics, 96(11):2209–2228, 2019.
  * [10] L. A. Grzelak and C. W. Oosterlee.  From arbitrage to arbitrage-free implied volatilities.  Journal of Computational Finance, 20(3):31–49, 2016.
  * [11] L. A. Grzelak, J. Witteveen, M. Suarez-Taboada, and C. W. Oosterlee.  The stochastic collocation Monte Carlo sampler: highly efficient sampling from “expensive” distributions.  Quantitative Finance, 19(2):339–356, 2019.
  * [12] I. Gyöngy.  Mimicking the one-dimensional marginal distribution of process having an Itô differential.  Probability Theory and Related Fields, 71:501–516, 1986.
  * [13] P. Hagan, D. Kumar, A. Leśniewski, and D. Woodward.  Managing smile risk.  Wilmott Magazine, pages 84–108, 2002.
  * [14] A. Langnau.  Introduction into “local correlation modelling”.  arXiv preprint arXiv:0909.3441, 2009.
  * [15] A. Langnau.  A dynamic model for correlation.  Risk, 23(4):74, 2010.
  * [16] A. Langnau and D. Cangemi.  Marking systemic portfolio risk with the Merton model.  Risk, 24(7):68, 2011.
  * [17] V. Lucic.  Correlation skew via product copula.  Global Derivatives Conference, 2013.
  * [18] R. Merton.  Option pricing when the underlying stocks are discontinuous.  J. Financ. Econ., 5:125–144, 1976.
  * [19] R. C. Merton.  Option pricing when underlying stock returns are discontinuous.  Journal of financial economics, 3(1-2):125–144, 1976.
  * [20] C. W. Oosterlee and L. A. Grzelak.  Mathematical modeling and computation in finance.  World Scientific, 2019.
  * [21] A. Reghai.  Breaking correlation breaks.  Risk, pages 90–95, 2010.
  * [22] H. C. Robert and T. M. William.  The orthogonal development of nonlinear functionals in series of Fourier-Hermite functionals.  Annals of Mathematics, 48(2):385–392, 1947.
  * [23] A. W. vd Stoep, L. A. Grzelak, and C. W. Oosterlee.  Collocating volatility: a competitive alternative to stochastic local volatility models.  IJTAF, 23(06), 2020.
  * [24] V. Zetocha.  Correlation skew via stochastic correlation and jumps.  Risk, 2015.

##  Appendix A Stochastic Collocation Method

Let’s assume two random variables, $Y$ and $X$, where the latter one is
cheaper to sample from (e.g., $X$ is a Gaussian random variable). These two
scalar random variables are connected, via,

$F_{Y}(Y)\stackrel{{\scriptstyle d}}{{=}}U\stackrel{{\scriptstyle
d}}{{=}}F_{X}(X),$ (A.1)

where $U\sim\mathcal{U}([0,1])$ is a uniformly distributed random variable,
$F_{Y}(\bar{y}):=\mathbb{P}[Y\leq\bar{y}]$ and
$F_{X}(\bar{x}):=\mathbb{P}[X\leq\bar{x}]$ are cumulative distribution
functions (CDF). Note that $F_{X}(X)$ and $F_{Y}(Y)$ are random variables
following the same uniform distribution. $F_{Y}(\bar{y}_{n})$ and
$F_{X}(\bar{x}_{n})$ are supposed to be strictly increasing functions, so that
the following inversion holds true,

$\bar{y}_{n}=F^{-1}_{Y}(F_{X}(\bar{x}_{n}))=:g(\bar{x}_{n}).$ (A.2)

where $\bar{y}_{n}$ and $\bar{x}_{n}$ are samples (numbers) from $Y$ and $X$,
respectively. The mapping function, $g(\cdot)=F^{-1}_{Y}(F_{X}(\cdot))$,
connects the two random variables and guarantees that $F_{X}(\bar{x}_{n})$
equals $F_{Y}(g(\bar{x}_{n}))$, in distributional sense and also element-wise.
The mapping function should be approximated, i.e., $g(\bar{x}_{n})\approx
g_{m}(\bar{x}_{n})$, by a function which is cheap. When function
$g_{m}(\cdot)$ is available, we may generate “expensive” samples,
$\bar{y}_{n}$ from $Y$, by using the cheaper random samples $\bar{x}_{n}$ from
$X$.

The Stochastic Collocation Monte Carlo method (SCMC) developed in [11] aims to
find an accurate mapping function $g(\cdot)$ in an efficient way. The basic
idea is to employ Equation (A.2) at specific collocation points and
approximate the function $g(\cdot)$ by a suitable monotonic interpolation
between these points. This procedure, see Algorithm I, reduces the number of
expensive inversions $F^{-1}_{Y}(\cdot)$ to obtain many samples from
$Y(\cdot)$.

Algorithm: SCMC Method Taking an interpolation function of degree $m-1$ (with
$m\geq 2$, as we need at least two collocation points), as an example, the
following steps need to be performed: 1. Calculate CDF $F_{X}(x_{j})$ on the
points $(x_{1},x_{2},...,x_{m})$, that are obtained, for example, from Gauss-
Hermite quadrature, giving $m$ pairs $({x}_{j},F_{X}({x}_{j}))$; 2. Invert the
target CDF ${y}_{j}=F^{-1}_{Y}(F_{X}({x}_{j}))$, $j=1,\ldots,m$, and form $m$
pairs $({x}_{j},{y}_{j})$; 3. Define the interpolation function, $y=g_{m}(x)$,
based on these $m$ point pairs $({x}_{j},{y}_{j})$; 4. Obtain sample $\hat{Y}$
by applying the mapping function $\hat{Y}=g_{m}(\hat{X})$, where sample
$\hat{X}$ is drawn from $X$.

The SCMC method parameterizes the distribution function by imposing
probability constraints at the given collocation points. Taking the Lagrange
interpolation as an example, we can expand function $g_{m}(\cdot)$ in the form
of polynomial chaos,

$Y\approx
g_{m}(X)=\sum_{j=0}^{m-1}\hat{\alpha}_{j}X^{j}=\hat{\alpha}_{0}+\hat{\alpha}_{1}X+...+\hat{\alpha}_{m-1}X^{m-1}.$
(A.3)

Monotonicity of interpolation is an important requirement, particularly when
dealing with peaked probability distributions.

The Cameron-Martin Theorem [22] states that any distribution can be
approximated by a polynomial chaos approximation based on the normal
distribution, but also other random variables may be used for $X$ (see, for
example, [11]).

##  Appendix B Proof of Lemma 3.1

###### Proof.

By substitution we find:

$\displaystyle\mathbb{E}[{\rm e}^{aX_{1}(T)+bX_{2}(T)}]={\rm
e}^{a\mu_{1}+b\mu_{2}}\mathbb{E}\left[\exp\left(a\sigma_{1}W_{1}(T)+b\sigma_{2}W_{2}(T)+\sum_{k=1}^{X_{\mathcal{P}}(T)}aJ_{k,{1}}+\sum_{k=1}^{X_{\mathcal{P}}(T_{j})}bJ_{k,2}\right)\right],$

with $X_{1}$, $X_{2}$ defined in (3.6), and where
$J_{k,i}\sim\mathcal{N}(\mu_{J},\sigma^{2}_{J})$. Due to the correlation
between Brownian motions we perform the following factorization,

$a\sigma_{1}W_{1}(T)+b\sigma_{2}W_{2}(T)\stackrel{{\scriptstyle\text{d}}}{{=}}\sqrt{a^{2}\sigma_{1}^{2}+b^{2}\sigma_{2}^{2}+2\rho_{1,2}ab\sigma_{1}\sigma_{2}}\widetilde{W}(T)=:\hat{\sigma}\widetilde{W}(T),$

with
$\hat{\sigma}^{2}=a^{2}\sigma_{1}^{2}+b^{2}\sigma_{2}^{2}+2\rho_{1,2}ab\sigma_{1}\sigma_{2}$
and where $\widetilde{W}(T)$ is an independent Brownian motion.

$\displaystyle\sum_{k=1}^{X_{\mathcal{P}}(T)}aJ_{k,{1}}+\sum_{k=1}^{X_{\mathcal{P}}(T)}bJ_{k,2}=\sum_{k=1}^{X_{\mathcal{P}}(T)}\left(aJ_{k,{1}}+bJ_{k,{2}}\right)=:\sum_{k=1}^{X_{\mathcal{P}}(T)}\hat{J}_{k},$

where
$\hat{J}\equiv\hat{J}_{k}=aJ_{k,{1}}+bJ_{k,{2}}\sim\mathcal{N}((a+b)\mu_{J},(a^{2}+b^{2})\sigma_{J}^{2})$;
thus we have:

$\displaystyle\mathbb{E}[{\rm
e}^{aX_{1}(T)+bX_{2}(T)}]=\mathbb{E}\left[\exp\left(\hat{\sigma}\widetilde{W}(T)+\sum_{k=1}^{X_{\mathcal{P}}(T)}\hat{J}_{k}\right)\right]\stackrel{{\scriptstyle\perp\\!\\!\\!\perp}}{{=}}{\rm
e}^{a\mu_{1}+b\mu_{2}}\mathbb{E}\left[{\rm
e}^{\hat{\sigma}\widetilde{W}(T)}\right]\mathbb{E}\left[{\rm
e}^{\sum_{k=1}^{X_{\mathcal{P}}(T)}\hat{J}_{k}}\right].$

The first expectation simply reads $\mathbb{E}\big{[}{\rm
e}^{\hat{\sigma}\widetilde{W}(T)}\big{]}={\rm
e}^{\frac{1}{2}\hat{\sigma}^{2}T},$ while the second expectation, thanks to
the tower property of expectations, gives:

$\displaystyle\mathbb{E}\left[\exp\left(\sum_{k=1}^{X_{\mathcal{P}}(T)}\hat{J}_{k}\right)\right]$
$\displaystyle=$ $\displaystyle\sum_{n\geq
0}\mathbb{E}\left[\exp\left(\sum_{k=1}^{n}\hat{J}_{k}\right)\right]\frac{{\rm
e}^{-\xi_{p}T(\xi_{p}T)^{n}}}{n!}$ $\displaystyle=$ $\displaystyle\sum_{n\geq
0}\frac{{\rm e}^{-\xi_{p}T(\xi_{p}T)^{n}}}{n!}\mathbb{E}^{n}\left[{\rm
e}^{\hat{J}}\right]$ $\displaystyle=$
$\displaystyle\exp\left(\xi_{p}T\mathbb{E}[{\rm e}^{\hat{J}}-1]\right),$

where the last expression can be recognized as a Taylor expansion of an
exponential function. Since

$\mathbb{E}[{\rm e}^{\hat{J}}]={\rm
e}^{(a+b)\mu_{J}+\frac{1}{2}(a^{2}+b^{2})\sigma_{J}^{2}},$

by substitutions the proof is complete. ∎

##  Appendix C Third Moment Derivations

For convenience we neglect the time argument $T_{j}$. The third moment of the
basket, by definition, is defined as:

$\displaystyle\mathbb{E}[B^{3}]=\sum_{i_{1}=1}^{N}\sum_{i_{2}=1}^{N}\sum_{i_{3}=1}^{N}\omega_{i_{1}}\omega_{i_{2}}\omega_{i_{3}}\mathbb{E}\left[S_{i_{1}}S_{i_{2}}S_{i_{3}}\right]=:\sum_{i_{1}=1}^{N}\sum_{i_{2}=1}^{N}\sum_{i_{3}=1}^{N}\omega_{i_{1}}\omega_{i_{2}}\omega_{i_{3}}\mathbb{E}\left[S_{i_{1},i_{2}}S_{i_{3}}\right],$

where $S_{i_{1},i_{2}}=S_{i_{1}}S_{i_{2}}$. Following the same procedure as
for the second moment we find:

$\displaystyle\mathbb{E}[B^{3}]$ $\displaystyle=$
$\displaystyle\sum_{i_{1}=1}^{N}\sum_{i_{2}=1}^{N}\sum_{i_{3}=1}^{N}\omega_{i_{1}}\omega_{i_{2}}\omega_{i_{3}}\left[\hat{\rho}({S_{i_{1},i_{2}},S_{i_{3}}})\sigma_{S_{i_{1},i_{2}}}\sigma_{S_{i_{2}}}+\mathbb{E}[S_{i_{1},i_{2}}]\mathbb{E}[S_{i_{3}}]\right].$

$\mathbb{E}[S_{i_{1},i_{2}}]:=\mathbb{E}[S_{i_{1}}S_{i_{2}}]$, can be computed
with (3.2), $\mathbb{E}[S_{i_{3}}]$ is given from the market data. Terms we
need to establish are: $\sigma_{S_{i_{1},i_{2}}}$ and
$\hat{\rho}({S_{i_{1},i_{2}},S_{i_{3}}})$. For the first term we have:

$\displaystyle\sigma_{S_{i_{1},i_{2}}}^{2}$ $\displaystyle=$
$\displaystyle\mathbb{E}[S_{i_{1}}S_{i_{2}}]-\mathbb{E}[S_{i_{1}}]\mathbb{E}[S_{i_{2}}],$

which again can be computed using Equation (3.2). For the correlation
coefficient we have:

$\displaystyle\hat{\rho}({S_{i_{1},i_{2}},S_{i_{3}}})$ $\displaystyle=$
$\displaystyle\rho\left({\rm
e}^{\hat{\alpha}_{i_{1},j,0}+\hat{\alpha}_{i_{1},j,1}X_{i_{1}}}{\rm
e}^{\hat{\alpha}_{i_{2},j,0}+\hat{\alpha}_{i_{1},j,1}X_{i_{2}}},{\rm
e}^{\hat{\alpha}_{i_{3},j,0}+\hat{\alpha}_{i_{3},j,1}X_{i_{3}}}\right)$
$\displaystyle=$ $\displaystyle\rho\left({\rm
e}^{\hat{\alpha}_{i_{1},j,1}X_{i_{1}}+\hat{\alpha}_{i_{2},j,1}X_{i_{2}}},{\rm
e}^{\hat{\alpha}_{i_{3},j,1}X_{i_{3}}}\right)=:\bar{\rho}({S_{i_{1},i_{2}},S_{i_{3}}}),$

which yields:

$\displaystyle\bar{\rho}({S_{i_{1},i_{2}},S_{i_{3}}})$
$\displaystyle\stackrel{{\scriptstyle\text{def}}}{{=}}$
$\displaystyle\frac{\mathbb{E}[{\rm
e}^{\hat{\alpha}_{i_{1},j,1}X_{i_{1}}+\hat{\alpha}_{i_{2},j,1}X_{i_{2}}+\hat{\alpha}_{i_{3},j,1}X_{i_{3}}}]-\mathbb{E}[{\rm
e}^{\hat{\alpha}_{i_{1},j,1}X_{i_{1}}+\hat{\alpha}_{i_{2},j,1}X_{i_{2}}}]\mathbb{E}[{\rm
e}^{\hat{\alpha}_{i_{3},j,1}X_{i_{3}}}]}{\sqrt{\mathbb{V}\text{ar}[{\rm
e}^{\hat{\alpha}_{i_{1},j,1}X_{i_{1}}+\hat{\alpha}_{i_{2},j,1}X_{i_{2}}}]\mathbb{V}\text{ar}[{\rm
e}^{\hat{\alpha}_{i_{3},j,1}X_{i_{3}}}]}}$ $\displaystyle=$
$\displaystyle\frac{\mathbb{E}[{\rm
e}^{\hat{\alpha}_{i_{1},j,1}X_{i_{1}}+\hat{\alpha}_{i_{2},j,1}X_{i_{2}}+\hat{\alpha}_{i_{3},j,1}X_{i_{3}}}]-\mathbb{E}[{\rm
e}^{\hat{\alpha}_{i_{1},j,1}X_{i_{1}}+\hat{\alpha}_{i_{2},j,1}X_{i_{2}}}]\phi_{X_{i_{3}}}(-i\hat{\alpha}_{i_{3},j,1})}{\sqrt{\mathbb{V}\text{ar}[{\rm
e}^{\hat{\alpha}_{i_{1},j,1}X_{i_{1}}+\hat{\alpha}_{i_{2},j,1}X_{i_{2}}}]\left(\phi_{X_{i_{3}}}(-2i\hat{\alpha}_{i_{3},j,1})-\phi_{X_{i_{3}}}^{2}(-i\hat{\alpha}_{i_{3},j,1})\right)}}.$

All the expression except for the first expectation can be computed using
(3.2):

$\displaystyle
a\sigma_{1}W_{1}(T)+b\sigma_{2}W_{2}(T)+c\sigma_{3}W_{3}(T)\stackrel{{\scriptstyle\text{d}}}{{=}}\hat{\sigma}\widetilde{W}(T),$
(C.1)

with

$\displaystyle\hat{\sigma}^{2}={a^{2}\sigma_{1}^{2}+b^{2}\sigma_{2}^{2}+c^{2}\sigma_{3}^{2}+2\rho_{1,2}ab\sigma_{1}\sigma_{2}}+2\rho_{1,3}ac\sigma_{1}\sigma_{3}+2\rho_{2,3}ab\sigma_{2}\sigma_{3},$
(C.2)

therefore

$\displaystyle\sum_{k=1}^{X_{\mathcal{P}}(T)}aJ_{k,{1}}+\sum_{k=1}^{X_{\mathcal{P}}(T)}bJ_{k,2}+\sum_{k=1}^{X_{\mathcal{P}}(T)}bJ_{k,3}=\sum_{k=1}^{X_{\mathcal{P}}(T)}\left(aJ_{k,{1}}+bJ_{k,{2}}+bJ_{k,{3}}\right)=:\sum_{k=1}^{X_{\mathcal{P}}(T)}\hat{J}_{k},$
(C.3)

where
$\hat{J}\equiv\hat{J}_{k}=aJ_{k,{1}}+bJ_{k,{2}}+cJ_{k,{3}}\sim\mathcal{N}((a+b+c)\mu_{J},(a^{2}+b^{2}+c^{2})\sigma_{J}^{2})$;
thus we have:

$\displaystyle\mathbb{E}[{\rm
e}^{aX_{1}+bX_{2}+cX_{3}}]\stackrel{{\scriptstyle\perp\\!\\!\\!\perp}}{{=}}{\rm
e}^{a\mu_{1}+b\mu_{2}+c\mu_{3}}\mathbb{E}\left[{\rm
e}^{\hat{\sigma}\widetilde{W}(T)}\right]\mathbb{E}\left[{\rm
e}^{\sum_{k=1}^{X_{\mathcal{P}}(T)}\hat{J}_{k}}\right],$ (C.4)

with $\hat{\sigma}\widetilde{W}(T)$ defined in (C.1) and $\hat{J}_{k}$ defined
in (C.3). Using the result from proof in B we find:

$\displaystyle\omega_{X}(a,b,c)$ $\displaystyle:=$
$\displaystyle\mathbb{E}[{\rm e}^{aX_{1}+bX_{2}+cX_{3}}]$ $\displaystyle=$
$\displaystyle\exp\left(a\mu_{1}+b\mu_{2}+c\mu_{3}+\frac{1}{2}\sigma
T+\xi_{p}T\mathbb{E}[{\rm e}^{\hat{J}}-1]\right)$ $\displaystyle=$
$\displaystyle\exp\left[a\mu_{1}+b\mu_{2}+c\mu_{3}+\frac{1}{2}\hat{\sigma}T+\xi_{p}T\left({\rm
e}^{(a+b+c)\mu_{J}+\frac{1}{2}(a^{2}+b^{2}+c^{2})\sigma_{J}^{2}}-1\right)\right],$

where $\sigma^{2}$ is defined in (C.2).

Combining these results with Lemma 3.1 we find the following expression for
the correlation: $\bar{\rho}({S_{i_{1},i_{2}},S_{i_{3}}})$,

$\displaystyle\bar{\rho}({S_{i_{1},i_{2}},S_{i_{3}}})=\frac{\omega_{X}(c_{1},c_{2},c_{3})-\omega_{X}(c_{1},c_{2})\phi_{X_{i_{3}}}(-ic_{3})}{\sqrt{\left(\omega_{X}(2c_{1},2c_{2})-\omega_{X}^{2}(c_{1},c_{2})\right)\left(\phi_{X_{i_{3}}}(-2ic_{3})-\phi_{X_{i_{3}}}^{2}(-ic_{3})\right)}},$

where: $c_{1}:=\hat{\alpha}_{i_{1},j,1}$, $c_{2}:=\hat{\alpha}_{i_{2},j,1}$,
$c_{3}:=\hat{\alpha}_{i_{3},j,1}$ and where $\omega_{X}(c_{1},c_{2},c_{3})$ is
defined in (C) and $\omega_{X}(c_{1},c_{2})$ in (3.7).