JMF  Vol.10 No.4 , November 2020
A Unified Stochastic Volatility—Stochastic Correlation Model
Abstract: This paper has two main contributions. First, we build a simple but rigorous stochastic volatility—stochastic correlation model. Mean-reverting and locally stochastic with dependent Brownian motions, our model proves to fit both marginal and joint distributions of the option market implied volatility and correlation. Second, asset correlations are currently modeled exogenously and then ad hoc assigned to an asset price process such as the Geometric Brownian Motion (GBM). This is conceptually and mathematically unsatisfying. We apply our approach to build a unified asset price—asset correlation model, which outperforms the standard GBM significantly.

1. Introduction and Motivation

Asset prices are typically modeled with the Geometric Brownian motion (GBM) of the form

d S t / S t = μ d t + σ d B t (1)

where St is the asset price, μ is the drift, σ is the volatility, and dBt is a standard Brownian motion.

Rapid developments in vanilla and exotic options markets, over the past several decades, have fundamentally challenged the static assumptions of μ and σ parameters in GBM. Merton [1] introduces jumps and shows that if the logarithm of the percentage jump is normally distributed, a closed form solution for European style options exists. Cox and Ross [2] create the constant elasticity of variance (CEV) model, where an exponential parameter α is added to the asset price. The value of α determines the dependence between asset price and volatility. In a pure jump extension, Madan et al. [3] follow a variance-gamma approach, which creates heavier tails and provides semi-analytic expressions for European style options.

Heston [4] in a seminal model correlates the asset return process with stochastic variance. Many Heston model extensions have since developed with Zhou [5], Hagan et al. [6], Brigo and Pallacinini [7], and Langnau [8] counted among many prominent examples.

The availability of market data of implied volatility surface1 has ushered in a new era of stochastic volatility modeling to reproduce a multitude of finer properties observed in the real world. Cont et al. [9] prove that 1-factor model is insufficient to represent the true dynamics seen in SPX volatility surfaces. Furthermore, the eigenfactors of the vol-surface time series are found not to be perfectly correlated to the underlying asset price movements—hence concluding—“Vega” risk cannot be reduced to “Delta” risk.

In a series of ground-breaking papers, Bergomi [10] and others set out to capture the term-structure of skew and vol-of-vol in a Forward Variance Swap framework, an analog of the Libor Market Model. A class of models has since multiplied aimed to price the “smiling” volatility derivatives consistently to vanilla options, therefore solving the “Joint S&P/VIX Smile Calibration Puzzle”. By adding simultaneous jumps to the Ornstein-Uhlenbeck process of the forward variance swap and the underlier GBM process, the flexible Lévy specification by Cont et al. [11] offers a greater analytical tractability and more efficient control of the vol-skew and the shifting correlation between spot and implied volatility.

Joint modeling of asset volatility and asset correlation should attract more research interest, as the need arises not only for asset allocation but also for derivatives pricing and risk management. By comparison to stochastic volatility, less progress has been made in the literature regarding asset correlation as a state-dependent market factor.

The summary statistics of realized stock correlations in Table 1 reveals that the level and variability of cross-sectional asset correlation both move higher conspicuously when market is in distress, from a distinctly lower regime observed during expansionary economic times. Meissner [12] further asserts that elevated correlations persist without necessarily accompanied by rising volatility during recessionary periods.

Table 1. Dow stock correlation level and Std. Dev. from Jan. 1972 to Mar. 20202.

Analogous to the irreducible “Vega” risk in the Implied Volatility (IV), the empirical evidence clearly supports the role of the Implied Correlation (IC) as a unique source of randomness to the financial system.

Hull et al. [13] model the asset process in a GBM setting and then sample the asset correlation from a beta distribution, with the distribution parameters exogenously derived and not time-varying. Ma [14] models stochastic volatility and stochastic correlation processes to price exchange rate options, and Ma [15] applies stochastic correlation to price and hedge multi-asset options, the volatility and correlation are however assumed independent from each other.

Emmerich [16] highlights the unique mathematical properties of a stochastic correlation process. Dϋllman et al. [17] model stochastic correlation with a Vasicek process. However, both papers do not combine the correlation process with the asset process.

Buraschi et al. [18] and Fonseca et al. [19] extend the Heston [4] model, by correlating a n-dimensional stochastic correlation process with a n-dimensional stochastic asset process.

The term of correlation is ubiquitous. It should be emphasized that the asset correlation in our study is restricted to the equity price returns. In contrast to Burtschell et al. [20], the local stochastic correlation model targets the default time correlation in a portfolio of Credit Default Swaps (CDS) instead.

The remaining paper is structured as follows: In Section 2 we introduce our unified stochastic asset volatility—stochastic asset correlation model. In Section 3 we show the real world fit of the model. Section 4 discusses the calibration of the eight parameters. Section 5 conducts the significance tests. Section 6 applies the model in an enhanced version of GBM.

2. The USVSC Model

Our proposed Unified Stochastic Volatility—Stochastic Correlation model (USVSC henceforth) consists of three stochastic differential equations:

{ d σ t = κ σ ( m σ σ t ) d t + ν σ σ t β d W t ( 2 ) d ρ t = κ ρ ( m ρ ρ t ) d t + υ ρ 1 ρ t 2 d Z t ( 3 ) d W t = ρ w d Z t + 1 ρ w 2 d Z t ( 4 )

where stochastic instantaneous volatility σ t and stochastic instantaneous correlation ρ t are set within a local volatility framework, mean-reverting with a rate of κ { σ , ρ } to a long-term mean of m { σ , ρ } . The diffusion coefficient is denoted as ν { σ , ρ } > 0 , and the positive skew in σ t is captured by a power parameter β . The modeled stochasticity is generated by three 1-dimensional standard Brownian motions { W t , Z t , Z t } , where Z t and Z t orthogonal, W t and Z t correlated with a coefficient ρ w .

Equation (2), referred to as CIRCEV hereafter, extends the Cox-Ingersoll-Ross (CIR) model with a β exponent parameter, otherwise known as Chan-Karolyi-Longstaff-Sanders (CKLS) model proposed by Chan et al. [21], or a variant of the Constant Elasticity of Variance (CEV) model. Evoking CIRCEV model specification is motivated by its non-negativity and mean-reversion properties fitting the implied volatility, and the presence of a positive skew found empirically and prominently displayed in Figure 1(c).

Equation (3) specifies a Jacobi process, otherwise known as Wright-Fisher diffusion process. Jacobi is a natural model choice for implied correlation, in which not only the mean-reverting empirical feature is preserved, the limiting behavior for a correlation is also obeyed naturally at −1 or +1. As discussed by Ahdida et al. [22], Emmerich [16] and Zetocha [23], Jacobi specification permits analytical solutions of the first two moments, as well as the average over an arbitrary length of time. Although, like Heston [4], a Feller-like inequality constraint should apply over model parameters to ensure the desired numerical stability.

Equation (4) gives a dependence structure between the implied volatility and the implied correlation, analogous to Heston [4] in treating the joint distribution between spot and local variance.

In summary, our USVSC model is formulated as a correlated CIRCEV and Jacobi processes in an attempt to capture the equity derivatives dynamics. The rational for disallowing the joint-dynamics between derivatives and underlying spot market is two-fold: 1) to keep the model at the lowest possible dimension in support of historical calibration; and 2) to expose the risk factors exclusively within derivatives, as opposed to a joint model traditionally done for efficient hedging purpose.

Our model consists of a total of eight parameters, namely, θ = { κ σ , m σ , ν σ , β , κ ρ , m ρ , ν ρ , ρ w } , forming the set of unknowns to be estimated empirically, with a goal to better understand the historically realized marginal and joint distributions of the risk neural option implied volatility and implied correlation through time.

3. Real World Fit

To measure the efficacy of our model, we let VIX represent IV (the option-implied volatility σ t ), and the CBOE Implied Correlation as IC (the option-implied correlation ρ t )3.

Figure 1. (a) Empirical relationship between IV and IC. (b) Time series plot of the empirical IV and IC. (c) Density of empirical IV. (d) Density of empirical IC.

Our dataset is comprised of the SPX, the VIX and the ICJ/JCJ indexes observed from 01/03/2007 to 08/01/2019, daily sampled and sourced from CBOE publishes both historical data series of the Implied Volatility (known as the VIX—the “fear gauge”), and the Implied Correlation, derived from standard 1-month options and dispersion trading strategies respectively, under the methodologies given by CBOE [24] and CBOE [25].

We first display the excellent real world fit of our model, which is followed by a detailed description of a 2-staged calibration and the correspondent statistical test results.

Figure 1 describes historical properties of the Implied Volatility (IV) and the Implied Correlation (IC).

Figure 2 displays the simulated properties of the Implied Volatility (IV) and the Implied Correlation (IC) using our calibrated model results.

With respect to Figure 2(a), the USVSC model is designed to produce a positive-“triangular” relationship between IV and IC as shown in Figure 1(a). A positivity in ρ w is expected to enable synchronized co-movements between σ t and ρ t . Furthermore, as the correlation approaches upper or lower bounds, i.e. +1 or −1, its stochasticity decreases due to the diffusion term υ ρ 1 ρ t 2 specified in Equation (3).

After outlining our calibration procedure in Section 3, we perform a series of statistical tests associated with Figures 2(a)-(d) in Section 4.

4. Model Calibration

The estimation scheme is described in this section. We approach the calibration in three Stages using the Maximum Likelihood Estimation (MLE) framework. We believe MLE is simplest to implement, numerically stable while allowing for the generation of the desired statistical inferences.

Figure 2. (a) Modeled relationship between IV and IC. (b) Time series plot of the modeled IV and IC. (c) Density of modeled IV. (d) Density of modeled IC.

4.1. Stage I—Marginal Distributions of IV and IC

By design of the proposed USVSC model, except for the correlated Brownian motions anchored by Equation (4), IV and IC have no interactions in the drift as well as the diffusion terms, divergent from that in Heston [4]. Under this specification, therefore, a decoupled and 1-dimensional historical calibration, separated for IV or IC process, is permissible.

To illustrate, we first suppose a generalized univariate stochastic differential equation (SDE) with drift and diffusion terms fully parameterized as follows:

d X t = f ( X t , θ ) d t + g ( X t , θ ) d W t (5)

where { X t } represents the time series of empirical observations, and θ the collection of model parameters. We assume the increments of observations d X t conditionally Gaussian N ( μ X , σ X 2 ) . This means, based on the Euler discretization scheme, we can derive explicitly the first two moments of the increments, similar to the approach taken by Ait-Sahalia [26]:

μ X = E ( Δ X t ) = f ( X t , θ ) Δ t , σ X 2 = E [ V a r ( Δ X t ) ] = g ( X t , θ ) 2 Δ t (6)

The transition density, under our assumption, can therefore be expressed in a closed-form:

p θ ( Δ t , x i | x i 1 ) = 1 2 π Δ t g ( x i 1 , θ ) 2 exp [ ( x i x i 1 f ( x i 1 , θ ) Δ t ) 2 2 Δ t g ( x i 1 , θ ) 2 ] (7)

The optimization objective function immediately follows, resulting in a maximization problem of the log-likelihood function to solve θ :

Θ = arg max θ log [ i = 1 N p θ ( Δ t , x i | x i 1 ) ] = arg max θ { log [ 2 π i = 1 N g ( x i 1 , θ ) 2 ] + i = 1 N ( x i x i 1 f ( x i 1 , θ ) Δ t ) 2 2 Δ t g ( x i 1 , θ ) 2 } (8)

With an unconstrained optimization routine4, we first report in Table 2 the results of USVSC model parameter estimation of θ , along with the corresponding estimation errors.

The calibrated parameters of IV, per Equation (2)—CIRCEV specification, including the long-term mean, the reversion rate (half-life), the magnitude in vol-of-vol, together with the CEV skew parameter β > 1 , are all of expected values with a good match to existing literature, such as Ait-Sahalia and Kimmel [27], and Bu et al. [28].

The calibrated parameters of IC, per Equation (3) Jacobi process, are shown to have met the inequality condition of the mean-reverting rate, κ ρ max ( υ ρ 2 1 + m ρ , υ ρ 2 1 m ρ ) , consistent to relevant discussions in Emmerich [16], Ma [14] and [15], as well as Zetocha [23].

To test the robustness, we produce 5k simulations to visualize in Figure 3. We confirm the fulfillment of correlation boundedness achieved by the calibrated parameters. In addition, the mean is fairly stable around the randomly assigned initial value of 0.2.

The red line in Figure 3 represents the time series mean of all simulated IC paths, which exhibits a relative stability as IC converges to a long-term mean prescribed by a Jacobi process as in Equation (3). Measured by each simulated single path shown in Figure 2(b), and relative to IV, however, a greater magnitude of variability in IC empirically observed in Figure 1(b) is largely captured, helped in part by a faster mean-reverting rate κ ρ ^ .

Figure 3. Time series of simulated implied correlation.

Table 2. Univariate calibration of IV and IC.

4.2. Stage II—Correlated Brownian Motions

Having estimated the seven parameter values in θ , { κ σ , μ σ , υ σ , β , κ ρ , μ ρ , υ ρ } , required by the model, we now turn our focus to the eighth parameter ρ w and the correlated Brownian motions in Equation (4).

We reconstruct the realized increments of Brownian motions, d W t and d Z t , making use of the estimated parameters from Stage I. Under the normality assumption, it follows that the asymptotic correlation coefficient estimation as specified in Equation (5) is just the expectation of the product of d W t and d Z t scaled by dt:

ρ w = E ( d W t × d Z t ) / d t (9)

The estimated correlation coefficient parameter ρ w ^ , together with the t-statistic and p-value are presented in Table 3. The zero-value null hypothesis is strongly rejected at 95+% confidence level. This is an important finding from our Stage II calibration.

The histogram of the product of the white noises is displayed in Figure 4 for illustration purpose.

5. Goodness-of-Fit Tests

In this section, we test if our calibrated model is statistically accurate in representing the empirical data. The marginal distributions of IV and IC resulted from calibration Stage I are first examined, which are followed by normality tests concerning the correlation parameter ρ w estimated in calibration Stage II.

Figure 4. Histogram of d W t × d Z t / d t .

Table 3. Correlation calibration result of pw of Equation (4).

H 0 : ρ w = 0 .

5.1. Testing for Stage I

To test the equality in aggregation distribution between the historical IV and IC densities in Figure 1(c) and Figure 1(d), and our model-simulated densities in Figure 2(c) and Figure 2(d), we proceed as follows.

Firstly, in Table 4, we perform non-parametric Kolmogorov-Smirnov (KS) test for the difference in univariate distributions. With the cumulative distribution functions (CDF) denoted as F(IV) and F(IC), respectively, between the empirical distributions and the modeled 1-path time series simulation, denoted as F ^ ( IV ) and F ^ ( IC ) , started with a pair of arbitrary points (IV0, IC0). Below 5% are the max distances between two CDFs measured by the D-Scores associated with IV and IC. Nearing 1% p-values found by the tests, however, suggest the null hypothesis should be rejected with 95% confidence level.

Secondly, in Table 5, we utilize a form of randomization test—the Permutation Test of Equality. With paired samples of univariate distributions, the high p-values registered for both IV and IC suggest a high degree of equality in distribution between the model generated IV or IC process and the correspondent empirical counterpart.

Graphical depiction rendered in Figure 5 shows the density functions of IV and IC, respectively, between the simulated series and the empirical observations. The blue color identifies the largest distances in density, along the x-axis continuous values of likely (IV Î [0, +1]) and feasible (IC Î [−1, +1]) ranges. This should shed additional light on the locality of the discrepancies in the interest of future studies.

Lastly, we study the sensitivity of inference test results to the numerical noises embedded in Monte Carlo simulations. When the number of simulation paths increases, the replicability of the distribution by our parameterized model is expected to significantly improve. This effect can be demonstrated via an unpaired 2-sample Wilcoxon test, also known as Mann-Whitney test.

Table 4. 1-Path Kolomogorov-Smirnov (KS) test for equivalence between empirical and the simulated 1-path time series with an arbitrary starting point of IV0 and IC0.

H 0 : F ( IV ) = F ^ ( IV ) , F ( IC ) = F ^ ( IC ) .

Table 5. 1-Path permutation test.

H 0 : F ( IV ) = F ^ ( IV ) , F ( IC ) = F ^ ( IC ) .

Figure 5. 1-Path IV and IC permutation test.

Judging by the large test statistics W and the high p-values (>>10%) reported in Table 6, resulted in by comparing 5K paths simulation against the single-path historical realizations, one can’t reject the null hypothesis of identical distribution for both IV and IC at 95% confidence level.

5.2. Testing for Stage II

Our SDE Equation (4) assumes a bivariate normal distribution generated by d W t and d Z t , two Wiener processes underlying IV and IC observables with an instantaneous and constant correlation parameter ρ w . We now examine the normality of the model derived Brownian increments.

First we present the density function plots with a normal curve overlay, labeled as Normality Test (1) in Figure 6. We observe a very mild case of fat-tail (leptokurtic) in distributions of d W t and d Z t .

Normality Test (2) is consisted of a pair of QQ-plots and ACF autocorrelation graphs. Figure 7 displays that a wide-ranged linearity and low autocorrelation have been tested in the time series of estimated residual, in support of the Gaussian distribution assumptions for d W t and d Z t .

With a reasonably high degree of confidence in normality, we conclude that our proposed Stage II calibration method for ρ w is consistent with empirical evidence.

5.3. Testing the Joint Distribution of IV and IC

Clearly Figure 2(a) makes a strong case in support of our USVSC model, as it successfully captures the very key feature of the real word implied volatility-correlation relationship as depicted by Figure 1(a).

To statistically measure the fit, we apply Fisher’s Z-transformation to test the Pearson correlation difference of two un-paired independent samples: one is the historical realization sample of IV and IC, and the other is estimated from the 5K paths of 2-dimensional simulation of IV and IC.

Table 6. 5K-Path Wilcoxon test.

H 0 : F ( IV ) = F ^ ( IV ) , F ( IC ) = F ^ ( IC ) .

Figure 6. d W t and d Z t normality tests (1).

Figure 7. d W t and d Z t normality tests (2).

With a zero correlation difference as the null hypothesis, we obtain a Z-Score of 1.59 reported in Table 7. This test statistics indicates that the Pearson correlation estimation, between simulated IV and IC, is indistinguishable in distribution from that of the empirical observed IV and IC, over the sampled time period, at 95% confidence level.

Table 7. IV and IC dependency test.

H 0 : ρ ( IV , IC ) = ρ ^ ( IV ^ , IC ^ ) .

We, therefore, conclude that the empirical joint distribution of IV and IC is replicated by our proposed USVSC model with sufficient accuracy.

6. Application of the Model—An Enhanced Geometric Brownian Motion

With the stochastic IV and IC fully parameterized and tested, we can now re-visit GBM model log-normal specification in Equation (1), the most foundational claim of the asset price returns in the finance literature.

Using GBM as the benchmark, Equation (1) parameters {μ, σ} are first calibrated through MLE. Table 8 gives the estimation results along with the corresponding standard errors.

Although the calibrated spot price volatility σ under GBM is robust and statistically significant, the same conclusion may not be drawn with respect to the parameter μ. The large estimation error suggests that the null hypothesis (μ = 0) is rather unrejectable with confidence, a subject to which we will return later.

Integrating our stochastic local volatility and correlation model, we proceed to replacing the constant diffusion coefficient σ with our time-varying and stochastic IV denoted by σt, as defined in Equation (2). The GBM Equation (1) thereby evolves into Equation (1a)—enhanced-GBM model, allowing the cash spot market infused with additional information sourced from the derivatives market:

d S t / S t = μ d t + σ t d B t (1a)

where d B t denotes a new Brownian motion for price returns, distinct from the original dBt from Equation (1).

To incorporate the proposed USVSC model, we further assume that the risk factors inherent to IV and IC markets are also integral to the dynamics of index trading, therefore can be thought of as priced state variables for asset pricing. This is not unreasonable, as the derivatives IV and IC have become standalone asset classes, albeit with a shorter history and a narrower investor base relative to the underlying cash market. VIX, for example, has evolved from being a latent risk parameter to a set of publicly investable trading product group largely helped by the creation of VIX Index and VIX futures commenced in 2004. Not unexpectedly, the stochasticity of IC is argued by Emmerich [16] as a fundamental source of risk, which is furthered by Driessen et al. [29] in a hunt for risk premium in the context of asset pricing. The roles of derivatives markets affecting index price discovering and the broader market have become more evident recently. Details and in-depth discussions on the market structure of VIX fixing, dispersion and worst-of option pricing, can be found in Osterrieder et al. [30], Li [31] and Zetocha [23].

Table 8. Equation (1) GBM model calibration.

Across three stylized state variables, i.e. the Brownian motion increments d B t together with dWt and dZt, a deterministic dependence structure is superimposed. This is manifested by Equation (10) and Equation (11) sequentially, with ρB and ρZ denoting the instantaneous correlation coefficients between spot-IV and spot-IC, respectively.

Since dWt and dZt are interconnected through ρw, as previously established by Equation (4), one can write:

d Z t = ρ w d W t + 1 ρ w 2 d W t (10)

with Wt W t as usual. Helped with the introduction of a third orthogonal factor W t , we can now decompose the Weiner process d B t in a 3-independent-factor representation:

d B t = ρ B d W t + ρ Z ρ w ρ B 1 ρ w 2 d W t + 1 ρ B 2 ( ρ Z ρ w ρ B ) 2 1 ρ w 2 d W t (11)

So that the recoveries of ρ B = E ( d B t × d W t ) / d t , ρ Z = E ( d B t × d Z t ) / d t and the identity property of E ( d B t × d B t ) / d t = 1 are satisfied. This particular setup should allow us to exploit the connections between index price dynamics with IV and IC observables5.

It should be noted that by now dWt and dZt have been attained from calibration Stage II. Therefore d W t and d W t are the two added sources of independent Gaussian white noise under the newly enhanced-GBM model Equation (1a).

Our task is therefore reduced to estimating the drift μ and correlation coefficients ρB and ρZ, and testing if the enhanced-GBM Equation (1a) is valid given our model setting. We proceed as follows.

Let xt = ln(St) to simplify the notation for asset price returns, after applying Ito’s lemma, Equation (1a) is expanded to a 3-factor model Equation (1b), admitting the information transformation from IV and IC markets rather explicitly:

d x t = ( μ 1 2 σ t 2 ) d t + σ t ( ρ B d W t + ρ Z ρ w ρ B 1 ρ W 2 d W t + 1 ρ B 2 ( ρ Z ρ w ρ B ) 2 1 ρ W 2 d W t ) (1b)

Taking the expectations on both sides, aided by the assumed Gaussianity of { d W t , d W t , d W t } ~ N ( 0 , d t ) , we arrive at the drift rate estimator

μ ^ = E [ d x t d t + 1 2 σ t 2 ] (12)

where all RHS variables are known, we obtain a calibration result and a statistical inference test in Table 9. At 95% confidence level, the null is once again unrejectable-same conclusion under the original GBM model. We thereby argue that the constant drift assumption under GBM (1) or (1a) is unsupported by the sampled historical data. That is, the linear drift term under GBM is not robust and likely misspecified, irrespective of the diffusion term being constant or locally stochastic.

It is clear now that the calibration of ρ B and ρ z is equivalent to locating the maximization of the log-likelihood (LL) function of a standard normal density over the full historical data sample:

d W t d t = d t σ t ( ρ z ρ B ρ w ) ( d x t + 1 2 σ t 2 d t μ d t ) 1 ρ w 2 ~ N ( 0 , 1 ) (13)

Figure 8 pictorializes the calibration approach, once again by utilizing the Maximum Likelihood Estimation approach. The RHS expression of Equation (13) containing 2 correlation parameters is assumed to follow a standard Gaussian distribution. Hence the inverse of a max likelihood function, i.e. LL = −Log(L), reaches the minimum while the likelihood function L is maximized.

An implementation for a non-linear optimization algorithm yields a solution pair of ρ B ^ and ρ z ^ , displayed in a matrix form in Table 10 together with the estimated ρ w ^ for completeness:

The estimated negative correlation ρ B ^ , between spot and IV, corroborates the leverage effect puzzle termed by Ait-Sahalia et al. [33], and converges with the results by Zetocha [23]. This is also a widely known phenomenon by retail investors, as described in CBOE [34].

The negative correlation ρ z ^ , between spot and IC, comes as no big surprises as stock prices tend to correlate more as the market is down than the reverse, consistent with a unique and well-known IC market feature—correlation skew by Reghai [35], Zetocha [23] and Delanoe [36]. In the academic literature, however, the interactivities of spot-IC have not been extensively studied relative to spot-IV.

Several normality tests undertaken for d B t d t estimated under GBM (1), as well as d B t d t under GBM (1a): Table 11 gives the numerical improvements in Jargue Bera statistics; Figure 9 charts the normalized density plots benchmarked against the standard normal; Figure 10 supplements QQ-plots. The substantial gains in performance, of our joint stochastic Implied Volatility and Implied Correlation model over the original GBM model, are demonstratively significant.

Our enhanced GBM model of Equation (1a) may also underlie a portfolio selection process. To avoid penalizing upside volatility potential, a semivariance volatility approach could be applied, see for example Estrada [37] or Chen et al [38].

7. Concluding Summary

The dynamics of asset prices is typically modeled by stochastic processes; however their correlations are exogenously modeled and ad hoc assigned to the price process. This is mathematically and conceptually unsatisfying.

Table 9. GBM (1a) calibration and inference test of drift μ.

H 0 : μ = 0 .

Table 10. Correlation matrix.

Table 11. GBM normality test.

H 0 : { d B t d t , d B t d t } ~ N ( 0 , 1 ) .

Figure 8. MLE solving ρ B and ρ z .

Figure 9. Density distribution comparison.

Figure 10. QQ-plot comparison.

We propose a simple and unified framework to model local stochastic volatility and local stochastic correlation, similar to the approach by Heston [4]. Powered by coupled CIRCEV and Jacobi processes, our structural USVSC model is parsimonious and capable of re-producing the non-trivial volatility and correlation skews observed in the real world. The calibration scheme developed is intuitive and tractable, with results proven significant. The goodness of fit is backed by robust inference results after undertaking multiple stringent statistical tests in both marginal and joint distributions of the implied volatility and implied correlation.

We show our well-calibrated model performs markedly better over the classical GBM approach, thus a suitable and practical challenger for asset pricing, especially under current market structure.

Future research of our USVSC model may entail analyzing the sensitivity and stability to historical data sampling periods, and testing its capabilities in fitting the volatility smiles, term structure and correlation skews, for the purpose of pricing basket options and risk management. Additionally, new trading strategies may be developed if the equity option-implied correlation skews have awareness of the CDO tranche-implied correlation skews, across firm capital structure and asset classes. Lastly, placing the implied correlation and joint dynamics of volatility-correlation in a wider context of risk premium study should add useful insights in the field of financial economics.


1The 3-D surface σ t B S ( K , T ) is implied from vanillaoption pricesvia Black-Scholes-Merton formula, with strike K and maturity T.

2Correlations estimated are the average of Pearson 30 × 30 correlation matrix of Dow stock price returns.

3Strictly speaking, the VIX only corresponds to the estimated implied volatility of 30-day SPX call/put contracts through a convergence framework between variance swap pricing and log-Price option contract. In order to model VIX Options and Futures over tenor T, consistently with SPX options, it is most common to express VIX as the result of the integration of an instantaneous forward variance rate σ ˜ τ , T 2 , as in Cont et al. [11] and Bergormi [10], rather than the instantaneous volatility stochastic factor σ t as we proposed in Equation (2): VIX τ , T = E ( 1 T τ τ + T σ ˜ τ , u 2 d u ) where T = 30 days. A similar expression exists for IC with respect to the quantity of ρ t . To simplify the calibration, without evoking the instantaneous forward variance curve nor the required instantaneous forward correlation curve, we implicitly assume σ t : = σ ˜ τ , τ + T and ρ t : = ρ ˜ τ , τ + T for illustration purpose. Essentially a 1-month forward term-structure is embedded in both stochastic variables. This slight abuse of notation is acceptable, since 1) T is fixed at 30 days everywhere, and 2) the emphasis of our study is structured around the systemic interactions between IC and IV, impacts from term-structure and optimal risk hedging strategies are left out for future studies.

4We thank Prof. Boukhetala for R package Sim. Diff. Proc. and the useful discussions.

5By this design, the risk transmission from IV into spot market is channeled through correlated Brownian motions, as the stochastic σ τ in Equation (1a) is set to historical observations directly. For additional insights of CKLS, readers are referred to Chan et al. [21] and Hu et al. [32].

Cite this paper: Lu, X. , Meissner, G. and Sherwin, H. (2020) A Unified Stochastic Volatility—Stochastic Correlation Model. Journal of Mathematical Finance, 10, 679-696. doi: 10.4236/jmf.2020.104039.

[1]   Merton, R.C. (1976) Option Pricing When the Underlying Stock Returns Are Discontinuous. Journal of Financial Economics, 3, 125-144.

[2]   Cox, C. and Ross, S. (1976) The Valuation of Option for Alternative Stochastic Processes. Journal of Financial Economics, 3, 145-166.

[3]   Madan, D., Carr, P. and Chang, E. (1998) The Variance Gamma Process and Option Pricing. Review of Finance, 2, 79-105.

[4]   Heston, S. (1993) A Closed-Form Solution for Options with Stochastic Volatility with Applications to Bond and Currency Options. The Review of Financial Studies, 6, 327-343.

[5]   Zhou, C. (2001) An Analysis of Default Correlations and Multiple Defaults. The Review of Financial Studies, 14, 555-576.

[6]   Hagan, P.S., Kumar, D., Lesniewski, A. and Woodward, D. (2002) Managing Smile Risk. Wilmott, 1, 84-108.

[7]   Brigo, C. and Pallacinini, A. (2008) Counterparty Risk and CCDSs under Correlation. Risk, 21, 84-88.

[8]   Langnau, A., Allianz, S.E. and Munch, L.M.U. (2009) Introduction into Local Correlation Modelling. Working Paper.

[9]   Cont, R. and Fonseca, J. (2002) Dynamics of Implied Volatility Surfaces. Quantitative Finance, 2, 45-60.

[10]   Bergomi, L. (2005) Smile Dynamics II. Risk Magazine.

[11]   Cont, R. and Kokholm, T. (2013) A Consistent Pricing Model for Index Options and Volatility Derivatives. Mathematical Finance, 23, 248-274.

[12]   Meissner, G. (2019) Correlation Risk Modeling and Management. RISK Books.

[13]   Hull, J., Predescu, M. and White, A. (2005) The Valuation of Correlation-Dependent Credit Derivatives Using a Structural Model. Journal of Credit Risk, 6, 99-132.

[14]   Ma, J. (2009) Pricing Foreign Equity Options with Stochastic Correlation and Stochastic Volatility. Annals of Economics and Finance, 10, 303-327.

[15]   Ma, J. (2009) A Stochastic Correlation Model with Mean Reversion for Pricing Multi-Asset Options. Asia-Pacific Financial Markets, 16, 97-109.

[16]   Emmerich, C. (2006) Modeling Correlation as a Stochastic Process. Bergische Universitaet Wuppertal Working Paper.

[17]   Dϋllmann, K., Kull, J. and Kunisch, M. (2008) Estimating Asset Correlations from Stock Prices or Default Rates—Which Method Is Superior? Deutsche Bundesbank Working Paper.

[18]   Buraschi, A., Porchia, P. and Trojani, F. (2010) Correlation Risk and Optimal Portfolio Choice. Journal of Finance, 65, 392-420.

[19]   Fonseca, J., Graselli, M. and Ielpo, F. (2008) Estimating the Wishart Affine Stochastic Correlation Model Using the Empirical Characteristic Function.

[20]   Burtschell, X., Gregory, J. and Laurent, J.-P. (2007) Beyond the Gaussian Copula: Stochastic and Local Correlation. Journal of Credit Risk, 3, 31-62.

[21]   Chan, K., Karolyi, A., Longstaff, F. and Sanders, A. (1992) An Empirical Comparison of Alternative Models of the Short-Term Interest Rate. Journal of Finance, 47, 1209-1227.

[22]   Ahdida, A. and Alfonsi, A. (2013) A Mean-Reverting SDE on Correlation Matrices. Stochastic Processes ad Their Applications, 123, 1472-1520.

[23]   Zetocha, V. (2014) Skewing Up Correlation: Understanding Correlation Skew in Equity Derivatives.

[24]   CBOE (2019) White Paper: CBOE Volatility Index.

[25]   CBOE (2009) CBOE S&P 500 Implied Correlation Index.

[26]   Ait-Sahalia, Y. (2002) Maximum Likelihood Estimation of Discretely Sampled Diffusions: A Closed-Form Approximation Approach. Econometrica, 70, 223-262.

[27]   Ait-Sahalia, Y. and Kimmel, R. (2007) Maximum Likelihood Estimation of Stochastic Volatility Models. Journal of Financial Economics, 83, 413-452.

[28]   Bu, R., Cheng, J. and Hadri, K. (2017) Specification Analysis in Regime-Switching Continuous-Time Diffusion Models for Market Volatility. Studies in Nonlinear Dynamics & Econometrics, 21, 65-80.

[29]   Driessen, J., Maenhout, P. and Vilkov, G. (2005) Option-Implied Correlations and the Price of Correlation Risk. Advanced Risk & Portfolio Management Paper.

[30]   Osterrieder, J., Roschli, K. and Vetter, L. (2019) The VIX Volatility Index—A Very Thorough Look at It.

[31]   Li, A. (2017) Improvements to the CBOE VIX Calculation Formula.

[32]   Hu, Y., Lan, G. and Zhang, C. (2014) The Explicit Solution and Precise Distribution of CKLS Model under Girsanov Transform. International Journal of Statistics and Probability, 4, 68.

[33]   Ait-Sahalia, Y., Fan, J. and Li, Y. (2013) The Leverage Effect Puzzle: Disentangling Sources of Bias at High Frequency. Journal of Financial Economics, 109, 224-249.

[34]   CBOE (2012) The Relationship of the SPX and the VIX® Index.

[35]   Reghai, A. (2010) Breaking Correlation Breaks. Risk Magazine, October, 90-95.

[36]   Delanoe, P. (2014) Local Correlation with Local Vol and Stochastic Vol.

[37]   Estrada, J. (2008) Mean-Semivariance Optimization—A Heuristic Approach. Journal of Applied Finance, 18.

[38]   Chen, L., Peng J., Zhang, B. and Rosyida, I. (2017) Diversified Models for Portfolio Selection Based on Uncertain Semivariance. International Journal of Systems Science, 48, 637-648.