Semi-analytic conditional expectations

A data-driven approach to computing expectations for the pricing and hedging of exotics

CLICK HERE TO VIEW THE PDF

Gaussian mixture model dynamically controlled kernel estimation (GMM-DCKE), a purely data-driven and model-agnostic method to compute conditional expectations, is introduced. Joerg Kienitz applies it to the pricing and hedging of (multi-dimensional) exotic Bermudan options and to calibration and pricing within stochastic local volatility models

Fast and accurate approximations of conditional expectations and their respective distributions are essential for many problems arising in quantitative finance. The least squares Monte Carlo (LSM) method (see Longstaff & Schwartz 2001) is the market standard for simulation-based estimation. It is well known that LSM requires a large number of paths and faces problems in estimating the tails. Tail estimates are crucial, eg, for risk management or pricing deep out-of-the-money options. Convergence depends on the choice of basis functions, but if LSM converges to the price, its derivative might not converge. Variance reduction, basis function selection and local regression have been considered for mitigation.

In this work, we built upon previous work by Geng et al (2022), but the proposed method is different from dynamically controlled kernel estimation (DCKE) since it does not rely on numerical kernel density estimation including (local) bandwidth selection and it does not need to apply Gaussian process regression (GPR) for interpolation or extrapolation and smoothing. The key difference is that it replaces the numerical methods with analytic expressions that lead simultaneously to smooth approximations.

GMM-DCKE allows us to apply further state variables in the sense of control variates (CVs). This improves upon the basic method and allows for analytical computation of proxy hedge sensitivities. Combining the GMM, CVs and properties of multivariate Gaussians, the contributions of this paper to the literature are the following:

• •

  the calculation of model-agnostic, data-driven conditional expectations on observed realisations at given time points even for multi-modal distributions with semi-analytic methods;

• •

  the GMM is fitted numerically, but all subsequent calculations are analytic, and thus they are fast; and

• •

  the method implies closed-form solutions for smooth proxy hedge sensitivities.

For illustration, the method is applied to pricing of (multi-dimensional) Bermudan derivatives and for pricing and calibration with stochastic local volatility models.¹¹ 1 A longer exposition and pedagogical illustrations can be found in Kienitz (2021) or at https://github.com/Lapsilago/GMM_DCKE.

To fix notation we take a discrete time grid $?=\{T_0,T_1,\mathrm{\dots },T_N\}$ with realisations of $d$ stochastic risk factors $S(T)={(S_1(T),\mathrm{\dots },S_d(T))}^{\mathrm{T}}$, $T\in ?$.

GMM

Let $K\in \mathbb{N}$ be a positive integer and let $Z_k$, $k=1,\mathrm{\dots },K$, be $d$-dimensional random variables with $Z_k\sim {?}_d({\mu }_k,{\mathrm{\Sigma }}_k)$ (denoting a $d$-dimensional Gaussian distribution, with ${\mu }_k\in {ℝ}^d$ and ${\mathrm{\Sigma }}_k\in {ℝ}^{d\times d}$ denoting the corresponding mean vectors and covariance matrices, respectively). Let $n_d(\cdot ;{\mu }_k,{\mathrm{\Sigma }}_k)$ be the corresponding probability density. We consider the $d$-variate $K$-component GMM $Z$ determined by the probability distribution $p_Z$:

$$p_Z(x)=\sum _{k=1}^K{\pi }_k{n}_d(x;{\mu }_k,{\mathrm{\Sigma }}_k),{\pi }_k>0,\sum _{k=1}^K{\pi }_k=1$$

${\pi }_k$ are called mixing weights and $\mathrm{GMM}(K)$ denotes a GMM with $K$ components. Once ${\pi }_k$, ${\mu }_k$ and ${\mathrm{\Sigma }}_k$, $k=1,\mathrm{\dots },K$, are fixed, the conditional distribution can be calculated analytically. Let:

$$Z={(Z_1,\mathrm{\dots },Z_m,Z_{m+1},\mathrm{\dots },Z_{m+n=d})}^{\mathrm{T}},Z\sim {?}_d({\mu }_Z,{\mathrm{\Sigma }}_Z)$$

We denote the realisations by $X$ and the labels by $Y$. We set $Z={(Y,X)}^{\mathrm{T}}$ with $m$-dimensional and $n$-dimensional Gaussians $Y={(Z_1,\mathrm{\dots },Z_m)}^{\mathrm{T}}$ and $X={(Z_{m+1},\mathrm{\dots },Z_{m+n})}^{\mathrm{T}}$, respectively; the conditionals, $Y\mid X$, and marginals for $Y$ and $X$ are available in closed form, with:

Taking $Y$, for instance, and denoting by $p_z(y,x)$ the joint density with regard to $Z$, we have:

	$p_Y(y)$	$={\displaystyle \int p_Z(y,x)dx}\sim {?}_m({\mu }_Y,{\mathrm{\Sigma }}_{YY})$
	$p_{Y\mid X=x}(y)$	$={\displaystyle \frac{p_Z(y,x)}{\int p_Z(y,x)dx}}\sim {?}_n({\mu }_{Y\mid X},{\mathrm{\Sigma }}_{Y\mid X})$

and the means and variances for the conditional distribution are given by:

	${\mu }_{Y\mid X}$	$={\mu }_Y-{\mathrm{\Sigma }}_{YX}{\mathrm{\Sigma }}_{XX}^{-1}(X-{\mu }_X)$		(1)
	${\mathrm{\Sigma }}_{Y\mid X}$	$={\mathrm{\Sigma }}_{YY}-{\mathrm{\Sigma }}_{YX}{\mathrm{\Sigma }}_{XX}^{-1}{\mathrm{\Sigma }}_{XY}$		(2)

respectively. Using these results for each single component $k=1,\mathrm{\dots },K$ of $\mathrm{GMM}(K)$ we have:

$$p_{k,Y\mid X=x}(y)=\frac{n_k(y,x)}{\int n_k(y,x)dx}\sim {?}_n({\mu }_{k,Y\mid X},{\mathrm{\Sigma }}_{k,Y\mid X})$$

and for the conditional distribution of $\mathrm{GMM}(K)$:

	$p_{Y\mid X}$	$={\displaystyle \sum _{k=1}^K}{\tilde{\pi }}_k{n}_{k,Y\mid X}$		(3)
	${\tilde{\pi }}_k$	$\sim {\displaystyle \frac{{\pi }_k{?}_m({\mu }_{k,X},{\mathrm{\Sigma }}_{k,XX})}{{\sum }_l{\pi }_l{?}_m({\mu }_{l,X},{\mathrm{\Sigma }}_{l,XX})}}$		(4)

For a given $\mathrm{GMM}(K)$, the conditional model is determined by (3) and (4). The key observation is that all conditional expectations, variances and covariances can be calculated analytically using (1) and (2).

The parameters

Let the set $?$ be of size $N$, and for $x\in ?$ we assume $x\in {ℝ}^d$. Then, $?$ is represented by an $N\times d$ matrix. Assume the number of components $K$ for $\mathrm{GMM}(K)$ is given. If $\theta$ denotes the collection of all parameters ${\pi }_k$, ${\mu }_k$ and ${\mathrm{\Sigma }}_k$, $k=1,\mathrm{\dots },K$, then $\theta$ is obtained by maximising the (log)likelihood:

$$l(\theta )=\sum _{n=1}^N\log \left(\sum _{k=1}^K{\pi }_k{n}_d(x_n,{\mu }_k,{\mathrm{\Sigma }}_k)\right)$$

by expectation maximisation (EM) (Dempster et al 1977).

GMM-DCKE algorithm

Input

Training set:

$?=(x_1,\mathrm{\dots },x_N)$, $x_n\in {ℝ}^d$, with labels $?=(y_1,\mathrm{\dots },y_N)$, $y_n\in ℝ$. $(?,?)$ are joint realisations of some random variables $X$, $Y$. The underlying risk factors are denoted by $X$, and $Y$ is a function of $X$, eg, a payoff at maturity $T$ for Bermudan derivatives or a variance value for a stochastic local volatility model.

Test set:

Let ${?}^{*}=(x_1^{*},\mathrm{\dots },x_M^{*})$, $x_n^{*}\in {ℝ}^d$, be a test set.

Output

The outputs are predictions $y^{*}=(y_1^{*},\mathrm{\dots },y_M^{*})$, $y_j^{*}\in ℝ$, $y_j^{*}\approx ?[Y\mid X=x_j^{*}]$. For our examples, these are the conditional expected values of a Bermudan derivative, or the expected variance values given the value of the risk factor in a stochastic volatility model. For in-sample performance we take ${?}^{*}=?$.

Parameters

The number of components is given by $K$ and, if we use a control variate, its realisations $?=(Z_1,\mathrm{\dots },Z_N)$. Including $?$ into the fitting of $\mathrm{GMM}(K)$, quantities such as the conditional mean, given by ${\mu }_Z(x):=?[Z\mid X=x]$, are calculated analytically.

Calculation

The first (numeric) step is to fit $\mathrm{GMM}(K)$ to approximate the joint distribution of $?$, $?$ and $?$ using EM. The next steps are analytic, using (1) and (2). To apply the CVs in order to stabilise the model we consider:

$Y^{*}:=Y\mid X+{\beta }_X(Z\mid X-{\mu }_{Z\mid X})$

(5)

with ${\beta }_X$ minimising the variance, (5):

${\beta }_{X=x}:=-{\displaystyle \frac{\mathrm{Cov}[Y,Z\mid X=x]}{\mathrm{Var}[Z\mid X=x]}}$

(6)

For $x\in ?$, with ${\mu }_{Z\mid X=x}$ a conditional mean, the ${\beta }_{X=x}\in {ℝ}^d$ are calculated analytically as a quotient of conditional covariance and conditional variance, respectively, for Gaussian distributions.

The quotient ${\beta }_X$ can be seen as a (conditional) trading strategy with respect to $Z$, given the values $X$ that minimise the variance for $Y^{*}$. If $Z$ is one of the risk factors, it is the (time-discrete) conditional minimal-variance delta, and thus a time-discrete, model-agnostic, data-driven hedge can be calculated analytically. Other choices, depending on the payoff, may be used to improve the estimate, eg, the maximum or minimum of the assets, when considering ‘best-of’ or ‘rainbow’ options.

Components and enhancements

The number $K$ of components can be determined by empirical results or by statistical methods, eg, a Bayesian Gaussian mixture, silhouette scores or minimising information criteria such as the Akaike and Bayesian information criteria:

$$\mathrm{AIC}:=2K-2\ln (\widehat{L});\mathrm{BIC}=K\ln (N)-2\ln (\widehat{L})$$

respectively, where $\widehat{L}:=p(x\mid \widehat{\theta },M)$, with $\widehat{\theta }$ the parameters maximising the likelihood, $x\in ?$ and $N$ the number data points. Taking the model where the AIC (or BIC) is smallest may lead to overfitting. Discrete gradients, ie, differences in AIC and BIC, for $\mathrm{GMM}(i)$ and $\mathrm{GMM}(i+1)$, respectively, can be considered. The number of components is chosen once this difference does not significantly increase further.

Silhouette scores exploit the compactness and separability of the clusters implied by the GMM. The Jensen-Shannon metric (JS) for $\mathrm{GMM}(i)$ and $\mathrm{GMM}(i+1)$ is calculated, and the number $i$ for which we observe a sudden jump in the JS value is chosen.

Known limitations are that the approximation is only valid for sample sizes $N$ that are much larger than the number of parameters in the model.

A method for enhancement is ‘bagging’, an ensemble method that consists of bootstrapping and aggregation (see Bishop 2006). In our setting we perform GMM-DCKE $k_{\mathrm{bag}}$ times with $N_{\mathrm{bag}}$ observations randomly sampled from $?$ with replacement and averaged over the number of experiments. Transformations can also be used. For instance, keeping values, especially of option prices or implied volatilities, positive is essential. (When applying Gaussian distributions, negative values are not ruled out.) If we have a set of strictly positive values $?$, we instead work with the set ${?}_{\epsilon ,\log }$ obtained by applying the transform $\log (a+\epsilon )$, $\epsilon >0$, for each $a\in ?$. Applying the GMM to this new training set and using the inverse transform $\exp (a-\epsilon )$ leads to strictly positive values.

Comparison with other approaches

Our work is related to that of Geng et al (2022), which introduced a kernel density estimator with a (local) bandwidth used to approximate the conditional expected values and Gaussian process regression used for smoothing and interpolation or extrapolation. Selecting a local bandwidth may be slow and depends strongly on the dimensionality. In contrast, our method is semi-analytic, stable and only depends on one parameter: the number of Gaussian regressions, which needs to be determined.

The work is also related to that of Halperin (2017), which presents a discrete-time option pricing model that is rooted in reinforcement learning. Our method essentially gives more direct access to this approach, circumventing the reinforcement learning part. Huge & Savine (2020) train a neural network, taking the hedge, ie, the sensitivities, into account to stabilise the calculations. Techniques such as kernel density estimation in multiple dimensions, reinforcement learning and neural networks are data and computationally intense when required to perform with high accuracy. GMM-DCKE derives the option price, conditional expectations and hedges accurately while applying sparse parameter sets and has significant advantages with regard to both the amount of data and the computation time, especially when considering high-dimensional problems.

Bermudan derivatives

For pricing Bermudan derivatives, for each consecutive pair $T_i$, $T_{i+1}\in ?$, we fit a GMM to the observed realisations and consider a discrete proxy hedging approach. Recall that if $S_t$ is the risky underlying, $r$ the risk-free rate and $B_t$ the bank account, then the self-financing condition of any hedging strategy $u_t$ implies that, for any time interval $[t,t+\mathrm{\Delta }]$:

$$u_t{S}_{t+\mathrm{\Delta }}+{\mathrm{e}}^{r\mathrm{\Delta }}{B}_t=u_{t+\mathrm{\Delta }}{S}_{t+\mathrm{\Delta }}+B_{t+\mathrm{\Delta }}$$

Denote by ${\mathrm{\Pi }}_t$ the replication portfolio of an option using a hedging strategy $u_t$ with respect to the available information, ie, the filtration ${ℱ}_t$. The hedging strategy $u_t^{*}$ that minimises the variance $\mathrm{Var}[{\mathrm{\Pi }}_T\mid {ℱ}_t]$ is given by:

$$\frac{\mathrm{Cov}[{\mathrm{\Pi }}_{t+1},\mathrm{\Delta }{S}_t\mid {ℱ}_t]}{\mathrm{Var}[\mathrm{\Delta }{S}_t\mid {ℱ}_t]}$$

The computation of the latter can be computationally intense. Notice that the minimal-variance hedging strategy $u_t^{*}$ is merely the minimal-variance delta. Taking this as a proxy, we include it in the calculation of conditional values to reduce the conditional variances. With GMM-DCKE the variance minimiser is available in closed form.

For pricing, we need to approximate the continuation value at each possible exercise date. This involves computing a conditional expectation. To illustrate our results we focus on a given time point $t$. For the full valuation this has to be applied in a backward algorithm. We choose the rough Bergomi and Bates models for illustration (Bates 1996; Bayer et al 2016).²² 2 We have actually applied the method to diffusion, jump-diffusion and pure jump processes. Let $W_1(t)$, $W_2(t)$ be Brownian motions with correlation coefficient $\rho \in [-1,1]$. Take $S_0=s_0$, $V_0=v_0$ as the initial values for the asset price and the variance, respectively. The stochastic differential equations (SDEs) for the evolution of the Bates model are:

(7)

The parameters in (7) are $r$, the risk-neutral rate; $\kappa$, the mean reversion of the variance; $\theta$, the long-term variance; and $\nu$, the volatility of variance. For a Poisson process $N$ with intensity $\lambda \ge 0$, $Y\sim {?}_1({\mu }_j,{\sigma }_j)$, ${\mu }_j\in ℝ$ and ${\sigma }_j\in {ℝ}^{+}$.

For the rough Bergomi model we take the logarithmic price process $X_t:=\log S_t$, $X_0=\log (s_0)=:x_0$. To model the instantaneous variance we consider the instantaneous forward variance ${\xi }_t^u$, $u\ge t$, for date $u$ observed at $t$. Let $\alpha =H-\frac{1}{2}\in (-\frac{1}{2},0]$ be a negative exponent depending on the Hurst exponent $H\in (0,\frac{1}{2}]$. For $\eta >0$, we have $X_0=x_0$ and ${\xi }_0^t=0$. The rough Bergomi model is given by:

(8)

with the variance obtained from (8) by setting:

$$V_t=\xi \exp (\eta Y-\frac{1}{2}{\eta }^2{t}^{2\alpha +1})$$

with $Y$ being the corresponding realisation of a Volterra process. We assume that a payoff $h$ is given and approximate the corresponding option price at time $t$.

We set the parameters as follows: $s_0=100$, $r=0.01$, $v_0={0.3}^2$, $\theta ={0.32}^2$, $\kappa =0.2$, $\nu =0.3$ and $\rho =-0.5$, $\lambda =0.2$, ${\mu }_j=0$ and ${\sigma }_j=0.25$ for the Bates model. For the rough Bergomi model we take $a=-0.43$, $\xi ={0.235}^2$, $\eta =0.9$ and $\rho =-0.5$.

Single underlying

We consider a call option and its delta for the Bates and rough Bergomi models. There is no obvious choice of $K$, since the number of components depends on the shape of the distribution, eg, skewness or kurtosis. Our numerical experiments have shown that the chosen number of Gaussian regressions ($3$ and $5$, respectively) lead to accurate results. In practice we recommend validating this assumption using cross validation or information criteria.

For numerical illustration we consider time $t=0.5$ and plot the conditional expectations for the calls and their deltas, together with quantiles $q\%$ and $(1-q)\%$, respectively ($q=0.5,1,5,10$), and compare them with both semi-analytical and nested Monte Carlo results (see figure 1). The absolute difference to these values, which serves as a benchmark, is only a few basis points and is much smaller than for the LSM. Some discrepancy is observed with respect to the proxy hedge, but this is due to the size of the time interval. The hedge still compares well to the instantaneously balanced hedge. Values of $t$ closer to $T$ lead to better results and improve upon the results stated in Geng et al (2022). Since we choose the underlying as the control variate, we derive the time-discrete minimal-variance delta analytically. All the simulations were done using 10,000 paths.

Figure 1: (a), (c) Conditional expectation and (b), (d) delta calculation for $t=\text{0.5}$, $T=\text{1}$ for (a), (b) the rough Bergomi model and (c), (d) the Bates model.

Let us summarise some of the findings. First, we observe that analytic smooth option prices and hedge sensitivities can be obtained using GMM-DCKE in a model-agnostic and data-driven way. Second, we observe the impact of the control variate; essentially, by including the proxy, we control the tail behaviour of the estimator. Thus, depending on the choice of underlying, we add linear tail behaviour, which means a slope of $0$ for the left tail and $1$ for the right tail, since our example is a call option. Compared with those in Geng et al (2022), the resulting hedges are smoother and more accurate, even for large time steps, eg, half-yearly ($t=0.5$). The difference from the displayed hedge sensitivities calculated using the model specifics is due to the numerical fit of the $\mathrm{GMM}(K)$ and because we consider time-discrete settings and not a continuous setting.

We are able to compute hedge sensitivities for the rough Bergomi model without nested Monte Carlo simulation.

For interest rate derivatives we consider the pricing of an at-the-money Bermudan swaption with the Cheyette model with and without stochastic volatility. Taking $x_0=y_0=0$ and $z_0=1$, the model is given by:

(9)

with $Z_i$, $i=1,2$, being uncorrelated Brownian motions. We take $\kappa =0.03$, $\beta =0.4$, $\theta =1$, $\epsilon =0.8$, ${\sigma }_0=0.2$, $m=0.15$ and $s_0=0.043011$. Using an alternating direction implicit (ADI) implementation with a grid of 200 points for the time component and 100 points for each space component. For the prices (yearly exercise), 5-year option maturity and 15-years underlying maturity, we obtain 572.4052 basis points (respectively, 552.2217bp) for the two-dimensional (2D) (respectively, three-dimensional (3D)) case. GMM-DCKE prices the options accordingly, with a maximum absolute difference of 1.284bp using 10,000 simulations. Both Python-based implementations were done on a Windows i9 machine with 16GB RAM and four cores using Python 3.6 with non-vectorised code. The runtimes are summarised in table A.

Table A: Runtimes (in seconds) pricing at-the-money Bermudan with ADI and GMM-DCKE for the Cheyette model.

	ADI		GMM-DCKE
Exercise	2D	3D	2D	3D
$Y$	0.28	0.94	0.37	0.83

GMM-DCKE can be considered as an alternative to other numerical methods for pricing interest rate derivatives, especially for more than two dimensions.

Multiple underlyings

We illustrate the performance of GMM-DCKE using a five-dimensional Heston model with asset correlation matrix $C_{\mathrm{aa}}$ and asset-variance correlation matrix $C_{\mathrm{av}}$ (compared to a full positive definite matrix $C$) defined, respectively, by:

	$C_{\mathrm{aa}}$
	$C_{\mathrm{av}}$	$=\mathrm{diag}{(-0.7,-0.8,-0.9,-0.8,-0.3)}^{\mathrm{T}}$

For the numerical experiments we show the performance of GMM-DCKE on the training set with 10,000 paths as well as using 500 paths for validation and we compare it to the DCKE method (Geng et al 2022). As an example, we consider a rainbow option with payoff given by:

$h_t^{\mathrm{rainbow}}$

$={\left(\max ({\displaystyle \frac{S_{1,t}}{S_{1,t_0}}},\mathrm{\dots },{\displaystyle \frac{S_{d,t}}{S_{d,t_0}}})-K\right)}^{+}$

(10)

We take $S_{u,0}=100$, $r=0.01$, $K=0.9$, $u=1,\mathrm{\dots },5$. For the $d$-dimensional setting we first derive $d$ $\mathrm{GMM}(K)$ representations for the random vectors $(S_{i,T},O_T,S_t)$, $i=1,\mathrm{\dots },d$, with $O_T$ being the option values at time $T$, $S_t\in {ℝ}^d$ being the asset prices at time $t$ and $S_{i,T}$ being the $i$th asset value, which serves as a control variate. We form the conditional distribution for each dimension on every $\mathrm{GMM}(K)$ component and calculate the conditional delta ${\delta }_{k,i}$ as well as conditional mean estimates ${\mu }_{k,i}$ for each $x_i$. The control variate is the weighted sum of the individual deltas estimated. For all samples $x_i$ we calculate the control variate using ${\sum }_{k=1}^K{\delta }_{k,i}(x_i-m_{k,i})$. As can be seen from figure 2, the results of GMM-DCKE and DCKE are close, and this is also observed for other payoffs, such as for basket options.

Figure 2: Rainbow option results for the five-dimensional Heston model with the (a) test and (b) validation sets for GMM-DCKE and DCKE.

Unlike DCKE, GMM-DCKE is based fully on the initial training set. To apply DCKE in multiple dimensions, the hedge sensitivities are computed by differentiating the payoff. Thus, DCKE depends on knowing the payoff function, and it is assumed that the corresponding derivatives can be calculated. With the hedge sensitivities calculated in this way, DCKE is not fully model-agnostic. With this additional adjustment it leads to results that are a little more accurate than those of GMM-DCKE. Fitting the full $d$-dimensional distribution can also be done, but it is computationally more intense.

The effect of using different control variates (ie, single assets, the max/min of all assets) is shown in figure 3. As expected, taking the $\max$ as a control variate leads to the best results.

Figure 3: Different control variates (single assets and the maximum and minimum) for GMM-DCKE.

Stochastic local volatility

For pricing and calibrating stochastic local volatility models, several techniques (ie, particle methods, binning and logarithmic normal approximations) have been considered (see Guyon & Henry-Labordère 2012; Muguruza 2020; van der Stoep et al 2014). We consider the dynamics of the forward, $S_t$, given by:

(11)

where $S_0=s_0$, $V_0=v_0$, $\mu$ and ${\sigma }_V$ are the drift and volatility functions, respectively, of the instantaneous variance, $\lambda :{ℝ}_{+}\times ℝ\to {ℝ}_{+}$ is the leverage function and $W$ and $Z$ are one-dimensional correlated Brownian motions. The function $\psi$ represents the stochastic volatility component. For the Heston model we have $\psi (X)=\sqrt{X}$. The SDE (11) can be reframed into a nonlinear equation in the sense of McKean and Vlasov, with the volatility depending on the probability distributions of $S$ and $V$. Showing the existence and uniqueness of the solutions to such equations is an involved mathematical problem, and for certain parameters, eg, for the stochastic volatility component, large values of the volatility of variance, there may not be any solution (see Guyon & Henry-Labordère 2012).

The market standard method for calibrating the model is to use a previously calibrated local volatility, ${\sigma }_D:{ℝ}^{+}\times {ℝ}^{+}\to {ℝ}^{+}$ (Dupire 1994), and a calibrated stochastic volatility model. The leverage function from (11) is derived by Markovian projection, and can be expressed as follows (see Guyon & Henry-Labordère 2012; van der Stoep et al 2014):

$${\lambda }^2(t,s)=\frac{{\sigma }_D{(t,s)}^2}{?[\psi {(V_t)}^2\mid S_t=s]}$$

(12)

The conditional expectation appearing in (12) needs to be calculated. This is approximated by GMM-DCKE. For calibration we take $T_1,T_2,\mathrm{\dots },T_N\in ?$ and consider:

• •

  $T_1$: no fit; use the initial value $\psi (V_t)$.

• •

  $T_2$: fit the $\mathrm{GMM}(K)$ model to the $T_1$ and $T_2$ values and store the $\mathrm{GMM}(K)$ parameters as ${\theta }_1$.

• •

  $T_n$: fit the $\mathrm{GMM}(K)$ model to the $T_{n-1}$ and $T_n$ values using ${\theta }_{n-1}$ as the initial values.

Using the most recently fitted parameters, we observed that the convergence can be improved. This relies (as do the other methods) on the input of the local volatility function ${\sigma }_D$. For widely fluctuating ${\sigma }_D$, it may be necessary to fully recalibrate at each step.

Figure 4: The conditional expectation $?[V_t\mid S_t=s]$ for the rough Bergomi model (a) with different numbers of Gaussians versus other market standard methods, (b) with bagging and (c) for stressed volatility of variance.

Figure 4 shows the accuracy of the method and also considers the version with bagging using only 20 runs with 10% (ie, 1,000) of samples, leading to speed-ups of 2–5 times and increased accuracy. Part (a) of the figure illustrates the application of GMM-DCKE to the rough Bergomi model for several choices of $K$. Part (b) shows the same model but with bagging applied. In our numerical experiments it turned out choosing $K=3$ is reasonable. This is also confirmed by taking the AIC and BIC. We also plotted the 1% and 99% quantiles for Gaussians versus binning, the 0.5%, 1%, 5% and 10% quantiles for the model with bagging and the 99.5%, 99%, 95% and 90% quantiles for the bagging case, respectively, to further illustrate the accuracy, even for the tails.

In terms of runtime we compared our method with the particle method, which is widely applied in financial institutions, and we observed our method performs slightly better when the parameter fitting is done using the previously calculated parameters and if we store the inverses of the matrices used for calculating the control variates (see (1), (2)). The particle method is often applied in conjunction with the Silverman rule, which allows a bandwidth to be picked for the kernel. Distributions that are very different from the Gaussian cannot be handled well with the assumptions made when applying this rule. A more sophisticated choice of the bandwidth (eg, a local bandwidth) may then be necessary. Determining the latter leads to optimisation problems that increase the runtime of the particle method.

For a Heston local volatility model we investigated a stressed parameter for a volatility of variance of 2.5 and observe that the particle and the method of Muguruza (2020) underestimate the smile. The binning method and GMM-DCKE, using $\mathrm{GMM}(3)$, lead to reasonable results (see figure 4(b)).

Finally, we have considered the same example as in Muguruza (2020). Since exactly the same data were not available, we calibrated to the data we observed and plotted a Heston smile that is close to the one shown in Muguruza (2020). We observe that applying GMM-DCKE for calculating the conditional expectations and determining the implied volatility leads to nearly identical results. The largest difference was less than 4.3bp for implied lognormal volatility, with an average error of 2.1bp. Furthermore, since GMM-DCKE does not assume a certain model, it can be applied to the more complex models considered in Muguruza (2020).

Conclusions and summary

GMM-DCKE is a purely data-driven and model-agnostic method to compute conditional expectations. We applied it to the pricing of (multi-dimensional) Bermudan options and the pricing and calibration of stochastic local volatility models. It also leads to hedging strategies with the realisations at given time points but without the need to know the underlying stochastic model. GMM-DCKE also can also be used for computing control variates for standard Monte Carlo simulation in multiple dimensions, where finding appropriate controls that can be computed fast is often difficult. In contrast to the method in Geng et al (2022), it is independent of the choice of (local) bandwidth for kernel density estimation. GMM-DCKE can be integrated into existing set-ups as an alternative to using other numerical methods. Due to its analytical nature and application of the GMM it only needs a decent amount of training data. Relevant practical examples show that its performance is at least comparable with existing numerical methods and can be enhanced by careful implementation of the method using vectorisation or efficient EM algorithms. Using methods such as bagging improves upon the suggested model. Further ideas – for instance, applying autoencoders with GMM-DCKE in low-dimensional latent space or using different distributions – are the subject of ongoing research.

Jörg Kienitz is an adjunct associate professor at the University of Cape Town and the University of Wuppertal. He owns finciraptor.de and is a partner in the quantitative services unit of Acadia. He is based in Bonn. He thanks Gordon Lee and Thomas McWalter for fruitful discussions, and two anonymous referees who helped to greatly improve the paper.

Email: [email protected]

References

• Bates D, 1996
  Jumps and stochastic volatility, exchange rate processes implicit in Deutsche mark options
  Review of Financial Studies 9, pages 69–107
• Bayer C, P Friz and J Gatheral, 2016
  Pricing under rough volatility
  Quantitative Finance 16(6), pages 887–904
• Bishop CM, 2006
  Pattern Recognition and Machine Learning
  Springer
• Dempster AP, NM Larid and DB Rubin, 1977
  Maximum likelihood from incomplete data via the EM algorithm
  Journal of the Royal Statistical Society B 39(1), pages 1–38
• Dupire B, 1994
  Pricing with a smile
  Risk July, pages 18–20
• Geng Q, J Kienitz, GT Lee and N Nowaczyk, 2022
  Dynamically controlled kernel estimation
  Risk February, pages 110–115, http://www.risk.net/7921186
• Guyon J and P Henry-Labordère, 2012
  Being particular about calibration
  Risk January, pages 88–93, http://www.risk.net/2135540
• Halperin I, 2017
  QLBS: Q-learner in the Black-Scholes (-Merton) worlds
  Preprint, SSRN, December (https://ssrn.com/abstract=3087076)
• Huge B and A Savine, 2020
  Differential machine learning: the shape of things to come
  Risk October, pages 76–81, http://www.risk.net/7688441
• Kienitz, J, 2021
  GMM-DCKE: semi-analytic conditional expectations
  Preprint, SSRN, August (https://dx.doi.org/10.2139/ssrn.3902490)
• Longstaff FA and ES Schwartz, 2001
  Valuing American options by simulation: a simple least-squares approach
  Review of Financial Studies 14(1), pages 113–147
• Muguruza A, 2020
  Not so particular about calibration: smile problem resolved
  Preprint, SSRN (https://ssrn.com/abstract=3461545)
• van der Stoep AW, LA Grzelak and CW Oosterlee, 2014
  The Heston stochastic-local volatility model: efficient Monte Carlo simulation
  International Journal of Theoretical and Applied Finance 17(7), article 1450045

• Tweet  
• Facebook  
• LinkedIn  
• Save this article
• Send to  
• Print this page  

• Topics
• Bermudan options
• Data
• Exotic options
• Gaussian model
• Hedging
• Pricing
• Stochastic local volatility (SLV)
• Cutting Edge
• Derivatives
• Banking
• Machine learning