The equivalent constant-elasticity-of-variance (CEV) volatility of the stochastic-alpha-beta-rho (SABR) model

Introduction

The stochastic-alpha-beta-rho (SABR) model proposed by Hagan et al. 2002 is one of the most popular stochastic volatility models adopted in the financial industry. Its commercial success is owing to a few factors. The model is intuitive and parsimonious. It provides a flexible choice of backbone, the trace of the at-the-money volatility against the spot price. Most importantly, Hagan et al. 2002 provide an analytic approximation of implied Black–Scholes (BS) volatility in closed form (hereafter, the HKLW formula), from which traders can readily convert to the option price using the BS formula.

The HKLW formula is an asymptotic expansion valid for small time to maturities (up to the first order in time) and near-the-money strike prices. Several authors have attempted to improve the HKLW formula. Based on the results of Berestycki et al. 2004, Obłój 2007 corrects the leading order term of the HKLW formula. Henry-Labordère 2005 derives the same leading order term from the heat kernel under hyperbolic geometry. Paulot 2015 further provides a second-order approximation which is accurate in a wider region of strike prices at the cost of numerical integrations for the second-order term. Further, Lorig et al. 2017 obtain implied BS volatility up to the third order in time, which unfortunately is valid only near the money. A more accurate solution of the SABR model, however, requires large-scale numerical methods such as the finite difference method (Hagan et al. 2014; Park 2014; von Sydow et al. 2019), continuous time Markov chain (Cui et al. 2018), multidimensional numerical integration (Henry-Labordère 2005; Islah 2009; Antonov et al. 2013; Korn and Tang 2013), or Monte-Carlo simulation (Chen et al. 2012; Cai et al. 2017; Choi et al. 2019).

By nature, analytic approximation methods suffer two important drawbacks when some parameters or strike price go beyond the comfort zones of the concerned methods: non-negligible price error from the true value and the occurrence of arbitrage opportunity. Nevertheless, these methods are still attractive to practitioners because they are fast and robust. Note that practitioners need to compute the prices and Greeks of thousands of European options (or swaptions) frequently during trading hours. The calibration of the model parameters to the observed volatility smile also requires fast option evaluation because the parameters must be found using iterative methods. The numerical methods mentioned above are computationally intensive and not fast enough to use for those purposes.

Fortunately, the errors in analytic approximations are not a significant issue for those who use the SABR model primarily to price and manage the risk of European options. Specifically, the model parameters should first be calibrated to the market prices of the options at several liquid strike prices near the money. Then, the calibrated model is used to price the options at other strike prices. In this sense, the SABR model serves as a tool to interpolate (as well as extrapolate) the volatility smile, meaning that the accuracy of the price formula is less a concern.

Arbitrage under the analytic approximations occurs because an absorbing boundary condition is actually not imposed at the origin, as the small-time asymptotics of the transition density does not feel the boundary. The SABR process has a non-zero probability of hitting zero for $0<\beta<1$, and an absorbing boundary condition should be explicitly imposed at the origin for $0<\beta<1/2$ for the price process to be a martingale and arbitrage-free. For this reason, the analytic approximations exhibit arbitrage opportunity in the low-strike region. The arbitrage outbreak is still an important concern to options market makers; savvy hedge funds can exploit them by purchasing a butterfly of options with a negative premium. To avoid such trades, market makers carefully keep track of the lower bound of the arbitrage-free region (and often patch a different arbitrage-free model below the bound). Therefore, the degree of arbitrage should be yet another performance measure for testing newly proposed analytic approximations, as in Obłój 2007, as much as the approximation error.

We contribute to the SABR model literature by proposing new analytic approximations that are more accurate and have a wider arbitrage-free strike region than existing studies. We derive the equivalent volatility under the constant-elasticity-of-variance (CEV) model, from which the option price is computed with the analytic CEV option price formula (Schroder 1989). We provide two formulas for the equivalent CEV volatility as spin-offs from existing studies. The first one (Theorem 1) is obtained by following the approximation method of Hagan et al. 2002. The second one (Theorem 2) is simplified from Paulot 2015’s original CEV volatility formula.

Our CEV-based approach is motivated by the simple intuition that the SABR model should converge to the CEV model when the volatility of volatility (vol-of-vol) approaches zero. Such motivation for using the CEV model is not novel in the SABR model literature. Yang et al. 2017 show that the CEV option price (with the CEV volatility being the initial SABR volatility) is a good approximation in certain parameter ranges and is, naturally, arbitrage-free. The practical use of the result, however, is limited because the parameters related to the volatility process (i.e., vol-of-vol and correlation) are ignored in the approximation and only one degree of freedom is left to fit the volatility smile. Our work extends Yang et al. 2017 as our CEV approximations have full dependency on the SABR parameters. Paulot 2015 also discusses the equivalent CEV volatility as an alternative to BS volatility and outlines the derivation. His emphasis, however, is placed on BS volatility and the CEV volatility approximation is not tested numerically. Further, there is no discussion about the implications such as the mass at zero. In short, the advantage of the CEV approach has not thus far been explored. We fill this research gap by advocating the use of CEV volatility for the SABR model.

The numerical results show that our CEV-based approximations are more accurate than the corresponding BS-based methods from which they stem. In particular, the presented CEV approaches are more accurate when the initial volatility is large. This finding complements Paulot 2015’s refinement, which makes the approximation more accurate for large vol-of-vol. Having both advantages, the second CEV approximation based on Paulot 2015 performs the best among all the methods over wide parameter ranges. In the numerical test for comparing the degree of arbitrage, the second CEV approximation also performs the best among all the methods in that negative implied probability density starts to appear at the lowest strike price (see Section 4.3).

Surprisingly, the projection of the SABR model to the CEV model has the effect of imposing an absorbing boundary condition at the origin because the CEV price formula assumes the same boundary condition. Our CEV approximations offer finite CEV volatility at a zero strike, making them capable of implying the probability of hitting the origin. The mass-at-zero approximation in a closed-form formula (Theorem 3) shows excellent agreement with the numerical results in small time and is consistent with the exponentially vanishing asymptotics in the limit (Chen and Yang 2019). To the best of our knowledge, our approximation is the first closed-form approximation that works for all $\rho$, as the existing estimation methods either work on the uncorrelated case only (Gulisashvili et al. 2018) or require numerical integration for the correlated case (Yang and Wan 2018). Even if the mass at zero from our CEV approximations becomes less accurate in large time or vov-of-vol, our CEV approximations remain internally consistent with the model-free smile shape determined by the (possibly incorrect) mass at zero (De Marco et al. 2017). This explains why our CEV approach has a wider arbitrage-free region.

The remainder of this paper is organized as follows. Section 2 reviews the SABR model and existing BS volatility approximations. Section 3 derives the equivalent CEV volatility and mass-at-zero approximation. Section 4 presents the numerical results. Finally, Section 5 concludes.

SABR model and analytic approximations

In this section, we introduce the SABR model and review various analytic approximation methods for the equivalent BS volatility in the order of increasing accuracy. Rather than simply repeating existing results, we reorganize the formulas in an insightful way, which should lead to our new approach in Section 3.

The stochastic differential equations (SDE) for the SABR volatility model (Hagan et al. 2002) are

where $F_{t}$ and $\sigma_{t}$ are the processes for the forward price and volatility, respectively, $\nu$ is the vol-of-vol, $\beta$ is the elasticity parameter, and $W_{t}$ and $Z_{t}$ are the standard Brownian motions (BM) correlated by $\rho$. Let $T$ be the time-to-maturity of the option; then, the SABR model is fully specified by the parameter set: $\{F_{0},\sigma_{0},\beta,\rho,\nu,T\}$. To simplify the notations, we also denote

In order for the process to have a unique solution, an explicit boundary condition has to be specified for $0<\beta<1/2$ and it has to be an absorbing boundary condition for the price process to be a martingale and arbitrage-free. The absorbing boundary condition holds naturally for $1/2\leq\beta<1$. Then, the SABR model has a probability mass at the origin for $0<\beta<1$ similar to the CEV model. For the avoidance of doubt, no boundary condition is necessary for $\beta=0$ (the normal SABR model) because the price can freely go negative.

We first standardize the SDE to not only simplify the notations but also help with the intuition and numerical implementation. In particular, we standardize the price, strike, and volatility by their typical scales:

The SDE with standardized variables become

Here, $\alpha$ is the initial volatility of the standardized price $f_{t}$ around $f_{0}=1$. This should not be confused with the notation of other studies (Hagan et al. 2002; Obłój 2007; Paulot 2015), where $\alpha$ is used for the initial volatility, $\alpha=\sigma_{0}$. Indeed, the formula for $\alpha$ is a quick way of converting CEV volatility $\sigma_{0}$ to BS volatility $\alpha$ from $\alpha F_{0}=\sigma_{0}F_{0}^{\beta}$. Since $f_{0}=1$, $\alpha$ also serves as an approximation of normal and CEV volatility. As it turns out, $\alpha$ is the 0-th order term of the equivalent volatility of the SABR model in all the base models we consider: normal, BS, and CEV. In the rest of this paper, we use standardized variables, $k$ and $\alpha$, in the volatility approximations. In particular, the equivalent volatility (and its error) is presented as a ratio to $\alpha$. If $\sigma(k)$ is the equivalent volatility in the standardized scale (i.e., strike price $k$ and $f_{0}=1$), volatility in the original scale, $\sigma^{\prime}(K)$, can be converted using

where $\beta^{\prime}$ is the elasticity parameter of the base model of the equivalent volatility (e.g., $\beta^{\prime}=1$ for the BS model, $\beta^{\prime}=0$ for the normal model, and $0<\beta^{\prime}<1$ for the CEV model). For BS volatility, $\sigma(k)$ and $\sigma^{\prime}(K)$ are same.

We can further reduce the dimension of the parameters by introducing a time scale. There are two choices of time scale: $T$ or $1/\nu^{2}$. Accordingly, the parameter set is reduced to $\{\alpha\sqrt{T},\nu\sqrt{T},\beta,\rho\}$ and $\{\nu/\alpha,\beta,\rho,\nu^{2}T\}$, respectively. While the two parameter sets are equivalent, the first set better suggests the patterns of the asymptotic expansion of $T$; for example, we expect each order of $T$ in a small-time expansion to be accompanied by $\alpha^{2}$, $\alpha\nu$, and $\nu^{2}$.

2 HKLW formula of Hagan et al. 2002

In the original paper, Hagan et al. 2002 derive the equivalent BS volatility using the singular perturbation in the limit of a small time-to-maturity and a near-the-money strike price. They first derive the equivalent normal volatility of the SABR model as (Hagan et al. 2002, Eq. (A.59))

and the function $H(z)$ is defined in a chain as follows The function $x(z)$ can be equivalently defined as $x(z)=-\log\left(\frac{V(z)-z-\rho}{1-\rho}\right)$. While the expression in Eq. (4) is used by Paulot 2015, the alternative expression, in the form of $-x(-z)$, is used by Hagan et al. 2002; Hagan et al. 2014 and Obłój 2007.:

We use the dummy variable $z$ here because the functions will also be used with arguments other than $\zeta$. Regarding the evaluation of $H(z)$, there are two special cases to comment. At $z=0$, $H(z)$ should be evaluated as 1 from Taylor’s expansion of $1/H(z)$ near $z=0$:

To obtain the equivalent BS volatility, Hagan et al. 2002 also derive the equivalent normal volatility of the BS model with volatility $\sigma_{\textsc{bs}}$ as a special case of the SABR model with $\alpha=\sigma_{\textsc{bs}}$, $\beta=1$, and $\nu=0$ (Hagan et al. 2002, Eq. (A.63)):

where $\log k$ appears as the limit of $(k^{\beta_{\ast}}-1)/\beta_{\ast}$ as $\beta_{\ast}\rightarrow 0$. Then, Eqs. (3) and (6) are equated to solve for $\sigma_{\textsc{bs}}$ up to $O(T)$:

Here, $\sigma_{\textsc{bs}}$ in Eq. (7b) is replaced by the leading order approximation of Eq. (7a) near $k=1$,

The approximation comes from the expansion (Hagan et al. 2002, Eq. (A.68b)) near $k=1$,

Further applying this expansion to $\beta_{\ast}\log k/(k^{\beta_{\ast}}-1)$ and $\zeta$, we finally arrive at the well-known HKLW formula (Hagan et al. 2002, Eq. (2.17)):

3 Corrected leading order term of Obłój 2007

Let us also separately define $q$ (and $z$) for the two special cases of $\beta=0$ and $1$:

Therefore, $z\rightarrow 0$ in the zero vol-of-vol limit ($\nu\downarrow 0$) and $z=q=0$ at the money ($k=1$).

Based on the results of Berestycki et al. 2004, Obłój 2007 corrects the leading order term $H(\zeta)$ of Hagan et al. 2002 by replacing $\zeta$ with $z$. After the correction, the equivalent normal and BS volatilities, Eqs. (3) and (7), respectively become

where $h_{\textsc{n}}$ and $h_{\textsc{bs}}$ are unchanged as defined in Eqs. (3b) and (7b), respectively. Both $q_{\textsc{n}}/q$ and $q_{\textsc{bs}}/q$ are numerically evaluated as 1 at $k=1$.

We briefly explain the variable $z$ in this correction. With the scaled time, $s=t\,\nu^{2}$, the SABR SDE are equivalently stated as

where $\hat{W}_{s}$ is a standard BM rescaled from $W_{t}$ by $\hat{W}_{s}=(1/\nu)W_{s\nu^{2}}$ (same for $\hat{Z}_{s}$). The variable $z$ is the Lamperti transformation,

evaluated with $f_{t}=k$. In fact, both $\zeta$ and $\zeta^{\prime}$ are approximations to $z$ near $k=1$.

We deliberately denote the new state variable by the same $z$ as the dummy variable in Eq. (4) because $z$ is the correct variable for the functions in Eq. (4). In the rest of this paper, we often omit the argument $z$ from the functions for the sake of conciseness (e.g., $H=H(z)$, $x=x(z)$, and $V=V(z)$), unless otherwise stated.

4 Improved normal volatility approximation of Hagan et al. 2014

While the HKLW formula, Eq (9), is widely used as the final outcome of Hagan et al. 2002, the normal volatility approximation, Eq. (3), is also popular among practitioners in fixed income trading, for which the SABR model was originally proposed. In fixed income trading, normal volatility is preferred to BS volatility for quoting and managing the risk of swaptions. Moreover, the normal volatility formula is considered to be more accurate because it avoids an extra step for approximating BS volatility.

In a follow-up paper, Hagan et al. 2014 present an equivalent normal volatility improved over Eq. (3):

This new approximation not only adopts the correction of Obłój 2007 in the leading order term but also further refines the first-order term, $h_{\textsc{n}}$. Rather than deriving the equivalent BS volatility, Hagan et al. 2014 promote using the normal model (Bachelier 1900) with this normal volatility to obtain the option price. If needed, implied BS volatility can be quickly inverted from the price using an accurate approximation such as Jäckel 2015.

In particular, the normal SABR ($\beta=0$) has been a popular model choice to admit negative interest rates (Antonov et al. 2015) when the interest rate hovered near zero since the global financial crisis in 2008. From both Eqs. (3) and (14), the normal volatility approximation of the normal SABR model is derived as

Here, several observations should be mentioned. For the normal SABR model, the equivalent normal volatility is a natural choice. Using the equivalent BS volatility would be non-sensical because it restricts the price to be non-negative, contradicting the motivation for choosing $\beta=0$. The leading order term is simply $H(z_{\textsc{n}})$, meaning that the zero vol-of-vol limit is correct; $\sigma_{\textsc{n}}\rightarrow\alpha$ as $\nu\downarrow 0$. Moreover, Eq. (15) does not have a divergence issue at $k=0$, unlike the general case of Eq. (14b). These observations are generalized to $0<\beta<1$ in our CEV approach in Section 3.

5 Refined first-order term of Paulot 2015

Using the heat kernel expansion of hyperbolic geometry, Paulot 2015 derives an equivalent BS volatility up to $O(T^{2})$. While Henry-Labordère 2005 also uses the heat kernel expansion to derive the correct leading order volatility, Paulot 2015 derives the $O(T)$ and $O(T^{2})$ terms without approximating the dependency on the strike price in each time order, thereby making the equivalent volatility valid for a wider region of strike prices. Although the $O(T^{2})$ term improves accuracy, we adopt the approximation only up to the order $O(T)$ because the $O(T^{2})$ term involves a numerical integration which defeats the purpose as an analytic approximation.

Paulot 2015 expresses the small-time asymptotics of option’s time value under a stochastic volatility model in the general form:

where $\sigma$ is the initial volatility of the model at $t=0$, $d$ is the geodesic distance between the initial ($f_{0}=1$) and final (i.e., $f_{T}=k$) points on the differential geometry characterized by the model, and $e$ is also to be determined from the model. Here, the time value is also understood as the out-of-the-money put option price, hence the notation $P(k)$, because we are concerned with the low-strike region ($k<1$). For example, the coefficients for the BS model with volatility $\sigma_{\textsc{bs}}$ are given by (Paulot 2015, Eq. (21))

The expansion of the option value under the SABR model in the short-time limit is given by (Paulot 2015, Eq. (32))

Here, $t_{1}$, $t_{2}$, and $G(t)$ in $A_{2}$ are defined by

We have rearranged the original expressions in Paulot 2015 into $A_{2}$ and $A_{3}$ (and $A_{1}$ later) to handle the limit $k\rightarrow 1$ later. Moreover, we further simplify the expression of $A_{2}$, which was originally given in terms of several layers of definitions. A provides the details of the simplification. The rearrangement and simplification not only facilitate the numerical implementation but also help explain the formula, as we discuss below.

Next, the two option value expansions are equated to solve for the equivalent BS model. Assuming an expansion in $T$, $\sigma_{\textsc{bs}}=\sigma_{\textsc{bs},0}(1+h_{\textsc{bs}}T)$, we obtain $\sigma_{\textsc{bs},0}$ and $h_{\textsc{bs}}$ sequentially as follows:

Finally, we arrive at the equivalent BS volatility formula of Paulot 2015:

where $H(z)$ is defined in Eq. (4); $q$, $z$, and $q_{\textsc{bs}}$ are defined in Eqs. (10)–(11); and $A_{1}$, $A_{2}$, and $A_{3}$ are defined in Eqs. (22) and (19).

It is worthwhile checking the circumstances in which the refinement of Paulot 2015 makes a difference to the approximations of Hagan et al. 2002 with Obłój 2007’s correction. To begin with, the leading order term, $(q_{\textsc{bs}}/q)H(z)$, is the same as that of Obłój 2007. Although the expressions of $h_{\textsc{bs}}$ in Eqs. (7) and (23) look different, they have the same value at the money ($k=1$) This is stated by Paulot 2015 without proof.. Indeed, we have purposely introduced $A_{1}$, $A_{2}$, and $A_{3}$ to decompose $h_{\textsc{bs}}$ in such a way that the three terms correspond to those in Eq. (7b), respectively:

The limits of $A_{1}$ and $A_{3}$ can be easily derived from the expansions, Eqs. (8) and (5), respectively. The limit of $A_{2}$ requires additional algebra, which is placed in B. Therefore, Paulot 2015 is distinguished from Hagan et al. 2002 and Obłój 2007 only for out-of-the-money strike prices ($k\neq 1$). Since the dependency on $k$ is manifested through $z=(\nu/\alpha)\,q$, we also expect Paulot 2015’s refinement to be more pronounced when the $\nu/\alpha$ ratio is large.

6 Low-strike smile and mass at zero of the BS-based approximations

We comment on the low-strike behavior of the BS volatility approximations. In general, analytic approximation does not take into account the boundary condition at the origin because the small-time asymptotics of the transition density ignores the boundary condition. Therefore, asymptotic approximation should not be trusted near $k=0$.

As $k$ approaches zero, approximation quality deteriorates and eventually causes negative implied probability density, inducing arbitrage. The occurrence of arbitrage can also be seen through the fact that the equivalent BS volatilities we have reviewed do not conform to the model-free volatility bound of Lee 2004 for the small strike: $\sqrt{2|\log k|/T}$ as $k\downarrow 0$ for all $T$. The corrected leading order of the equivalent BS volatilities scale as $(q_{\textsc{bs}}/q)H(z)\sim O(|\log k|)$ as $k\downarrow 0$ for $0<\beta<1$, and this breaches the upper bound of Lee 2004.

where $N^{-1}(\cdot)$ is the inverse cumulative distribution function (CDF) of the normal distribution. The first term in Eq. (24) corresponds to Lee 2004’s upper bound.

Similar to the CEV model, the SABR model also has the probability mass at the origin for $0<\beta<1$. Therefore, the small-strike smile under the SABR model is subject to Eq. (24). Several researches estimate the mass at zero under the SABR model and, eventually, obtain the small-strike smile by taking advantage of Eq. (24). Gulisashvili et al. 2018 derive the approximations for $M_{T}$ for the uncorrelated case ($\rho=0$) in small- and large-time limits. Yang and Wan 2018 derive $M_{T}$ for the correlated case in the small-time limit. The practical use of the formulas, however, is limited because Eq. (24) is valid for a small strike (diverges at $k=1$) and it is non-trivial to merge the small-strike smile into the approximations reviewed earlier that are valid near the money.

CEV-based approximation

Since we advocate the use of implied CEV volatility, we briefly review the CEV model. The standardized CEV model with volatility $\sigma_{\textsc{cev}}$ is given by

The standardized prices of the call and put options with strike price $k$ and time-to-maturity $T$ are respectively (Schroder 1989)

where $F_{\chi^{2}}(x;r,x_{0})$ and $\bar{F}_{\chi^{2}}(x;r,x_{0})$ are respectively the CDF and complementary CDF of the non-central chi-squared distribution with degrees of freedom $r$ and non-centrality parameter $x_{0}$. The Greeks, namely, the sensitivity of the price with respect to the parameters, are also analytically available; see Larguinho et al. 2013. The functions related to the non-central chi-squared distribution are included in many standard numerical libraries. Several approximation methods (Sankaran 1963; Fraser et al. 1998) are also available to speed up the evaluation; see Larguinho et al. 2013 for further details.

Under the CEV model, the absorbing boundary condition has to be imposed at the origin in order to make the price process a martingale and arbitrage-free. The option formulas above are indeed derived with the absorbing boundary condition. Eqs. (25)–(26) imply $C_{\textsc{cev}}(0)=1$ and $P_{\textsc{cev}}(0)=0$, satisfying the put–call parity at $k=0$. The CDF of the price distribution is given by

In particular, the mass at zero for $0<\beta<1$ is analytically available as

where $\bar{\Gamma}(x;a)$ is the complementary CDF of the gamma distribution This is equivalent to the upper incomplete gamma function normalized by the gamma function $\Gamma(a)$. with the shape parameter $a>0$. The definition of $\bar{\Gamma}(x;a)$ and its asymptotic expansions for large $x$ (Abramowitz and Stegun 1972, (6.5.32)) are respectively given by

For later use, we summarize the small-strike asymptotics of the put option price under the CEV model. If a price distribution in general is defined to be non-negative, the CDF and put option price, $P(k)$, satisfy

because the term in the middle is the price of the put spread struck at $0$ and $k$ with $P(0)=0$. Therefore, if $M_{T}>0$,

In the CEV model context, for all $T$ and $0<\beta<1$, we have

It can also be directly shown from Eq. (26) using the fact that the second term vanishes faster than $O(k)$:

Finally, the small-strike and small-time asymptotics are

2 Observations and insights

Before we proceed to the explicit derivation, it is possible to postulate the form of the equivalent CEV volatility based on insights and observations. The equivalent BS volatility approximations in Section 2, after adopting the correction of Obłój 2007, can be cast into the following generic form:

In Eq. (7), for example, the breakdown of $h_{\textsc{bs}}$ is obvious. In Eq. (23), $A_{1}$, $A_{3}$, and $A_{2}$ are recognized as the terms corresponding to $O(\alpha^{2})$, $O(\alpha\nu)$, and $O(\nu^{2})$, respectively, based on their values in the limit $k\rightarrow 1$.

To understand the roles of the parts of Eq. (34) better, let us consider the limit $\nu\downarrow 0$. At this limit, the SABR model converges to the CEV model with $\sigma_{\textsc{cev}}=\alpha$ and, therefore, Eq. (34) plays the role of converting CEV volatility, $\sigma_{\textsc{cev}}=\alpha$, to BS volatility $\sigma_{\textsc{bs}}$:

Indeed, this form is close to the well-known local volatility conversion of Hagan and Woodward 1999. However, it is only an approximation; the converted $\sigma_{\textsc{bs}}$ does not exactly reproduce the CEV option value for volatility $\sigma_{\textsc{cev}}$. Since the next order term is $O(\alpha^{4})\,T^{2}$, approximation quality is expected to be poor for large initial volatility, $\alpha\sqrt{T}>1$.

To preserve the CEV model limit in the volatility approximation, we naturally consider the equivalent CEV volatility. As we recognize that $q_{\textsc{bs}}/q$ and $O(\alpha^{2})$ are purely involved in the conversion between the CEV and BS models, we do not expect these two terms to appear in the expression of CEV volatility. Dividing Eq. (34) by Eq. (35) on each side of the equation, we can factor out these terms and obtain the expected CEV volatility up to $O(T)$ in the form of

This form yields the desired limit, $\sigma_{\textsc{cev}}\rightarrow\alpha$ as $\nu\downarrow 0$, owing to the absence of the $O(\alpha^{2})$ term and the limit $H\rightarrow 1$ as $z\rightarrow 0$. Since the multiplication of Eqs. (35) and (36) up to $O(T)$ conversely yields Eq. (34), the BS volatility approximation is understood as a two-step conversion: (i) from the SABR model to CEV volatility and (ii) then to BS volatility. Therefore, the CEV volatility approximation in the form above is expected to be more accurate because the second approximation step is omitted. In particular, the advantage of the CEV volatility approach over the BS volatility approach becomes clear for $\alpha\sqrt{T}>1$ (when the quality of the second approximation is poor).

Indeed, we can already appreciate the postulated form for the two special cases of $\beta=0$ and 1. The equivalent normal volatility for the normal SABR model ($\beta=0$), Eq. (15), follows the form. When $\beta=1$, the equivalent BS volatilities, Eqs. (7), (9), and (23), all turn into the expected form since $q=q_{\textsc{bs}}$ and $\beta_{\ast}=0$. In the next two subsections, we generalize to $0<\beta<1$ by providing explicit derivations.

3 Equivalent CEV volatility based on Hagan’s approach

We present our first CEV volatility approximation. While we similarly follow the derivation of Hagan et al. 2002, we use the improved normal volatility, Eq. (14), instead of Eq. (3).

The equivalent CEV volatility of the SABR model up to $O(T)$ is given by

We derive the equivalent CEV volatility in a manner similar to obtaining the equivalent BS volatility in Section 2.2. We equate the normal volatilities for the SABR and CEV models from Eq. (14) as follows:

Then, we solve for $\sigma_{\textsc{cev}}$ up to $O(T)$ and obtain the result:

where we use $\sigma_{\textsc{cev}}\approx\alpha$ from $H=1$ when $k=1$. $\square$

This approximation indeed meets our expectation, Eq. (36). Consequently, the SABR model converges to the CEV model; $\sigma_{\textsc{cev}}\rightarrow\alpha$ at all $k$ as $\nu\downarrow 0$.

4 Equivalent CEV volatility based on Paulot’s approach

The equivalent CEV volatility is discussed in Paulot 2015 as an alternative to the equivalent BS volatility. For the CEV model with volatility $\sigma_{\textsc{cev}}$, the coefficients of the option value expansion, Eq. (16), are given by (Paulot 2015, p 52)

This is a special case of Eq. (18) with $\nu\downarrow 0$ and $\alpha=\sigma_{\textsc{cev}}$. In turn, Eq. (17) is a special case of Eq. (39) for $\beta\rightarrow 1$ ($\beta_{\ast}\rightarrow 0$). Based on this option value expansion, our second CEV volatility approximation is given below.

The equivalent CEV volatility from Paulot 2015 is simplified to

where $H(z)$ and $z$ are defined in Eqs. (4) and (10), respectively and $A_{2}$ and $A_{3}$ are given by Eq. (19).

Using the expansion in $T$, $\sigma_{\textsc{cev}}=\sigma_{\textsc{cev},0}(1+h_{\textsc{cev}}T)$, we sequentially obtain the two terms as

The second CEV approximation is compared with the first one in Theorem 1 in the same way that the BS volatility of Paulot 2015 is compared with that of Hagan et al. 2002. The two CEV volatility approximations have the same leading order term, $H$, and their first-order terms, $h_{\textsc{cev}}$, also converge to the same value at the money:

Similarly, they differ at out-of-the-money strike prices and this difference is expected to be pronounced if $\nu/\alpha$ is large. The second CEV volatility is more accurate and is valid for a wider region of strike prices than the first one because the first-order term, $h_{\textsc{cev}}$, is obtained without approximating the dependency on the strike price.

5 Probability mass at zero implied from the CEV volatility approximations

We examine the mass at zero and low-strike smile of the two CEV approximations in Theorems 1 and 2. We show that, unlike the BS-based approximations, the CEV-based approximations can properly imply the probability mass at the origin. We start with the observation that the equivalent CEV volatility has a finite value at $k=0$.

For $0<\beta<1$, the CEV approximations in Theorems 1 and 2 have finite CEV volatilities at a zero strike:

where the first-order term at $k=0$, $h_{\textsc{cev}}^{k=0}$, is finite. If $\rho=0$, the expression is further simplified to

Here, the CEV volatilities at $k=0$ are evaluated with $z|_{k=0}=-\xi$. In Theorem 1, $h_{\textsc{cev}}^{k=0}$ is trivially given by

For $0<\beta<1$, the mass at zero implied from the CEV approximations (Theorems 1 and 2) is given by

where $\bar{\Gamma}(x;a)$ is defined in Eq. (28).

Given that the equivalent CEV volatility is bounded near $k=0$ and the asymptotic behavior in Eq. (32), the small-strike put price from the CEV volatility approximation is similar to that of the CEV model Eq. (31):

The mass at zero follows from $M_{T}=\lim_{k\downarrow 0}P(k)/k$ in Eq. (30). $\square$

The numerical experiments in Section 4.2 will show that the above approximation to the mass at zero is accurate. Therefore, Theorem 3 serves as an alternative to existing methods. Unlike Gulisashvili et al. 2018, our method can be used for correlated cases. While Yang and Wan 2018 require a two-dimensional numerical integration for correlated cases, our method is in a closed-form.

Given that $0<M_{T}<1$ under the CEV approximation, we expect that the small-strike smile, when converted to BS volatility, is consistent with Eq. (24). Although the result of De Marco et al. 2017 is model independent, they use the CEV model as a benchmark for the numerical test because the CEV model is a rare case that renders the exact arbitrage-free option prices and mass at zero. As argued earlier, our CEV approximations are expected to have a lower boundary of arbitrage in the small-strike region.

The mass at zero in Theorem 3 is also consistent with the findings of Chen and Yang 2019. With the absorbing boundary condition at the origin, Chen and Yang 2019 proves that a positive constant $T_{0}$ (depending on the SABR parameters) exists such that

thereby characterizing the asymptotics of the not-feeling-boundary principle as $M_{T}=O(e^{-T_{0}/T})$ as $T\downarrow 0$. Theorem 3 not only is consistent with the asymptotics, but also provides a closed-form expression for $T_{0}$.

In the mass at zero in Theorem 3, the time scale of the exponential decay, $T_{0}$, is given by a closed form:

From Theorem 3 and the leading order asymptotics, $\bar{\Gamma}(x;a)\sim x^{a-1}e^{-x}/\Gamma(a)$,

The expression for $T_{0}$ naturally follows. Because $T_{0}$ involves only the leading order CEV volatility at $k=0$, it is independent of the choice of the CEV volatility, presented in either Theorem 1 or Theorem 2. $\square$

It should be noted that the results in Theorem 3 and Corollary 3.1 are implied from our CEV volatility approximations rather than explicitly derived from the SABR dynamics. As such, further mathematical (dis)proof is required to show that the results conform to the true behavior of the SABR model. Below we prove that is the case for $\rho=0$. The proof for the general case, however, is beyond the scope of this paper. The proof for $\rho=0$ is based on the literature on the exponential functional of BM. The exponential functional of BM, defined by

is a heavily studied topic in stochastic analysis. See Matsumoto and Yor 2005a; Matsumoto and Yor 2005b for an extensive review. Originally, it was inspired by the pricing of continuously monitored Asian options. In the context of the SABR model, however, the integrated variance is related to the functional by

When $\rho=0$, conditional on the path of $\hat{\sigma}_{t}$ over $0\leq t\leq T$, the forward price $f_{T}$ under the SABR model is distributed according to the CEV model with variance $\sigma_{\textsc{cev}}^{2}T$ replaced by the integrated variance. Therefore, the option price can be expressed as the expectation of the CEV option price over stochastic variance $(\alpha/\nu)^{2}A^{[-1/2]}_{\nu^{2}T}$ (Islah 2009). The mass at zero for $\rho=0$ is similarly expressed as (Gulisashvili et al. 2018)

When $\rho=0$ and $0<\beta<1$, the decay time scale $T_{0}$ from Corollary 3.1, simplified to

is consistent with the true SABR dynamics.

We begin the proof with the Laplace transform of $1/A_{t}$ from Matsumoto and Yor 2005a:

By applying Laplace’s method, we show that

because the minimum of $(\xi^{2}/2)e^{-x}+\cosh(x)$ occurs when $e^{x}=\sqrt{1+\xi^{2}}$. It can be also derived from an alternative expression of the Laplace transform (Antonov et al. 2019, Eq. (3.108)). With Girsanov’s theorem and the asymptotic expansion of $\bar{\Gamma}(x;a)$ in Eq. (29), $M_{T}$ can be expressed as

where $\mu=-1/2$, $a=1/(2\beta_{\ast})$, and $t=\nu^{2}T$. Since $A_{t}=O(t)$ as $t\downarrow 0$, the first term, $-\xi^{2}/(2A_{t})$, of the exponent asymptotically dominates the rest as $t\downarrow 0$. The small-time asymptotics of $M_{T}$ should have the same leading order term of Eq. (44). Therefore, the result follows as

Note that two special cases of the expectation in similar forms are exactly known for $t>0$ (Matsumoto and Yor 2005a, Corollary 4.6):

We comment on (in)consistencies between the approximation methods for the mass at zero. In the $\nu\downarrow 0$ limit, the mass at zero from the three methods, Yang and Wan 2018, Gulisashvili et al. 2018, and Theorem 3, all converges to that of the CEV model with $\sigma_{\textsc{cev}}=\alpha$, which is consisent with the intution. In Theorem 3, this is the case because $\sigma_{\textsc{cev}}^{k=0}\rightarrow\alpha$ as $\nu\downarrow 0$. We note, however, the difference in decay time scales (i.e., the $T\downarrow 0$ limit) between Theorem 4 and Yang and Wan 2018. When $\rho=0$, the mass at zero from Yang and Wan 2018 reduces to that of the CEV model even though $\nu$ is not small. That implies $-T\log M_{T}\rightarrow 1/(2\beta_{\ast}^{2}\alpha^{2})$ as $T\downarrow 0$. However, this is significantly different from $T_{0}$ in Theorem 4 if $\xi=\nu/(\beta_{\ast}\alpha)$ is in the order of 1. The numerical experiment in Section 4.2 (Figure 2) is in favor of our asymptotic time scale.

The success of our CEV approximations in estimating the mass at zero is surprising considering that the boundary condition has not been explicitly considered in the derivation of the CEV approximations. Below are our insights on how this happens. We argue that the CEV-based approach is effective because it works as a control variate method. Recall that the CEV volatility based on Paulot 2015 in Theorem 2 is derived by matching the small-time option values between the CEV and SABR models:

where the left-hand side from Eqs. (16) and (39) is for the CEV model and the right-hand side from Eqs. (16) and (18) is for the SABR model. The asymptotics on each side are originally intended to hold in the small-time limit near $k=1$ in general. As such, they are not accurate near $k=0$. The correct CEV price asymptotics in Eq. (33) are different from the CEV value on the left-hand side in terms of the powers of $k$ and $T$ in the prefactor. However, the incorrect term $T^{3/2}k^{\beta/2}$ also arises on the right-hand side for the SABR model and they just cancel out. The choice of the CEV model also makes the exponents of the two sides closer (e.g., $q^{2}$ versus $(q/H)^{2}$). As a result, the CEV approach prevents the equivalent volatility from diverging to infinity as $k\downarrow 0$. Finally, the SABR price asymptotics, Eq. (41), obtained via the equivalent CEV volatility are similar to those of the CEV model, Eq. (31), which correctly carries the absorbing boundary condition. Notably, the output, Eq. (41), is correct even though the input, Eq. (45), is incorrect because the CEV approach works as a control variate method.

Numerical results

In this section, we numerically test the two CEV-based approximations in comparison to existing BS-based approximations. The Python code used in this study can be found at https://github.com/PyFE/PyfengForPapers. For easier reference, the methods are labeled as follows:

BS-A: Eq. (9), the HKLW formula of Hagan et al. 2002.

BS-B: Eq. (23), Paulot 2015’s equivalent BS volatility of order $O(T)$.

BS-C: Lorig et al. 2017’s equivalent BS volatility up to $O(T^{3})$

DMHJ: Eq. (24), the small-strike BS volatility smile determined by the mass at zero (De Marco et al. 2017)

CEV-A: Eq. (37) in Theorem 1, the equivalent CEV volatility based on Hagan et al. 2002.

CEV-B: Eq. (40) in Theorem 2, the equivalent CEV volatility based on Paulot 2015.

The three BS-based approximations are good representatives of existing methods. Let us make several comments on the selection of the methods. In BS-A, we do not adopt Obłój 2007’s leading order correction since the HKLW formula is widely used already. In BS-B, we do not use the $O(T^{2})$ term because it is not in a closed form. The naive CEV approximation of Yang et al. 2017 is not included because it has no dependency on $\nu$ and $\rho$. Similarly, the price approximation of Jordan and Tier 2011 is not included as it works only for $\rho=0$. Although it is not reviewed in Section 2, we include the $O(T^{3})$ approximation of Lorig et al. 2017, which is labeled BS-C. In Lorig et al. 2017, the higher-order result is achieved at the expense of the lower accuracy for off-the-money strike; the approximation is valid only for near-the-money strike prices because the dependency on $k$ is expressed in terms of the expansions of $\log k$ in all orders of time. This is in contrast to BS-A and BS-B, in which the accuracy over all $k$ is maintained to a certain degree through $H(z)$ in the leading order term The leading and first-order terms of BS-C coincide with those of BS-A and BS-B at $k=1$.. For this reason, the approximation quality of BS-C may deteriorate faster than that of BS-A or BS-B as $k$ moves away from one. Therefore, it is of additional interest to numerically compare the three BS-based methods. The DMHJ asymptotics is not self-contained, as it needs the externally estimated mass at zero and are only valid for $k<1$. However, they provide a useful reference as a model-free volatility smile at small strikes. We thus evaluate DMHJ with the two values of the mass at zero; one estimated from CEV-A or CEV-B and the other the true value.

We also compare the following methods to estimate the mass at zero:

GHJ-LN: The moment-matched lognormal approximation of Gulisashvili et al. 2018 for $\rho=0$ in small-time limit.

MC: Eq. (43) with Monte-Carlo simulated $A^{[-1/2]}_{\nu^{2}T}$ for $\rho=0$.

CEV-A: Theorem 3 with the zero strike CEV volatility from Theorem 1.

CEV-B: Theorem 3 with the zero strike CEV volatility from Theorem 2.

In GHJ-LN, the integrated variance is approximated by a lognormal distribution by matching to the first and second moments. Thus, the expectation in Eq. (43) is computed as a one-dimensional integration. While the reference does not specify the numerical method, we use the Gauss–Hermite quadrature for efficient integration with fast convergence (Choi and Wu 2021). In YW, the mass at zero is expressed as an expansion consist of the $O(1)$ and $O(\nu\sqrt{T})$ terms. As mentioned earlier, when $\rho=0$, the first order term vanishes and the leading order term is reduced to that of the CEV model. For accuracy benchmarking, we also implement the Monte-Carlo (MC) method. We draw the random values of $A^{[-1/2]}_{\nu^{2}T}$ by simulating $\hat{\sigma}_{t}$ on a discretized time grid and integrating $\hat{\sigma}_{t}^{2}$ by Simpson’s rule for quadrature. We generate 40,000 paths with 20 time steps between $t=0$ and $T$.

Table 1 shows the parameter sets used for numerical test. Sets 1–3 are for the volatility approximation (Section 4.1) and the arbitrage boundary (Section 4.3), while Sets 1, 4, and 5 are for the mass at zero (Section 4.2).

First, we test the accuracy of the analytic approximations using Sets 1–3. The three parameter sets are carefully chosen for comparing the relative strength of various methods. Those sets have been used in past studies (von Sydow et al. 2019; Cai et al. 2017; Antonov et al. 2019) and the exact option prices are available through numerical means (e.g., the finite difference method and the Monte-Carlo simulation). These exact values are additionally verified within reasonable accuracy using the continuous-time Markov chain codes https://github.com/jkirkby3/PROJ_Option_Pricing_Matlab adopted by Cui et al. 2018; Cui et al. 2019 The result for Set 2 is also reported in Cui et al. 2018 and Cui et al. 2019.. The table also displays $\alpha\sqrt{T}$ and $\nu\sqrt{T}$ as references to the performance of each method. The accuracy of the approximations tends to deteriorate as the two variables become larger. The ratio $\nu/\alpha$ also indicates how much Paulot’s refinement is noticeable.

Tables 2–4 respectively compare the approximation accuracy for the three parameter sets. Each table shows the standardized BS volatility error, $(\sigma_{\textsc{bs}}-\sigma_{\textsc{bs}}^{\text{exact}})/\alpha$, where $\sigma_{\textsc{bs}}$ is the implied BS volatility from the analytic approximation methods In the BS-based methods, $\sigma_{\textsc{bs}}$ is directly obtained from the volatility formulas, whereas in the CEV-based methods, $\sigma_{\textsc{bs}}$ is converted from the CEV option prices in Eqs. (25)–(26). and $\sigma_{\textsc{bs}}^{\text{exact}}$ is the exact BS volatility reported in the literature. The $\sigma_{\textsc{bs}}^{\text{exact}}$ value and corresponding call option price (in the original scale) are also provided in the tables for reference.

Table 2 shows that all the methods are accurate for Set 1 with an error below 0.02. This is expected from the values of $\nu\sqrt{T}$ and $\alpha\sqrt{T}$ not exceeding one. Among the methods, BS-C is the most accurate because it contains higher-order terms and the tested strikes are near the money. The two CEV-based methods are more accurate than the first two BS-based methods by a small margin.

Set 2 is often used in the literature (Cai et al. 2017; Cui et al. 2018) to reveal the failure of the HKLW formula. Table 3 shows that the CEV-based approximations are superior to their BS-based counterparts because of the large value of $\alpha\sqrt{T}=3.257$. Owing to the absence of the $O(\alpha^{2})\,T$ term, the CEV-based approximations have a much smaller error when $\alpha\sqrt{T}>1$, whereas both BS-A and BS-B largely over-predict volatility. The difference between BS-A and CEV-A (and between BS-B and CEV-B) is believed to be as the error generated from converting CEV volatility to BS volatility (see Section 3.2), and such error can be avoided by adopting CEV volatility instead. In addition, the numerical results are close between BS-A and BS-B and between CEV-A and CEV-B because the small $\nu/\alpha$ ratio makes $z$ small for fixed $k$ and, therefore, makes the type–B methods equipped with Paulot 2015’s refinement indistinguishable from the type–A methods. BS-C shows an interesting result. While it is the most accurate at the money ($k=1$) as expected, the error rises in the out-of-the-money region. Notably, the volatility skew (i.e., the slope of $\sigma_{\textsc{bs}}$) by BS-C is in the opposite direction to those from the other methods as well as the exact skew.

Set 3 is perhaps the most challenging case because of the large values of $\nu\sqrt{T}$ and $\alpha\sqrt{T}$, which are both higher than one. Table 4 shows the results. Although the errors are large overall as expected, CEV-B shows higher accuracy over the other methods except BS-C. In general, the B-type methods perform better than the A-type methods. This grouping of the results is different from the two previous test sets because $z$ is amplified by the large $\nu/\alpha$ ratio and Paulot 2015’s refinement dominates the CEV volatility effect. BS-C is the most accurate among the tested methods, although not by a significant amount.

2 Mass at zero and its small-time asymptotics

To demonstrate the accuracy of the mass-at-zero approximation in Theorem 3, we test with Set 4. Chen and Yang 2019 report the mass-at-zero values computed with the finite difference method for Set 4 and several modifications from the base values. Figure 1 shows the comparison. To efficiently compare the small-time asymptotics of $M_{T}$, we plot the ratio of decay time scale, $-T\log M_{T}/T_{0}$, as a function of $T$, where $T_{0}$ is from Corollary 3.1. If $M_{T}$ follows Corollary 3.1, $-T\log M_{T}/T_{0}\rightarrow 1$ as $T\downarrow 0$. Therefore, Corollary 3.1 is verified from the plot. CEV-A and CEV-B from Theorem 3 show good agreement with the exact values overall even though it is a correlated case ($\rho=-0.5$). In particular, they are accurate for small $T$, supporting Corollary 3.1. Among the parameter variations, the increased vol-of-vol in subplot (e) shows the largest error for the same level of $T$. This is consistent with fact that the analytic approximation deteriorates as $\nu\sqrt{T}$ rises. Unexpectedly, CEV-A shows slightly higher accuracy than CEV-B, which may just be a coincidence for this parameter set.

Table 5 shows the result for Set 5 and the variations with respect to $\beta$ and $\rho$. The parameter has been tested in Yang and Wan 2018 and the values for finite difference (all $\rho$) and YW ($\rho\neq 0$) are from the reference. The finite difference for $\rho\neq 0$ and MC for $\rho=0$ are considered the most accurate method. Although GHJ-LN is only available for $\rho=0$, it again shows the best agreement with MC. Overall, both CEV-B The CEV-A results are very close to those of CEV-B. and YW agree with the benchmark values for all cases.

3 Small-strike smile and arbitrage boundary

Using Sets 1–4, we examine the implication of the mass at zero embedded in the CEV approximations on the volatility smile, and also investigate the arbitrage boundary where negative implied density starts to occur. Figure 3 plots the small-strike BS volatility smile for the four parameter sets. We compare BS-B, CEV-B, and DMHJ along with the exact volatilities from Tables 2–4. The DMHJ asymptotics uses the mass at zero from CEV-B volatility at the origin. In all the sets, BS-B diverges to infinity ($\sim O(|\log k|)$) faster than DMHJ ($\sim O(\sqrt{|\log k|})$) as $k\downarrow 0$, indicating arbitrage. By contrast, CEV-B converges to DMHJ, indicating a lower degree of arbitrage. In near-the-money strike region, CEV-B merges to BS-B instead, while DMHJ diverges at $k=1$. Therefore, CEV-B bridges DMHJ in a low strike and BS-B near the money.

Although the mass at zero implied from Theorem 3 may not necessarily be accurate, CEV-B is still consistent with DMHJ. To demonstrate this point, we also plot the DMHJ asymptotic smile with a more accurately estimated $M_{T}$ (DMHJ-2). For Sets 1 and 2, $M_{T}$ is estimated from MC. For Set 4, $M_{T}$ is obtained from Chen and Yang 2019. For Set 3, we use an upper bound of $M_{T}$ obtained from Eq. (30) and the option value at $k=0.1$. We know that this upper bound is tight because the DMHJ-2 smile passes the exact volatilities at $k=0.1$ and 0.4. Because the mass-at-zero obtained from CEV-B may deviate significantly from the true value, DMHJ and DMHJ-2 can be quite different as shown in Figure 3(c) for Set 3. Nevertheless, the CEV-B smile is consistent with that of DMHJ up to the possibly incorrect mass at zero. Therefore, we expect CEV-B to exhibit less arbitrage than BS-B.

To elaborate this point further, we examine the occurrence of arbitrage under the analytic approximations by explicitly detecting the negative probability density function (PDF). The PDF at $k$ is implied from the volatility approximations by the second-order difference equation with small $h$:

where $C(k)$ is the price of the option struck at $k$, computed using the equivalent volatility at $k$. We compare the degree of arbitrage by locating the arbitrage boundary. Given that the analytic approximations exhibit arbitrage at a low strike, the location of the arbitrage boundary is defined as the first $k$ value having a negative implied PDF when $k$ decreases from one to zero. Since the difference equation can also be understood as the premium of the option butterfly, options market makers can trade options above the boundary without being arbitraged against if they use the analytic approximations of the SABR model. Therefore, the lower the arbitrage boundary is, the better the analytic approximation is.

Figure 4 depicts the PDF for Set 3 implied from the five methods as well as the exact prices. The PDFs from the approximation methods deviate from the exact value and eventually fall below zero as $k$ approaches zero, indicating arbitrage opportunities at low strikes. The PDFs from BS-B and CEV-B, however, have the lowest arbitrage boundary ($k=0.19$). The PDF deviation of BS-C is the most severe; the PDF falls below zero at the highest strike price ($k=0.30$) and then behaves erratically.

We further compare the arbitrage location in wide parameter ranges. In this experiment, we vary $\nu$, $\alpha$, and $\beta$ from the parameter values of Set 3. The upper panel of Figure 5 shows the result for varying $\nu$. While the location of the arbitrage boundary generally increases as $\nu$ increases, the arbitrage location is much lower for BS-B and CEV-B than for BS-A and CEV-A. Again, this is because Paulot 2015’s refinement becomes pronounced as the $\nu/\alpha$ ratio increases. In the middle panel where $\alpha$ is varied, the grouping pattern is different; the two CEV-based methods exhibit lower arbitrage boundary than the first two BS-based methods do because the CEV effect dominates Paulot 2015’s refinement as $\nu/\alpha$ decreases. The lower panel shows the result for varying $\beta$. Across all the methods, the arbitrage boundary is higher for $\beta$ closer to zero, reaffirming that the arbitrage is related to the absorbing boundary condition. BS-B and CEV-B have an arbitrage boundary consistently lower by about 0.1. Overall, CEV-B has the lowest arbitrage boundary among all the methods. CEV-B performs the best in all three parameter variations because the method is equipped with the two enhancements that complement each other: Paulot 2015’s refinement and the CEV-based approximation. The arbitrage boundary of BS-C behaves erratically, probably due to the side effect of the higher-order terms in $k$.

Conclusion

The SABR model is the dominant stochastic volatility model in the financial industry. Since finding an accurate solution is computationally burdensome, robust analytic approximations of the equivalent volatility are favored in practice despite their imperfection. We show that the quality of the approximation can be significantly improved by deriving the equivalent CEV volatility instead of the BS volatility. Projecting the SABR model on the CEV model takes advantage of an absorbing boundary condition at zero and thus exhibits less arbitrage for small strikes than projecting on the BS model. With this article, we hope to see that the equivalent CEV volatility receives more attention in both research and application.

Acknowledgements

Jaehyuk Choi was supported by the 2019 Bridge Trust Asset Management Research Fund. Lixin Wu was supported by Grant #16306717 of the Research Grants Council of Hong Kong. The authors thank Nan Chen and Nian Yang for providing the data points of Chen and Yang 2019.

Appendix A Simplification of Paulot’s formula

We simplify $t_{1}$, $t_{2}$, and $G(t)$ from the original definition of Paulot 2015. The intermediate variables in Paulot 2015 used to define $t_{1}$, $t_{2}$, and $G(t)$ are first simplified to

Using the expressions above, we simplify $t_{1}$ and $t_{2}$, respectively, to

The function $G(t)$ is also simplified to the form presented in Eq. (21), using

Appendix B Convergence of Paulot’s approximation

First, we handle the special case of $\beta=1$:

For all three cases of (21), the first and second derivatives are the same as

Evaluated at $t=t_{2}=(1+\rho)/\rho_{\ast}$,

The first derivative of $G(t)$ at $t=t_{2}$ is 0, as

Putting all terms together, we have, as $k\rightarrow 1$ ($z\rightarrow 0$),

References