Most-likely-path in Asian option pricing under local volatility models

TL;DR

This paper develops a novel approximation method for Asian option prices under local volatility models using a most-likely path approach, combining path integrals, asymptotic analysis, and large deviation theory.

Contribution

It introduces a new framework for approximating Asian option prices via the most-likely path, integrating path integrals and asymptotics in local volatility models.

Findings

Derived a path-integral expression for Asian option prices.

Established the existence of the most-likely path using large deviation theory.

Provided asymptotic formulas for small-time option pricing and implied volatility.

Abstract

This article addresses the problem of approximating the price of options on discrete and continuous arithmetic average of the underlying, i.e. discretely and continuously monitored Asian options, in local volatility models. A path-integral-type expression for option prices is obtained using a Brownian bridge representation for the transition density between consecutive sampling times and a Laplace asymptotic formula. In the limit where the sampling time window approaches zero, the option price is found to be approximated by a constrained variational problem on paths in time-price space. We refer to the optimizing path as the most-likely path (MLP). Approximation for the implied normal volatility follows accordingly. The small-time asymptotics and the existence of the MLP are also recovered rigorously using large deviation theory.

Tables

Table 1. Table 1. Prices of Asian calls in the Black-Scholes model in the benchmark scenarios. The last four columns correspond to the approximate price from ( 4.15 ) (ALW), the third order approximation from [ 18 ] (FPP3), the precise evaluation in [ 27 ] (Linetsky), and the relative discrepancy of ALW to Linetsky.

Case	$S_0$	$K$	$r$	$\sigma$	$T$	ALW	FPP3	Linetsky	rel. discrp.
1	2	2	0.02	0.1	1	0.056042	0.055986	0.055986	0.10%
2	2	2	0.18	0.3	1	0.219607	0.218387	0.218387	0.56%
3	2	2	0.0125	0.25	2	0.172939	0.172267	0.172269	0.39%
4	1.9	2	0.05	0.5	1	0.195034	0.193164	0.193174	0.96%
5	2	2	0.05	0.5	1	0.248277	0.246406	0.246416	0.73%
6	2.1	2	0.05	0.5	1	0.308029	0.306210	0.306220	0.59%
7	2	2	0.05	0.5	2	0.355167	0.350040	0.350095	1.45%

Table 2. Table 2. Prices of Asian calls in the CIR model in the benchmark scenarios considered in [ 12 ] . The last four columns correspond to the approximate price from ( 4.15 ) (ALW), from [ 12 ] (DN), the third order approximation from [ 18 ] (FPP3), and the relative discrepancy of ALW to FPP3.

Case	$S_0$	$K$	$r$	$\sigma$	$T$	ALW	DN	FPP3	rel. discrp.
1	2	2	0.02	0.14	1	0.055591	0.0197	0.055562	0.05%
2	2	2	0.18	0.42	1	0.218521	0.2189	0.217874	0.30%
3	2	2	0.0125	0.35	2	0.171331	0.1725	0.170926	0.24%
4	1.9	2	0.05	0.69	1	0.191950	0.1902	0.190834	0.58%
5	2	2	0.05	0.72	1	0.252333	NA	0.251121	0.48%
6	2.1	2	0.05	0.72	1	0.309864	0.3098	0.308715	0.37%
7	2	2	0.05	0.71	2	0.356411	0.3339	0.353197	0.91%

Equations

$$d{S}_t=S_t{\sigma }_{\ell }(S_t,t)d{W}_t=a(S_t,t)d{W}_t{S}_0=s_0>0.$$

$$0\le a(s,t)\le C(1+|s|)\text{ for all t∈[0,T] and for all s∈R; }$$

$$|a(s,t)-a(s^{\prime },t)|\le C_R|s-s^{\prime }|.$$

$$x=\phi (s,t)={\int }_{s_0}^s\frac{d\xi }{a(\xi ,t)}\text{mbox}fors>0.$$

$$p(T,s_T|t,s_t)=\frac{g(T-t,\phi (s_T,T)-\phi (s_t,t))}{a(s_T,T)}{\tilde{\mathbb{E}}}_{\phi (s_t,t),\phi (s_T,T)}\left[e^{{\int }_t^{T}h(X_s,s)d{X}_s-\frac{1}{2}{\int }_t^T{h}^2(X_s,s)ds}\right]$$

$$p(T,s_T|t,s_t)$$

$$\Psi (s_t,s_T)={\tilde{\mathbb{E}}}_{\phi (s_t,t),\phi (s_T,T)}\left[e^{{\int }_t^{T}h(X_s,s)d{X}_s-\frac{1}{2}{\int }_t^T{h}^2(X_s,s)ds}\right].$$

$$\Psi (s_t,s_T)$$

$${\tilde{\mathbb{E}}}_{\phi (s_t,t),\phi (s_T,T)}\left[e^{-\frac{1}{2}{\int }_t^T{h}^2(X_s,s)+h_x(X_s,s)+2H_t(X_s,s)ds}\right].$$

$${\left(\frac{1}{n}\sum _{i=1}^n{S}_{t_i}-K\right)}^{+},$$

$$p(t_1,s_{t_1}|t_0,s_{t_0})p(t_2,s_{t_2}|t_1,s_{t_1})\cdots p(t_n,s_{t_n}|t_{n-1},s_{t_{n-1}})$$

$$=$$

$$C$$

$$=\iint {\left(\frac{1}{n}\sum _{i=1}^n{s}_{t_i}-K\right)}^{+}\prod _{i=1}^{n}p(t_i,s_{t_i}|t_{i-1},s_{t_{i-1}})d{s}_{t_1}\cdots d{s}_{t_n}$$

$$=\frac{1}{{\left(2\pi \Delta t\right)}^{\frac{n}{2}}}\iint {\left(\frac{1}{n}\sum _{i=1}^n{s}_{t_i}-K\right)}^{+}{e}^{-\frac{D(\mathbf{s},\mathbf{t})}{\Delta t}}W(\mathbf{s},\mathbf{t})d\mathbf{s},$$

$$D(\mathbf{s},\mathbf{t})$$

$${\int }_R{e}^{-\frac{\theta (x)}{\tau }}f(x)dx$$

$$\underset{\mathbf{s}}{\min }\frac{1}{2}\sum _{i=1}^n{\left|\phi (s_{t_i},t_i)-\phi (s_{t_{i-1}},t_{i-1})\right|}^2$$

$$\frac{1}{n}\sum _{i=1}^n{s}_{t_i}\ge K.$$

$$C$$

$$\mathbb{E}\left[{\left(\frac{1}{T}{\int }_0^T{S}_{t}dt-K\right)}^{+}\right]=\underset{n\to \infty }{\lim }\mathbb{E}\left[{\left(\frac{1}{n}\sum _{i=1}^n{S}_{t_i}-K\right)}^{+}\right]$$

$$\log C$$

$$+\log \left[\frac{\nabla W({\mathbf{s}}^{\ast },\mathbf{t})\cdot \nabla D({\mathbf{s}}^{\ast },\mathbf{t})}{\sqrt{\det {\partial }_{\mathbf{t}}^{2}D({\mathbf{s}}^{\ast },\mathbf{t})}|\nabla D({\mathbf{s}}^{\ast },\mathbf{t}){|}^2}+\frac{1}{2}\mathrm{tr}\left\{{\partial }_{\mathbf{t}}^{2}W({\mathbf{s}}^{\ast },\mathbf{t}){\left[{\partial }_{\mathbf{t}}^{2}D({\mathbf{s}}^{\ast },\mathbf{t})\right]}^{-1}\right\}+\mathcal{O}\left(\Delta t\right)\right],$$

$$\frac{1}{2\Delta t}\sum _{k=1}^n{\left|\phi (s_{t_i},t_i)-\phi (s_{t_{i-1}},t_{i-1})\right|}^2.$$

$$\phi (s_{t_i},t_i)-\phi (s_{t_{i-1}},t_{i-1})={\phi }_s(s_{t_{i-1}},t_{i-1})\Delta s_{t_i}+{\phi }_t(s_{t_{i-1}},t_{i-1})\Delta t+o(\Delta s_{t_i},\Delta t),$$

$$\underset{\Delta t\to 0}{\lim }\frac{1}{2\Delta t}\sum _{k=1}^n{\left|\phi (s_{t_i},t_i)-\phi (s_{t_{i-1}},t_{i-1})\right|}^2$$

$$=\underset{\Delta t\to 0}{\lim }\frac{1}{2}\sum _{k=1}^n{\left|\frac{\frac{\Delta s_i}{\Delta t}}{a(s_{i-1},t_{i-1})}\right|}^2\Delta t$$

$$=\frac{1}{2}{\int }_0^T{\left[\frac{\dot{s}(t)}{a(s(t),t)}\right]}^{2}dt.$$

$$\underset{s:t\mapsto s(t)}{\min }\frac{1}{2}{\int }_0^T{\left[\frac{\dot{s}(t)}{a(s(t),t)}\right]}^{2}dt$$

$$\frac{1}{T}{\int }_0^{T}s(t)dt=K,s(0)=s_0.$$

$$\frac{d}{dt}\left(\frac{\dot{s}}{a}\right)-\frac{a_t}{a^2}\dot{s}+\frac{\lambda }{T}a=0$$

$$s(0)=s_0,\dot{s}(T)=0,$$

$$L(s,\lambda )=\frac{1}{2}{\int }_0^T{\left[\frac{\dot{s}(t)}{a(s(t),t)}\right]}^{2}dt-\lambda \left(\frac{1}{T}{\int }_0^{T}s(t)dt-K\right),$$

$$0$$

Full text

1. Introduction

Asian options, also known as average price options, are among the most liquidly traded exotic options in commodities such as agriculture, energy, and fixed incomes markets. Nowadays average price options represent a high percentage of options on oil; some are directly on the futures contracts of oil, while others on spreads between two types of oil futures. Asian options are also commonly used as a risk management vehicle for, owing to its averaging feature, a) the underlying average price is more difficult to manipulate; b) the average price is less sensitive to abrupt shocks; and c) such options are cheaper than similar vanilla options.

There is a rich literature on the problem of pricing Asian options, in part due to the difficulty of finding analytical solutions even for simpler cases like the Black-Scholes model. We mention only the following few and refer the interested readers to the references therein. To our knowledge, [25] was the first published work tackling the problem of Asian option pricing in Black-Scholes model. As numerically pricing Asian options by Monte Carlo simulations is concerned, [25] also introduced a variance reduction technique by using the price of geometric average option, whose analytic form is readily obtained, as a control variate. See also [11] for further analysis and extension to portfolio options on the technique of conditioning. Further development and improvement of Monte Carlo schemes since then were followed up by [6] (for a direct method of estimating the Greeks of an Asian option by Monte Carlo simulation, see Section 4.2 on P.275), [36] (combining control variate and change of measure/likelihood ratio), and [23] (importance sampling), etc. Attempts to find closed or semi-closed form expressions for Asian option pricing in Black-Scholes models first appeared in the seminal work of [22]. Among other interesting results in [22], the Laplace transform of the Asian option price with respect to time to maturity is derived and has been known as the celebrated Geman-Yor formula. An extension of the Geman-Yor approach to a jump diffusion model can be found in [7]. However, numerical inversion of the Geman-Yor formula was shown slow and needed to be handled with care, see for example the discussions in [15] and [17]. Analytical approximation of the risk neutral density for the average price, which in turn yields an approximation of the Asian option, under Black-Scholes model dates back to the work of [35] and [31]. Both papers applied the Edgeworth expansion to the density of discretely monitored average price around lognormal distribution and obtained Black-Scholes type formula (up to correction terms) for the price of Asian option. Such analytical approximations are more appealing in practice than the Monte Carlo solutions because explicit expressions for the Greeks are readily accessible. Finally, PDE method is pursued by [32] and [37] (see also [38] for a more recent development) after an ingeniously chosen change of variable.

Literatures on the pricing of Asian options under more general dynamics for the underlying security such as local volatility or stochastic volatility models, as opposed to those in the Black-Scholes model, are comparatively little. Less ambitious approaches to the pricing of Asian options include arbitrage free bounds as in [1] and [32], and also approximative and asymptotic solutions such as [12], [18], and [28]. Approximations resorting to asymptotic expansions are mostly Itô-Taylor type expansion based as originated from the work of [26], see also [8] for a more recent development along this line. Despite being straightforward but tedious, such expansions usually require calculations up to third or fourth order in order to achieve satisfactory accuracy. Finally, though not directly related to the current paper, Asian option pricing under stochastic volatility models is discussed in [16] and [34].

In the current paper, we address the problem of approximating the price of options on the discrete arithmetic average, and its continuous-time limit, with the underlying following a local volatility model. For the discretely monitored Asian option, we assume that the average is over a set of equally spaced discrete time samples before expiry. The application of the Brownian bridge representation for the transition density (see Theorem 2.1 below) obtained in [40] between consecutive sampling time points leads to a “path-integral” type expression for the Asian option price, see (2.12). A direct application of a Laplace asymptotic formula (in this case high-dimensional, see Lemma 2.1) yields an approximation of the option price (see Theorem 2.2). In the limit where the sampling time window approaches zero, the leading order term (in small time to expiry) can be expressed as a constrained variational problem of finding an optimal path, referred to as the most-likely-path (MLP), in the time-price space. An approximation of the continuously monitored Asian option price is obtained once the variational problem is solved, see Definition 2.1 and Theorem 2.3. The MLP approximation coincides with a rigorous derivation of the leading order based on a recent extension of the Freidlin-Wentzell theorem for a large class of models, see Theorem 3.1. As for implied volatility, we opt to use the Bachelier model as benchmark rather than the Black-Scholes model because of the lack of analytical expression for Asian options in the Black-Scholes model. Such defined implied volatility in the European option case is sometimes referred to as the implied normal volatility in practice. By comparing corresponding expansions from the benchmark Bachelier model and from the local volatility model, we obtain, as the main result of the paper, the lowest order approximation of the implied normal volatility for Asian option in Theorem 2.5.

The paper is organized as follows. Section 2 lays out the model and provides derivations of the Laplace type approximation for discretely monitored Asian calls and the most-likely-path approximations for continuously monitored Asian call options. Section 3 is devoted to a rigorous derivation of the asymptotic behavior obtained in Section 2 based on large deviation theory. It was drawn to our attention that a rigorous proof of the most-likely-path approximation for the price of continuously monitored Asian option by the theory of large deviation has been done independently in [29]. Finally, we conclude by presenting numerical tests of the most-likely-path approximation.

Throughout the text, $(W_{t},t\geq 0)$ denotes the standard Brownian motion defined on the filtered probability space $(\Omega,(\mathcal{F}_{t})_{t\geq 0},\mathbb{P})$ satisfying the usual conditions. Dot will always refer to the partial derivative with respect to the time variable and prime to the space variable $s$.

2. Asian option pricing in local volatility models

In this section, we derive asymptotic expansions of the prices of both discretely and continuously monitored Asian calls in local volatility models. First, the most-likely-path approximation in the discrete case is obtained in Theorem 2.2 using a Brownian bridge representation for the transition density obtained in [40] (see Theorem 2.1 below) and a high-dimensional Laplace asymptotic formula (Lemma 2.1). Second, an expression for the most-likely-path approximation in the continuous case is derived in Theorem 2.3 by taking a formal limit of the approximation in the discrete case. This provides a heuristic derivation for the leading order of the most-likely-path approximation for continuously monitored Asian calls.

We assume that the underlying asset $S$ follows the local volatility model

[TABLE]

We will suppose that the diffusion function $a(s,t)$ is strictly positive (except possibly at $s=0$ where it could be [math]), and grows at most linearly in $s$: there exists $C>0$ such that

[TABLE]

and is locally Lipschitz: for every $R>0$, there exists $C_{R}$ such that for all $t\in[0,T]$ and for every $x,y\in\mathbb{R}$ with $|s|,|s^{\prime}|<R$

[TABLE]

These assumptions are sufficient to imply the existence and uniqueness of a strong solution to the SDE (2.1). Such solutions will also satisfy a large deviation principle based on the work in [9] as explained in Section 3.

Let $p(T,s_{T}|t,s_{t})$, $t<T$, be the transition probability density from $(t,s_{t})$ to $(T,s_{T})$ for the local volatility model (2.1). Consider the Lamperti transformation from $s$ to $x$

[TABLE]

Since $a$ is assumed locally Lipschitz and is strictly positive except at $s=0$, the transformation is well-defined for all $s>0$, except possibly at $s=0$. The following representation for the transition density $p$ is derived in [40].

Theorem 2.1**.**

Let $S=(S_{t},t\geq 0)$ be the diffusion process given by (2.1). Define the function $h$ by $h(x,t)=\varphi_{t}(s,t)-a_{s}(s,t)/2$, with $s=\varphi^{-1}(x,t)$, where $\varphi$ is the Lamperti transformation (2.4) and subindices refer to corresponding partial derivatives. Then the transition density $p$ of $S$ from $(t,s_{t})$ to $(T,s_{T})$ has the representation:

[TABLE]

where $g$ denotes the centered Gaussian density with variance $t$: $g(t,\xi)=\exp(-\xi^{2}/2t)/\sqrt{2\pi t}$ and $\tilde{\mathbb{E}}_{x,y}[\cdot]$ is the expectation under the Brownian bridge measure from $x$ to $y$.

Equivalently, if $H$ is an antiderivative of $h$ with respect to $x$, namely, $\partial_{x}H(x,t)=h(x,t)$, for all $x$ and $t$, then

[TABLE]

For notational simplicity, hereafter we shall denote the expectation term in (2.5) as

[TABLE]

With this notation, we also have for the term in (2.6)

[TABLE]

2.1. Small time asymptotic for discretely monitored Asian call

Assume an (arithmetic) Asian call is sampled discretely at the time points $t_{1}<t_{2}<\cdots<t_{n}$ with $t_{0}=0$ and $t_{n}=T$. In other words, the payoff of such an Asian call is

[TABLE]

where the time interval between the sampling points is assumed equal, i.e., $t_{i}-t_{i-1}=\Delta t=T/n$ for $i=1,\cdots,n$. By using the Brownian bridge representation (2.5) or (2.6), the joint density for $S_{t_{1}},\cdots,S_{t_{n}}$ can be written as

[TABLE]

In this notation, the price $C=C(s_{0},0;K,T)$ of a discretely monitored Asian call struck at $K$ can be written as

[TABLE]

where

[TABLE]

for $\mathbf{t}=(t_{1},\cdots,t_{n})$, $\mathbf{s}=(s_{t_{1}},\cdots,s_{t_{n}})$ and $d\mathbf{s}=ds_{t_{1}}\cdots ds_{t_{n}}$.

Hence, the evaluation of the Asian call price $C$ becomes the computation of the multidimensional integral (2.12). The following Laplace asymptotic formula will be applied to estimate the $n$-dimensional integral in (2.12) for small $\Delta t$.

Lemma 2.1**.**

*(Laplace asymptotic formula)

Let $R$ be a closed set in $\mathbb{R}^{n}$ with nonempty and smooth boundary $\partial R$. Suppose $\theta$ is a continuous function in $R$ and attains its minimum uniquely at $x^{*}\in\partial R$ and, given any $\epsilon>0$, there exists $\delta>0$ such that $\theta(x)\geq\theta(x^{*})+\delta$ for all $x\in R\setminus B_{\epsilon}(x^{*})$, where $B_{\epsilon}(x^{*})=\{x:|x-x^{*}|<\epsilon\}$ is the open ball of radius $\epsilon$ centered at $x^{*}$. Assume that $f$ is integrable in $R$, i.e., $\int_{R}|f(x)|dx<\infty$ and that $f$ vanishes identically in $R^{c}$ and on the boundary $\partial R$ but the inward normal directional derivative of $f$ at $x^{*}$ is nonzero. Then we have the asymptotic expansion as $\tau\to 0^{+}$*

[TABLE]

where $\partial_{\mathbf{t}}^{2}f(x^{*})$ and $\partial_{\mathbf{t}}^{2}\theta(x^{*})$ are the Hessian matrices of $f$ and $\theta$ respectively in the tangential direction to $R$ at $x^{*}$.

The proof of the lemma is standard and straightforward. See for instance Section 8.3 in [4].

Now we apply the Laplace asymptotic formula (2.1) to the multidimensional integral (2.12) by taking $\theta=D(\mathbf{s},\mathbf{t})$, $f=\left(\frac{1}{n}\sum_{i=1}^{n}s_{t_{i}}-K\right)W(\mathbf{s},\mathbf{t})$, and $R$ as the half space $R:=\left\{\mathbf{s}:1/n\sum_{i=1}^{n}s_{t_{i}}\geq K\right\}$. The crucial step in applying the Laplace asymptotic formula (2.1) is the determination of the minimum point of $D$ in the half space $R$, which boils down to solving the constrained optimization problem:

[TABLE]

subject to

[TABLE]

Remark 2.1**.**

We shall assume $s_{0}<K$ in the following for if $s_{0}\geq K$, the value of the constrained optimization problem (2.15):(2.16) is 0 since one can simply take $s_{t_{1}}=s_{t_{2}}=\cdots=s_{t_{n}}=s_{0}$. Then $s_{1}+\cdots+s_{n}=ns_{0}\geq nK$ and $\varphi(s_{t_{i}},t_{i})=\varphi(s_{t_{0}},t_{i})=0$ for all $1\leq i\leq n$. Thus the objective function in (2.15) attains its global minimum [math]. **

Remark 2.2**.**

We show in Section 6 Appendix I that, for $\Delta t$ small enough, the objective function in (2.15) is in fact convex, which in turn implies that the minimizer, if there exists any, is unique since the constraint is a linear inequality. **

We summarize the result for the price of a discretely monitored Asian call in Theorem 2.2 whose proof in the time homogeneous case is simply a direct application of the Laplace asymptotic formula (2.1) since the functions $D$ and $W$ are independent of $\mathbf{t}$. We remark that, modulo the exponential term, the result suggests the next order term is of order $3/2$ in $\Delta t$, which coincides with the order in the case for European options, see for example Theorem 2.3 in [19].

Theorem 2.2**.**

*(Discrete monitored Asian option)

The price $C=C(s_{0},0;K,T)$ of a discretely monitored Asian call struck at $K$ with $K>s_{0}$ and expiry time $T$ has the following asymptotic expansion as $T\to 0^{+}$, for fixed $n$,*

[TABLE]

where $\mathbf{s}^{*}=(s_{t_{1}}^{*},\cdots,s_{t_{n}}^{*})$ is the minimizer of the minimization problem (2.15) subject to the constraint $1/n\sum_{i=1}^{n}s_{t_{i}}\geq K$.

2.2. Continuously monitored Asian call and the most-likely-path

The approximate price of an Asian call obtained in Theorem 2.2 is subject to solving a high-dimensional constrained optimization problem which is daunting in general. However, in the limit as $\Delta t$ approaches zero, the optimization problem converges to a variational problem to which, in certain cases such as Black-Scholes and CIR, the associated Euler-Lagrange equations have closed form solution. The heuristic computation is given in this section. The rigorous derivation for the leading order using large deviation is done in Section 3.

Let $\{0=t_{0}<t_{1}<\cdots<t_{n}=T\}$ be a partition of the interval $[0,T]$ with $t_{i}-t_{i-1}=\Delta t:=T/n$, for $i=1,\cdots,n$. Then, the price of a continuously monitored Asian call can be written as the limit of the prices of discretely monitored Asian calls as $n\to\infty$. Precisely,

[TABLE]

by applying the Lebesgue dominated convergence theorem. Hence, to the lowest order, it is natural to approximate the price of a continuously monitored Asian call by taking the limit of the approximate price of discretely monitored Asian call in (2.17) as $\Delta t\to 0$. To be specific, rewrite the logarithm of (2.17) as

[TABLE]

where we recall that $\mathbf{s}^{*}$ is the minimizer of (2.15). Notice that the first term in the last expression is dominant as $\Delta t\to 0$.

To determine the limit as $\Delta t\to 0$ of the dominant term, we bring $\Delta t$ back to the objective function (2.15) as

[TABLE]

Note that, since

[TABLE]

where $\Delta s_{t_{i}}=s_{t_{i}}-s_{t_{i-1}}$, we have

[TABLE]

Therefore, in the limit as $\Delta t$ approaches zero, the optimization problem (2.15) becomes the following variational problem111An equivalent formulation of the same problem in the Black-Scholes case was obtained in [10] (see (1.4)) for the analysis of fundamental solution. on the space of paths $s:t\mapsto s(t)$

[TABLE]

subject to

[TABLE]

Definition 2.1**.**

---

*The optimal path of the variational problem (2.23):(2.24) is referred to as the most-likely-path (MLP) for the Asian call struck at $K$. **

In view of the above, one expects that, at the heuristic level and for the leading order, the logarithm of the price $C$ of an out-of-the-money Asian call is approximately given by the solution to the constrained variational problem (2.23):(2.24), which is equivalent to determining the most-likely-path. To determine the MLP, the Euler-Lagrange equation associated with the constrained variational problem (2.23):(2.24) along with proper boundary conditions is derived as follows.

Lemma 2.2**.**

The optimal path of the constrained variational problem (2.23):(2.24) satisfies the Euler-Lagrange equation

[TABLE]

with boundary conditions

[TABLE]

where $\lambda$ is chosen such that $1/T\int_{0}^{T}s(t)dt=K$.

Proof.

We first rewrite the constrained variational problem (2.23):(2.24) in the Lagrangian form

[TABLE]

where $\lambda$ is the Lagrange multiplier. Let $f:[0,T]\mapsto\mathbb{R}$ be a perturbation around the optimal path $s(t)$ with $f(0)=0$. The first order criterion of optimality yields

[TABLE]

By applying integration by parts to the second integral and noting that $f(0)=0$, the last equality becomes

[TABLE]

Finally, since $f$ is arbitrary and $a(s,t)>0$, we obtain the Euler-Lagrange equation

[TABLE]

which simplifies to

[TABLE]

with boundary conditions $s(0)=s_{0},\dot{s}(T)=0$. ∎

We summarize the final result in the following theorem. A rigorous proof based on large deviation theory is postponed to Theorem 3.1 in Section 3.

Theorem 2.3**.**

*(Log price of a continuously monitored Asian call)

The price $C(S_{t},t;K,T)$ at time $t$ of a continuously monitored out-of-the-money Asian call struck at $K>s_{0}$ and expiry time $T$, is approximately equal to*

[TABLE]

where $s:t\mapsto s(t)$ is the solution to the constrained variational problem (2.23):(2.24).

We conclude the section by deriving closed form expressions for the most-likely-path in the Bachelier, the Black-Scholes, and the Cox-Ingersoll-Ross (CIR) models by solving their corresponding boundary value problems (2.25):(2.26). However, for more general cases we will have to resort to an iteration scheme for numerical computations of the most-likely-path. See Section 4 for a numerical scheme.

Example 2.1**.**

*(Bachelier model)

*In the Bachelier model, $a(s,t)=\sigma$, a constant. The Euler-Lagrange equation (2.25) reduces to

[TABLE]

whose general solution is $s(t)=-\lambda\sigma^{2}t^{2}/(2T)+c_{1}t+c_{0}$. From the boundary conditions (2.26) together with $1/T\int_{0}^{T}s(t)dt=K$, we conclude that $\lambda=3(K-s_{0})/(\sigma^{2}T^{2})$, $c_{1}=3(K-s_{0})/T$, and $c_{0}=s_{0}$. Therefore, the most likely path for Asian call in the Bachelier model is a downward parabola in the $t$-$s$ plane given by

[TABLE]

We remark that the most-likely-path in this case does not depend on $\sigma$. Also, (2.35) for the Bachelier model becomes

[TABLE]

Thus it captures the exponential decay of an out-of-the-money Asian call in the Bachelier model, see (2.64) below. **

Example 2.2**.**

*(Black-Scholes model)

*In the Black-Scholes model, $a(s,t)=\sigma s$ with $\sigma$ being a constant. The Euler-Lagrange equation (2.25) reads

[TABLE]

to which the general solutions are

[TABLE]

and

[TABLE]

where $c_{1}$ and $c_{2}$ are (to be determined) constants. The condition $\dot{s}(T)=0$ implies that $c_{2}=T$ in either case and the conditions $s(0)=s_{0}$, $\int_{0}^{T}s(t)dt=TK$ imply respectively that

[TABLE]

and

[TABLE]

By dividing the two equations in (2.42) and rearranging terms, $c_{1}$ is given by the solution to the equation

[TABLE]

or $c_{1}=f^{-1}\left(s_{0}/K\right)/T$, where $f(x)=x/\sinh x$. On the other hand, from (2.43) we obtain that

[TABLE]

or $c_{1}=g^{-1}\left(K/s_{0}\right)/T$, where $g(x)=\sin x/x$ for $x\in[0,\pi]$.

Hence, the most-likely-path $s(t)$ for $0\leq t\leq T$ in the Black-Scholes model after simplification becomes

[TABLE]

where $c_{1}=f^{-1}\left(s_{0}/K\right)/T$ for $s_{0}<K$; whereas for $s_{0}>K$

[TABLE]

and $c_{1}=g^{-1}\left(K/s_{0}\right)/T$. Moreover, for $s_{0}<K$,

[TABLE]

and for $s_{0}>K$

[TABLE]

Example 2.3**.**

*(Cox-Ingersoll-Ross model)

*For the CIR model, $a(s,t)=\sigma\sqrt{s}$, where $\sigma$ is a constant. The Euler-Lagrange equation (2.25) becomes

[TABLE]

for some constant $\phi$ with general solution given by

[TABLE]

if $\phi>0$ and

[TABLE]

if $\phi<0$. The boundary conditions (2.26) imply respectively that

[TABLE]

Thus

[TABLE]

The parameter $\phi$ is determined by the solution to the equation

[TABLE]

if $s_{0}<K$ and by

[TABLE]

if $s_{0}>K$. Finally, we have, subject to the determination of $\phi$, that

[TABLE]

for $s_{0}<K$ and

[TABLE]

if $s_{0}>K$. **

2.3. Implied normal volatility for Asian options

As far as implied volatility is concerned, we opt to use the Bachelier model as benchmark rather than the conventional Black-Scholes model because of the lack of analytical expression for Asian options in the Black-Scholes model. Such defined implied volatility in the European counterpart is sometimes used and referred to as the implied normal volatility in practice. For European calls, this approximation is good whenever $\sigma\sqrt{T}$ is small, at least for at-the-money options, see for example [33]. We expect the same should hold for Asian options, though the lack of analyticity makes it hard to check. In particular, note that the approximation might be problematic for small strike prices as the ratio $S/K$ will be large.

Recall that, by straightforward calculations, under the Bachelier model

[TABLE]

where $\sigma_{b}$ is a constant, the price of a continuously monitored Asian call option struck at $K$ with expiry $T$ has the closed form expression

[TABLE]

since the average price $1/T\int_{0}^{T}S_{t}dt$ is normally distributed with mean $s_{0}$ and variance $\sigma_{b}^{2}T/3$. In (2.60), $N(\cdot)$ denotes the cumulative distribution function of standard normal distribution. On the other hand, for a discretely monitored Asian call struck at $K$ with expiry $T$ in the Bachelier model, since $S_{t_{1}}+\cdots+S_{t_{n}}$ is normally distributed with

[TABLE]

its price $C_{b}^{d}$ is given by

[TABLE]

where $A_{n}:=6n^{2}/(n+1)(2n+1)$. Apparently, $A_{n}\to 3$ as $n\to\infty$.

Regarding the small time expansion, we remark that by using the asymptotic expansion $N(x)=N^{\prime}(x)[-1/x+1/x^{3}+\mathcal{O}(1/x^{5})]$ as $x\to-\infty$, we have

[TABLE]

as $t\to T^{-}$. Similarly, in the discrete case as $t\to T$,

[TABLE]

We remark that once the small time asymptotic for the price of an out-of-the-money call is established on the model side, it is a common practice to derive the small time asymptotic of implied volatility thereby. To that end, the following expansion for the implied normal volatility in terms of the call price given in [24] will be helpful.

Proposition 2.1**.**

For a fixed strike out-of-the-money call with time to expiry $T$, let $C=C(T)$ be the price of the European call regarded as a function of $T$. Then as time to expiry $T\to 0$, the implied normal volatility $\sigma_{N}$ has the asymptotic

[TABLE]

as $T\to 0$.

2.3.1. Implied normal volatility for discretely monitored Asian option

The implied normal volatility for a discretely monitored Asian call is defined by solving the following equation for $\sigma_{b}$

[TABLE]

where $A_{n}:=6n^{2}/(n+1)(2n+1)$ and $C^{d}$ is the price of a discretely monitored Asian call obtained by the model or from the market. Obviously, among other parameters, such defined $\sigma_{b}$ depends on $K$ and $T$.

To derive an asymptotic expansion for the implied normal volatility $\sigma_{b}$ defined in (2.67) in small time, the idea, as in [19], is to compare the corresponding terms in the expansions on both side of (2.67). The lowest order term is thus obtained by matching the exponential terms on both side of (2.67). Precisely, recall the small time asymptotic of the price of a discrete Asian call from Theorem 2.2

[TABLE]

where $\Delta t=T/n$. By matching the exponential terms in the above expression and in (2.65), we have

[TABLE]

Hence, the lowest (zeroth) order approximation of the implied normal volatility is given by

[TABLE]

where recall that $(s_{t_{1}}^{*},\cdots,s_{t_{n}}^{*})$ is the solution to the $n$-dimensional constrained optimization problem (2.15):(2.16). We summarize the result in the following theorem.

Theorem 2.4**.**

*(Implied normal volatility asymptotic for discrete Asian call)

For a discretely monitored out-of-the-money Asian call struck at $K$, i.e., $s_{0}<K$, in which the underlying is driven by the local volatility model (2.1), the implied normal volatility $\sigma_{b}^{d}$ defined as in (2.67) has the asymptotic expansion as $T\to 0^{+}$*

[TABLE]

where

[TABLE]

with $A_{n}=6n^{2}/(n+1)(2n+1)$ and $(s_{t_{1}}^{*},\cdots,s_{t_{n}}^{*})$ the solution to the constrained optimization problem (2.15):(2.16).

2.3.2. Implied normal volatility for continuously monitored Asian option

For continuously monitored Asian calls, the implied normal volatility is thus defined by solving the following equation for $\sigma_{b}$

[TABLE]

where $C$ is the price of an Asian call obtained by the model or from the market. Notice that the right hand side of (2.72) can be regarded as the price of a European call option in the Bachelier model but with one third of $\sigma_{b}^{2}$ as the variance (volatility squared) parameter. In other words, in the Bachelier world, the price of an Asian call is equal to the price of its European counterpart but with only one third of variance. By combining the asymptotics (2.35) for the call price and (2.66) for the implied normal volatility, a small time asymptotic for the implied normal volatility of an out-of-the-money Asian call is established. Alternatively, we may also obtain the same approximation by straightforwardly taking the limit as $n\to\infty$ of $\sigma^{d}_{b,0}$ in Theorem 2.4. We summarize the result in the following theorem.

Theorem 2.5**.**

*(Implied normal volatility asymptotic for continuous Asian call)

For a continuously monitored out-of-the-money Asian call struck at $K$, i.e., $s_{0}<K$, in which the underlying is driven by the local volatility model (2.1), the implied normal volatility $\sigma_{b}$ defined in (2.72) has the asymptotic expansion as $T\to 0^{+}$*

[TABLE]

where

[TABLE]

and $\tilde{s}(t)$ is the most-likely-path for Asian option determined by solving the variational problem (2.23):(2.24).

Proof.

Recall from (2.66) that the implied normal volatility $\sigma$ of a European option has the asymptotic

[TABLE]

Thus, by substituting $\log C(T)$ with (2.35) and using the fact that Asian variance equals one third of its European counterpart in Bachelier world, we obtain for implied normal volatility $\sigma_{b}$ of an Asian call

[TABLE]

where we used (2.35) in the last equality. Finally, the result is obtained by rearranging terms. ∎

Example 2.4**.**

*(Time-dependent Bachelier model) *Finally, we consider a time-dependent Bachelier model in which

[TABLE]

Note that in this case $a(s,t)=\sigma\theta(t)$, the Euler-Lagrange equation (2.25) reduces to

[TABLE]

Integrating the last equation subject to the condition $\dot{s}(T)=0$ gives

[TABLE]

for some constant $\lambda$ set using the condition

[TABLE]

This gives

[TABLE]

Also, from (2.82) one easily obtain

[TABLE]

This gives

[TABLE]

On the other hand, we have $\bar{S}_{T}=1/T\,\int_{0}^{T}S_{t}\,dt$ where

[TABLE]

A straightforward computation then yields

[TABLE]

so the exact Bachelier implied volatility in the time-dependent Bachelier case is given by

[TABLE]

Thus, the most-likely-path approximation is exact in the time-dependent Bachelier case. **

2.4. Approximation of Greeks

The approximation of implied volatility in Theorem 2.5 is also applicable for approximations of Greeks. For example, we may calculate the delta as follows. For notational simplicity, denote by $v_{b}:=\sigma_{b}^{2}T$. Suppressing and holding the other parameters fixed, since the Bachelier implied volatility is defined through

[TABLE]

the delta $\Delta$ satisfies

[TABLE]

where $C_{b}$ is the function defined in (2.72). Note that by straightforward calculations we have

[TABLE]

Thus, with $v_{b}\approx v_{b,0}:=\sigma_{b,0}^{2}T$, it follows that

[TABLE]

The expressions on the right hand side of (2.93) can be calculated easily except the last term which can be calculated as follows. Recall that

[TABLE]

we have

[TABLE]

The integral in (2.95) in the general case needs to be evaluated numerically. Similarly, since the gamma $\Gamma$ satisfies

[TABLE]

an approximation of $\Gamma$ is given by substituting $v_{b}$ with $v_{b,0}$, subject to numerically evaluations of the derivatives $\partial_{s}v_{b,0}$ and $\partial_{s}^{2}v_{b,0}$. However, we remark that in the Black-Scholes (Example 2.2) and the CIR (Example 2.3) models, since the $\sigma_{b,0}$’s have closed form expressions, closed expressions for Greeks are available. We refer the reader to [30] for more detailed discussions on Greeks of Asian options in the Black-Scholes model.

3. Large deviation principle

In this section, we prove Theorem 3.1 which is a large deviation reformulation of Theorem 2.3 for continuously monitored options. Finer tools would be needed to go beyond the leading order of the asymptotic expansion in Theorem 2.2 suggested by the discrete approximation. However, in particular cases, these terms are given by a Girsanov change of measure. We will illustrate this using the Bachelier model at the end of the section. A more refined expansion for more general processes does not appear tractable with current methods.

Theorem 3.1** (Large deviation for the log-price of continuous monitored Asian call).**

The price $C(s_{0},0;K,T)$ at time $t=0$ of a continuously monitored out-of-the-money Asian call struck at $K>s_{0}$ with expiry time $T$ for the price process $(S_{t},t\geq 0)$ of (2.1) admits the following expansion in $T$

[TABLE]

where

[TABLE]

Moreover, if $f(t)=x+\int_{0}^{t}a(f(s),s)g(s)ds$ for some $g\in L^{2}([0,T])$ then

[TABLE]

and $I(f)$ is $\infty$ otherwise.

The idea to prove the theorem is to treat $\varepsilon=T$ as a small parameter and expand around $\varepsilon=0$. We do a simple time-change $t=uT$ and write for simplicity

[TABLE]

Note that with this notation:

[TABLE]

By the scaling property of Brownian motion, the process $(S^{\varepsilon}_{u},u\in[0,1])$ satisfies the SDE

[TABLE]

where $a_{\varepsilon}(x,u):=a(x,u\varepsilon)$. Since $a(x,\cdot)$ is assumed to be continuous uniformly in $x$, we can write

[TABLE]

A large deviation principle for the process means that there exists a rate function $I:\mathcal{C}([0,T])\to[0,\infty]$ that is lower semi-continuous and has compact level sets such that for any Borel subset $A$ of paths in $\mathcal{C}([0,1])$,

[TABLE]

where int($A$) denotes the interior of $A$ and cl($A$) is its closure. Roughly speaking, a large deviation principle quantifies the probability of atypical path at the exponential scale with the help of the rate function. The proof of Theorem 3.1 is based on an extension of the Freidlin-Wentzell theorem, see e.g. [14]. The standard statement of the theorem holds for time-homogeneous drift and volatility under Lipschitz and boundedness assumptions. Here we will use a recent result of [9]. Note that the locally Lipschitz condition (2.3) is not satisfied in the CIR model. However, in this case, the volatility is time-homogeneous and a weaker assumption is needed for a large deviation to hold. This is the content of Theorem 4 in [9].

Theorem 3.2** (Theorem 2, Theorem 4 and Example 1 in [9]).**

The family of diffusions (3.6) satisfy a large deviation principle with rate function

[TABLE]

whenever the set $\{g\in L^{2}([0,1]):f(t)=x+\int_{0}^{t}a(f(s),0)g(s)ds\}$ is non-empty and $I(f)$ is $\infty$ otherwise.

Formally, it is good to think of the set $\{g\in L^{2}([0,T]):f(t)=x+\int_{0}^{t}a(f(s),0)g(s)ds\}$ as the set of “white-noise paths” $t\mapsto\dot{W}_{t}$. The map $g\mapsto f$ where $f$ is the solution of $f(t)=x+\int_{0}^{t}\bar{a}(f(s),s)g(s)ds$ can then be thought of as the map sending an underlying Brownian path to the corresponding diffusion path. In the case where the volatility is non-zero for the path $f$, the map can be inverted and the rate function reduces to the simplest case

[TABLE]

This is certainly the case for geometric Brownian motion and the CIR model. The reader is referred to Theorem 1 in [9] and Proposition 3.11 in [3] for general sufficient conditions for (3.10) to hold.

Proof of Theorem 3.1.

Observe that for any random variable $X$ and $K>0$, we have the identity

[TABLE]

For simplicity, write $T$ for the functional $T:\mathcal{C}([0,1])\to\mathbb{R}$ with $T(f)=\int_{0}^{1}f(u)du$. Note that $T$ is continuous on $\mathcal{C}([0,1])$ equipped with the topology of uniform convergence. By the contraction principle (see e.g. [14]) and Theorem 3.2, the family of random variables $(T(S^{\varepsilon}),\varepsilon>0)$ satisfies a large deviation principle with rate function $J:\mathbb{R}\to[0,\infty]$ with

[TABLE]

By (3.11), we have

[TABLE]

On one hand, we have for any $M>K$

[TABLE]

For $M$ large enough (possibly dependent on $\varepsilon$), the term in the second parenthesis is smaller than $2$. Pick $M=\varepsilon^{-1}$. For $\varepsilon$ small enough we thus have

[TABLE]

This proves that

[TABLE]

since $(T(S^{\varepsilon}),\varepsilon>0)$ satisfies a large deviation with rate function $J$.

On the other hand, for $\delta>0$

[TABLE]

Therefore, for any $\delta>0$, the above with the use of (3.8) becomes

[TABLE]

It remains to show that $\lim_{\delta\to 0}J(K+\delta)=J(K)$. First, notice that $\liminf_{\delta\to 0}J(K+\delta)\geq J(K)$ since $J$ is lower semi-continuous (it is a rate function). So it suffices to show that $\limsup_{\delta\to 0}J(K+\delta)\leq J(K)$. By definition, $J(K+\delta)=\inf\{I(x):x\in\mathcal{C}([0,1])\text{ and }\int_{0}^{1}x(u)du=K+\delta\}$. Pick a sequence $(y_{n})\in\mathcal{C}([0,1])$ such that $I(y_{n})\to J(K)$ and $\int_{0}^{1}y_{n}(u)du=K$. By definition of the infimum, this sequence can be picked such that $I(y_{n})<J(K)-1/n$. Pick $z$ a differentiable function on $[0,1]$ such that $\int z=1$ and $z(0)=0$. Then $\int(y_{n}(u)+\delta z(u))du=K+\delta$. Moreover

[TABLE]

It is easy to check that for fixed $n$, $\lim_{\delta\to 0}I(y_{n}+\delta z)=I(y_{n})$. Therefore

[TABLE]

Since $n$ is arbitrary, this proves the claim. ∎

The asymptotic expansion (2.35) suggests that the lower order corrections to (3.1) are much smaller than the dominant term in $T^{-1}$. We expect more precisely that

[TABLE]

This can be verified rigorously for the Bachelier model (Example 1) using a Girsanov change of measure designed to tilt towards the most-likely path minimizing the rate function $J$ for a given $K$. The correction is expected to be of the same order for other models as well as for other types of options.

For the Bachelier model, the diffusion on $[0,1]$ is $dS_{u}^{\varepsilon}=\sqrt{\varepsilon}dW_{u}$ and the most-likely path (2.37) when written as a path on $[0,1]$ is

[TABLE]

By (3.11), the price of the call option becomes

[TABLE]

Write $\mathbb{Q}$ for the measure with $d\mathbb{Q}/d\mathbb{P}=\exp\left(\varepsilon^{-1/2}\int_{0}^{1}\dot{s}(u)dW_{u}-\frac{\varepsilon^{-1}}{2}\int_{0}^{1}\dot{s}(u)^{2}du\right)$. Under $\mathbb{Q}$, the process $(W_{u},u\in[0,1])$ is a Brownian motion with drift $\varepsilon^{-1/2}\dot{s}(u)$. With this notation, (3.23) becomes

[TABLE]

where $(\widetilde{W}_{u},u\in[0,1])$ is a standard Brownian motion under $\mathbb{Q}$ and we use the fact that $\int_{0}^{1}s(u)du=K$ by definition of the most likely path. By doing the change of variable $y=x-K$, this reduces to

[TABLE]

The first term is $e^{-\varepsilon^{-1}J(K)}$ and gives the first order. To evaluate the second term, it is convenient to first integrate $x$ to get

[TABLE]

Note that $\dot{s}(u)=3(K-s_{0})(1-u)$, therefore

[TABLE]

Write $X$ for the random variable $\int_{0}^{1}\widetilde{W}_{u}du$ which is Gaussian with mean [math] and variance $1/3$. We have

[TABLE]

The integral is of order $1$. We conclude that

[TABLE]

4. Numerical tests

From Theorem 2.5, we have the following approximate formula for Asian Bachelier implied volatility:

[TABLE]

where $\tilde{s}(t)$ is the most-likely-path for an Asian option determined by solving the variational problem

[TABLE]

subject to

[TABLE]

We now proceed to test numerically the approximate implied volatility formula (4.1) for various definitions of the local volatility function $a(s,t)=s\,\sigma_{\ell}(s,t)$ in the local volatility model (2.1). In each case, to evaluate (4.1), we need to compute the most-likely-path $\tilde{s}(t)$. To this end, we exploit the following iteration scheme.

Lemma 4.1**.**

The most-likely path $\tilde{s}(t)$ satisfies the recursive formula

[TABLE]

where

[TABLE]

Proof.

From Lemma 2.2, $\tilde{s}(t)$ satisfies the Euler-Lagrange equation (2.25) which we reiterate in the following for convenience

[TABLE]

with boundary conditions

[TABLE]

where $\lambda$ is chosen such that $1/T\int_{0}^{T}s(t)dt=K$.

For ease of notation, define $\tilde{a}(t):=a\left(\tilde{s}(t),t\right)$ and $f(s,t)=a_{t}/a=\partial_{t}\log a$. Also, $\tilde{f}(t)=f(\tilde{s}(t),t)$. Applying the integrating factor $\exp\left(-\int_{0}^{t}\tilde{f}(u)du\right)$ and integrating (4.5) with the boundary condition $\dot{s}(T)=0$ gives

[TABLE]

It follows that

[TABLE]

Rearranging and integrating again gives

[TABLE]

where $I(t)=\int_{0}^{t}\int_{r}^{T}\,\tilde{a}(r)\tilde{a}(u)\exp\left(-\int_{t}^{u}\tilde{f}(v)dv\right)\,dudr$. Now apply the boundary condition $1/T\int_{0}^{T}s(t)dt=K$ to get

[TABLE]

and the result follows. ∎

Lemma 4.1 leads to an efficient fixed-point algorithm for solving for the most-likely-path. The natural choice of first guess is the Bachelier most-likely-path (2.37):

[TABLE]

The resulting algorithm typically converges sufficiently after three or four iterations.

4.1. A time-dependent CIR model

We consider the model

[TABLE]

with $S_{0}=1$, $\sigma=0.2$ and $\lambda=1$. Thus

[TABLE]

Though we computed a quasi-closed-form for the most-likely-path for the time-homogeneous case in Example 2.3, in this time-inhomogeneous case, we choose to compare our approximate implied volatility formula (4.1) evaluated using the fixed-point iteration algorithm against Monte Carlo simulations generated with 1000 time steps and 2 million sample paths. The results are shown in Figure 4.1. We remark that the most-likely-path approximation in this example slightly underestimates the normal implied volatility for Asian option inferred from simulation.

4.2. Time-dependent quadratic local volatility

Next we consider the following quadratic local volatility model:

[TABLE]

with $\sigma=0.2$, $\psi=-0.5$, $\gamma=0.1$, and $\lambda=1$. We remark that though in this example the function $a$ grows quadratically to infinity as $|x|\to\infty$ which violates the linear growth condition (2.2) required for the theoretical argument, we did the numerical experiment for testing the applicability of the most-likely-path methodology.

Though a closed-form solution for European options with these parameters is given in [2], we again resort to Monte Carlo simulation to estimate the value of Asian options in this model. Likewise, simulations are generated with 1000 time steps and 2 million sample paths. The results are shown in Figure 4.2. Similarly, the most-likely-path approximation in this example also underestimates the normal implied volatility for Asian option inferred from simulation.

4.3. Benchmark scenarios in Black-Scholes and CIR models

In this subsection, for the Black-Scholes and CIR models, we compare the most-likely-path approximation to a few existing approximations with the benchmark scenarios proposed in [21] and [17] that were commonly used in the literatures on Asian option pricing, see for instance [12], [15], [18], [27], [39].

With the approximation of Bachelier implied volatility $\sigma_{b,0}$ given in Theorem 2.5, we may calculate the approximate price of the Asian call struck at $K$ and expired at $T$ by

[TABLE]

with

[TABLE]

and $\mu=r-q$. Note that (4.15) is indeed the price of the Asian call struck at $K$ and expired at $T$ under the Bachelier model

[TABLE]

with risk free rate $r$ and dividend rate $q$, both assumed constant. We remark that, since many of the benchmark scenarios are of ATM, the ATM approximate implied volatility is obtained by taking the limit of $\sigma_{b,0}$ as $K$ approaches $S_{0}$. Explicitly, the limits in the Black-Scholes and CIR cases are given by

[TABLE]

Table 1 exhibits the numerical results for the asymptotic approximation for the Asian options obtained from (4.15) for the scenarios considered in [17]. Based on the results in [27], the relative discrepancies of the approximate prices are less than 1.5% in all the seven benchmark scenarios and within 1% for options expiring in a year. Table 2 shows the numerical tests for Asian option pricing in the CIR model for the seven scenarios proposed in [12]. Data is quoted from Table 5 in [18]. Based on the results in [18], the relative discrepancies of the approximate prices are less than 1% in all the seven benchmark scenarios and within 0.6% for options expiring in a year. We may thus arguably conclude that the short expiry approximation results of the current paper provide a reasonable approximation for Asian option prices in these models.

5. Conclusion

We have derived a small time asymptotic of the price of discretely monitored Asian options up to first order. The most-likely-path approximation for continuously monitored Asian calls has been derived heuristically by taking the limit of the corresponding term from discretely monitored case and proved rigorously using the theory of large deviation. Numerical experiments in both time-dependent CIR and time-dependent quadratic local volatility models showed that there is room to improve the performance of the most-likely-path approximation. One possibility is to include higher order terms. Large deviation theory generically provides no insights beyond the term of exponential decay. On the other hand, it is conceivable that taking the limit of corresponding terms from the discretely monitored case may yield tractable expressions for numerical evaluations. The calculations of higher order terms are considerably more involved and henceforth were left to further work.

Acknowledgement

We thank the anonymous referees for their careful reading and valuable comments. We would like to express our gratitude to Jim Gatheral for numerical results and Dan Pirjol for valuable and interesting discussions, and sending us an early version of their papers. L.-P. A. is supported by NSF CAREER 1653602, NSF grant DMS-1513441, and a Eugene M. Lang Junior Faculty Research Fellowship. N.-L. L. is supported by JSPS KAKENHI Grant Numbers 25780213. T.-H. W. is partially supported by the Natural Science Foundation of China grant 11601018.

6. Appendix - Convexity of the constrained optimization problem

We analyze the convexity of the constrained optimization problem (2.15):(2.16) in this appendix. Recall the Lamperti transformation $\varphi(s,t):=\int_{s_{0}}^{s}1/a(\xi,t)d\xi$ and the objective function $D$

[TABLE]

For notational simplicity, we shall write the function $a(\cdot,t_{i})$ as $a_{i}(\cdot)$ and similarly $\varphi(\cdot,t_{i})$ as $\varphi_{i}(\cdot)$. By straightforward calculations, the second partial derivatives of $D$ are given by

[TABLE]

We decompose the Hessian matrix $H=\left[\partial_{s_{i}}\partial_{s_{j}}D\right]$ as $H=H^{1}+H^{2}$, where $H^{1}$ is the symmetric tridiagonal matrix with diagonal entries given by

[TABLE]

and off-diagonal entries by

[TABLE]

$H^{2}$ is a diagonal matrix with diagonal entries given by

[TABLE]

We claim that $H^{1}$ is positively definite. Let $\hat{H}^{1}$ be the square submatrix of $H^{1}$ consisting of the first $n-1$ rows and columns of $H^{1}$. We partition $H^{1}$ as

[TABLE]

where $c=\left[0\,\cdots\,0\;H^{1}_{n-1,n}\right]$ is an $(n-1)$ row vector. Note that, by induction, one can show that the principle minors of $\hat{H}^{1}$ are $(k+1)\prod_{j=1}^{k}a_{j}^{-2}(s_{j})>0$, for $k=1,\cdots,n-1$. By applying the identity for determinant, if $A$ is invertible,

[TABLE]

we calculate the determinant of $H^{1}$ as

[TABLE]

where $\hat{\hat{H}}^{1}$ is the square submatrix of $H^{1}$ by deleting the last two rows and columns of $H^{1}$. Thus, by Sylvester’s criterion $H^{1}$ is positive definite.

Finally, since

[TABLE]

where $\Delta s_{k}=s_{k}-s_{k-1}$, the entries of $H^{2}$ are small if the $s_{i}$’s are close to each other and $\Delta t$ is small. In this case, $H$ is positive definite. Hence, the objective function $D$ is convex.