1 Introduction↩︎

Option pricing is one of the main tasks of a derivatives desk. Not only to obtain the derivatives prices and trade them but also to calibrate more complex models or manage all the associated risks. In practice, most institutions consider stochastic volatility models to better capture market dynamics. However, including this feature adds more complexity to the calculations, making the processes more complicated and requiring more time to execute.

Over the last decades, a research line has consisted of finding different methods to approximate the price accurately. On the one hand, there are approximations focused on finding an expansion of the implied volatility. For example, using perturbation methods as in [@HaganWoodward99], [@Hagan02] and [@Jacquier14], using heat kernel expansions as in [@Gatheral12] or [@Labordere08], or asymptotic expansion method as in [@KaramiShiraya]. On the other hand, others focused on prices. For example, the Alòs decomposition formula in [@Alos12] expanded later in [@Alos20] and [@GLMV19a], adjoint expansions as in [@Pagliarani], the cos methods as in [@Fang08] or Watanabe’s expansion in [@Kunitomo03], [@Takahashi] and [@Nagami21].

Despite the vast literature on option approximations, there are few written references on quadratic payoffs. These products are not commonly traded but are useful to calculate the convexity’s correction for CMS. The main references are [@Hagan19] and [@HaganWoodward20]. The value of a quadratic call, put, and swaps are given by \[\begin{align} \label{quadratic95payoff95definition} V^{QC}(T, F, K) :&=& N_{T_p} \mathbb{E}^{\mathbb{Q}^{N}}\left[\left(F_{T}- K\right)^{2}_{+}\right],\\ V^{QP}(T, F, K) :&=& N_{T_p} \mathbb{E}^{\mathbb{Q}^{N}}\left[\left(K-F_{T}\right)^{2}_{+}\right],\\ V^{QS}(T, F, K) :&=& N_{T_p} \mathbb{E}^{\mathbb{Q}^{N}}\left[\left(F_{T}- K\right)^{2}\right]. \end{align}\tag{1}\] where \(T\) is the expiration date, \(T_p\) the payment date and \(F_T\) is the forward rate process associated to the derivative. On the other hand, \(\mathbb{E}^{\mathbb{Q}^{N}}\left[\cdot\right]\) is the expected value associated to the numeraire \(N_t\). Although there is a closed formula for normal and log-normal when volatility is constant, it is not good enough since these products depend more on market skews and smiles than vanilla options. Many banks use replication to price these structures. Note that \[\begin{align} \left(F_{T}- K\right)^{2}_{+} = 2 \int^{\infty}_{K} \left(F_{T}- K\right)_{+} dK', \end{align}\] therefore, it is possible to approximate the value as the sum of options with different strikes. That is the main advantage of this approach, as it is consistent with vanilla option pricing, easy to implement, and automatically satisfies the put-call parity. Unfortunately, from a practical standpoint, it is necessary a dense option grid with different strikes. Only a few options are traded per expiration and strike, so it is essential to calculate the grid, making it computationally intensive. In addition, the option prices with higher strikes, which are required for the calculation, are illiquid and have a ‘fictitious’ market price. To overcome this, an alternative methodology based on finding an approximation formula under the SABR model has been proposed by [@Hagan19] and [@HaganWoodward20].

The paper aims to obtain an alternative approach for pricing CMS derivatives. Firstly, we will establish a connection between the price of a CMS derivative and a quadratic payoff using a bit of Mallivin calculus. This relationship is obtained without relying on any specific model. Then, we will extend the results of Watanabe’s expansions to the case of quadratic payoffs over local and stochastic local volatility. The approximations are generic. To evaluate their accuracy, we will compare the approximations numerically under the normal SABR model against the market standards: Hagan’s approximation, and a Monte Carlo simulation.

The structure of the paper is as follows. In Section 2, we will introduce the main concepts of the CMS market to make the paper self-contained. In Section 3, we explain how to derive a general expression for the convexity adjustment for CMS caps, CMS floors, and CMS swaps using Malliavin calculus. This expression will depend on the different quadratic payoffs explained above. In Section 4, we provide a brief explanation of the Watanabe calculus. We will also give an example of the normal SABR model contrasting the accuracy of the approximation. In Section 5, we prove how to use Watanabe calculus to obtain the price of a quadratic call, quadratic put, and quadratic swaps under different models. The approximation is applied to local and stochastic local volatility models. In particular, the accuracy is tested for the normal SABR model by comparing the results against Monte Carlo and Hagan’s approximation. Finally, in Section 7, we will present our conclusions.

We provide the code for running all the examples on the papper at https://github.com/Dagalon/PyStochasticVolatility.

2 A brief introduction to the CMS market↩︎

One of the principal uses of quadratic payoffs is the convexity adjustment calculation, especially for CMS products. Therefore, to make the article self-contained, in this section, we will do a brief introduction to the CMS market.

We will consider the market standard of two curves: one to discount the cash flows and another to estimate the floating coupon. These two curves are equivalent to the discount and repo curves in the equity market.
The following notation will be used throughout the paper:

• We will denote the discount curve by \(d\) and the estimated curve by \(e\).

• \(P^{i}(T)\) is the spot discount factor to T using the curve i=d, e.

• \(P^{i}(t,T)= \frac{P^{i}(T)}{P^{i}(t)}\) is the forward discount factor associated to the curve i=d, e from t to T.

• \(\Delta(T_a,T_b)\) is the fraction’s year between \(T_a\) and \(T_b\).

• \(L^{i}(t,T_a,T_b) = \frac{1}{\Delta(T_a,T_b)}\Big(\frac{P^{i}(t,T_a)}{P^{i}(t,T_b)} - 1\Big)\) is the floating rate associated with the curve i=d, e.

The swap is the most commonly traded derivative in the interest rate market. While there are various swaps types, the market standard is to exchange a fixed coupon for a floating one until the end of the contract. Fixed and floating rate payments can be scheduled on different dates and periods. But, to simplify the notation, we assume that payments are made on the same dates and in the same fractions of the year. Today is denoted as \(T_{0}\), and the future payment dates are scheduled by \[\begin{align} \mathcal{T}=\left\{T_{a}, T_{a+1}, \ldots, T_{b}\right\}. \end{align}\]

The standard swap price is \[\begin{align} \label{swap95payoff} V^{Swap}(T_0,T_{a}, T_{b}, C) &= w * \Biggl(\sum_{i=a}^{b} \Delta(T_{i-1},T_{i})P^{d}(T_0,T_{i}) C \nonumber \\ & \qquad- \sum_{i=a}^{b} \Delta(T_{i-1},T_{i})P^{d}(T_0,T_{i}) \mathbb{E}_{T_0}^{T_i}\bigg[L^{e}(T_{i-1},T_{i-1},T_{i})\bigg]\Biggl). \end{align}\tag{2}\]

When we receive a fixed coupon, \(w=1\), we are in a receiver swap. On the other hand, if we pay the fixed coupon, \(w=-1\), we are in a payer swap. The measure \(\mathbb{Q}^{T_i}\) used on the floating leg is associate to the numeraire \(N^{T_p}_t:= P^{d}(t, T_p)\).

There are two key concepts when working with swaps: the swap rate and the annuity factor. The swap rate is the fixed coupon that makes the structure fair. It is defined by \[\begin{align} S_{a,b}(T_0) &= \frac{\sum_{i=a}^{b} \Delta(T_{i-1},T_{i})P^{d}(T_0,T_{i}) \mathbb{E}_{T_0}^{T_i}\bigg[L^{e}(T_{i-1},T_{i-1},T_{i})\bigg]}{\sum_{i=a}^{b} \Delta(T_{i-1},T_{i})P^{d}(T_0,T_{i})}. \end{align}\] The annuity factor is defined by \[\begin{align} 01(T_0, T_a, T_b) &= \sum_{i=a}^{b} \Delta(T_{i-1},T_{i})P^{d}(T_0,T_{i}). \end{align}\]

A constant maturity swap (CMS) is a fixed-for-floating swap, where the floating leg pays, or receives, periodically a swap rate with a fixed time to maturity. It has the same structure as the standard swap, but changing the floating payment from \(\mathbb{E}^{T_i}\big[L^{e}(T_{i-1},T_{i-1},T_{i})\big]\) to \(\mathbb{E}_{T_0}^{T_p}\big[S_{i-1,i-1+b}(T_{i-1})\big]\) with \(T_{i-1} \leq T_{p_i} < T_i\).
The CMS swap price is given by \[\begin{align} \label{cms95swap95payoff} V^{Swap}_{CMS}(T_0,T_{a}, T_{b}, C) &= w * \Biggl(\sum_{i=a}^{b} \Delta(T_{i-1},T_{i})P^{d}(T_0,T_{i}) C \nonumber \\ &\qquad - \sum_{i=a}^{b} \Delta(T_{i-1},T_{i})P^{d}(T_0,T_{p_i}) \mathbb{E}_{T_0}^{T_{p_i}}\bigg[S_{i-1,i-1+b}(T_{i-1})\bigg] \Biggl). \end{align}\tag{3}\]

There is a significant difference between valuing a standard swap and a CMS swap. The floating coupon of the standard swap is a martingale under the forward measure \(N^{T_p}_t\). In other words, \[\mathbb{E}_{T_0}^{T_i}\bigg[L^{e}(T_{i-1},T_{i-1},T_{i})\bigg] = L^{e}(T_0,T_{i-1},T_{i}).\] But the swap rate is not, it implies that \[\mathbb{E}_{T_0}^{T_p}\bigg[S_{i-1,i-1+b}(T_{i-1})\bigg] \neq S_{i-1,i-1+b}(T_0).\] The swap rate is a martingale respect to the measure induced by the numeraire \(N^{a,b}_t:=01(t, T_a, T_b)\). Since the floating rate and the swap rate are martingales with respect to different numeraries, a convexity adjustment is necessary for CMS products.

In fact, we have that \[\mathbb{E}_{T_0}^{T_{p_i}}\bigg[S_{i-1,i-1+b}(T_{i-1})\bigg] = \frac{N^{i-1,i-1+b}_{T_{0}}}{P^{d}(T_{0}, T_{p_i})} \mathbb{E}_{T_0}^{T_{p_i}}\bigg[S_{i-1,i-1+b}(T_{i-1}) \frac{N^{i-1,i-1+b}_{T_{i-1}}}{P^{d}(T_{i-1}, T_{p_i})} \bigg].\]

The essence of CMS modeling lies on the definition of a mapping function \(M(t,a,b)\). In [@AndreasenPiterbargIII], a comprehensive discussion is presented regarding the selection of the mapping function and the requirements it must satisfy. We use the following mapping function: \[M(t,a,b,p)=\frac{P^{d}(t,T_p)}{N_t^{a,b}}.\] Then, we have that \[\begin{align} \label{convexity95adjustment95cms95swap} \mathbb{E}_{T_0}^{T_{p_i}}\bigg[S_{i-1,i-1+b}(T_{i-1})\bigg] &= \mathbb{E}_{T_0}^{N^{i-1,i-1+b}}\bigg[S_{i-1,i-1+b}(T_{i-1}) \frac{M(T_{i-1},i-1,i-1+b, p_i)}{M(T_{0},i-1,i-1+b, p_i)} \bigg] \nonumber \\ &= S_{i-1,i-1+b}(T_{0}) \nonumber \\ &+ \mathbb{E}_{T_0}^{N^{i-1,i-1+b}}\bigg[S_{i-1,i-1+b}(T_{i-1}) \left(\frac{M(T_{i-1},i-1,i-1+b, p_i)}{M(T_{0},i-1,i-1+b, p_i)} - 1 \right) \bigg]. \end{align}\tag{4}\] The second term in the summand of the right side is known as convexity adjustment for the CMS rate, see [@Hagan03] for more details.
Concluding this section, we introduce CMS Caps and Floors, the most common CMS derivative products. The paper aims to develop different price expansions under various models for calculating their convexity adjustments. A CMS Cap and CMS Floor consist of a series of calls and puts, named caplets or floorlets, with the underlying \(S_{a,b}(t)\). The payoff is given by \[\begin{align} \label{cms95cap95floor95payoff} V^{Cap}_{CMS}(T_0,T_a,T_b, K) &= \left(\sum_{i=a}^{b} \Delta(T_{i-1},T_{i})P^{d}(T_0,T_{p_i}) \mathbb{E}_{T_0}^{T_{p_i}}\bigg[S_{i-1,i-1+b}(T_{i-1})-K)^{+}\bigg] \right), \nonumber \\ V^{Floor}_{CMS}(T_0,T_a,T_b, K) &= \left(\sum_{i=a}^{b} \Delta(T_{i-1},T_{i})P^{d}(T_0,T_{p_i}) \mathbb{E}_{T_0}^{T_{p_i}}\bigg[K-S_{i-1,i-1+b}(T_{i-1}))^{+}\bigg] \right). \nonumber \end{align}\tag{5}\] Now, we will focus on calculating the value of a single caplet, that is \[P^{d}(0,T_p) \mathbb{E}^{T_p}\bigg[ \Big(S_{a,b}(T_a) - K\Big)^{+}\bigg] = P^{d}(0,T_a) \mathbb{E}^{T_a}\bigg[ \Big(S_{a,b}(T_a) - K\Big)^{+} P^{d}(T_a,T_p)\bigg].\] Using a change of numeraire, we have that \[\begin{align} N^{a,b}_0 \mathbb{E}^{N^{a,b}}\bigg[ \Big(S_{a,b}(T_a) - K\Big)^{+} M(T_a,a,b,p)\bigg] &= N^{a,b}_0 \mathbb{E}^{N^{a,b}}\bigg[ \Big(S_{a,b}(T_a) - K\Big)^{+} \bigg]\nonumber \\ &+ N^{a,b}_0 \mathbb{E}^{N^{a,b}}\bigg[ \Big(S_{a,b}(T_a) - K\Big)^{+} \Big(M(T_a,a,b,p) - 1\Big) \bigg] \end{align}\] where \(N^{a,b}_0 = 01(0,a,b)\) and \(M(t,a,b,p)=\frac{P(t,T_p)}{01(t,a,b)}.\) The second term is the convexity adjustment for a CMS Cap. To calculate this term, it is usually assumed that \[M(t,a,b,p)= M_{p}\big(S_{a,b}(t)\big),\] and use the following replication formula \[\bigg(g\Big(S_{a,b}(t)\Big)\bigg)^{+}=g^{\prime}(K)\Big(S_{a,b}(t) - K\Big)^{+} + \int_{K}^{\infty} \Big(S_{a,b}(t) - u\Big)^{+} g^{\prime \prime}(u) du\] with \[g\Big(S_{a,b}(t)\Big)= \Big(S_{a,b}(t) - k\Big) \Big(M_{p}(S_{a,b}(t))-1\Big).\]

The primary concern with this formula is that it requires knowledge of the swaption price on the strike grid where the right side tends to \(\infty\). Consequently, the convexity adjustment heavily relies on the extrapolation assumption.

In recent papers such as [@Hagan19] and [@HaganWoodward20], the authors have proposed a new methodology that expresses the convexity adjustment for CMS Caps, CMS Floors, and CMS Swaps as a quadratic payoff. Using the singular perturbation technique, they obtain an expansion of the implied volatility for the quadratic payoff when the dynamics of the underlying follow a SABR-type model.

In the following sections, we will demonstrate how to use the Watanabe calculus to compute the price of quadratic payoff under a general stochastic local volatility model.

3 The connection between Quadratic Payoffs and the CMS Cap-Floor convexity↩︎

In this section, we will illustrate the connection between the quadratic payoffs and the convexity adjustment for CMS Cap, Floor, and Swap instruments. We will use Malliavin’s techniques to derive a generic representation for this adjustment. For more detailed content on Malliavin calculus see [@AlosLorite] or [@Nualart].

We will start with the convexity adjustment for a CMS Cap, as we saw in Section 2. We have that \[\begin{align} \label{ca95cms95cap} CA_{CMS}^{cap}(0,a,b)&= N^{a,b}_0 \mathbb{E}^{N^{a,b}}\Bigg[ \Big(S_{a,b}(T_a) - K\Big)^{+} \bigg(M(T_a,a,b,p) - 1\bigg) \Bigg] \nonumber \\ &= N^{a,b}_0 \mathbb{E}^{N^{a,b}}\Bigg[ \Big(S_{a,b}(T_a) - K\Big)^{+} \bigg( \mathbb{E}^{N^{a,b}}\Big[M(T_a,a,b,p)|S_{a,b}(T_a)\Big] - 1\bigg) \Bigg] \nonumber \\ &= N^{a,b}_0 \mathbb{E}^{N^{a,b}}\Bigg[ \Big(S_{a,b}(T_a) - K\Big)^{+} \bigg(M(S_{a,b}(T_a)) - 1\bigg) \Bigg] \nonumber \\ \end{align}\tag{6}\] where we will define \[M\big(S_{a,b}(T_a)\big):= \mathbb{E}^{N^{a,b}}\Big[M(T_a,a,b,p)|S_{a,b}(T_a)\Big].\] Here, the function \(M\big(S_{a,b}(T_a)\big)\) is called the mapping function. A popular choice for this function is the linear model. In this model, it is assumed that \[M\big(S_{a,b}(T_a)\big)= \mathcolor{highlightNEW}{\alpha}(T_a) S_{a,b}(T_a) + \beta(T_a).\] The parameters \(\mathcolor{highlightNEW}{\alpha}(T_a)\) and \(\beta(T_a)\) are chosen according to the properties of \(S_{a,b}(T_a)\) under the measure associated with the numeraire \(N_t^{a,b}\). Using the Clark-Ocone formula with \(F_t = S_{a,b}(t)\), we have \[M_{S_{a,b}(T_a)}= M_{S_{a,b}(0)} + \int_{0}^{T_a} \mathbb{E}_s^{N^{a,b}}\bigg[ \partial_{S} M\big(S_{a,b}(T_{a})\big) D^{W}_s S_{a,b}(T_{a}) \bigg] dW_s.\] If we freeze \(\partial_{S} M\big(S_{a,b}(T_{a})\big)\) in \(S_{a,b}(0)\), we get the following approximation \[\begin{align} M_{S_{a,b}(T_a)} &\approx M_{S_{a,b}(0)} + \partial_{S} M\big(S_{a,b} (0)\big) \int_{0}^{T_a} \mathbb{E}_s^{N^{a,b}}\bigg[ D^{W}_s S_{a,b}(T_{a}) \bigg] dW_s \\ &= M\big(S_{a,b}(0)\big) + \partial_{S} M\big(S_{a,b}(0)\big) \Big(S_{a,b}(T_a) - S_{a,b}(0)\Big). \end{align}\] Therefore, if we substitute the above equality in 6 , we get that \[CA_{CMS}^{cap}(0,a,b) \approx \Big(M\big(S_{a,b}(0)\big) - 1\Big) V^{C}(T_a,S_{a,b},K) + \Big(\partial_{S} M\big(S_{a,b}(0)\big)-1\Big) V^{QC}(T_a,S_{a,b},K).\] Likewise, in the CMS Floor case, we obtain that \[CA_{CMS}^{floor}(0,a,b) \approx \Big(M\big(S_{a,b}(0)\big) - 1\Big) V^{P}(T_a,S_{a,b},K) + \Big(\partial_{S} M\big(S_{a,b}(0)\big)-1\Big) V^{QP}(T_a,S_{a,b},K)\] and, in the swap case, we find that \[CA_{CMS}^{swap}(0,a,b) \approx \Big(M\big(S_{a,b}(0)\big) - 1\Big) \Big(S_{a,b}(0) - K\Big) + \Big(\partial_{S} M\big(S_{a,b}(0)\big)-1\Big) V^{QS}(T_a,S_{a,b},K).\]

4 Basic introduction to Watanabe’s exapansion↩︎

In this section, we will give a brief introduction to Watanabe’s expansion in the field of quantitative finance, see [@Kahl] for more details. This method supposes that we can find an expansion of the underlying \(F_t\) in terms of a small parameter \(\epsilon << 1\) \[\begin{align} \label{base95expansion} F_{\epsilon}(T) &= F_{0}+ \epsilon \cdot g_{1}(T) + \epsilon^{2} \cdot g_{2}(T) + \epsilon^{3} \cdot g_{3}(T) + \cdots \nonumber \\ &= F_{0} + \epsilon \sigma_{0}\sqrt{T} \Big(\hat{g}_{1}(T) + \epsilon \hat{g}_{2}(T) + \epsilon^{2} \hat{g}_{3}(T) + \cdots \Big) \end{align}\tag{7}\] where \[\begin{align} \hat{g}_{i}(T) = \frac{g_{i}(T)}{\sigma_{0}\sqrt{T}}. \end{align}\] The expansion is constructed such that \(\hat{g}_{1} \sim \mathcal{N}(0,1)\). We refer to \(\sigma_{0}\) as the initial volatility and \(T\) as the time to maturity. We use the auxiliary term \(\epsilon\) to derive the asymptotic option pricing formula. To calculate \[\begin{align} V^{C}(T,F_T,K):=\mathbb{E}^{\mathbb{Q}_N}\Big[(F_T - K)^{+}\Big] \end{align}\] we will use the forward expansion 7 . It is easy to show that \[\begin{align} \label{expansion95base95generic95payoff} \frac{V^{C}(T,F,K)}{N_T} &= \mathbb{E}^{\mathbb{Q}_N}\Bigg[\bigg(x_{0} - K + \epsilon \sigma_{0} \sqrt{T} \sum_{j=1}^{\infty} \epsilon^{j-1} \hat{g}_{j}(T)\bigg)^{+} \Bigg] \nonumber \\ &= \epsilon \sigma_{0} \sqrt{T} \mathbb{E}^{\mathbb{Q}_N}\Bigg[\bigg(\hat{g}_{1}(T) - y + \epsilon{}\sum_{j=1}^{\infty} \epsilon^{j-1}\hat{g}_{j+1}(T) \bigg)^{+} \Bigg] \end{align}\tag{8}\] where \[\begin{align} \label{label95y} y = \frac{K - F_{0}}{\epsilon \sigma_{0} \sqrt{T}} \end{align}\tag{9}\] and \[\begin{align} V^{C}_{(j)} := \partial^{(j)}_{F} (F - K)^{+}_{|F=\hat{g}_{1}(T),K=y}. \end{align}\] Then, using a Taylor’s expansion on \(V^{C}(T,T,T,T,\cdot,K)\) around \(\hat{g}_{1}(T) - y\), we get the following price expansion \[\begin{align} \label{expansion95generic95payoff95order955} \frac{V^{C}(T,F,K)}{N_T} &= \epsilon \sigma_{0} \sqrt{T} \mathbb{E}^{\mathbb{Q}^{N}}\Bigg[\Big(\hat{g}_{1}(T) - y\Big)^{+} + \epsilon V^{C}_{(1)} \hat{g}_{2}(T) \nonumber \\ & \left.+ \epsilon^{2} \bigg(V^{C}_{(1)} \hat{g}_{3}(T) + \frac{1}{2} V^{C}_{(2)} \hat{g}_{2}^{2} \bigg) \right. \nonumber \\ & \left. +\epsilon^{3} \bigg( V^{C}_{(1)}\hat{g}_{4} +V^{C}_{(2)}\hat{g}_{2}(T) \hat{g}_{3}(T) + \frac{1}{6} V^{C}_{(3)}\hat{g}^{2}_{2}(T) \bigg) \right. \nonumber \\ &+ \epsilon^{4} \bigg(V^{C}_{(1)} \hat{g}_{5}(T) + V^{C}_{(2)}\Big(\hat{g}_{2}(T)\hat{g}_{4}(T) + \frac{1}{2}\hat{g}^{2}_{3}(T) \Big) \nonumber \\ &+ \frac{1}{2} V^{C}_{(3)} \hat{g}^{2}_2(T) \hat{g}_3(T) + \frac{1}{24} V^{C}_{(4)} \hat{g}^{2}_{2}(T) \bigg) + O(\epsilon^{5})\Bigg]. \end{align}\tag{10}\]

Note that to use the price expansion given in 10 , we need to know the terms \(\hat{g}_i(T)\) for \(i=1, \cdots, n\). These terms are expressed as \[\label{g95i95definition} \hat{g}_i(T) = \beta_i(T) I_{(\underbrace{0,\cdots,0}_\text{n-times},\underbrace{1,\cdots,1}_\text{m_B-times},\underbrace{1,\cdots,1}_\text{m_W-times})}(T)\tag{11}\] where \(\beta_i(T)\) is a non-stochastic function and \[\begin{align} I_{(\underbrace{0,\cdots,0}_\text{n-times},\underbrace{1,\cdots,1}_\text{m_B-times},\underbrace{1,\cdots,1}_\text{m_W-times})}(T) \end{align}\] is the next iterated integral

\[\label{g95i95integral95representation} \int_{0}^{T}\;\cdots \int_{0}^{t_{m_B}} \int_{0}^{t^{\prime}_1} \cdots \int_{0}^{t^{\prime}_{m_W}} \int_{0}^{t^{\prime \prime}_1} \cdots \int_{0}^{t^{\prime \prime}_{n}} du^{\prime \prime}_{n + 1} \cdots du^{\prime \prime}_1 dB_{u^{\prime}_{m_B + 1}} \cdots dB_{u^{\prime}_1} dW_{u_{m_B+ 1}} \cdots dW_{u_1}\tag{12}\] with \(W_t\) and \(B_t\) correlated Brownian motions of the underlying and volatility processes, respectively. The integral (12 ) may seem complex, but using the Brownian bridge and conditioning on \(\hat{g}_1(T)\), the expectation of the integral (11 ) is simplified to

\[\mathbb{E}^{\mathbb{Q}_N}\Bigg[ I_{(\underbrace{0,\cdots,0}_\text{n-times},\underbrace{1,\cdots,1}_\text{m_B-times},\underbrace{1,\cdots,1}_\text{m_W-times})}(T) \Bigg] = \int_{0}^{T} g_{(\underbrace{0,\cdots,0}_\text{n-times},\underbrace{1,\cdots,1}_\text{m_B-times},\underbrace{1,\cdots,1}_\text{m_W-times})}(u) \hat{\phi}_1(u) du\] where \(\hat{\phi}_1(u)\) is the pdf of \(\hat{g}_1(T)\) and \[\begin{align} g_{(\underbrace{0,\cdots,0}_\text{n-times},\underbrace{1,\cdots,1}_\text{m_B-times},\underbrace{1,\cdots,1}_\text{m_W-times})}(u) \end{align}\] is a function to be defined.
To illustrate, we will consider the normal SABR model and we will obtain the order \(O(\epsilon^{2})\) approximation as described above.

5 Quadratic option price↩︎

In this section, following the methodology presented in Section 4, we will derive a Watanabe expansion for quadratic payoff derivatives when the underlying follows a stochastic volatility model, a local volatility model, and a local stochastic volatility. It is worth noting that the quadratic payoff \(V^{QC}(F, K)\) is more influenced by the tail of the distribution than the call option price \(V^{C}(T, F, K)\). In the case of local volatility, we will explain why this happens.
Consider that the forward price dynamics are given by \[\label{forward95local95model} dF_t = \sigma(F_t) dW^{N}_t\tag{14}\] where \(\sigma(\cdot)\) is Lipschitz function to ensure existence and uniqueness of the solution. Now, if we use Itô-Tanaka formula to \(\Big((F_t - K)^{+}\Big)^{2}\), we get that \[\label{qc95lv95representation} \frac{V^{QC}(T,F,K)}{N_T} = \bigg((F_0 - K)^{+}\bigg)^{2} + \int_{0}^{T} \mathbb{E}^{\mathbb{Q}_N}\bigg[I\Big(F_u > K\Big) \sigma^{2}(F_u)\bigg] du.\tag{15}\] It is noteworthy that the value of \(\frac{V^{QC}(T,F,K)}{N_T}\) depends on the value of \(\sigma(F)\) when \(F \in (K, \infty)\). For this reason, the quadratic call price is highly sensitive to our assumptions about the behavior of the local volatility function at extreme strikes.
As seen in Section 4, we will find the Watanabe’s price expansion for the quadratic call, put, and swap.

5.1 Quadratic option under local volatility↩︎

When the volatility is constant, that is \(\sigma(\cdot)=\sigma\), the quadratic call option price is \[\label{bachelier95quadratic95option} \frac{V^{QC}(T,F,K)}{N_T} = \bigg(\Big(F_0 - K\Big)^{2} + \sigma^2 T\bigg) \Phi\bigg( \frac{K - F_0}{\sigma \sqrt{T}}\bigg) + \Big(F_0 - K\Big) \sigma \sqrt{T} \phi\bigg(\frac{K - F_0}{\sigma \sqrt{T}}\bigg).\tag{16}\] where \(\Phi(\cdot)\) and \(\phi(\cdot)\) are, respectively, the density and cumulative function for a standard normal distribution. We are interested in building an option price expansion with base (16 ). We will consider the small diffusion model \[\label{local95vol95model95small} dF_{\epsilon,t} = \epsilon \sigma(F_{\epsilon,t}) dW^{N}_t\tag{17}\] where \(\sigma(\cdot)\) is a positive function such that belongs to \(\mathcal{C}^{3}_b(\mathbb{R})\).

As previously stated, to obtain the Watanabe’s price expansion, we need to find an expression for the \(\hat{g}_{i}\) terms. Comparing Lemma 2 with the price expansion 7 , we have that

\[\begin{align} \label{coeff95lv95expansion} \hat{g}_{1}(T) &= \frac{W_T}{\sqrt{T}} \nonumber \\ \hat{g}_{2}(T) &= \partial_F \sigma(F_0) \frac{I_{(1,1)(T)}}{\sqrt{T}} \nonumber \\ \hat{g}_{3}(T) &= \bigg( \frac{1}{2} \Big(\partial^{2}_F \sigma(F_0)\Big) \sigma(F_0) + \Big(\partial_F \sigma(F_0)\Big)^{2}\bigg) \frac{I_{(1,1,1)}(T)}{\sqrt{T}} + \frac{1}{2} \partial^2_F \sigma(F_0) \sigma(F_0) \frac{I_{(0,1)}(T)}{\sqrt{T}}. \end{align}\tag{18}\]

5.2 Quadratic option under stochastic local volatility↩︎

SLV-type models are widely used in the industry; see, for example, [@Bang2019] and [@Labordere]. These models allow for control of the backbone through a layer of local volatility and its movement over time. In the interest rate market, a SABR dynamic is the market standard. We will consider a general SABR dynamics given by: \[\begin{align} \label{general95SLV95mode} dF_t &= C(F_t) \sigma_t dW_t, \nonumber \\ d\sigma_t &= \nu \sigma_t dB_t. \end{align}\tag{19}\] Here, \(<dW_t,dB_t>= \rho dt\) and \(C(\cdot)\) is a positive function that belongs to \(\mathcal{C}^{3}_b(\mathbb{R})\). To achieve an expansion for \(F_t\) similar to ?? , we will consider the next parametrization of the above dynamic \[\begin{align} \label{general95SLV95mode95small95e} dF_{\epsilon, t} &= \epsilon C(F_{\epsilon, t}) \sigma_{\epsilon, t} dW_t, \nonumber \\ d\sigma_{\epsilon, t} &= \epsilon \nu \sigma_{\epsilon, t} dB_t . \end{align}\tag{20}\]

Then, as we have seen before, we can imply the \(\hat{g}_{i}\) terms. \[\begin{align} \label{g95s95slv} \hat{g}_1(T) &= \frac{W_T}{\sqrt{T}} \nonumber \\ \hat{g}_2(T) &= \partial_F C(F_0) \frac{\mathcolor{highlightNEW}{\alpha}}{\sqrt{T}}I_{(1,1)}(T) + \frac{\nu }{\sqrt{T}}\int_{0}^{T}\int_{0}^{u_1} dB_{u_2} dW_{u_2} \nonumber \\ \hat{g}_3(T) &= \frac{\nu^2}{\sqrt{T}}\int_0^T \int_0^{u_1}\int_0^{u_2} dB_s dB_{u_2} dW_{u_1} + \frac{\mathcolor{highlightNEW}{\alpha}\nu}{\sqrt{T}} \partial_F C(F_0) \int_{0}^{T} B_u W_u dW_u \nonumber \\ &+ \frac{\mathcolor{highlightNEW}{\alpha}^2 }{2\sqrt{T}}\Big(\partial^2_F C(F_0)\Big)C(F_0) \int_{0}^{T} W^2_u dW_u + \frac{\mathcolor{highlightNEW}{\alpha}\nu}{\sqrt{T}}\partial_F C(F_0) \int_0^{T} \int_0^{u_1} \int_0^{u_2} dB_s W_{u_2} dW_{u_1} \nonumber \\ &+ \mathcolor{highlightNEW}{\alpha}^{2} \partial_{F}C(F_0) \frac{I_{(1,1,1)}(T)}{\sqrt{T}}. \end{align}\tag{21}\]

6 Future works↩︎

In a recent paper [@Skoufis2024], the author proposes a closed formula for the convexity adjustment in the case of average RFR indexes when the underlying follow a time-dependent SABR model. This problem was addressed in [@DavidRaul2023] in the case of the Cheyette model with local volatility using Malliavin calculus. To make the paper self-contained, let us introduce the basic concepts. Let us define \[\begin{align} I(T_a,T_b)&=\mathbb{E}^{\mathbb{Q}}\left[\frac{e^{\int_{T_a}^{T_b} r_u du} - 1}{T_b -T_a} \right],\\ I_{avg}(T_a,T_b)&=\frac{\int_{T_a}^{T_b} r_u du}{T_b - T_a}, \end{align}\] and \(R(t)=\mathbb{E}^{\mathbb{Q}}\left[I(T_a,T_b) \right]\). In addition, we will suppose that \[\begin{align} dR_{t}(T_a,T_b) &= h_t C(R_{t}(T_a,T_b)) \sigma_{t} dW_t, \nonumber \\ d\sigma_{t} &= \nu \sigma_{t} dB_t \nonumber \\ <dW_t,dB_t> &= \rho dt \end{align}\] where \(h_t = \min \left\{1, \left(\frac{T_b - t}{T_b - T_a} \right)^{q} \right\}\) with \(0 < q < 1\). From the definition of \(I_{avg}(T_a,T_b)\), we have that \[I_{avg}(T_a,T_a,T_b) = \frac{\ln\Big(1 + (T_b - T_a) I(T_a,T_b)\Big)}{T_b - T_a} = \frac{\ln\Big(1 + (T_b - T_a) R(T_a, T_a,T_b)\Big)}{T_b - T_a}.\] Then, the issue is to compute \(\mathbb{E}\left[ I_{avg}(T_a,T_b) \right]\) when \(R_{t,T_a,T_b}\) follow the dynamic [sabr95time95dependent95rfr]. Using Itô formula, we obtain that \[\begin{align} \mathbb{E}\left[ I_{avg}(T_a,T_b) \right] &= I(0,T_a,T_b) + \frac{1}{2} \mathbb{E}^{\mathbb{Q}} \left[\int_{0}^{T_b} G^{\prime \prime}(R(u,T_a,T_b)) h_u^{2} \sigma^{2}_u C^{2}(R(u,T_a,T_b)) ds \right] \\ & \approx I(0,T_a,T_b) + \frac{1}{2} \mathbb{E}^{\mathbb{Q}} G^{\prime \prime}(R(0,T_a,T_b)) \mathbb{E}^{\mathbb{Q}}\bigg[ \Big(R(T_b,T_a,T_b) - R(0,T_a,T_b)\Big)^2 \bigg] \end{align}\] where \[\begin{align} G(R)=\frac{\ln(1 + (T_b - T_a) R)}{T_b - T_a}. \end{align}\] Let us suppose that we have an expansion similar to 7 . Then, it is easy to show that \[\begin{align} \mathbb{E}^{\mathbb{Q}}&\bigg[ \Big(R(T_b,T_a,T_b) - R(0,T_a,T_b)\Big)^2 \bigg] \\ &= \mathcolor{highlightNEW}{\alpha}^2 C^2(R(0,T_a,T_b)) \mathbb{E}^{\mathbb{Q}}\bigg[ \hat{g}^2_{1}(T_b) + \hat{g}^2_{1}(T_b) + 2\hat{g}_{1}(T_b) \hat{g}_{3}(T_b) + O(\epsilon^3)\bigg] T. \end{align}\] The challenge for future papers is to obtain the expansion of \(F_{\epsilon, t}\) under a general SLV dynamic, this way could compute the convexity adjustment of CMS or average RFR, even the pricing of generic terminal payoff for a general SLV model.

7 Conclusions↩︎

In conclusion, we have proposed a new approach for pricing CMS derivatives. Our method uses Mallaivin’s calculus and Watanabe’s expansions. Firstly, we establish a model-independent connection between the price of a CMS derivative and the quadratic payoffs using Mallaivin calculus. Then, we extend the results of Watanabe’s expansions to the quadratic payoff case under local and stochastic local volatility. The approximations obtained are generic. To evaluate its accuracy, we compare the approximations numerically under the SABR model against the industry standards: Hagan’s approximation, and a Monte Carlo simulation. Our approximations are as precise as Hagan’s approximation, but they cover a wider range of models.

Appendix↩︎

We will consider \(W_{t}\), \(W_{\perp,t}\) two independent Brownian motions.

8 Normal SABR↩︎

8.1 Calculation of \(\mathbb{E}^{\mathbb{Q}^{N}}\left[ \hat{g}_{2}(T) |\hat{g}_1(T)\right]\)↩︎

We have that \[\mathbb{E}^{\mathbb{Q}^{N}}\Bigg[ \hat{g}_{2}(T) \Big|\hat{g}_1(T)\Bigg]=\mathbb{E}^{\mathbb{Q}^{N}}\Bigg[\frac{\nu}{\sqrt{T}} \bigg(\rho \int_{0}^{T} W_{u} dW_{u} + \bar{\rho} \int_{0}^{T} W_{\perp,u} dW_{u}\bigg) \Big| \hat{g}_{1}(T)\Bigg].\] Therefore, \[\mathbb{E}^{\mathbb{Q}^{N}}\Bigg[ \hat{g}_{2}(T) \Big|\hat{g}_1(T)\Bigg]= \frac{\rho \nu}{\sqrt{T}} \mathbb{E}^{\mathbb{Q}^{N}}\Bigg[ \int_{0}^{T} W_{u} dW_{u} \Big| \hat{g}_{1}(T) \Bigg] = \frac{1}{2}\rho \nu\sqrt{T}\Bigg(\hat{g}^{2}_{1}(T) - 1\Bigg).\]

8.2 Calculation of \(\mathbb{E}^{\mathbb{Q}^{N}}\left[ \hat{g}_{3}(T) |\hat{g}_1(T)\right]\)↩︎

We will follow the same approach as in the previous section. Therefore, we have \[\begin{align} \hat{g}_{3}(T) &=& \frac{\nu^{2}}{2\sqrt{T}} \int_{0}^{T} B_{u}^{2} - u dW_{u} \\ &=& \frac{\nu^{2}}{2\sqrt{T}} \int_{0}^{T} B_{u}^{2} dW_{u} - \frac{\nu^{2}}{2\sqrt{T}} \int_{0}^{T} u dW_{u}. \end{align}\] In one hand, we obtain that \[\begin{align} \mathbb{E}^{\mathbb{Q}^{N}}\Bigg[B(T)\Big|W_{T} = \hat{g}_{1}(T) \sqrt{T} \Bigg] =\frac{1}{4}\nu^{2} \hat{g}_{1}(T) T. \end{align}\] On other hand, we have that \[\begin{align} \frac{\rho^{2} \nu^{2}}{2\sqrt{T}} \mathbb{E}^{\mathbb{Q}^{N}}\Bigg[\int_{0}^{T} W_{u} dW_{u} \Big|W_{T} = \hat{g}_{1}(T) \sqrt{T} \Bigg] &= \frac{1}{12}\rho^{2} \nu^{2} \Bigg( 2 \hat{g}^{3}_{1}(T) - \hat{g}_{1}(T) \Bigg) T. \\ \frac{\hat{\rho}^{2} \nu^{2}}{2\sqrt{T}} \mathbb{E}^{\mathbb{Q}^{N}}\Bigg[\int_{0}^{T} W^{2}_{u,\perp} dW_{u} \Big|W_{T} = \hat{g}_{1}(T) \sqrt{T} \Bigg] &= \frac{1}{4}\hat{\rho}^{2} \nu^{2}\hat{g}_{1}(T)T. \end{align}\] Therefore, when combine it, we get that \[\label{g395g1} \mathbb{E}^{\mathbb{Q}^{N}}\Bigg[\hat{g}_{3}(T) \Big| \hat{g}_{1}(T) \Bigg] = \frac{1}{12} \rho^{2} \nu^{2} \Bigg( 2 \hat{g}^{3}_{1}(T) - \hat{g}_{1}(T) \Bigg) T + \frac{1}{4}\hat{\rho}^{2} \nu^{2}\hat{g}_{1}(T) T - \frac{1}{4}\nu^{2} \hat{g}_{1}(T) T.\tag{22}\]

8.3 Calculation of \(\mathbb{E}^{\mathbb{Q}^{N}}\left( \hat{g}^{2}_{2}(T) |\hat{g}_1(T)\right)\)↩︎

From the definition of \(\hat{g}_{2}(T)\), we have that \[\begin{align} \frac{\nu^{2}}{T} \Bigg(\underbrace{\rho^{2} \bigg(\int_{0}^{T} W_{u} dW_u\bigg)^{2}}_{A(T)} + \underbrace{\bar{\rho}^{2} \bigg(\int_{0}^{T} W_{\perp, u} dW_u\bigg)^{2}}_{B(T)} + \underbrace{2 \rho \bar{\rho} \bigg(\int_{0}^{T} W_{\perp, u} dW_u\bigg) \bigg(\int_{0}^{T} W_{u} dW_u\bigg)}_{C(T)} \Bigg) \end{align}\] Now, we will compute \(\mathbb{E}^{\mathbb{Q}^{N}}\Big[\cdot | \hat{g}_1(T)\Big]\) for each term. The easiest one is \(\mathbb{E}^{\mathbb{Q}^{N}}\Big[C(T)|\hat{g}_{1}(T)\Big]\), due to the independence of \(W_{T}\) and \(W_{\perp,t}\). We get \[\mathbb{E}^{\mathbb{Q}^{N}}\Big[C(T)|\hat{g}_{1}(T)\Big] = 0.\] Then, we will compute the rest of the terms \[\label{g95295295first95term} \frac{\nu^{2} \rho^{2}}{T} \mathbb{E}^{\mathbb{Q}^{N}} \Bigg[\bigg(\int_{0}^{T} W_{u} dW_{u}\bigg)^{2} \Big|W_{T} = \hat{g}_{1}(T) \sqrt{T}\Bigg] = \frac{1 }{4}\nu^{2} \rho^{2} T\Bigg(\hat{g}^{2}_{1}(T) - 1 \Bigg)^{2}.\tag{23}\]

For the last term, we must remember the next identity for a Brownian bridge \(X_s\) with \(X_{0} = 0\) and \(X_{T} = \hat{g}_{1}(T)\sqrt{T}\) \[\label{brownian95bridge95dw95dw} \mathbb{E}^{\mathbb{Q}^{N}}\Bigg[dX_u dX_s \Bigg] = \frac{\hat{g}^{2}_{1}(T)}{T} + \frac{\min(u,s)}{t(t -\min(u,s))} + \delta(u -s) - \frac{1}{(t -\min(u,s))} du ds.\tag{24}\] Therefore, we have that \[\begin{align} \label{g95295295second95term} &\frac{\nu^{2} \hat{\rho}^{2}}{T} \mathbb{E}^{\mathbb{Q}^{N}}\Bigg[ \left(\int_{0}^{T} W_{u,\perp} dW_{u} \right)^{2} \Big|W_{T} = \hat{g}_{1}(T) \sqrt{T} \Bigg] \nonumber \\ &= \frac{\nu^{2} \hat{\rho}^{2}}{T} \left( \int_{0}^{T}\int_{0}^{T} \min(u_{1},u_{2}) \mathbb{E}^{\mathbb{Q}^{N}}\Bigg[dW_{u_{1}} dW_{u_{2}} \Big|W_{T} = \hat{g}_{1}(T) \sqrt{T} \Bigg]\right) \nonumber \\ &= \frac{\nu^{2} \hat{\rho}^{2}}{T} \int_{0}^{T}\int_{0}^{T} \min(u_{1},u_{2}) \left( \frac{\hat{g}^{2}_{1}(T)t}{t^{2}} + \frac{\min(u_{1},u_{2})}{t(t - \min(u_{1},u_{2}))} + \delta(u_1-u_2) - \frac{1}{t- \min(u_{1},u_{2})} \right) du_{1}du_{2} \nonumber \\ &= \frac{1}{6}\nu^{2}\hat{\rho}^{2}\bigg(2 \hat{g}^{2}_{1}(T)+1\bigg)T. \end{align}\tag{25}\] Now, if we join 23 and 25 , we obtain that \[\label{sec:a13:g295power952} \mathbb{E}^{\mathbb{Q}^{N}}\Bigg[\hat{g}^{2}_{1}(T)\Big|W_{T} = \hat{g}_{1}(T) \sqrt{T}\Bigg] =\frac{1}{12}\nu^{2}\Bigg( 3\rho^{2} \bigg(\hat{g}^{2}_{1}(T) - 1 \bigg)^{2} + 4\hat{\rho}^{2} \hat{g}^{2}_{1}(T)+ 2\hat{\rho}^{2}\Bigg) T.\tag{26}\]

9 SABR↩︎

9.1 Calculation of \(\mathbb{E}^{\mathbb{Q}^{N}}\left[ \hat{g}_{2}(T) |\hat{g}_1(T)\right]\)↩︎

From 21 , we have that \[\mathbb{E}^{N}\Bigg[ \hat{g}_{2}(T) \Big| \hat{g}_{1}(T) \Bigg] = \mathcolor{highlightNEW}{\alpha}\partial_F C(F_0) \mathbb{E}^{N}\Bigg[ \frac{I_{(1,1)}(T)}{\sqrt{T}} \Big| \hat{g}_{1}(T) \Bigg] + \nu \rho \mathbb{E}^{N}\Bigg[ \frac{I_{(1,1)}(T)}{\sqrt{T}} \Big| \hat{g}_{1}(T) \Bigg]\] Now, if we use 24 , we have that \[\mathbb{E}^{N}\Bigg[\hat{g}_{2}(T) \Big| \hat{g}_{1}(T) \Bigg] = \frac{\sqrt{T}}{2} \bigg(\hat{g}^{2}_1(T) - 1\bigg) \bigg(\nu \rho + \frac{1}{2}\mathcolor{highlightNEW}{\alpha}\partial_F C(F_0) \bigg).\]

9.2 Calculation of \(\mathbb{E}^{\mathbb{Q}^{N}}\big[ \hat{g}_{3}(T) |\hat{g}_1(T)\big]\)↩︎

From the expression of \(\hat{g}_3(T)\) in 21 , we obtain that \[\mathbb{E}^{N}\Bigg[\hat{g}_{3}(T) \Big| \hat{g}_{1}(T) \Bigg] = A(T) + B(T) + C(T)\] where \[\begin{align} A(T) &= \frac{\nu^2}{\sqrt{T}}\int_0^T \int_0^{u_1}\int_0^{u_2} dB_s dB_{u_2} dW_{u_1} + \mathcolor{highlightNEW}{\alpha}\nu \Big(\partial_F C(F_0)\Big) \frac{1}{\sqrt{T}}\int_{0}^{T} B_u W_u dW_u, \\ B(T) &= \frac{1}{2} \mathcolor{highlightNEW}{\alpha}^2 \Big(\partial^2_F C(F_0)\Big)C(F_0)\frac{1}{\sqrt{T}}\int_{0}^{T} W^2_u dW_u + \mathcolor{highlightNEW}{\alpha}\nu \Big(\partial_F C(F_0)\Big) \frac{1}{\sqrt{T}}\int_0^{T} \int_0^{u_1} \int_0^{u_2} dB_s W_{u_2} dW_{u_1}, \\ C(T) &= \mathcolor{highlightNEW}{\alpha}^{2} \Big(\partial_{F}C(F_0)\Big)^2 \frac{1}{\sqrt{T}}I_{(1,1,1)}(T). \end{align}\] Through a direct computation and employing the representation of the distinct iterated integrals, we get that \[\begin{align} \mathbb{E}^{N}\Bigg[A(T)\Big|\hat{g}_1(T)\Bigg]&= \frac{1}{6}\nu^2 \rho^2 \bigg(\hat{g}^3_1(T) - 3\hat{g}_1(T)\bigg) T + \frac{1}{6}\mathcolor{highlightNEW}{\alpha}\nu \Big(\partial_F C(F_0)\Big) \bigg(2\hat{g}^3_1(T)- 3\hat{g}_1(T)\bigg) T \\ \mathbb{E}^{N}\Bigg[B(T)\Big|\hat{g}_1(T)\Bigg]&= \frac{1}{12}\mathcolor{highlightNEW}{\alpha}^2 \Big(\partial^2_F C(F_0)\Big)C(F_0)\bigg(2\hat{g}^3_1(T) - 3\hat{g}_1(T)\bigg)T + \frac{1}{3}\mathcolor{highlightNEW}{\alpha}\rho \nu \Big(\partial_F C(F_0)\Big) \bigg(\hat{g}^3_1(T) - 3\hat{g}_1(T)\bigg)T \\ \mathbb{E}^{N}\Bigg[C(T)\Big|\hat{g}_1(T)\Bigg]&= \frac{1}{3}\mathcolor{highlightNEW}{\alpha}^{2} \Big(\partial_{F}C(F_0)\Big)^2 \bigg(\hat{g}^3_1(T) - 3\hat{g}_1(T) \bigg)T \end{align}\] Rearranging the previous equalities, we obtain that \[\begin{align} \mathbb{E}^{N}\Bigg[ \hat{g}_{3}(T) | \hat{g}_{1}(T) \Bigg] &=\frac{1}{3}\Bigg(\frac{1}{2}\nu^2 \rho^2 + \mathcolor{highlightNEW}{\alpha}\rho \nu \partial_F C(F_0) + \mathcolor{highlightNEW}{\alpha}^{2} \Big(\partial_{F}C(F_0)\Big)^2 \Bigg) \Bigg(\hat{g}^3_1(T) - 3\hat{g}_1(T) \Bigg)T \\ &+\frac{1}{6}\Bigg(6\mathcolor{highlightNEW}{\alpha}\nu \partial_F C(F_0) + 3\mathcolor{highlightNEW}{\alpha}^2 \Big(\partial^2_F C(F_0)\Big)C(F_0) \Bigg) \Bigg(2\hat{g}^3_1(T) - 3\hat{g}_1(T)\Bigg)T. \end{align}\]

9.3 Calculation of \(\mathbb{E}^{\mathbb{Q}^{N}}\big[ \hat{g}^{2}_{2}(T) |\hat{g}_1(T)\big]\)↩︎

From 9.1 and 26 , we have that \[\mathbb{E}^{N}\Bigg[ \hat{g}^{2}_{2}(T) \Big| \hat{g}_{1}(T) \Bigg] = A(T) + B(T)\] where \[\begin{align} A(T)&= \frac{1}{4} \mathcolor{highlightNEW}{\alpha}^2 \Big(\partial_F C(F_0)\Big)^2 \Big(\hat{g}^{2}_1(T) - 1\Big)^2 T ,\\ B(T) &= \frac{\nu^{2}}{12}\bigg(3 \rho^{2}\Big(\hat{g}^{2}_{1}(T) - 1 \Big)^{2} + 4\hat{\rho}^{2} \hat{g}^{2}_{1}(T) +2\hat{\rho}^{2}\bigg)T. \end{align}\] Then, if we combine both expressions, we get that \[\mathbb{E}^{N}\Bigg[ \hat{g}^{2}_{2}(T) | \hat{g}_{1}(T) \Bigg] = \frac{1}{12}\Bigg(3\mathcolor{highlightNEW}{\alpha}^2 \Big(\partial_F C(F_0)\Big)^2+ \nu^{2} \rho^{2} \bigg) \Big(\hat{g}^{2}_{1}(T) - 1 \Big)^{2} +4\nu^{2}\hat{\rho}^{2} \hat{g}^{2}_{1}(T) + 2\nu^{2}\hat{\rho}^{2}\Bigg)T.\]

10 Call option expansion normal SABR↩︎

10.1 \(\epsilon^0\)-coefficient call option expansion↩︎

In this case, we must to compute \(\mathbb{E}^{\mathbb{Q}^{N}}\bigg[\Big(\hat{g}_1(T) - y\Big)^{+}\bigg]\). From the expression of \(\hat{g}_1\), see ?? , it is easy to show that

\[\label{order95095price} \mathbb{E}^{\mathbb{Q}^{N}}\bigg[\Big(\hat{g}_1(T) - y\Big)^{+}\bigg] = \phi(y) - y \bar{\Phi}(y)\tag{27}\] where \(\phi(\cdot)\) is the PDF for a standard normal and \(\bar{\Phi}(y) = 1 - \Phi(y)\) with \(\Phi(\cdot)\) the CDF for a standard normal.

10.2 \(\epsilon\)-coefficient call option expansion↩︎

We have to compute \(\mathbb{E}^{\mathbb{Q}^{N}}\bigg[I\Big(\hat{g}_{1}(T) > y\Big) \hat{g}_{2}(T) \bigg]\). From 8.1, we have that \[\begin{align} \label{order95195price} \mathbb{E}^{\mathbb{Q}^{N}}\Bigg[ I\Big(\hat{g}_{1}(T) > y\Big) \hat{g}_2(T)\Bigg]&=\mathbb{E}^{\mathbb{Q}^{N}}\Bigg[ I\Big(\hat{g}_{1}(T) > y\Big) \mathbb{E}^{\mathbb{Q}^{N}}\Big[\hat{g}_2(T) | \hat{g}_1(T) \Big] \Bigg] \nonumber\\ &= \frac{\rho \nu \sqrt{T}}{2} y \phi(y). \end{align}\tag{28}\]

10.3 \(\epsilon^2\)-coefficient call option expansion↩︎

As we did before, we have that \[\begin{align} \label{order95395term951} \mathbb{E}^{\mathbb{Q}^{N}}\Bigg[ I\Big(\hat{g}_1(T) > y\Big) \hat{g}_3(T) \Bigg] &= \mathbb{E}^{\mathbb{Q}^{N}}\Bigg[I\Big(\hat{g}_1(T) > y\Big) \frac{\nu^{2}}{12}\bigg(\rho^{2} \Big( 2 \hat{g}^{3}_{1}(T) - \hat{g}_{1}(T)\Big) + 3\hat{\rho}^{2} \hat{g}_{1}(T) - 3 \hat{g}_{1}(T) \bigg)T\Bigg] \nonumber \\ &= \frac{1}{12}\nu^{2}\phi(y) \bigg(\rho^{2}\left(2y^{2} + 1\right) + 3\hat{\rho}^{2} - 3 \bigg) T, \end{align}\tag{29}\] and \[\begin{align} \label{order95395term952} \mathbb{E}^{\mathbb{Q}^{N}}\Bigg[ \frac{1}{2}\delta\Big(\hat{g}_1(T) - y\Big) \hat{g}^{2}_{2}(T)\Bigg] &= \frac{1}{24}\nu^{2}\phi(y)\bigg( 3 \rho^{2}\Big(y^{2} - 1 \Big)^{2} + 4\hat{\rho}^{2} y + 2\hat{\rho}^{2}\bigg) T. \end{align}\tag{30}\]

11 Quadratic call option local volatility↩︎

11.1 \(\epsilon^0\)-coefficient quadratic call option expansion↩︎

We must to compute \(\mathbb{E}^{\mathbb{Q}^{N}}\Bigg[\Big((\hat{g}_1 - y)^{+}\Big)^{2}\Bigg]\). From the expression of \(\hat{g}_1\), see ?? , it is easy to show that \[\label{qc95order95095price95lv} \mathbb{E}^{\mathbb{Q}^{N}}\Bigg[\Big((\hat{g}_1(T) - y)^{+}\Big)^{2}\Bigg] = \Big(1+y^{2}\Big) \bar{\Phi}(y) - y \phi(y).\tag{31}\]

11.2 \(\epsilon\)-coefficient quadratic call option expansion↩︎

To calculate this term, we have to calculate \(\mathbb{E}^{\mathbb{Q}^{N}}\Bigg[2 \Big(\hat{g}_1(T) - y\Big)^{+} \hat{g}_2(T)\Bigg]\). It is easy to show using integration by parts that

\[\label{qc95order95195price95lv} \mathbb{E}^{\mathbb{Q}^{N}}\Bigg[2\Big(\hat{g}_1(T) - y\Big)^{+} \hat{g}_{2}(T) \Bigg] = \Big(\partial_{F}\sigma(F_0)\Big) \sqrt{T} \phi(y).\tag{32}\]

11.3 \(\epsilon^2\)-coefficient quadratic call option expansion↩︎

The last term that we will calculate is

\[\mathbb{E}^{\mathbb{Q}^{N}}\Bigg[2 \Big(\hat{g}_1(T) - y\Big)^{+} \hat{g}_3(T) + I\Big(\hat{g}_1(T) > y\Big) \hat{g}^{2}_2(T)\Bigg] = \mathbb{E}^{\mathbb{Q}^{N}}\Bigg[A(T) + B(T)\Bigg]\] where \[\begin{align} A(T) &= 2 \Big(\hat{g}_1(T) - y\Big)^{+} \hat{g}_3(T), \\ B(T) &= I\Big(\hat{g}_1(T) > y\Big) \hat{g}^{2}_2(T). \end{align}\] From 18 , we have that \[B(T) = \frac{1}{4}I\Big(\hat{g}_1(T) > y\Big) \Big(\partial_F\sigma(F_0)\Big)^2 \Big(\hat{g}^{2}_1(T) - 1\Big)^2 T.\] Also, we get that \[\begin{align} A(T) = \Big(\hat{g}_1(T) - y\Big)^{+} \Bigg(\frac{1}{3}\bigg( &\Big(\partial^{2}_F \sigma(F_0)\Big) \sigma(F_0) + \Big(\partial_F \sigma(F_0)\Big)^{2}\bigg) \bigg(\hat{g}^{3}_1(T) - 3\hat{g}_1(T)\bigg)T\\ & + \Big(\partial^{2}_F \sigma(F_0)\Big) \sigma(F_0)\bigg(\hat{g}_1(T)T - \frac{1}{\sqrt{T}}\int_{0}^{T} W_u du \bigg) \Bigg). \end{align}\] Therefore, we have that \[\begin{align} \mathbb{E}^{\mathbb{Q}^{N}}\Bigg[2 \Big(\hat{g}_1(T) - y\Big)^{+} \hat{g}_3(T)\Bigg] &=\frac{1}{6}\Bigg( 2\bigg( \Big(\partial^{2}_F \sigma(F_0)\Big) \sigma(F_0) + \Big(\partial_F \sigma(F_0)\Big)^{2}\bigg) y \phi(y) + 3\Big(\partial^{2}_F\sigma(F_0)\Big)\sigma(F_0)\bar{\Phi}(y) \Bigg) T , \nonumber \\ \mathbb{E}^{\mathbb{Q}^{N}}\Bigg[ I\Big(\hat{g}_1(T) > y\Big) \hat{g}^{2}_2(T) \Bigg] &= \frac{1}{4}\Big(\partial_F \sigma(F_0)\Big)^{2} \bigg(\Big(y^3+y\Big)\phi(y)+ 2\bar{\Phi}(y)\bigg)T. \end{align}\]

12 Quadratic call option SABR↩︎

12.1 \(\epsilon^0\)-coefficient quadratic call option expansion↩︎

We do the same calculation that in 31 .

12.2 \(\epsilon\)-coefficient quadratic call option expansion↩︎

We have to compute \(\mathbb{E}^{\mathbb{Q}^{N}}\Bigg[2\Big(\hat{g}_{1}(T) - y\Big)^{+} \hat{g}_{2}(T) \Bigg]\). Then, \[\mathbb{E}^{\mathbb{Q}^{N}}\Bigg[2\Big(\hat{g}_{1}(T) - y\Big)^{+} \hat{g}_{2}(T) \Bigg] = \mathbb{E}^{\mathbb{Q}^{N}}\Bigg[2\Big(\hat{g}_{1}(T) - y\Big)^{+} \mathbb{E}^{\mathbb{Q}^{N}}\Big[ \hat{g}_{2}(T)| \hat{g}_{1}(T)\Big]\Bigg].\] Using 9.1, we have that \[\begin{align} \mathbb{E}^{\mathbb{Q}^{N}}\Bigg[2\Big(\hat{g}_{1}(T) - y\Big)^{+} \mathbb{E}^{\mathbb{Q}^{N}}\Big[ \hat{g}_{2}(T)| \hat{g}_{1}(T)\Big]\Bigg] &= \sqrt{T} \bigg(\nu \rho + \frac{1}{2}\mathcolor{highlightNEW}{\alpha}\Big(\partial_F C(F_0)\Big) \bigg) \mathbb{E}^{\mathbb{Q}^{N}}\bigg[\Big(\hat{g}_{1}(T) - y\Big)^{+} \Big(\hat{g}^{2}_{1}(T) - 1\Big) \bigg] \\ &= \sqrt{T} \bigg(\nu \rho + \frac{1}{2}\mathcolor{highlightNEW}{\alpha}\Big(\partial_F C(F_0)\Big) \bigg) \int_{y}^{\infty} (s - 1)(s^2 - 1) \phi(s) ds \nonumber \\ & = \left(\nu \rho + \frac{1}{2}\mathcolor{highlightNEW}{\alpha}\Big(\partial_F C(F_0)\Big) \right) \phi(y) \sqrt{T} \end{align}\]

12.3 \(\epsilon^{2}\)-coefficient quadratic call option expansion↩︎

As we did in 11.3, we must calculate \[\mathbb{E}^{\mathbb{Q}^{N}}\Bigg[ 2 \Big(\hat{g}_1(T) - y\Big)^{+} \hat{g}_3(T) + I\Big(\hat{g}_1(T) > y\Big) \hat{g}^{2}_2(T) \Bigg]= A(T) + B(T)\] where \[\begin{align} A(T) &= 2 \mathbb{E}^{\mathbb{Q}^{N}}\Bigg[\Big(\hat{g}_1(T) - y\Big)^{+} \hat{g}_3(T) \Bigg],\\ B(T)&= \mathbb{E}^{\mathbb{Q}^{N}}\Bigg[I\Big(\hat{g}_1(T) > y\Big) \hat{g}^{2}_2(T) \Bigg]. \end{align}\] From 21 and using 9.1,9.2 and 9.3, we obtain that \[\begin{align} A(T) &= \frac{2}{3}\bigg(\frac{1}{2}\nu^2 \rho^2 + \mathcolor{highlightNEW}{\alpha}\rho \nu \Big(\partial_F C(F_0)\Big) + \mathcolor{highlightNEW}{\alpha}^{2} \Big(\partial_{F}C(F_0)\Big)^2 \bigg)y\phi(y)T \\ &= \bigg(\mathcolor{highlightNEW}{\alpha}\nu \Big(\partial_F C(F_0)\Big) + \frac{1}{2} \mathcolor{highlightNEW}{\alpha}^2 \Big(\partial^2_F C(F_0)\Big)C(F_0)\bigg) \bigg( \frac{2}{3}y\phi(y) + \bar{\Phi}(y) \bigg)T, \end{align}\] and \[B(T)= \frac{1}{12}\Bigg(\bigg(3\mathcolor{highlightNEW}{\alpha}^2 \Big(\partial_F C(F_0)\Big)^2+ \nu^{2} \rho^{2} \bigg) \bigg(2\bar{\Phi}(y)+\Big(y^3+y\Big)\phi(y) \bigg) + 4 \nu^2 \hat{\sigma}^2y \phi(y) + 6 \nu^2 \hat{\sigma}^2\bar{\Phi}(y) \Bigg)T .\]

---

1. Corresponding author, david.garcia.lorite@gmail.com↩︎