Pricing Exotic Derivatives for Cryptocurrency Assets—A Monte Carlo Perspective

Mesias Alfeus¹, Shiam Kannan²

¹Department of Statistics and Actuarial Science, Stellenbosch University, Stellenbosch, South Africa.
²Dyson School of Applied Economics and Management, Cornell University, Ithaca, USA.

DOI: 10.4236/jmf.2021.114033                

Abstract

In the current paper, we develop a methodology to price lookback options for cryptocurrencies. We propose a discreetly monitored window average lookback option, whose monitoring frequencies are randomly selected within the time to maturity, and whose monitoring price is the average asset price in a specified window surrounding the instant. We price these options whose underlying asset is the CCI30 index of various Cryptocurrencies, as opposed to a single cryptocurrency, with the intention of reducing volatility, and thus, the option price. We employ the Normal Inverse Gaussian (NIG) and Rough Fractional Stochastic Volatility (RFSV) models to the cryptocurrency market, and using the Black-Scholes as the benchmark model. In doing so, we intend to capture the extreme characteristics such as jumps and volatility roughness for cryptocurrency price fluctuations. Since there is no availability of a closed-form solution for lookback option prices under these models, we utilize the Monte Carlo simulation for pricing, and augment it using the antithetic method for variance reduction. Finally, we present the simulation results for the lookback options, and compare the prices resulting from using the NIG model, RFSV model with those from the Black-Scholes model. We find that the option price is indeed lower for our proposed window average lookback option, than for a traditional lookback option. We found the Hurst parameter to be H=0.09 which confirms that the cryptocurrency market is indeed rough.

Keywords

Cryptocurrency, Monte Carlo Simulation, Lookback Options, Rough Volatili-ty, CCI30 Index

1. Introduction

Lookback options are exotic options that enable the holder to buy or sell an underlying asset at any price the underlying asset took in a lookback time window within the time to maturity [1] [2]. For example, the payoff for a call option of this type is defined as the difference between the maximum asset price within the lookback window and the fixed strike price. For floating strike lookback call options, the payoff is the difference between the asset price at maturity and the minimum asset price within the lookback window. Lookback options belong to a class of path dependent options. Using these options, a trader is able to exploit or hedge extreme movements in the underlying asset price. Lookback options also help investors in reducing exposure to high volatility arising from extreme events, such as geopolitical developments, epidemics such as COVID-19, or natural disasters.

Due to their flexibility in reducing exposure to market risk, lookback options are well-suited for application in the cryptocurrency market, whose movements are marked by unpredictable swings with high volatility (roughly over 40% in the past few years). However, the major drawback we face is that such high volatility results in similarly high option prices, which is one of the main concerns we address in this paper by using rough fractional stochastic volatility model. See also similar studies by [3] and [4].

Our novel approach to reducing cryptocurrency option prices involves choosing an entire index, the CCI30 index in our case, as the underlying asset, as opposed to a single cryptocurrency such as Bitcoin or Ethereum. We also monitor the underlying asset across a few dates, similar to the method proposed in [5], instead of discreetly monitoring it over the entire time period until maturity. However, the problem in this method is that it enables a trader to potentially inflate or deflate the closing prices on these selected dates to produce a more favourable outcome for him. To combat the efficacy of such actions, we propose an arithmetic or geometric average of the underlying asset price over a time window centred at the discrete monitor date, as opposed to basing the option on merely the closing price on the monitor date.

For pricing, we apply the Monte Carlo simulation to price lookback options under both a Black-Scholes, a pure jump process Normal Inverse Gaussian (NIG) and in the case for rough fractional stochastic volatility model. When using Black-Scholes, we assume that the distribution of CCI30 log returns are normally distributed. However, this assumption is unlikely considering the volatility and extreme fluctuations of cryptocurrency prices. Thus, following the method proposed in [6], we employ the Normal Inverse Gaussian (NIG) model, to capture extreme variations and jumps in the CCI30 returns. Using MATLAB as a computational tool, we develop a technology to price these options and compare the prices of lookback options priced under the Black-Scholes model, the NIG model, rough volatility. We also compare prices between traditional lookback options, specifically fixed strike European lookback options, and our window-average options. We could further apply our scheme to lookback put and call options with floating strikes.

The remainder of the paper is structured as follows: Section 2 highlights the key aspects of the Monte Carlo pricing method. Section 3 discusses types of lookback options and the methodology of incorporating averaging techniques. Section 4 introduces the selected models that will be used to price lookback options, namely Black-Scholes (a benchmark model), Normal Inverse Gaussian (for jumps) and Rough Fractional Stochastic Volatility (for rough volatility). In Section 5, we discuss how model parameters are estimated from the cryptocurrency data, the CCI30 Index. Section 6 discusses numerical implementation results. Section 7 concludes the paper.

2. The Monte Carlo Pricing Method

The Monte Carlo simulation is a versatile method for pricing exotic derivatives, where analytical close-formed solutions are unavailable. Suppose that you want to compute the expectation $\mathbb{E}\left[X\right]$ of a random variable X, with distribution F. The Monte Carlo basic idea is to construct N independent draws $X_1\mathrm{,}{X}_2\mathrm{,}\cdots \mathrm{,}{X}_N$ from the distribution of X. The Strong Law of Large Numbers guarantees that if

${\bar{X}}_N=\frac{1}{N}{\displaystyle \underset{l=1}{\overset{N}{\sum }}}{X}_l,$

Then

${\bar{X}}_N\to \mathbb{E}\left[X\right]\text{almost}\text{surely}\text{as}N\to \infty \mathrm{.}$

From the Central Limit Theorem we have that

$\frac{1}{{\sigma }_X\sqrt{N}}{\displaystyle \underset{l=1}{\overset{N}{\sum }}}\left(X_l-\mathbb{E}\left[X\right]\right)\to \mathcal{N}\left(\mathrm{0,1}\right)\text{in}\text{distribution}\text{as}N\to \infty \mathrm{,}$

where ${\sigma }_X$ is the standard deviation of X¹. It follows that, for large N, the error

${\bar{X}}_N-\mathbb{E}\left[X\right]\sim \mathcal{N}\left(\mathrm{0,}\frac{{\sigma }_X}{\sqrt{N}}\right)\mathrm{.}$

In other words, for a large number of simulations, the Monte Carlo method computes the expectation accurately. One can also encapsulate methods for variance reduction in Monte Carlo Method². In this paper, we will employ the method of antithetic variates. This method works in the same way, by introducing negative dependence runs of the Monte Carlo simulation. For example, if $Z\sim \mathcal{N}\left(\mathrm{0,1}\right)$, then so is $-Z$. Pricing financial derivatives by Monte Carlo involves simulating a driving process. For example, for a model driven by a Brownian motion, each path requires the generation of many independent normal random variables $Z_1\mathrm{,}{Z}_2\mathrm{,}\cdots \mathrm{,}{Z}_N$. The use of the independent normal random variables $-Z_1\mathrm{,}-Z_2\mathrm{,}\cdots \mathrm{,}-Z_N$ generates a path which is a mirror image of the original. Thus for the price of one set of normal random variables $Z_1\mathrm{,}{Z}_2\mathrm{,}\cdots \mathrm{,}{Z}_N$, one can generate two Brownian sample paths which are negatively correlated. For more about Monte Carlo method refer to [7].

The fundamental concept for pricing an European style lookback call option is presented below, see also [4] [8].

Given the number of paths (Monte-Carlo simulations) $N_S$, number of time steps per path $N_T$, and time to maturity T: Discretize the time period $\left[\mathrm{0,}T\right]$ into sub-intervals such that $T_i=i\Delta t$ ( $\Delta t=\frac{T}{N_T}$, for $i=0,1,\cdots ,N_T$ ).

Do the following steps $N_S$ times:

· Using risk-neutral asset dynamical equation for crypto index (Black Scholes or NIG), simulate evolution of asset price.

· Perform random simulation of asset price movement each time step for a total of $N_T$ time steps.

· Compute pay-off at time to maturity.

· Discount pay-off at maturity at risk-free interest rate.

Average the discounted the pay-off for all $N_S$ paths. This average is the price of the option.

3. Lookback Options

There are many types of lookback options that can be employed for risk management or capitalizing from extreme movements in the underlying asset. Let $S_0\mathrm{,}{S}_1\mathrm{,}{S}_2\mathrm{,}\cdots \mathrm{,}{S}_T$ be the values the underlying asset takes in time instances $t=0,1,2,\cdots ,T$. The following Table 1 depicts the payoffs for different lookback options:

For example, the value of a European fixed strike lookback call option at current time 0 is given by

$c_{\text{fix}}={\text{e}}^{-rT}\mathbb{E}\left[\max \left(M_0^T-K\mathrm{,0}\right)\right]={\text{e}}^{-rT}{\displaystyle {\int }_0^{\infty }}\mathbb{P}\left(M_0^T-K\ge x\right)\text{d}x$ (1)

Averaging Methods

In this paper, we explore discrete lookback options for which the asset price is

Table 1. Lookback option type.

Where, $M_0^T$ is the Maximum of $\left\{S_0\mathrm{,}{S}_1\mathrm{,}{S}_2\mathrm{,}\cdots \mathrm{,}{S}_T\right\}$ ; $m_0^T$ is the Minumum of $\left\{S_0\mathrm{,}{S}_1\mathrm{,}{S}_2\mathrm{,}\cdots \mathrm{,}{S}_T\right\}$ ; $S_T$ is the Stock price at the expiry date; K is the Strike price; $S_0$ is today’s stock price; $x^{+}=\max \left\{x,0\right\}$.

not monitored daily, but rather on a fixed subset of days. We measure the impact of averaging the asset price over a time period surrounding the monitoring points on option price. A similar approach has been explored by [5], but his approach does not account for window averaging as we propose below. We use arithmetic averaging to perform averaging over the window period, but other methods, such as geometric or harmonic averaging can also be used. Notably, the harmonic mean is the smallest of the aforementioned averages: $A_h\le A_g\le A_a$, where $A_a$ is the arithmetic mean, $A_g$ is the geometric mean, and $A_h$ is the harmonic mean. This indicates that the arithmetic mean yields a higher option price than other averaging techniques.

In the above Figure 1, we have three monitor instants, $T_1\mathrm{,}{T}_2$ and $T_3$, within the interval $\left[\mathrm{0,}T\right]$. The window periods for averaging are $W_1\mathrm{,}{W}_2$ and $W_3$, respectively. The width of the windows are not constant, and neither is the spacing between the monitor time instants.

If we use harmonic averaging over the window period for a lookback call option, the option price will be lower than if we were to use geometric or arithmetic averaging. Similarly, these averaging methods could be used for floating strike options as well to reduce the option price.

Let $S_1\mathrm{,}{S}_2\mathrm{,}\cdots \mathrm{,}{S}_L$ be the asset price within a time window of size L. Then, the averages are calculated as follows:

$A_a=\frac{{\displaystyle {\sum }_{i=1}^L{S}_i}}{L}$ (2)

$A_g=\sqrt[L]{{\displaystyle \underset{i=1}{\overset{L}{\prod }}}{S}_i}$ (3)

$A_h=\frac{L}{{\displaystyle {\sum }_{i=1}^L\frac{1}{S_i}}}$ (4)

Option value in Equation (1) can now be given as

$c_{\text{fix}}={\text{e}}^{-rT}\mathbb{E}\left[\max \left({\bar{M}}_0^T-K\mathrm{,0}\right)\right]\mathrm{,}$ (5)

where

${\bar{M}}_0^T=A_{\cdot }{\mathbb{I}}_{\left(A_{\cdot }>M_0^T\right)}+M_0^T{\mathbb{I}}_{\left(A_{\cdot }\le M_0^T\right)}\mathrm{,}$

and where the subscript ∙ of A represents the averaging technique, $\mathbb{I}$ is the indicator function.

Figure 1. Monitoring window averaging technique for lookback options.

4. Pricing Model

4.1. Black-Scholes Model

The Black-Scholes model [9] is the most widely used model for financial derivative pricing. It states that the log return of a stock (or any underlying asset) can be specified as follows:

$\log S_T-\log S_0\sim \mathcal{N}\left[\left(\mu -\frac{{\sigma }^2}{2}\right)T\mathrm{,}{\sigma }^{2}T\right]$ (6)

where $\mathcal{N}\left(\mu \mathrm{,}\sigma \right)$ is a random normal value with mean $\mu$ and variance $\sigma$, $\mu$ is the expected mean of the daily returns, and $\sigma$ is the standard deviation of the daily returns. The above Equation (6) can additionally be written as follows:

$\log \left(\frac{S_T}{S_0}\right)\sim \mathcal{N}\left[\left(\mu -\frac{{\sigma }^2}{2}\right)T\mathrm{,}{\sigma }^{2}T\right]$ (7)

In order to derive a dynamical model of the underlying asset for the Monte Carlo simulation, we use the following risk-neutral representation of the evolution of asset price:

$S_t=S_0\exp \left(\left(r-\frac{{\sigma }^2}{2}\right)\Delta t+\sigma B_t\right)$ (8)

where r is the annual risk-free rate of return and $B_t$ is the the Wiener process.

From Equation (8), a discretized asset dynamical equation for generating Monte-Carlo simulation paths can be expressed as follows

$S_{t+\Delta t}=S_t\exp \left(\left(r-\frac{{\sigma }^2}{2}\right)\Delta t+\sigma \sqrt{\Delta t}\epsilon \right)$ (9)

Given that $B_{\Delta t}=\sqrt{\Delta t}\epsilon$ and $\epsilon \sim \mathcal{N}\left(\mathrm{0,1}\right)$. In Equation (9)

· $\Delta t$ is the time interval per time step

· $S_{t+\Delta t}$ is the stock price at the next time step

· $S_t$ is the stock price at the current time step

· r is the risk-free rate of return

· $\sigma$ is the annual volatility

· $\epsilon$ is the standard normal random variable

As shown in Equation (9), generating the asset evolution over time only requires one parameter to be estimated from the crypto-index logarithm of the return time series: $\sigma$. The other parameters are user-entered; for our purposes, they are treated as constants and specified in advance.

To generate a vector of paths for each time step, we use the vectorized version of the above equation, which in MATLAB vector notation can be expressed as follows:

$S\left(t+\Delta t\right)=S\left(t\right)\exp \left(\left(r-\frac{{\sigma }^2}{2}\right)\Delta t+\sigma \sqrt{\Delta t}\epsilon \right)$ (10)

In this case, for each time step, a vector of random normal variables is generated for all $N_S$ paths. $S\left(t+\Delta t\right)$ and $S\left(t\right)$ are vectors of stock prices, and $\epsilon$ is a random value of size $N_S\times 1$.

For each time step, the vectorized Monte Carlo simulation is run as follows:

1) Generate $N_S$ randomized asset prices for the current time step using Equation (9); repeat $N_T$ times

2) Once time to maturity T is reached, calculate payoffs for each trajectory, for $N_S$ trajectories

3) Discount payoffs to the present day for each trajectory; average discounted payoffs to compute option price

In addition, for our Antithetic Monte Carlo simulation, the vector equations are modified

$S^{+}\left(t+\Delta t\right)=S^{+}\left(t\right)\exp \left(\left(r-\frac{{\sigma }^2}{2}\right)\Delta t+\sigma \sqrt{\Delta t}\epsilon \right)$ (11)

$S^{-}\left(t+\Delta t\right)=S^{-}\left(t\right)\exp \left(\left(r-\frac{{\sigma }^2}{2}\right)\Delta t-\sigma \sqrt{\Delta t}\epsilon \right)$ (12)

The evolutions of both $S^{+}$ and $S^{-}$ are used for pay-off computation at maturity to improve Monte Carlo convergence.

4.2. Normal Inverse Gaussian (NIG) Model

The Normal Inverse Gaussian process is a four-parameter distribution that captures skewness and kurtosis better than the Normal distribution. It was introduced to financial application by As the cryptocurrency market returns are highly volatile and characterized by extreme events, we propose to apply the NIG distribution to model the log returns (see also [8] ). If X is a random number generated from the NIG distribution, it can be represented as

$X\sim \mathcal{N}\mathcal{I}\mathcal{G}\left(\alpha \mathrm{,}\beta \mathrm{,}\mu \mathrm{,}\delta \right)$

where

· $\alpha$ is the tail parameter, controlling tail behavior (large values of $\alpha$ implies light tails and small values of $\alpha$ indicate heavier tails)

· $\beta$ is the skewness parameter (a negative $\beta$ implies left-skewedness and a positive $\beta$ implies right-skewedness);

· $\mu$ is the location parameter;

· $\delta$ is the scale parameter (determines the spread of returns).

The asset dynamical equation used in the Monte Carlo simulation for the NIG model is:

$S_{t+\Delta t}=S_t\exp \left(m\Delta t+\omega \Delta t+X_{\Delta t}^{\text{NIG}}\right)$ (13)

where m is the martingale-correcting term, given by

$m=r-d$

where r is the risk-free interest rate and d is the dividend rate, $\omega$ is given by:

$\omega =\delta \left(\sqrt{{\alpha }^2-{\left(1+\beta \right)}^2}-\sqrt{{\alpha }^2-{\beta }^2}\right)$

and where

$X_{\Delta t}^{\text{NIG}}\sim \mathcal{N}\mathcal{I}\mathcal{G}\left(\alpha \mathrm{,}\beta \mathrm{,}\mu \mathrm{,}\delta \right)$

To generate the asset evolution over time, four parameters, $\left[\alpha \mathrm{,}\beta \mathrm{,}\mu \mathrm{,}\delta \right]$, need to be estimated from the log return time series.

For fast simulations, we use the vector form of the above equation, given by:

$S_{t+\Delta t}=S_t\exp \left(m\Delta t+\omega \Delta t+X_{\Delta t}^{\text{NIG}}\right)$ (14)

where $S\left(t+\Delta t\right)$ and $S\left(t\right)$ are vectors of stock prices, and $X_{\Delta t}^{\text{NIG}}$ is the vector of NIG random values of size $N_S\times 1$.

4.3. Rough Fractional Stochastic Volatility (RFSV) Model

The dynamics of the Rough Fractional Stochastic Volatility (RFSV) model is given by:

$\begin{array}{l}\text{d}{S}_t=r_t{S}_t\text{d}t+{\sigma }_t{S}_t\text{d}{W}_t\\ {\sigma }_t=\exp \left\{X_t\right\}\\ \text{d}{X}_t=-\alpha \left(X_t-m\right)\text{d}t+\nu \text{d}{W}_t^H\end{array}$ (15)

where $m\in ℝ$, and $\nu$ and $\alpha$ are positive parameters, with a classical Brownian motion W and a fraction Brownian motion (fBM) $W^H$ ³. Here we assume that the two driving processes are uncorrelated. $H\in \left(\mathrm{0,1}\right)$ is the Hurst parameter that characterize the increment of the driving process. In particular, when

$H=\frac{1}{2}$, fractional Brownian motion boils down to just Brownian motion.

· If $H>1/2$, increments are positively correlated.

· If $H<1/2$, increments are negatively correlated and the measured smoothness of the volatility.

While standard Brownian motion has independent increments, fBM displays auto-correlation, i.e., it does not exhibit independent increments. The covariance between an fBM process at times t and s is

$\mathbb{E}\left[W_t^H{W}_s^H\right]=\frac{1}{2}\left({\left|t\right|}^{2H}+{\left|s\right|}^{2H}-{\left|t-s\right|}^{2H}\right)\mathrm{.}$ (16)

This covariance implies an auto-correlation function that decays slower at $H\ne \frac{1}{2}$. The auto-correlation function has the slowest decay when $\frac{1}{2}<H<1$.

This implies that an fBM process will exhibit longer memory as H gets closer to 1. This gave rise to the model of [10] using fBM to model volatility as a long

memory process, hence their restriction of $H\in \left(\frac{1}{2}\mathrm{,1}\right)$. Recently, however, [11] argues that volatility exhibits shorter memory, hence their restriction of $H\in \left(\mathrm{0,}\frac{1}{2}\right)$. [3] showed that the Hurst paremeter estimated from realised volatility of Bitcoin is closer to zero.

To simulate the model in Equation (15), we need to be able to simulate the fractional Brownian motion. In particular, considering the method which completely capture the covariance structure of fBM, as opposed to approximate methods which aim to reduce computation times by approximating the covariance structure.

In this paper we use Hosking Method introduced by Hosking [12] to simulate the sample paths of the fBM. This method is concerned with simulating fractional Gaussian noise (fGn). A sample fBM path can then be recovered by using a cumulative sum on the generated fGn sequence. The sequence ${\left(X_n\right)}_{n\in \mathbb{N}}$ of fractional Gaussian noise is computed recursively by computing the conditional distribution of $X_{n+1}$ given $X_n\mathrm{,}\cdots \mathrm{,}{X}_0$. The required sample is found by generating a standard normal random variable $X_0$ and calculating the remaining $X_{n+1}$ recursively. The conditional distribution is itself Gaussian with mean

${\mu }_n=c{\left(n\right)}^{\prime }\Sigma {\left(n\right)}^{-1}\left[\begin{array}{c}{X}_n\\ ⋮\\ X_1\\ X_0\end{array}\right]\mathrm{,}$ (17)

and variance

${\sigma }_n^2=1-c{\left(n\right)}^{\prime }\Sigma {\left(n\right)}^{-1}c\left(n\right)\mathrm{,}$ (18)

where $c\left(n\right)$ is an $\left(n+1\right)$ -column vector with elements $c{\left(n\right)}_k=\sigma \left(m\mathrm{,}m+k+1\right)$, for $k=0,\cdots ,n$. The algorithm presented by Hosking computes $\Sigma {\left(n\right)}^{-1}c\left(n\right)$ recursively to ensure greater efficiency.

5. Model Parameter Estimation

We employ Crypto Currency Index 30 (CCI30) daily data⁴ covering the period from 2015/01/01-2021/03/17. These crypto currency trade 365 days in a year and in total we have 2268 data point in our sample. Figure 2 shows the index prices as well the the return series. In Figure 3, we show 10-day realised volatility of the index⁵. A sharp increase in the index price is seen during the COVID-19 pandemic.

5.1. Black-Scholes Model Parameter

In order to generate asset paths based on the Black-Scholes model, we need the volatility parameter from the log returns, which are estimated as:

Figure 2. CCI30 index prices and returns. (a) Index price trajectory; (b) Index return.

Figure 3. 10-days realised volatility of crypto currency index 30.

$r_t=\log \left(\frac{P_t}{P_{t-1}}\right)$ (19)

where $P_t$ is the closing price of the CCI30 index on day t, and $P_{t-1}$ is the closing price of CCI30 on day $t-1$. After computing the log return time series, we estimate the variance.

For our purposes, we used price data from 1/1/2015 to 1/19/2020 to compute the return series and its variance, see also [13] [14] [15].

The annual volatility parameter $\sigma$ that is needed for the Black-Scholes asset dynamical model is computed as:

$\sigma =\sqrt{{\hat{\sigma }}^2\times 365},$

where ${\hat{\sigma }}^2$ is the variance estimated from the return series.

In our case, the number of trading days in a year is 365 days for the annualized volatility calculation, as cryptocurrency markets do not close any day of the year, as opposed to the stock market, for example. Table 2 shows the estemates for Black-Scholes model parematers.

5.2. NIG Model Parameter

In order to generate asset paths based on the NIG model, we need to estimate four parameters, $\left[\alpha \mathrm{,}\beta \mathrm{,}\mu \mathrm{,}\delta \right]$, from the CCI30 log returns. The log returns for the CCI30 index are computed exactly in the same way as in Equation (18).

After computing the log return time series, we estimate its mean, variance, skewness and kurtosis.

These four parameters are then used as input to the function from the MATLAB built-in NIG toolbox function which then returns the parameters $\left[\alpha \mathrm{,}\beta \mathrm{,}\mu \mathrm{,}\delta \right]$. The estimates for these parameters are given in Table 3.

5.3. RFSV Model Parameter

Recall that fBM is a Gaussian process with the property that

$\mathbb{E}\left[{\left|W_{t+\Delta }^H-W_t^H\right|}^q\right]=K_q{\Delta }^{qH}\mathrm{.}$

verify that the empirical distributions of log-volatility are approximately Gaussian for various time lags.

Thus, to estimate the smoothness of the volatility process, that is, H, [11] use

Table 2. Black-Scholes parameter estimates.

Table 3. NIG parameter estimates.

the following approach. Suppose that we have access to N discrete observations of the volatility process ${\sigma }_k$ on $\left[\mathrm{0,}T\right]$. Calculate

$m\left(q\mathrm{,}\Delta \right)=\frac{1}{N-\Delta }{\displaystyle \underset{k=1}{\overset{N-\Delta }{\sum }}}{\left|\log \left({\sigma }_{k+\Delta }\right)-\log \left({\sigma }_k\right)\right|}^q\mathrm{,}$

where $\Delta \in \mathbb{N}$ is the lag.

Now, assuming that the log-volatility process has stationary increments, then $m\left(q\mathrm{,}\Delta \right)$ can be seen as an estimate of

$\mathbb{E}\left[{\left|\log \left({\sigma }_{\Delta }\right)-\log \left({\sigma }_0\right)\right|}^q\right]=K_q{\Delta }^{qH}\mathrm{.}$

Taking logs, we get

$\log \mathbb{E}\left[{\left|\log \left({\sigma }_{\Delta }\right)-\log \left({\sigma }_0\right)\right|}^q\right]=\log K_q+qH\log \Delta \mathrm{.}$

We can then compute $m\left(q\mathrm{,}\Delta \right)$ for different values of $\Delta$ for each q and regress $\log m\left(q\mathrm{,}\Delta \right)$ against $\log \Delta$. The slope of each line of best fit is then an estimate of qH and the intercept of the best fit is $\nu$.

We then make use of 10-day realised volatility to estimate H and $\nu$ as in [11]. Having estimated these two parameters from 10-day realised volatility, we are left with three other RFSV model parameters to be specified, that is $X_0\mathrm{,}\alpha$ and m. To get these model parameters, we use Black-Scholes call prices for Lookback option and calibrate RFSV model to these prices. We set up our objective function that we minimise to infer the remaining RFSV model parameters.

We define calibration error based on the relative mean squared error (RMSE) of the objective function for parameter calibration as follows:

$\frac{1}{N_O}{\displaystyle \underset{j=1}{\overset{N_O}{\sum }}}{\left(P_j^{\text{model}}\left(\Theta \right)-P_j^{\text{BS}}\right)}^2$ (19)

where $N_O$ is the number of options used in the calibration, $\Theta$ is the model parameter set, $P_j^{\text{model}}\left(\Theta \right)$ is the jth theoretical model price dependent on $\Theta$, $P_j^{\text{BS}}$ is the jth Black-Scholes option price. Calibration results are given in Table 4.

6. Numerical Results

In this section we discuss numerical implementation results. In Tables 5-7, we show discretely monitoring lookback option prices consider all time steps for all three models. Option prices across various strikes are almost close to each others for all models. We also report the standard errors from the Monte Carlo and the

Table 4. RFSV paremeter estimates.

confident interval. Table 8 shows option prices using selected time spacing interval. In all models, we see that the larger the time spacing, the smaller the option prices. In Table 9 we show lookback options for all models using arithemtic, geometric and harmonic averaging methods with varying window size. Harmonic averaging method has a slighly low prices as compared to other methods. Figure 4 shows prices under each model when time spacing for window size increase. Rough volatility model has lower price compare to the other models. This agrees with finding by [16] that rough volatility model prices are lower in general. Figure 5 shows implied volatility for all three models. We find that BS lies in-between NIG and RFSV.

6.1. Discrete Monitoring Across all Timesteps

In this section, we provide numerical results for lookback options under discrete monitoring cross all available time-steps for all models.

(a) (b)

Table 5. Black-Scholes Lookback option prices. (a) Plain Monte Carlo method; (b) Antithetic Monte Carlo method.

(a) (b)

Table 6. NIG Lookback option prices. (a) Plain Monte Carlo method; (b) Antithetic Monte Carlo method.

(a) (b)

Table 7. RFSV Lookback option prices. (a) Plain Monte Carlo method; (b) Antithetic Monte Carlo method.

6.2. Discrete Monitoring at Select Timesteps

In this section, we provide numerical results for lookback options under discrete monitoring cross selected time-steps for all models.

(a) (b) (c)

Table 8. Discrete Monitoring at Select Timesteps. (a) Black-Scholes Monte Carlo method; (b) NIG Monte Carlo method; (c) RFSV Monte Carlo method.

6.3. Discrete Monitoring Spaced with Window Averaging

In this section, we provide numerical results for lookback options under discrete monitoring with window averaging all models.

(a) (b) (c)

Table 9. Discrete Monitored Spaced Lookback with Window Averaging. (a) BS Monte Carlo method; (b) RFSV Monte Carlo method; (c) NIG Monte Carlo method

7. Conclusions

Upon running our simulations, we observe many noticeable differences between the methods we use to compute the price of the options, as well as between different types of options. Firstly, when comparing option prices between the Monte Carlo and Antithetic Monte Carlo simulations for the Black-Scholes Model, we find that employing the latter significantly reduces standard error.

Regarding select lookback options, we notice a clear decrease in price as time spacing between the monitor instances increases for both Black-Scholes and NIG implementations, a result of a lower number of asset return values being considered for pricing the option. With respect to our proposed select window average lookback options, we observe a similar phenomenon for both the Black-Scholes and NIG implementations. We find that as the size of the window increases, the option price decreases, regardless of the type of averaging that is

Figure 4. Changes in Option Price as Time Spacing Increases for Black-Scholes, NIG and RFSV model.

Figure 5. Model implied volatility.

used. Moreover, our prediction that the harmonic mean yields the lowest option price holds true across all window sizes, the arithmetic mean results in the highest option price, followed by the geometric and harmonic means, respectively. It must be noted, however, that we employ the same random trajectory in the simulations for comparing the different averages, so as to fairly compare their performances. Finally, for all simulations, we observe that the option prices resulting from using the NIG model are generally slightly lower than from using the Black-Scholes model, irrespective of any other considerations (e.g. time spacing and type of average). Lookback options prices obtained via rough volatility volatility model are lower than the benchmark model.

Further research direction related to this work will be to consider window average lookback options by incorporating model regime switching features into our implementation.

NOTES

¹ ${\sigma }_X$ may be unknown but we can always find its unbiased estimate: $s_N=\sqrt{\frac{1}{N-1}{\displaystyle {\sum }_{l=1}^N}{\left(X_l-{\bar{X}}_N\right)}^2}$.

²Using the fact that a portfolio of two negatively correlated asset is less risky than the risk inherent in each asset alone.

³A fractional Brownian motion ${\left(W_t^H\right)}_{t\in ℝ}$ is a centred self-similar Gaussian process with stationary increments such that for all $t\in ℝ$,

$\mathbb{E}\left[{\left|W_{t+\Delta }^H-W_t^H\right|}^q\right]=K_q{\Delta }^{qH}\mathrm{,}$

where $\Delta \ge 0$ is a small time interval, $K_q$ is the moment of order $q>0$ of the absolute value of a standard Gaussian random variable and $H\in \left(\mathrm{0,1}\right)$ is the Hurst parameter which defines the fBM.

⁴This data is freely available online: https://www.cci30.com.

⁵The formula of Realised volatility is the square root of realized variance (RV). We first compute the return.

$r_t=\log \frac{P_t}{P_{t-1}},$

then

${\text{RV}}_t={\displaystyle \underset{i=t-10}{\overset{t-1}{\sum }}}{r}_i^2$

Conflicts of Interest

The authors declare no conflicts of interest.

References

[1]	Barndimarte, P. (2002) Numerical Methods in Finance: A Matlab-Based Introduction. John Wiley & Sons, London.
[2]	Clewlow, C. and Strickland, L. (2003) Implementing Derivative Models. John Wiley & Sons, London.
[3]	Takaishi, T. (2020) Rough Volatility of Bitcoin. Finance Research Letters, 32, Article ID: 101379. https://doi.org/10.1016/j.frl.2019.101379
[4]	Umeorah, N. (2017) Pricing Barrier and Look-Back Options. North-West University Thesis.
[5]	Chang, H.-C.H. and Li, K. (2018) The Amnesiac Lookback Option: Selectively Monitored Lookback Options and Cryptocurrencies. Frontiers in Applied Mathematics and Statistics, 4, 10. https://doi.org/10.3389/fams.2018.00010
[6]	Alfeus, M. (2013) Fast Pricing of Barrier Options. BSc Hons Thesis, Stellenbosch University, Stellenbosch.
[7]	Glasserman, P. (2004) Monte Carlo Methods in Financial Engineering. Springer, New York.
[8]	Saebo, K.K. (2009) Pricing Exotic Options with NIG Model Using Path Integration. Thesis Report, Norwegian University of Science & Technology, Trondheim.
[9]	Black, F. and Scholes, M. (1973) The Pricing of Options and Corporate Liabilities. Journal of Political Economy, 81, 637-654. https://doi.org/10.1086/260062
[10]	Comte, F. and Renault, E. (1998) Long Memory in Continuous-Time Stochastic Volatility Models. Mathematical Finance, 8, 291-323. https://doi.org/10.1111/1467-9965.00057
[11]	Gatheral, J.J., Jaisson, T. and Rosenbaum, M. (2018) Volatility Is Rough. Quantitative Finance, 18, 933-949. https://doi.org/10.1080/14697688.2017.1393551
[12]	Hosking, J.R.M. (1984) Modeling Persistence in Hydrological Time Series Using Fractional Differencing. Water Resources Research, 20, 1898-1908. https://doi.org/10.1029/WR020i012p01898
[13]	Chen, L., Gao, R., Bian, Y. and Di, H. (2021) Elliptic Entropy of Uncertain Random Variables with Application to Portfolio Selection. Soft Computing, 25, 1925-1937. https://doi.org/10.1007/s00500-020-05266-z
[14]	Chen, L., Peng, J., Zhang, B. and Rosyida, I. (2017) Diversified Models for Portfolio Selection Based on Uncertain Semivariance. International Journal of Systems Science, 48, 637-648. https://doi.org/10.1080/00207721.2016.1206985
[15]	Gross, P. (2006) Parameter Estimation for Black-Scholes Equation. Technical Report.
[16]	Alfeus, M., Overbeck, L. and Schlogl, E. (2019) Regime Switching Rough Heston Model. The Journal of Futures Markets, 39, 538-552. https://doi.org/10.1002/fut.21993