Forecasting semi-stationary processes and statistical arbitrage

ABSTRACT

If a financial derivative can be traded consecutively and its terminal payoffs can be adjusted as the sum of a bounded process and a stationary process, then we can use the moving average of the historical payoffs to forecast and the corresponding errors form a generalised mean reversion process. Thus we can price the financial derivatives by its moving average. One can even possibly get statistical arbitrage from certain derivative pricing. We particularly discuss the example of European call options. We show that there is a possibility to get statistical arbitrage from Black–Scholes's option price.

Keywords:

• Stationary process
• statistical arbitrage
• Black–Scholes model
• mean reversion

1. Introduction

In economics, arbitrage is the act of exploiting price differences between two or more markets: building a portfolio of deals upon the imbalance, expecting to get the spread profit. In academia, arbitrage refers to the cost of the possibility obtaining a risk-free profit after the transaction is completed. For example, when there is an arbitrage opportunity, there will be a low buy and sell high immediately.

The practice of ‘statistical arbitrage’ actually came together with gambling. Suppose that a dice is rolled continuously and a gambler may lose 1 dollar each time when the result is ‘1’ or wins 1 dollar otherwise. Denote by $S_n$ the gambler's gain after n tries, then $S_n/n$ is approximately $\frac{2}{3}$ for large n. That is the law of large numbers in probability theory on which based the theory of mathematical statistics. Furthermore the strong law of large numbers tells us that when $n\to \mathrm{\infty }$, $S_n/n$ approaches to $\frac{2}{3}$ with probability 1. However, that is not an arbitrage opportunity for this gambler as he still has a chance to lose no matter how large n is. Morgan Stanley started in the early 1980s (see Gregory van Kipnis' foreword to Pole, 2007) to apply ‘statistical arbitrage’ to get profit from the stock market. Nevertheless, the returns in the stock market are not independent as rolling dices, so many mathematical models were created, of which most would be possibly never published due to profitability. In Pole (2007), Andrew Pole showed many important historical examples, ‘rules’, and structural models of statistical arbitrage. In academic use, an arbitrage is risk-free; in common use, as in statistical arbitrage, it may refer to expected profit, though losses may occur, and in practice, there are always risks in arbitrage, some minor (such as fluctuation of prices decreasing profit margins), some major (such as devaluation of a currency or derivative). In academic use, an arbitrage involves taking advantage of differences in price of a single asset or identical cash flows; in common use, it is also used to refer to differences between similar assets (relative value or convergence trades), as in merger arbitrage.

There have been a few more general academic definitions of statistical arbitrage since that time on. Hogan, Jarrow, and Warachka (2002) defined ‘Statistical Arbitrage’ by four conditions that the discounted cumulative value $v\left(t\right)$ should satisfy: (1) $v\left(0\right)=0;$ (2) $\underset{t\to \mathrm{\infty }}{lim}E\left[v\right(t\left)\right]>0;$ (3) $\underset{t\to \mathrm{\infty }}{lim}\mathrm{V}\mathrm{a}\mathrm{r}\left(v\right(t\left)\right)/t=0;$ (4) $\underset{t\to \mathrm{\infty }}{lim}P\left[v\right(t)<0]=0$. Their Condition (3) excluded the cases where $v\left(t\right)$ is the sum of independent identically distributed outcomes, and thus excluded our dice game and most of casino games. According to Lo (2010), ‘Statistical Arbitrage’ refers to highly technical short-term mean-reversion strategies involving large numbers of securities (hundreds to thousands, depending on the amount of risk capital), very short holding periods (measured in days to seconds), and substantial computational, trading and information technology (IT) infrastructure. Wang and Zheng (2014) simplified the definition to repeatedly trading a basket of assets according to the same algorithm and get an accumulated profit with statistically stable positive rate. Moreover, the concept of the ergodic theorem of stationary process was applied into the practice of statistical arbitrage in Wang and Zheng (2014). We will further this discussion to derivative pricing in this paper. We will discuss the time series with bounded trend component, of which the increments have some nice property as mean reversion.

The ‘mean reversion’ phenomena is a special case of stationarity, which suggests that prices and returns eventually move back towards the mean or average. This mean or average can be the historical average of the price or return or another relevant average such as the growth in the economy or the average return of an industry. A typical example of mean reversion phenomena is the increments of the Ornstein–Uhlenbeck process (see, e.g., Karatzas & Shreve, 1987, p. 358). However, the Ornstein–Uhlenbeck process comes together with Gaussian law which could not be justified by some real market data. Thus we will define a generalised discrete time mean reversion process as the increments of a semi-stationary process which is the sum of a bounded process and a stationary process. The time average of a generalised mean reversion process vanishes faster than that of an ordinary stationary processes with null mean. Therefore a generalised mean reversion process is more stable than an ordinary stationary one. We show in the next section that if we use the moving average to estimate a time series with bounded trend component, then the time average of the error vanishes in inverse proportion to the time length. As its application, we give an example of statistical arbitrage from Black–Scholes' option price based on estimate of payoffs with generalised mean reversion errors.

We introduce the concept of generalised mean reversion process as the increment process of a semi-stationary process in Section 2. We use Section 3 to explain some basic facts of option pricing. We use Section 4 to study the options of the ETFs which track the main stock indices in the US market. Finally we show that one may get statistical arbitrage from Black–Scholes model in the last section.

2. Forecasting with MA

Given a discrete time stochastic process $\left\{X\right(t){\}}_{t=0,1,2,\dots },$ it is very important to give an estimate $\tilde{X}(t+1)$ for $X(t+1)$ based on the observed information $\left\{X\right(0),X(1),\dots ,X(t\left)\right\}$. Certainly, the selection of $\tilde{X}(t+1)$ depends on the selection of the norm of the error $X(t+1)-\tilde{X}(t+1).$ When we need to estimate continuously for $t=m+1,m+2,\dots ,N$, then the most popular ones are the sum of $l_p-$norms ($p\ge 1$) of errors $\sum _{t=m+1}^N\left\{E\right[|X\left(t\right)-\tilde{X}\left(t\right){|}^p]{\}}^{1/p}$ or the similar ones. However, in many financial applications, an investor is more interested in the $l_p$-norm of the accumulated sum: (1) $\begin{array}{rl}& {\left\{E\left[|\sum _{t=m+1}^N\right(X\left(t\right)-\tilde{X}\left(t\right)\left){|}^p\right]\right\}}^{1/p}\\ & ={\left\{E\left[|\sum _{t=m+1}^{T}X\right(t)-\sum _{t=m+1}^N\tilde{X}(t\left){|}^p\right]\right\}}^{1/p}.\end{array}$(1) The economic reason is very simple. If $X(t+1)$ is the value of certain financial derivative at time t + 1 and $X(t+1)-\tilde{X}(t+1)$ is the loss caused by estimate $\tilde{X}(t+1),$ then an investor is more interested in the accumulated lose $\sum _{t=m+1}^N\left(X\right(t)-\tilde{X}(t\left)\right)$ rather than the sum $\sum _{t=m+1}^N|X\left(t\right)-\tilde{X}\left(t\right)|$. For the error norm (1(1) $\begin{array}{rl}& {\left\{E\left[|\sum _{t=m+1}^N\right(X\left(t\right)-\tilde{X}\left(t\right)\left){|}^p\right]\right\}}^{1/p}\\ & ={\left\{E\left[|\sum _{t=m+1}^{T}X\right(t)-\sum _{t=m+1}^N\tilde{X}(t\left){|}^p\right]\right\}}^{1/p}.\end{array}$(1) ), the simplest estimate is to use the moving average (ab. MA) (2) $\begin{array}{r}{M}^{\left(m\right)}\left(t\right)=\frac{X\left(t\right)+X(t-1)+\cdots +X(t-m+1)}{m}\end{array}$(2) to forecast $X(t+1),$ which is called the MA model in time series analysis. We will impose a condition ‘semi-stationary’ on $\left\{X\right(t\left)\right\}$ under which (1(1) $\begin{array}{rl}& {\left\{E\left[|\sum _{t=m+1}^N\right(X\left(t\right)-\tilde{X}\left(t\right)\left){|}^p\right]\right\}}^{1/p}\\ & ={\left\{E\left[|\sum _{t=m+1}^{T}X\right(t)-\sum _{t=m+1}^N\tilde{X}(t\left){|}^p\right]\right\}}^{1/p}.\end{array}$(1) ) will be bounded in T for some $p\ge 0.$

We say that $\left\{X\right(t\left)\right\}$ is a ‘strongly stationary’ process, if for each positive a, the processes $\left\{X\right(t+a\left)\right\}$ and $\left\{X\right(t\left)\right\}$ obey the same probability law. Thus a sequence of independent identically distributed random variables is a strongly stationary process. In the previous example of continuously rolling dice, if we denote $X\left(t\right)=-1$ when the tth outcome is ‘1’, and $X\left(t\right)=1$ otherwise, then $\left\{X\right(t\left)\right\}$ is a strongly stationary process. The weak form is the weakly stationary process. $\left\{X^{\prime }\right(t){\}}_t$ is a ‘weakly stationary’ process if: (i) $E\left[X^{\prime }\right(t\left)\right]$ is a constant; (ii) $\mathrm{C}\mathrm{o}\mathrm{v}\left(X^{\prime }\right(t),X^{\prime }(t+a\left)\right)=\mathrm{C}\mathrm{o}\mathrm{v}\left(X^{\prime }\right(0),X^{\prime }(a\left)\right)$ for each t. Certainly, if a strongly stationary process $\left\{X\right(t\left)\right\}$ has its second moments, then $\left\{X\right(t\left)\right\}$ is weakly stationary. The most important property of a strongly stationary process is the Birkhoff's theorem (Loeve, 1977, p. 76). The time-average $\left(X\left(1\right)+X\left(2\right)+\cdots +X\left(N\right)\right)/N$ of a strongly stationary process $\left\{X\right(t\left)\right\}$ converges with probability 1. Furthermore, when $\left\{X\right(t\left)\right\}$ has the so-called ‘ergodic’ property, then the limit is just its mean.

Definition 2.1

We call $\left\{X\right(t\left)\right\}$ a strongly (or weakly) semi-stationary process if there is a bounded process $\left\{A\right(t\left)\right\}$ such that $\left\{X\right(t)-A(t\left)\right\}$ is a strongly (or weakly, respectively) stationary process.

It is easy to see that for a semi-stationary process $\left\{X\right(t\left)\right\}$, $X\left(T\right)/T\to 0$ with the rate $1/T$ and the time average $\frac{\mathrm{\Delta }X\left(1\right)+\mathrm{\Delta }X\left(2\right)+\cdots +\mathrm{\Delta }X\left(N\right)}{N}$ of the increments $\mathrm{\Delta }X\left(t\right)=X(t+1)-X\left(t\right)$ of a semi-stationary process $\left\{X\right(t\left)\right\}$ vanishes at the rate $1/N$ which is faster than that of an ordinary sequence of independent identically distributed random variables, which is $1/\sqrt{N}$.

A time series $\left\{x\right(t\left)\right\}$ is often to be considered as a sequence of data which can be written as $x\left(t\right)=m\left(t\right)+S\left(t\right)+ϵ\left(t\right),$ where $m\left(t\right)$ is the deterministic trend component, $S\left(t\right)$ is the seasonal component and $ϵ\left(t\right)$ is a weakly stationary process (see Wang & Zheng, 2014, p. 101, for example). Since $S\left(t\right)$ is periodic, it is bounded. If $m\left(t\right)$ is also bounded, then $\left\{x\right(t\left)\right\}$ is a weakly semi-stationary process. We are going to use MA to forecast a weakly semi-stationary time series and the cumulated error has mean vanishes at the rate $1/T$.

The increment of a semi-stationary process has the properties which are very similar to the so-called mean reversion process. The Mean Reversion phenomena in stock prices have been studied for more than three decades. This theory suggests that prices and returns eventually move back towards the mean or average. This mean or average can be the historical average of the price or return or another relevant average such as the growth in the economy or the average return of an industry (see, e.g., Ansley, Spivey, & Wrobleski, 1977; Fama & French, 1988; Mukherji, 2011; Poterba & Summers, 1988). Its typical mathematical model is the Ornstein–Uhlenbeck process of which the increments can be written as (3) $\mathrm{d}X\left(t\right)=\theta (\mu -X(t\left)\right)\mathrm{d}t+\sigma \mathrm{d}W\left(t\right),$(3) where μ, $\theta (>0)$ and $\sigma (>0)$ are constants. $\mathrm{d}X\left(t\right)$ and $\mathrm{d}W\left(t\right)$ are the increments of $X\left(t\right)$ and Brownian motion $W\left(t\right)$ respectively. Since the above equation is not taught in ordinary text books of probability theory, we briefly explain its meaning here. One considers $\mathrm{d}W\left(t\right)$ as the noise and σ as its magnitude. $\left(X\right(t)-\mu )$ is the distance from the mean μ. So this process has a drift $\theta (\mu -X(t\left)\right)$ push its path back to the mean. If we denote $X\left(t\right)-\mu =Y\left(t\right),$ then $\mathrm{d}Y\left(t\right)=-\theta Y\left(t\right)\mathrm{d}t+\sigma \mathrm{d}W\left(t\right),$ which has its integral form (see Karatzas & Shreve, 1987, p. 358). $Y\left(t\right)=Y\left(0\right)\exp \{-\theta t\}+\sigma {\int }_0^t\exp \left\{\theta \right(s-t\left)\right\}\mathrm{d}W\left(s\right).$ Thus $\begin{array}{rl}X\left(t\right)& =\left(X\right(0)-\mu )\exp \{-\theta t\}\\ & +\sigma {\int }_0^t\exp \left\{\theta \right(s-t\left)\right\}\mathrm{d}W\left(s\right)+\mu ,\end{array}$ which is strongly stationary if $X\left(0\right)$ has Gaussian distribution with mean μ and variance ${\sigma }^2/2\theta$. However, it is known that this process does not precisely match a lot of data experimentally. Therefore we need a more general definition. If we fix $\left\{Y\right(t\left)\right\}$ and set $X\left(t\right)=Y\left(t\right)+\mu$, then $X\left(t\right)$ satisfies (3(3) $\mathrm{d}X\left(t\right)=\theta (\mu -X(t\left)\right)\mathrm{d}t+\sigma \mathrm{d}W\left(t\right),$(3) ) and they have same increments: $\begin{array}{rl}& X(t+1)-X\left(t\right)\\ & =Y(t+1)-Y\left(t\right),(t=0,1,2,\dots ).\end{array}$

Hence, from the uniqueness in law of the Ornstein–Uhlenbeck process (for fixed initial stationary distribution, μ, θ and σ), we get the following.

Lemma 2.2

For fixed θ and σ, the probability law of the increments $\left\{X\right(t+1)-X(t){\}}_{t=0,1,2,\dots }$ of stationary Ornstein–Uhlenbeck process do not depend on the mean μ.

The above Lemma inspired us to introduce the following.

Definition 2.3

Let $\left\{U\right(t){\}}_{t=1,2,\dots }$ be a stochastic process. If there is a strongly (or weakly) semi-stationary process $\left\{X\right(t){\}}_{t=0,1,2,\dots }$ such that $U\left(t\right)=\mathrm{\Delta }X\left(t\right),(t=1,2,\dots ),$ then $\left\{U\right(t\left)\right\}$ is called a generalised strongly (or weakly, respectively) Mean Reversion process and $\left\{X\right(t\left)\right\}$ is its integrated process. Furthermore, if the integrated process is strongly (or weakly) stationary, then $\left\{U\right(t\left)\right\}$ is called strongly (or weakly, respectively) Mean Reversion process.

Thus the time average $\frac{U\left(1\right)+\cdots +U\left(N\right)}{N}=\frac{X\left(N\right)-X\left(0\right)}{N}\to 0$ with rate $1/N$.

We have easily that $\begin{array}{rl}& \sum _{t=1}^m\frac{t}{m}X\left(t\right)+\sum _{t=m+1}^N\left[X\right(t)-M^{\left(m\right)}(t-1\left)\right]\\ & =\sum _{t=1}^m\frac{t}{m}X\left(t\right)+\sum _{t=m+1}^N\left[X\left(t\right)-\frac{1}{m}\sum _{i=1}^{m}X(t-i)\right]\\ & =\sum _{t=1}^m\frac{t}{m}X\left(t\right)+\sum _{t=m+1}^{N}X\left(t\right)\\ & -\frac{1}{m}\sum _{i=1}^m\sum _{t=m+1}^{N}X(t-i)\\ & =\sum _{i=1}^m\frac{m+1-i}{m}X(N+1-i)\\ & =\sum _{t=1}^m\frac{t}{m}X(N-m+t).\end{array}$ Thus we have

Lemma 2.4

$\begin{array}{rl}& \sum _{t=m+1}^N\left[X\right(t)-M^{\left(m\right)}(t-1\left)\right]\\ & =\sum _{i=1}^m\frac{i}{m}X(N-m+i)-\sum _{i=1}^m\frac{i}{m}X\left(i\right).\end{array}$

Since $\sum _{i=1}^m\left(i/m\right)X(t-m+i)$ is a finite linear combination of the values at different time of a strongly (or weakly) stationary processes, which is also strongly (weakly, respectively) stationary, we get

Theorem 2.5

When $\left\{X\right(t\left)\right\}$ is a strongly (or weakly) semi-stationary process, the error process $\left\{X\right(t+1)-M^{\left(m\right)}(t\left)\right\}$ is a generalised strongly (or weakly, respectively) mean reversion process and the accumulate error is strongly (or weakly, respectively) semi-stationary such that $\frac{\sum _{t=m+1}^N\left(X\right(t)-M^{\left(m\right)}(t-1\left)\right)}{N}\to 0$ with rate $1/N$.

Therefore, if $\left\{X\right(t\left)\right\}$ is a sequence of financial derivatives in the market, which is semi-stationary, then we can use $\left\{M^{\left(m\right)}\right(t-1\left)\right\}$ to forecast its values with errors forming a generalised mean reversion process. Hence the errors have mean nearly 0 in long run. In particular, we will consider option prices in the next section.

The ergodic theorem for stationary process has its weak point. It cannot be applied to an arbitrary subsequence of a stationary process. For example, when X is a standard Gaussian random variable, then $\{X,-X,X,-X,\dots \}$ is a strongly stationary process. However, one can easily choose a subsequence which has no convergent time average. Indeed, if we choose $\{X,-X,-X,-X,X,X,X,X,X,X,X,X,X,X,-X,-X,\dots \},$ then its time average will take values $X/2$ and $-\left(X/2\right)$ repeatedly. Therefore, we need the following.

Lemma 2.6

Suppose that $\left\{X_j\right\}$ is a weakly stationary sequence. If $\left\{X_{j_i}\right\}$ is a subsequence and $\underset{i}{min}\{|j_i-j_{i+1}|\}\ge J,$ then for any $ϵ>0,$ $\begin{array}{rl}& P\left[\left|\frac{1}{N}\sum _i^N{X}_{j_i}-E\left[X_1\right]\right|\ge ϵ\right]\\ & \le \frac{1}{{ϵ}^2}\left[\frac{\mathrm{V}\mathrm{a}\mathrm{r}\left(X_1\right)}{N}+\underset{i\ge J}{sup}\left|\mathrm{C}\mathrm{o}\mathrm{v}(X_1,X_{1+i})\right|\right].\end{array}$

Proof. $\begin{array}{rl}& E{\left|\frac{1}{N}\sum _i^N{X}_{j_i}-E\left[X_1\right]\right|}^2\\ & \le \frac{1}{N^2}\sum _{i=1}^N\mathrm{V}\mathrm{a}\mathrm{r}\left(X_{j_i}\right)+\frac{2}{N^2}\sum _{i<k}^N\mathrm{C}\mathrm{o}\mathrm{v}(X_{j_i},X_{j_k})\end{array}$ Thus we get the result by simplification and the classical Chebyshev's inequality in probability theory.

3. Price of option

The basic idea of Black–Scholes–Merton's theory (Black & Scholes, 1973; Merton, 1973) is that the option price of an asset depends only on the current price $S\left(t\right)$ of the asset, mature time, strike price, volatility and risk-free interest rate. It is also known that the real traded prices of various options in the market are quite different from those theoretical ones, so that the concepts of ‘implied volatility’ and ‘stochastic volatility’ were introduced (Canina & Figlewski, 1993; Chrisensen & Prabhala, 1998; Dumas, Fleming, & Whaley, 1998; Heston, 1993; Poon & Granger, 2003) as attempts to fill the gap between the theoretical ones and the real market ones. In the last two decades, many mathematicians and statisticians introduced a lot of more sophisticated stochastic models in order to describe more accurately the movement of stock price and its option prices. Among them, the volatility study always attracts the main interests.

In the basic Black–Scholes's model, stock share price is assumed to be a geometric Brownian motion. That is, $S\left(t\right)=S\left(0\right)\exp \left\{\sigma W\left(t\right)-\left(\frac{{\sigma }^2}{2}-\mu \right)t\right\}.$ The main mathematical tool in Black–Scholes's theory is Cameron–Martin–Girsanov's theorem (Karatzas & Shreve, 1987, p. 190) which states that up to a bounded time T, the induced probability measures corresponding to different μ are equivalent. Therefore one can choose the probability measure which makes the geometric Brownian motion a martingale ($\mu =0$) to get the call option price c through taking the mathematical expectation of $\left(S\right(T)-K{)}^{+}=c+{\int }_0^{T}H(t\left)\mathrm{d}S\right(t),$ where $H\left(t\right)$ is the hedge, K is the strike price and T is the mature time. Therefore the call option price is just the mathematical expectation of $\left(S\right(T)-K{)}^{+}$ with respect to the risk neutral ($\mu =0$) probability measure. However, some limit property cannot be kept when $T\to \mathrm{\infty }$. For example, $\log S\left(T\right)/T\to -({\sigma }^2/2-\mu )$ depending on μ. Therefore in long run, some results deduced from real market price might be different to those deduced from the risk neutral one.

Let us fix a positive integer T. The payoff of $a_i$ shares of European call option with mature time at $(i+1)T$ and strike price $K_i$ ($i=0,1,2,3,\dots$) is $a_i{V}_{i+1}=a_i\left(S\right((i+1)T)-K_i{)}^{+}.$ Denote by $Q_i$ the option price per share paid at time iT, then a buyer's profit at time $(i+1)T$ will be $\begin{array}{rl}{a}_i{C}_{i+1}& =a_i{V}_{i+1}-a_i{Q}_i\\ & =a_i\left(S\right((i+1)T)-K_i{)}^{+}-a_i{Q}_i.\end{array}$ If we can find a sequence of $\left\{\right(K_i,a_i){\}}_i$ such that $\{a_i{V}_{i+1}{\}}_i$ form a semi-stationary sequence, then we may apply Theorem 2.5 to use its moving averages to estimate the payoffs $\{a_i{V}_{i+1}{\}}_i$ based on the market data and taken as $\left\{a_i{Q}_i\right\}$ such that the error terms $\left\{a_i{C}_i\right\}$ are generalised mean reversion.

In order to get such a stationary sequence, we may choose (at time iT) $K_i=kS\left(iT\right)$ and $a_i=1/S\left(iT\right)$ where k is a prefixed positive percentage constant. Denote by $q_i$ the cost (at time iT) for purchasing $a_i$ shares of option, then the terminal payoff and the buyer's profit will be $\begin{array}{rl}{v}_{i+1}& ={\left(\frac{S\left(\right(i+1\left)T\right)}{S\left(iT\right)}-k\right)}^{+}\mathrm{a}\mathrm{n}\mathrm{d}\\ c_{i+1}& ={\left(\frac{S\left(\right(i+1\left)T\right)}{S\left(iT\right)}-k\right)}^{+}-q_i,\end{array}$ which can also be considered the corresponding values percentage to the original asset prices. In practice, in order to get profit percentage continuously, we may fix a large constant C and just trade $C/S\left(jT\right)$ shares of options with strike price $kS\left(jT\right)$ (both rounded out to the nearest adequate digits). In Black–Scholes's model, $\{S\left(\right(i+1\left)T\right)/S\left(iT\right){\}}_i$ are just the exponential functions of the increments of Brownian motion, which are independent and identically distributed Gaussian random variables. So the price $q_i$ is the mean of $(S\left(\right(i+1\left)T\right)/S\left(iT\right)-k{)}^{+}$ under the risk-neutral measure, which is a constant. Therefore $\{v_i{\}}_i$ and $\{c_i{\}}_i$ are both strongly stationary processes in Black–Scholes' theory. Furthermore, the error terms $\{c_i{\}}_i$ are independent identically distributed in Black–Scholes' formula, of which the time average converges slower than that of generalised mean reversion process.

By Lemma 2.6, if $\mathrm{C}\mathrm{o}\mathrm{v}(v_1,v_j)\to 0,$ and $\{v_{j_k}{\}}_k$ is a subsequence such that $\underset{k}{min}\{j_{k+1}-j_k\}$ is sufficiently large, then $1/N\sum _i^N{v}_{j_i}$ will be close to $E\left[v_1\right]$. Its special case is the logarithmic return of Heston model (Heston, 1993), which is strongly stationary and ergodic, $\mathrm{C}\mathrm{o}\mathrm{v}(c_1,c_N)\to 0$ (when $N\to \mathrm{\infty }$).

In Black–Scholes–Merton's theory, the only undetermined factor is the volatility σ. In practice, the investors quite often use the near term sample variance to estimate the volatility at time t. That is, (4) $\begin{array}{rl}{\hat{\sigma }}^2\left(t\right)& =\frac{1}{n-1}\sum _{j=t-n}^{t-1}\left[\phantom{\sum _{i=t-n}^{t-1}}(\log S(j+1)-\log S(j\left)\right)\right.\\ & -{\left.\frac{1}{n}\sum _{i=t-n}^{t-1}(\log S(i+1)-\log S(i\left)\right)\right]}^2\end{array}$(4) for large n. However, if we use the above $\hat{\sigma }\left(n\right)$ with large n for a unified volatility in Black–Scholes formula, the result would be not good. So very often the investors use smaller n to estimate the volatility of $S\left(t\right)$ near t, which is just Heston's price with $\rho =0$ (Heston, 1993; Poon & Granger, 2003).

From the above discussion, we have four different candidates for option prices: (1) the percentage prices $\left\{q_i\right\}$ are the moving averages $\left\{M\right(iT\left)\right\}$; (2) $\left\{q_i\right\}$ are the historical call option prices in the real market; (3) $\left\{q_i\right\}$ are Black–Scholes' prices with unified volatility $\hat{\sigma }$; (4) $\left\{q_i\right\}$ are Black–Scholes' prices with moving volatility $\hat{\sigma }\left(t\right)$. We will use Market data to show that the cumulated error is stationary in Case (1). We will also illustrate the cumulated errors in the other three cases for comparison.

4. Statistics with real market data

The strong stationarity of logarithmic return is implied in the most popular mathematical models for asset prices. If we choose the trading volume inversely proportional to the current price, then the consecutive payoffs will be strongly stationary as discussed in the previous section. Thus we can use the moving averages $\left\{M\right(iT\left)\right\}$ to estimate a this sequence and the errors will be mean reversion according to Theorem 2.5. In this section, we use the data of ETFs DIA, QQQ, SPY for main stock indices in the US market to illustrate the cumulated buyer's profit percentage, which is also the cumulated percentage errors, to compare with that of Black–Scholes'.

We consider the buyers successive profits $c_{i+1}={\left(\frac{S\left(\right(i+1\left)T\right)}{S\left(iT\right)}-k\right)}^{+}-q_i$ for the prices $\left\{q_i\right\}$ defined in the four cases listed at the end of the previous section, which are also the cumulated errors of the corresponding option pricing.

The first thing is to test the stationarity. We have checked by the ADF test with Matlab software package and the hypothesis that ‘$\left\{c_{i+1}\right\}$ has a unit root’ is rejected with 95% of confidence, for the trading data of SPY (the ETF tracing S&P500), QQQ (the ETF tracing NASDAQ) and DIA (the ETF tracing DOW) under various k and T. In other words, $\left\{c_{i+1}\right\}$ can be statistically considered as a stationary process. For example, when we take $S\left(t\right)$ as the price of SPY, take k = 0.98 and T = 7, Figure 1 shows the successive weekly trading profit sequence $\left\{c_{i+1}\right\}$ when $q_i=M\left(iT\right)$ (m = 6) from 6th January 1995 through the end of July 2015.

Figure 1. The successive profit sequence $\left\{X\right(j\left)\right\}$ from 6th January 1995.

We also found that if $q_i$ is the weekly expired percentage option prices in the real market (with $kS\left(jT\right)$ rounded out to the nearest strike price in the market), $\left\{c_{i+1}\right\}$ still pass the ADF test for SPY, QQQ and DIA data with various k and T.

The parameter k may assume any value. For simplicity, we just list the results for k = 0.98, 1 and 1.02 here. We discuss two time scales for T: (a) weekly expired options; (b) monthly expired options. For weekly expired options, we consider to trade call options between the closing time of neighbouring Fridays. Thus if we purchase $1/S\left(jT\right)(j=1,2,\dots ,n)$ unit of call option at time jT with mature time $(j+1)T$ and strike price $K=kS\left(jT\right)$, the profit by time $(j+1)T$ will be $c_{j+1}$, and the cumulated profit by time $(j+1)T$ is $\sum _{i=1}^j{c}_{i+1}$. In that case, we can take T = 7 calendar days except for very few holidays. Similarly, when we consider the monthly expired options, we trade the options at the 3rd Friday of neighbouring months. In that case, actually T = 30 or 31 (except for February). We compare the cumulated profits of the following four prices:

1. $\left\{q_i\right\}$ are the moving average $\left\{M\right(iT\left)\right\}$;

2. $\left\{q_i\right\}$ are the historical call option price in the real market;

3. $\left\{q_i\right\}$ are Black–Scholes' price with unified volatility $\hat{\sigma }$;

4. $\left\{q_i\right\}$ are Black–Scholes' price with moving volatility $\left\{\hat{\sigma }\right(iT\left)\right\}$.

We should mention some technical points here: (a) when there was no trade at the closing time, we use the mean of bid and ask price as our last price; (b) the strike prices $\left\{kS\right(jT\left)\right\}$ are rounded out to the nearest available strike prices in the market; (c) all our data are percentage priced; (d) the parameter n in (4(4) $\begin{array}{rl}{\hat{\sigma }}^2\left(t\right)& =\frac{1}{n-1}\sum _{j=t-n}^{t-1}\left[\phantom{\sum _{i=t-n}^{t-1}}(\log S(j+1)-\log S(j\left)\right)\right.\\ & -{\left.\frac{1}{n}\sum _{i=t-n}^{t-1}(\log S(i+1)-\log S(i\left)\right)\right]}^2\end{array}$(4) ) is chosen according to the usual method applied in the real market (Poon & Granger, 2003) and the selection of m in (2(2) $\begin{array}{r}{M}^{\left(m\right)}\left(t\right)=\frac{X\left(t\right)+X(t-1)+\cdots +X(t-m+1)}{m}\end{array}$(2) ) is less crucial. Actually we found that there is no much differences when m is chosen between 2 and 12.

4.1. SPY analysis

Figure 2 shows the daily price of SPY from 3rd January 1995 through 2nd July 2015.

Figure 2. The daily price of SPY from 3rd January 1995 through 2nd July 2015.

We have SPY trading data of monthly expired options from 21st January 2005 through 20th March 2015, and its trading data of weekly expired options from 1st July 2011 through 27th March 2015. Thus we can show our cumulative profits of corresponding options in those two periods respectively.

4.1.1. SPY results under different k values

Figures 3–5 show the cumulated profits of both monthly and weekly expired options of SPY under different k values.

Figure 3. The cumulated profits of SPY options when k = 0.98: (a) monthly expired options and (b) weekly expired options.

Figure 4. The cumulated profits of SPY options when k = 1: (a) monthly expired options and (b) weekly expired options.

Figure 5. The cumulated profits of SPY options when k = 1.02: (a) monthly expired options and (b) weekly expired options.

4.2. QQQ analysis

Figure 6 shows the daily price of QQQ from 10th March 1999 through 2nd July 2015.

Figure 6. The daily price of QQQ from 10th March 1999 through 2nd July 2015.

We have QQQ trading data of monthly expired options from 15th February 2002 through 17th April 2015, and trading data of weekly expired options from 1st July 2011 through 27th March 2015.

4.2.1. QQQ results under different k values

Figures 7–9 show the cumulated profits of both monthly and weekly expired options of SPY under different k values.

Figure 7. The cumulated profits of QQQ options when k = 0.98: (a) monthly expired options and (b) weekly expired options.

Figure 8. The cumulated profits of QQQ options when k = 1: (a) monthly expired options and (b) weekly expired options.

Figure 9. The cumulated profits of QQQ options when k = 1.02: (a) monthly expired options and (b) weekly expired options.

4.3. DIA analysis

Figure 10 shows the daily price of DIA from 20th January 1998 through 2nd July 2015.

Figure 10. The daily price of DIA from 20th January 1998 through 2nd July 2015.

We have DIA trading data of monthly expired options from 21st June 2002 through 17th April 2015, and trading data of weekly expired options from 10th August 2012 through 27th March 2015.

4.3.1. DIA results under different k values

Figures 11–13 show the cumulated profits of both monthly and weekly expired options of SPY under different k values.

Figure 11. The cumulated profits of DIA options when k = 0.98: (a) monthly expired options and (b) weekly expired options.

Figure 12. The cumulated profits of DIA options when k = 1: (a) monthly expired options and (b) weekly expired options.

Figure 13. The cumulated profits of DIA options when k = 1.02: (a) monthly expired options and (b) weekly expired options.

5. Possibility of statistical arbitrage

From the above data analysis, we can easily find that the errors of our moving average pricing of option are mean reversion and the cumulated error is stationary. When the strike price is sufficiently low, then the asset price will be always above the strike price and the call option will be always executed. So that the buyer at Black–Scholes' option price can always take statistical arbitrage in long run. Let us take T = 7. If one investor continuously purchase at time iT call option for $1/S\left(iT\right)$ share expired at $(i+1)T$ according to Black–Scholes' price, then Figure 14 shows his cumulated profit in the last 22 years (01/29/1993–07/02/2015) when k = 0 and k = 0.9.

Figure 14. Cumulated gain against Black–Scholes price: (a) when k = 0 and (b) when k = 0.9.

However, there are still many problems left. For example, what happens if only some of the profits $\left\{c_{i+1}\right\}$ occur? That is, the partial sum of a strongly stationary sequence may not have a convergent mean if its terms are not uncorrelated as shown in Section 2. However, we can apply Lemma 2.6 to construct statistical arbitrages if the covariance function of $\left\{c_{i+1}\right\}$ tends to 0 quick enough.

We show here another more practical example to get statistical arbitrage from Black–Scholes' call option. The following graph shows that when $q_i=M\left(iT\right),$ m = 6 and T = 7, the sample autocorrelation function of $c_{i+1}$ tends to 0 (Figure 15).

Figure 15. The autocorrelation function.

Thus if $c_{i_k}$ is a subsequence such that $\underset{k}{min}\{i_{k+1}-i_k\}$ is sufficiently large, we can apply Lemma 2.1 to get their mean sufficiently close to 0. Denote by $Q\left(iT\right)$ the Black–Scholes' price. Our strategy is: (a) when $M\left(iT\right)-Q\left(iT\right)>0.01$ and the previous trade was made at least 9 weeks ago, buy $1/S\left(iT\right)$ of call option at $Q\left(iT\right)$; (b) when $Q\left(iT\right)-M\left(iT\right)>0.0005$ and the previous trade was made at least 9 weeks ago, sell $1/S\left(iT\right)$ call option at $Q\left(iT\right)$. Figure 16 shows our cumulated profit and Figure 17 gives our profit against Black–Scholes price with moving volatility. Thus we can get statistical arbitrage from the mean, which is similar to the case of high-frequency trading (see Wang & Zheng, 2014).

Figure 16. Cumulated gain against Black–Scholes price.

Figure 17. Cumulated gain against Black–Scholes price with stochastic volatility.

Disclosure statement

No potential conflict of interest was reported by the authors.

Additional information

Notes on contributors

Si Bao

Si Bao is now working in Xiangcai security Co., LTD. She studied for her Ph.D. in School of Statistics from ECNU.

Shi Chen

Shi Chen is a data scientist of PayPal Holding Inc. He received his Ph.D in statistics from ECNU in 2017.

Wei An Zheng

Wei An Zheng is Professor of ECNU and Professor Emeritus of University of California, Irvine, USA.

Yu Zhou

Yu Zhou works in Guotai Junan Securities. He received his Ph.D in statistics from ECNU in 2016.