# A hidden Markov model for statistical arbitrage
in international crude oil futures markets

Viviana Fanelli  Department of Economics, Management and Business Law,
University of Bari Aldo Moro, Italy<EMAIL_ADDRESS>,  Claudio
Fontana  Department of Mathematics “Tullio Levi - Civita”, University of
Padova, Italy<EMAIL_ADDRESS>and  Francesco Rotondi  Department of
Finance, Bocconi University, Italy<EMAIL_ADDRESS>
###### Abstract.

In this work, we study statistical arbitrage strategies in international crude
oil futures markets. We analyse strategies that extend classical pairs trading
strategies, considering the two benchmark crude oil futures (Brent and WTI)
together with the newly introduced Shanghai crude oil futures. We document
that the time series of these three futures prices are cointegrated and we
model the resulting cointegration spread by a mean-reverting regime-switching
process modulated by a hidden Markov chain. By relying on our stochastic model
and applying online filter-based parameter estimators, we implement and test a
number of statistical arbitrage strategies. Our analysis reveals that
statistical arbitrage strategies involving the Shanghai crude oil futures are
profitable even under conservative levels of transaction costs and over
different time periods. On the contrary, statistical arbitrage strategies
involving the three traditional crude oil futures (Brent, WTI, Dubai) do not
yield profitable investment opportunities. Our findings suggest that the
Shanghai futures, which has already become the benchmark for the Chinese
domestic crude oil market, can be a valuable asset for international
investors.

###### Key words and phrases:

Pairs trading; crude oil futures; Shanghai crude oil futures; cointegration;
spread process; mean-reverting process; regime switching; stochastic
filtering.

JEL classification. C51, C58, G11, G15.
Financial support from the University of Padova (research programme STARS StG
PRISMA - “Probabilistic Methods for Information in Security Markets”) is
gratefully acknowledged.
Corresponding author: F. Rotondi.

##  1\. Introduction

Pairs trading strategies represent a well-known instance of statistical
arbitrage strategies, exploiting temporary deviations of the prices of similar
securities from their long-term equilibrium in order to achieve profits when
convergence to the equilibrium is reached (see, e.g., [Vidyamurthy, 2004],
[Elliott et al., 2005], [Gatev et al., 2006]). Since one cannot determine ex
ante when prices are going to realign, such strategies do not constitute pure
arbitrage strategies but rather statistical arbitrage strategies, that are
expected to deliver a profit over a sufficiently long time horizon. The
profitability of this type of strategies presumes a strong long-term
relationship, so that mispricings are temporary and likely to quickly revert
back. This is typically captured through the existence of a cointegration
relation among the considered security prices.

In this work, we analyse statistical arbitrage strategies in international
crude oil futures markets, inspired by pairs trading strategies. The futures
contracts we consider are two standard benchmark crude oil futures, the Brent
and the West Texas Intermediate (WTI), and the Shanghai crude oil futures
introduced in March 2018. It is well known that the Brent and the WTI futures
prices are strongly cointegrated (see, e.g., [Kristoufek and Vosvrda, 2014],
[Cerqueti et al., 2019], [Cerqueti and Fanelli, 2021], [Cotter et al., 2022]),
while preliminary evidence (see [Yang and Zhou, 2020], [Huang and Huang,
2020], [Ji et al., 2022]) shows that the newly introduced Shanghai crude oil
futures cointegrate with them as well. These empirical findings suggest the
possibility of achieving statistical arbitrage opportunities when these three
crude oil futures are jointly traded. This represents the starting point and
the motivation of this work.

The Shanghai International Energy Exchange (INE) introduced crude oil futures
on March 26, 2018. These contracts are known as the Shanghai crude oil futures
and are traded on the Shanghai Futures Exchange (SHFE). They are based on the
INE’s own benchmark, which is a blend of various Middle Eastern and Asian
crude oil grades that reflects the specific needs of China’s oil imports.
Shanghai crude oil futures are denominated in Chinese yuan, differently from
the widely adopted crude oil futures denominated in US dollars. While the
worldwide relevance of the Brent and WTI futures is well established, it is
interesting to investigate how the Shanghai futures, which are nowadays traded
mainly by local investors, can become a profitable asset for international
portfolios. While recent studies (see [Wang et al., 2022] and [Shao and Hua,
2022]) show that Shanghai futures have become the benchmark for Chinese
domestic crude oil markets, their role in global financial markets is still to
be addressed. Indeed, probably for the first time, Shanghai crude oil futures
provide international investors with direct access to China’s oil market,
allowing them to participate in the pricing and trading of crude oil in China.

In this work, we investigate whether the Shanghai crude oil futures can be
efficiently integrated into international portfolios of crude oil futures
contracts in order to achieve profitable investments in the form of
statistical arbitrage strategies. Pairs trading on traditional crude oil
benchmarks has already been studied (see, e.g., [Dunis et al., 2006], [Cummins
and Bucca, 2012], [Baviera and Baldi, 2019]), while [Niu et al., 2023] test
pairs trading strategies on the different Shanghai futures contracts traded at
the INE. To the best of our knowledge, this work is the first to consider
statistical arbitrage strategies involving simultaneously traditional crude
oil futures together with the newly introduced Shanghai crude oil futures.

From a methodological viewpoint, we improve on traditional approaches to pairs
trading by considering a stochastic model for the cointegration spread driven
by a mean-reverting process with regime-switching modulated by a hidden Markov
chain. This generalizes the approach of [Tenyakov and Mamon, 2017] and
[Elliott and Bradrania, 2018] and enables us to capture different
cointegration regimes. As the Markov chain is unobserved, stochastic filtering
techniques are employed to estimate the current regime and the model
parameters, similarly as in [Erlwein and Mamon, 2009] and [Erlwein et al.,
2010]. Parameter estimation is done by means of a filter-based version of the
Expectation Maximization (EM) algorithm, as introduced in [Elliott et al.,
1995] (see also [Fontana and Runggaldier, 2010] for an application in credit
risk). The filtering approach enables us to estimate the most likely regime
and the model parameter in a dynamic way, thereby ensuring that the model is
constantly tuned to the actual market evolution. In turn, this will enable us
to construct statistical arbitrage strategies that are constantly updated as
new observations of the cointegration spread become available.

After estimating the model over a training sample, we move out-of-sample and
analyse several types of statistical arbitrage strategies. Our empirical
analysis shows that statistical arbitrage strategies that make use of our
hidden Markov model for the cointegration spread are remarkably profitable
when compared to more traditional strategies based only on the time series of
futures prices. Moreover, we find that including the Shanghai futures in the
portfolio, rather than a standard crude oil benchmark as the Dubai crude oil
futures, delivers a significant investment performance even when large
transaction costs are accounted for. Our findings indicate that the greater
profitability of statistical arbitrage strategies involving the Shanghai
futures can be explained by the much larger speed of adjustment of the
Shanghai futures time series compared to the Brent and to the WTI. This
implies that the Shanghai futures prices tend to revert back quickly to the
long-run relationship, making easier the exploitation of temporary
mispricings. This highlights the potential attractiveness of including the
Shanghai futures in a portfolio of international crude oil futures contracts.

The contribution of the paper is threefold. First, we corroborate the evidence
of cointegration among the relatively new Shanghai crude oil futures and two
long-standing crude oil futures benchmarks, the Brent and the WTI. To the best
of our knowledge, this is the first time that cointegration is assessed
considering the three contracts together, rather than on a pairwise basis.
Second, we propose a mean-reverting stochastic model for the cointegration
spread that allows for regime switching by means of a hidden Markov chain
determining the parameters of the process. We apply filtering techniques to
dynamically estimate the most likely regime and the model parameters. Third,
we study different statistical arbitrage strategies involving three crude oil
futures contracts. Our evidence suggests that strategies including the
Shanghai crude oil futures are remarkably profitable and robust with respect
to several different tests.

The paper is structured as follows. Section 2 describes the cointegration
structure and the stochastic model for the resulting long-run relationship
among the futures prices. Section 3 describes the statistical arbitrage
strategies we implement and how their performance is assessed. Section 4
contains the results of our empirical analysis. In particular, we first
consider a fixed time window and afterwards we assess the robustness of our
results over different time periods and with respect to different levels of
transaction costs. Section 5 concludes. The results of an additional empirical
analysis based on daily data are reported as Supplementary Material.

##  2\. Cointegration analysis and spread modelling

In this section, we introduce the stochastic model for the long-run
relationship among the three crude oil future prices (Section 2.1) and
describe how the model parameters can be estimated by means of filtering
techniques combined with the EM algorithm (Section 2.2).

###  2.1. A mean-reverting hidden Markov model

Let $F^{i}=(F^{i}_{t})_{t\geq 0}$ denote the futures price process, for
$i\in\\{B,S,W\\}$, of the Brent, the Shanghai and the WTI futures,
respectively. The processes $F^{B}$ and $F^{W}$ are typically found to be
cointegrated, as documented by [Lanza et al., 2005], [Hammoudeh et al., 2008],
[Kristoufek and Vosvrda, 2014] among others. Similarly, [Yang and Zhou, 2020]
and [Huang and Huang, 2020] show that $F^{S}$ is cointegrated with $F^{B}$ and
with $F^{W}$. Therefore, it is reasonable to guess that the three futures
price processes jointly considered are cointegrated, similarly to the case of
the Brent, the Dubai and the WTI futures (see [Galay, 2019]). To the best of
our knowledge, cointegration among $F^{B}$, $F^{S}$, $F^{W}$ has not been
tested in previous works and will be shown to hold in the following. In the
presence of cointegration, there exists a linear combination of $F^{B}$,
$F^{S}$, $F^{W}$ that follows a stationary and possibly mean-reverting
process. As usual in the context of cointegration among commodity prices, we
call this linear combination the _spread process_ and denote it by
$S=(S_{t})_{t\geq 0}$. More specifically,

(2.1) $S_{t}:=\lambda^{0}+\sum_{i=B,S,W}\lambda^{i}F_{t}^{i},$

for suitable coefficients $\lambda^{0},\lambda^{B},\lambda^{S},\lambda^{W}$
that constitute the _cointegration vector_.

Using the cointegration vector estimated in-sample an investor can evaluate
day-by-day the spread out-of-sample and set up statistical arbitrage
strategies based on it, taking as portfolio weights the elements of the
cointegration vector. These strategies turn out to be profitable as long as
the cointegrating relationship among the assets is stable, namely as long as
the spread $S$ remains stationary. In a standard cointegration analysis, the
cointegration vector is assumed to be constant. However, previous works on the
Brent and on the WTI document the existence of structural breaks in the
cointegration vector (see, e.g., [Caporin et al., 2019]), showing also that
the cointegration vector depends significantly on the time window over which
the relationship is estimated. As already noticed by [Lee and Papanicolaou,
2016], this may impact pairs trading strategies based on cointegration. For
this reason, it seems appropriate to model the spread as a stochastic process
with regime switching, updating the model parameters dynamically by
incorporating new information as soon as it becomes available. In this work,
we shall adopt this approach, starting from the explicit modelling of the
spread process.

Let $(\Omega,\mathcal{F},\mathbb{F}=(\mathcal{F}_{t})_{t\geq 0},\mathbb{P})$
be a filtered probability space supporting a Brownian motion $W=(W_{t})_{t\geq
0}$ and a Markov chain $\textbf{X}=(\textbf{X}_{t})_{t\geq 0}$ with $N$ states
and transition matrix $\boldsymbol{\Pi}$. The states of the Markov chain X
represent different market regimes, thereby accounting for the time-varying
nature of the cointegrating relationship. The Markov chain X admits a
semimartingale representation of the form

$\textbf{X}_{t}=\textbf{X}_{0}+\int_{0}^{t}\boldsymbol{\Pi}\textbf{X}_{s}\mathrm{d}s+\textbf{M}_{t},$

where $\textbf{M}=(\textbf{M}_{t})_{t\geq 0}$ is a martingale with
$\textbf{M}_{0}=0$ (see [Elliott et al., 1995], Chap. 7).

For the spread process $S$, we assume mean-reverting dynamics with regime-
switching:

(2.2)
$\mathrm{d}S_{t}=a(\textbf{X}_{t})\left(\beta(\textbf{X}_{t})-S_{t}\right)\mathrm{d}t+\xi(\textbf{X}_{t})\mathrm{d}W_{t},\qquad\text{
with }S_{0}=s_{0}\in\mathbb{R},$

where $\alpha$, $\beta$, $\xi$ are real-valued functions on $\mathbb{R}^{N}$.
Without loss of generality and similarly as in [Elliott et al., 1995], we can
assume that the state space of the Markov chain X is the canonical basis of
$\mathbb{R}^{N}$, denoted by $\\{\textbf{e}_{1},\ldots,\textbf{e}_{N}\\}$, so
that $a(\textbf{X}_{t})=\textbf{a}^{\top}\textbf{X}_{t}$ with
$\textbf{a}\in\mathbb{R}^{N}$. For each $i=1,\dots,N$, the component $a_{i}$
represents the speed of mean-reversion of $S$ in the $i^{\text{th}}$ regime.
Analogously, $\beta_{i}$ (resp. $\xi_{i}$) represents the long-run mean (resp.
volatility) of $S$ in the $i^{\text{th}}$ regime. The resulting model for $S$
is a _Ornstein–Uhlenbeck process driven by a hidden Markov model_ (OU-HMM
henceforth) which has been previously used in the literature for the modelling
of interest rates ([Erlwein and Mamon, 2009], [Grimm et al., 2020]), commodity
spot prices ([Erlwein et al., 2010]) and spreads between equity price
processes ([Elliott and Bradrania, 2018]).

Since the current market regime is generally unknown by market participants,
we assume that the Markov chain X is unobservable. As a consequence, the
current state of X has to be filtered from the observations of the spread.
Since we want to implement statistical arbitrage strategies that exploit the
probabilistic structure of the model (2.2) for $S$, we need to estimate
$\boldsymbol{a}$, $\boldsymbol{\beta}$, $\boldsymbol{\xi}$ and the transition
matrix $\boldsymbol{\Pi}$. In this incomplete observation setup, the
estimation of $\boldsymbol{a}$, $\boldsymbol{\beta}$, $\boldsymbol{\xi}$ and
$\boldsymbol{\Pi}$ can be done by a _filter-based expectation maximization
algorithm_ that we now describe.

###  2.2. The filter-based expectation maximization algorithm

We now describe the filter-based expectation maximization (EM) algorithm for
the estimation of the OU-HMM in (2.2). The following results are adapted from
the general techniques described in [Elliott et al., 1995] and also applied to
an Ornstein-Uhlenbeck process in [Erlwein and Mamon, 2009].

First of all, it is necessary to discretize the continuous-time model for the
spread in (2.2). This discretization is needed not only because only discrete
observations of $S$ are available in practice, but also because in a
continuous-time setting the EM algorithm does not allow to estimate
$\boldsymbol{\xi}$, since the EM algorithm is based on an equivalent change of
probability.

Consider the spread process $S$ in (2.2) over a time interval $[t,t+1]$ and
let $\Delta$ denote the time step. If $\mathbf{X}$ is constant over this
interval, we can explicitly compute the solution of (2.2) as

$S_{t+1}=e^{-a(\textbf{X}_{t})\Delta}S_{t}+\beta(\textbf{X}_{t})(1-e^{-a(\textbf{X}_{t})\Delta})+\xi(\textbf{X}_{t})\int_{t}^{t+1}e^{-a(\textbf{X}_{t})(t+1-s)}\mathrm{d}W_{s}.$

From this equation we can derive the discrete-time version
$y=(y_{t})_{t\in\mathbb{N}}$ of the spread:

(2.3)
$y_{t+1}=\gamma(\textbf{X}_{t})+\alpha(\textbf{X}_{t})y_{t}+\eta(\textbf{X}_{t})z_{t},$

where $(z_{t})_{t\in\mathbb{N}}$ is a sequence of i.i.d. standard normal
random variables and

(2.4)
$\begin{array}[]{l}\alpha(\textbf{X}_{t})=e^{-a(\textbf{X}_{t})\Delta},\\\
\gamma(\textbf{X}_{t})=\beta(\textbf{X}_{t})\left(1-e^{-a(\textbf{X}_{t})\Delta}\right),\\\
\eta(\textbf{X}_{t})=\xi(\textbf{X}_{t})\sqrt{\frac{1-e^{-2a(\textbf{X}_{t})\Delta}}{2a(\textbf{X}_{t})}}.\end{array}$

Analogously, we discretize the Markov chain $\mathbf{X}$ as

$\mathbf{X}_{t+1}=\boldsymbol{\Pi}\mathbf{X}_{t}+\mathbf{v}_{t+1},$

with $\mathbf{v}_{t+1}$ denoting a martingale increment. We also introduce the
following auxiliary quantities that will be needed for the recursive
estimation of the model parameters:

  * •

$J^{ij}_{t}:=\sum_{n=1}^{t}\langle\mathbf{X}_{n-1},\mathbf{e}_{i}\rangle\langle\mathbf{X}_{n},\mathbf{e}_{j}\rangle$,
representing the cumulative number of jumps of $\mathbf{X}$ from state
$\mathbf{e}_{i}$ to state $\mathbf{e}_{j}$ until time $t$. Notice that
$J^{ij}_{t}=J^{ij}_{t-1}+\langle\mathbf{X}_{t-1},\mathbf{e}_{i}\rangle\langle\mathbf{X}_{t},\mathbf{e}_{j}\rangle$.

  * •

$O^{i}_{t}:=\sum_{n=1}^{t}\langle\mathbf{X}_{n},\mathbf{e}_{i}\rangle$,
representing the occupation time of $\mathbf{X}$ in state $\mathbf{e}_{i}$
until time $t$. Notice that
$O^{i}_{t}=O^{i}_{t-1}+\langle\mathbf{X}_{t},\mathbf{e}_{i}\rangle$.

  * •

$T^{i}_{t}(f)=\sum_{n=1}^{t}\langle\mathbf{X}_{n-1},\mathbf{e}_{i}\rangle
f_{n}$, where $f_{n}$ is a generic function of the observations of $y$ up to
time $n$ (in our case, the function $f_{n}$ will be given by $f_{n}=y_{n}$,
$f_{n}=y_{n}^{2}$ or $f_{n}=y_{n}y_{n-1}$). Notice that
$T^{i}_{t}(f)=T^{i}_{t-1}(f)+\langle\mathbf{X}_{t-1},\mathbf{e}_{i}\rangle
f_{t}$.

We denote by $\mathbb{F}^{y}=(\mathcal{F}_{t}^{y})_{t\in\mathbb{N}}$, with
$\mathcal{F}^{y}_{t}:=\sigma\left\\{y_{0},\dots,y_{t}\right\\}$ for all
$t\in\mathbb{N}$, the filtration generated by the process $y$. We assume that
$\mathbb{F}^{y}$ represents the information available to the investor, who
cannot observe the Markov chain X and can only access discrete-time
observations of the spread.

We now state the recursive equations for the filters of the unobserved Markov
chain $\mathbf{X}$ and of the quantities $J$, $O$, $T$, as well as the
recursive equations for the EM estimators of the parameters
$\boldsymbol{\alpha}$, $\boldsymbol{\gamma}$, $\boldsymbol{\eta}$ and
$\boldsymbol{\Pi}$. For $t\in\mathbb{N}$, we denote by
$\hat{\mathbf{X}}_{t}:=\mathbb{E}\left[\mathbf{X}_{t}|\mathcal{F}^{y}_{t}\right]$
the filtered estimate of the latent Markov chain at time $t$, obtained using
only the information contained in the observation process $y$ up to time $t$.
In an analogous way we define the filtered estimates $\hat{J}_{t}^{ij}$,
$\hat{O}_{t}^{i}$ and $\hat{T}_{t}^{i}(f)$ of the quantities introduced above.
Moreover, we denote by
$\hat{\boldsymbol{\alpha}}^{(t)}=(\hat{\alpha}_{i}^{(t)})_{i=1,\dots,N}$ the
estimate of the parameter $\boldsymbol{\alpha}$ on the basis of the
information generated by the observation process $y$ up to time $t$. The
estimators $\hat{\boldsymbol{\gamma}}^{(t)}$, $\hat{\boldsymbol{\eta}}^{(t)}$,
$\hat{\boldsymbol{\Pi}}^{(t)}$ are defined in an analogous way.

As explained in [Elliott et al., 1995, Chapter 8], in order to filter the
unobserved Markov chain $\mathbf{X}$ one can resort to a change of measure
passing from $\mathbb{P}$ to an equivalent probability measure
$\tilde{\mathbb{P}}$ under which the observation process $y$ is independent of
the Markov chain $\mathbf{X}$. For each $T\in\mathbb{N}$, the probability
measures $\mathbb{P}$ are $\tilde{\mathbb{P}}$ are related by the Radon-
Nikodym derivative

$\Lambda_{T}=\left.\frac{\mathrm{d}\mathbb{P}}{\mathrm{d}\tilde{\mathbb{P}}}\right|_{\mathcal{F}_{T}}=\prod_{t=1}^{T}\lambda_{t},$

where

(2.5)
$\lambda_{t}:=\frac{1}{\eta(\mathbf{X}_{t-1})}\exp\left(-\frac{1}{2}\left(\frac{\left(y_{t}-y_{t-1}\alpha(\mathbf{X}_{t-1})-\gamma(\mathbf{X}_{t-1})\right)^{2}}{\eta(\mathbf{X}_{t-1})^{2}}-y_{t}^{2}\right)\right),$

with $\Lambda_{0}=1$. For convenience of notation, it is useful to represent
the possible values of $\lambda_{t}$ in (2.5) associated to each of the $N$
states of $\mathbf{X}$ by a diagonal matrix
$\mathbf{D}_{t}=[d^{ij}_{t}]_{i,j=1,...,N}$ with

(2.6)
$d^{ii}_{t}=\displaystyle\frac{1}{\eta_{i}}\exp\left(-\frac{1}{2}\left(\frac{\left(y_{t}-y_{t-1}\alpha_{i}-\gamma_{i}\right)^{2}}{\eta_{i}^{2}}-y_{t}^{2}\right)\right),\qquad\text{
for }i=1,\ldots,N.$

The starting point of the parameter estimation algorithm is a set of initial
guesses for the quantities $\mathbf{X}_{0}$ (that delivers also the initial
values of the $O^{i}_{0}$’s) and for $\boldsymbol{\Pi}$,
$\boldsymbol{\gamma},\boldsymbol{\alpha}$, $\boldsymbol{\eta}$.

By performing similar computations as in [Erlwein and Mamon, 2009, Section 4],
the recursive filtering equations for $\mathbf{X}$ and for the auxiliary
quantities $J$, $O$, $T$ are given as follows, where we denote by
$\boldsymbol{1}$ the unit vector in $\mathbb{R}^{N}$:

(2.7)
$\begin{array}[]{rl}\hat{\mathbf{X}}_{t+1}&=\boldsymbol{\Pi}\mathbf{D}_{t+1}\hat{\mathbf{X}}_{t},\\\
\hat{J}^{i,j}_{t+1}&=\langle\boldsymbol{1},\boldsymbol{\Pi}\mathbf{D}_{t+1}\hat{J}^{i,j}_{t}+\langle\hat{\mathbf{X}}_{t},\mathbf{e}_{i}\rangle\langle\mathbf{D}_{t+1}\mathbf{e}_{i},\mathbf{e}_{i}\rangle\pi_{ij}\mathbf{e}_{j}\rangle,\\\
\hat{O}^{i}_{t+1}&=\langle\boldsymbol{1},\boldsymbol{\Pi}\mathbf{D}_{t+1}\hat{O}^{i}_{t}+\langle\hat{\mathbf{X}}_{t},\mathbf{e}_{i}\rangle\langle\mathbf{D}_{t+1}\mathbf{e}_{i},\mathbf{e}_{i}\rangle\boldsymbol{\Pi}\mathbf{e}_{i}\rangle,\\\
\hat{T}_{t+1}^{i}(f)&=\langle\boldsymbol{1},\boldsymbol{\Pi}\mathbf{D}_{t+1}\hat{T}_{t}^{i}(f)+\langle\hat{\mathbf{X}}_{t},\mathbf{e}_{i}\rangle\langle\mathbf{D}_{t+1}\mathbf{e}_{i},\mathbf{e}_{i}\rangle
f_{t}\boldsymbol{\Pi}\mathbf{e}_{i}\rangle,\end{array}$

for all $t\in\mathbb{N}$. Observe that the filters in (2.7) depend on
$\boldsymbol{\Pi}$ and, through the matrix $\mathbf{D}_{t+1}$, also on
$\boldsymbol{\gamma}$, $\boldsymbol{\alpha}$, $\boldsymbol{\eta}$. Therefore,
when implementing the filtering procedure in practice these quantities must be
replaced by their estimates $\hat{\boldsymbol{\Pi}}^{(t)}$,
$\hat{\gamma}_{i}^{(t)}$, $\hat{\alpha}_{i}^{(t)}$, $\hat{\eta}_{i}^{(t)}$
obtained at the previous time point $t$. In analogy to [Erlwein and Mamon,
2009, Section 5], the equations for the EM estimates of $\boldsymbol{\Pi}$,
$\boldsymbol{\gamma},\boldsymbol{\alpha}$, $\boldsymbol{\eta}$ at time $t$ are
given by

(2.8)
$\begin{array}[]{rl}\hat{\pi}_{ij}^{(t+1)}&=\displaystyle\frac{\hat{J}^{i,j}_{t+1}}{\hat{O}^{j}_{t+1}},\\\
\hat{\gamma}_{i}^{(t+1)}&=\displaystyle\frac{\hat{T}_{t+1}^{i}(y_{t+1})-\hat{\alpha}_{i}^{(t)}\hat{T}_{t}^{i}(y_{t})}{\hat{O}^{i}_{t+1}},\\\
\hat{\alpha}_{i}^{(t+1)}&=\displaystyle\frac{\hat{T}_{t}^{i}(y_{t+1},y_{t})-\hat{\gamma}_{i}^{(t+1)}\hat{T}_{t}^{i}(y_{t})}{\hat{T}_{t}^{i}(y_{t}^{2})},\\\
\hat{\eta}_{i}^{(t+1)}&=\displaystyle\frac{\hat{T}_{t+1}^{i}(y_{t+1}^{2})+(\hat{\alpha}_{i}^{(t+1)})^{2}\hat{T}_{t}^{i}(y_{t}^{2})+(\hat{\gamma}_{i}^{(t+1)})^{2}\hat{O}^{i}_{t+1}}{\hat{T}_{t}^{i}(y_{t}^{2})}+\\\
&\quad+\displaystyle\frac{\hat{\gamma}_{i}^{(t+1)}\hat{T}_{t+1}^{i}(y_{t+1})-2\hat{\alpha}_{i}^{(t+1)}\hat{T}_{t}^{i}(y_{t+1},y_{t})-2\hat{\alpha}_{i}^{(t+1)}\hat{\gamma}_{i}^{(t+1)}\hat{T}_{t}^{i}(y_{t})}{\hat{O}^{i}_{t+1}}.\\\
\end{array}$

It is worth noticing that every time a new observation becomes available all
the filters/estimates in (2.7) and (2.8) are updated without re-running the
whole filtering algorithm, due to its recursive nature. This is an advantage
of the online filter-based EM algorithm, which allows the model to be
constantly tuned to the actual market information.

######  Remark 2.1 (Numerical aspects).

In this remark we address a number of practical numerical issues that arise
when implementing the above filtering algorithm:

  1. (1)

As far as the initialization of the algorithm is concerned, the guess for the
initial state of the $N$-state Markov chain $\mathbf{X}_{0}$ determines also
$O^{i}_{0}$. In our empirical analysis we set by convention
$\hat{\mathbf{X}}_{0}=\mathbf{e}_{1}$. This implies $\hat{O}^{1}_{0}=1$ and
$\hat{O}^{i}_{0}=0$ for $i=2,\dots,N$. With respect to the initial guess for
$\boldsymbol{\Pi}$, if $N=2$ we set $\hat{\pi}_{11}^{(0)}=0.6$ and
$\hat{\pi}_{22}^{(0)}=0.5$ since the perfectly symmetric initial guess
$\hat{\pi}_{11}^{(0)}=\hat{\pi}_{22}^{(0)}=0.5$ tends to create numerical
instabilities as there is no difference across states in the transition
matrix. If $N=3$ we set $\hat{\pi}_{11}^{(0)}=0.5,\hat{\pi}_{12}^{(0)}=0.25$,
$\hat{\pi}_{21}^{(0)}=0.3,\hat{\pi}_{22}^{(0)}=0.4$ and
$\hat{\pi}_{31}^{(0)}=0.2,\hat{\pi}_{33}^{(0)}=0.6$. The initial guesses for
$\boldsymbol{\gamma},\boldsymbol{\alpha}$ and $\boldsymbol{\eta}$ can be
derived from the data. More specifically, if the Markov chain $\mathbf{X}$ was
constant, the process $y$ would follow a simple AR(1) process with scalar
parameters satisfying

(2.9) $y_{t+1}=\gamma+\alpha y_{t}+\varepsilon_{t+1},$

where $\varepsilon=(\varepsilon_{t})_{t\in\mathbb{N}}$ is a Gaussian white
noise process with variance $\eta^{2}$. In this case, the parameters $\gamma$,
$\alpha$, $\eta$ can be easily estimated by ordinary least squares. In our
empirical analysis we do so using the first $20$ datapoints of $y$,
corresponding to one month of daily observations (this number is twice the
_batch parameter_ $m$ that we introduce below and that we set equal to $10$).
For the case $N=1$ we therefore set $\hat{\gamma}^{(0)}=\hat{\gamma}_{OLS}$,
$\hat{\alpha}^{(0)}=\hat{\alpha}_{OLS}$ and
$\hat{\eta}^{(0)}=\sqrt{\hat{\eta}^{2}_{OLS}}$, where $\hat{\gamma}_{OLS}$,
$\hat{\alpha}_{OLS}$ and $\hat{\eta}^{2}_{OLS}$ denote the ordinary least
square estimates of the parameters of (2.9). For the case $N=2$ we let the
initial guess for the two states for all the parameters be equally spaced with
respect to the OLS estimates determined as above. Therefore, we set
$\hat{\gamma}_{1}^{(0)}=1.3\hat{\gamma}_{OLS}$ and
$\hat{\gamma}_{2}^{(0)}=0.7\hat{\gamma}_{OLS}$ and we do the same for the
other two parameters. For the case $N=3$ we follow the same reasoning and set
$\hat{\gamma}_{1}^{(0)}=1.3\hat{\gamma}_{OLS}$,
$\hat{\gamma}_{2}^{(0)}=\hat{\gamma}_{OLS}$ and
$\hat{\gamma}_{3}^{(0)}=0.7\hat{\gamma}_{OLS}$, and we do the same for the
other two parameters. We want to point out, however, that the EM algorithm is
robust with respect to the specification of the initial guesses.

  2. (2)

During the first iterations of the algorithm it might happen that the quantity
$\hat{O}^{i}_{t+1}$, which appears at the denominator in
$\hat{\gamma}_{i}^{(t+1)}$ and $\hat{\eta}_{i}^{(t+1)}$ in (2.8), stays equal
to zero. This happens if the Markov chain has never visited state $i$ before
time $t+1$. If that is the case, there is no way to update
$\hat{\gamma}_{i}^{(t+1)}$ and $\hat{\eta}_{i}^{(t+1)}$, which are left
unchanged throughout the current iteration of the algorithm.

  3. (3)

Quantities like $\hat{T}_{t}^{i}(y_{t-1}^{2})$, which appears at the
denominator in $\hat{\alpha}_{i}^{(t+1)}$ and $\hat{\eta}_{i}^{(t+1)}$ in
(2.8), might take very small values (especially when the $y$’s are close to
zero) and this might create numerical instabilities. The same kind of
numerical instabilities might arise when computing the diagonal elements of
$\mathbf{D}$ if the $\hat{\eta}_{i}^{(t+1)}$’s, which appear at the
denominator of the $d_{t}^{ii}$’s in (2.6), are too close to zero. If not
controlled for, these instabilities propagate, preventing the algorithm to
reach convergence. Therefore, for all denominators in (2.8) we control that
the new estimate of the quantity of interest is not smaller/larger than ten
times the previous estimate and, in case, we truncate the new estimate to that
level.

  4. (4)

The two previous issues can be mitigated by adopting a _batchwise_ approach,
as pointed out in [Grimm et al., 2020]. According to this approach, the
parameters in (2.8) are not updated at every time point $t$ as the filters in
(2.7) but rather every $m$ steps. This technique stabilizes the estimates of
the filtered quantities in (2.7) and, as a consequence, the ones of the
parameters in (2.8). In our analysis we set $m=10$, which in practice implies
that the model parameters are updated every two weeks.

##  3\. Statistical arbitrage strategies

In this section we present the statistical arbitrage strategies we test on
data (Section 3.1) and we describe the performance measures we use to assess
their profitability (Section 3.2).

###  3.1. Statistical arbitrage strategies

As explained in [Vidyamurthy, 2004, Chapter 8], the basic intuition behind
pairs trading is that whenever the spread deviates “sufficiently” from its
equilibrium value, the investor should open a trading position, appropriately
investing in the underlying assets. This position is then closed when this
divergence corrects itself and the spread goes back to its equilibrium value.
The explicit specification of what we mean above by “sufficiently” determines
the precise trading strategy, which can be summarized by a simple rule based
on _opening_ and _closing_ signals.

Regardless the way in which opening/closing signals are determined, another
relevant point is how to build the trading position once the opening signal is
received. In our analysis we adopt a standard approach and we assume that the
trading position is determined by the cointegration vector
$\lambda=(\lambda^{0},\lambda^{B},\lambda^{S},\lambda^{W})$ in (2.1), where
$\lambda^{0}$ determines the borrowing/lending of money from a riskless money
market account (earning zero interest rate, for simplicity of presentation).
Abstracting from transaction costs, if an opening signal is received at time
$t$ when $S_{t}>0$ (resp. $S_{t}<0$), the investor has to go short (resp.
long) on portfolio $\lambda$. This generates an inflow of money at time $t$
equal to $|S_{t}|$. Closing this position when the spread reverts back to zero
will deliver no cashflow. Therefore, these strategies are designed to deliver
positive payoffs at their opening, with zero cashflows when they are closed,
similarly to standard arbitrage opportunities (of the second kind). However,
since it is not known ex ante if and when $S_{t}$ will revert back to zero,
these strategies are not pure arbitrage opportunities but rather statistical
arbitrage strategies (see [Bondarenko, 2003] and the recent generalizations in
[Rein et al., 2021] for a mathematical formalization of the concept of
statistical arbitrage strategies).

In our analysis we shall consider five different statistical arbitrage
strategies. The first one is a simple benchmark of pairs trading and requires
no explicit modelling of the spread process. The second and the fourth rules
rely only on sample estimates of the moments of the spread. On the contrary,
the third and the fifth rules exploit our stochastic model for the spread, as
they are based on the one-step ahead forecast of the spread and on its
forecast interval. While the timing of opening signals differs across the five
different rules, the closing one is assumed to be the same for simplicity of
comparison. In particular, if there was an open position at time $t-1$, the
position is closed at time $t$ if the spread changes sign (in discrete time,
this corresponds to the first passage time of the spread at zero, coherently
with the above description).

Strategy 1: plain vanilla (PV). This first strategy is a benchmark in the
pairs trading literature and dates back to [Burgess, 1999]. According to this
rule, the investor should open the position as soon as the spread differs from
zero. Following this strategy, the investor exploits every minimal deviation
of the assets from their long-run relationship and might be highly profitable.
However, this rule involves a lot of trading (since the spread oscillates
frequently around zero) and might not be profitable when transaction costs are
taken into account.

Formally, assuming that there is no open position at $t-1$, a position is
opened at $t$ if $S_{t}\neq 0$.

Strategy 2: probability interval (ProbI). If the investor wants to avoid too
much trading and save on transaction costs, she should open a position only
when the spread deviates significantly from zero. One common choice to
quantify this deviation (see, e.g., [Endres and Stübinger, 2019]) is to trade
when the spread exceeds the so-called Bollinger band ([Bollinger, 2001]),
which coincides with the 95% probability interval under a normality
assumption, where the two moments of the distribution are estimated
dynamically over a rolling window of $n$ days. This improves on the related
strategy analyzed in the seminal work of [Gatev et al., 2006], where the bands
are just assumed to be constant and equal to twice the sample standard
deviation over the entire test sample. The ProbI strategy generalizes the
strategy analyzed in [Bock and Mestel, 2009], who consider two different bands
associated to two alternating market regimes.

Formally, assuming that there is no open position at $t-1$, a position is
opened at $t$ if $S_{t}\notin\left(\hat{\mu}_{t-1-n:t-1}\pm
q_{\alpha}\hat{\sigma}_{t-1-n:t-1}\right)$, where $\hat{\mu}_{n_{1}:n_{2}}$
(resp. $\hat{\sigma}_{n_{1}:n_{2}}$) is the sample estimate of the expected
value (resp. standard deviation) estimated using the datapoints from $n_{1}$
to $n_{2}$ and $q_{\alpha}$ is the $\alpha$-quantile of a standard normal
random variable. In our empirical analysis we use the $0.975^{\text{th}}$
quantile. Notice that this strategy is equivalent to the one based on the
z-score (considered for instance in [Avellaneda and Lee, 2010]), defined as
$z_{t}:=\frac{S_{t}-\hat{\mu}_{t-1-n:t-1}}{\hat{\sigma}_{t-1-n:t-1}}$, that
prescribes to open the position at $t$ if
$z_{t}\notin(-q_{\alpha},+q_{\alpha})$.

Strategy 3: prediction interval (PredI). By exploiting our stochastic model
for the spread process, we can improve on the previous strategy by adopting a
forward-looking approach, inspired by the approach of [Elliott et al., 2005].
According to the PredI rule, the band the investor should consider is not the
backward-looking confidence interval considered in strategy 2 above, but
rather the prediction interval where the first two moments of the distribution
are computed according to the probabilistic model described in Section 2.1.

Formally, assuming that there is no open position at $t-1$, a position is
opened at $t$ if $S_{t}\notin(\mathbb{E}[S_{t}|\mathcal{F}^{y}_{t-1}]\pm
q_{\alpha}\sqrt{\mathbb{V}ar[S_{t}|\mathcal{F}^{y}_{t-1}]})$, where
$\mathbb{E}[S_{t}|\mathcal{F}^{y}_{t-1}]$ is the one step-ahead forecast of
$S_{t}$ and $\mathbb{V}ar[S_{t}|\mathcal{F}^{y}_{t-1}]$ its variance, both
computed by relying on the filtering algorithm described in Section 2.2 on the
basis of the information gneerated by the observation of the discrete-time
spread process $y$ itself. More specifically, using the model of Section 2, we
have that

$\mathbb{E}[S_{t}|\mathcal{F}^{y}_{t-1}]=\hat{\gamma}^{(t-1)}(\hat{\mathbf{X}}_{t-1})+\hat{\alpha}^{(t-1)}(\hat{\mathbf{X}}_{t-1})S_{t-1}\qquad\text{and}\qquad\mathbb{V}ar[S_{t}|\mathcal{F}^{y}_{t-1}]=(\hat{\eta}^{(t-1)}(\hat{\mathbf{X}}_{t-1}))^{2}.$

Observe that these quantities depend on the filtered estimate of the latent
Markov chain X. Similarly to the ProbI strategy described above, we use the
$0.975^{\text{th}}$ quantile.

Strategy 4: realized increment (RI). As suggested by [Dunis et al., 2006],
another way to avoid too much trading and save on transaction costs is to look
for increments of the spread that are significantly “larger than usual”.
According to this strategy, the investor has to first compute the time series
of the spread increments, that we label by $x=(x_{t})_{t\in\mathbb{N}}$ with
$x_{t}:=S_{t}/S_{t-1}-1$. Then, she should open a position at $t$ if the
increment $x_{t}$ is significantly larger than the previous ones. One simple
way to formalize this is to rely on two historical quantiles.

Formally, assuming that there is no open position at $t-1$, a position is
opened at $t$ if
$x_{t}\notin\left(q_{\underline{\alpha},x},q_{\overline{\alpha},x}\right)$
with $\underline{\alpha}<\overline{\alpha}$ and where $q_{\alpha,x}$ is the
empirical quantile of order $\alpha$ of variable $x$. In our analysis we
consider the $\underline{\alpha}=0.025^{\text{th}}$ and the
$\overline{\alpha}=0.975^{\text{th}}$ quantiles.

Strategy 5: predicted increment (PI). By exploiting the stochastic model for
the spread process described in Section 2, we can improve on strategy RI
described above. More specifically, the investor can compute the series
$\hat{x}=(\hat{x}_{t})_{t\in\mathbb{N}}$ of predicted spread increments, with
$\hat{x}_{t}:=\mathbb{E}[S_{t+1}|\mathcal{F}^{y}_{t}]/S_{t}-1$, and trade
whenever the predicted increment lies outside a given band.

Formally, assuming that there is no open position at $t-1$, a position is
opened at $t$ if
$\hat{x}_{t}\notin\left(q_{\underline{\alpha},\hat{x}},q_{\overline{\alpha},\hat{x}}\right)$.
In our analysis we consider the $\underline{\alpha}=0.025^{\text{th}}$ and the
$\overline{\alpha}=0.975^{\text{th}}$ quantiles and we compute the one step-
ahead forecast $\mathbb{E}[S_{t+1}|\mathcal{F}^{y}_{t}]$ as described in
Strategy 3 above.

######  Remark 3.1 (Alternative statistical arbitrage strategies).

Other statistical arbitrage strategies inspired by pairs trading proposed in
the literature involve: moving-average trading strategies based on the
difference between a short and a long moving average for the spread (see
[Alizadeh and Nomikos, 2008]); strategies based on first-passage times when
the spread is modelled as a standard Ornstein-Uhlenbeck process (see [Bertram,
2010] and [Cummins and Bucca, 2012]); strategies constructed by relying on a
machine learning approach that analyses possible strategies within a training
sample (see [Sarmento and Horta, 2021]).

Finally, we shall also consider a passive strategy that does not incur into
transaction costs: a buy-and-hold position on the S&P500 index that will serve
as a benchmark (S&P strategy).

###  3.2. Performance measures

As pointed out in [Gatev et al., 2006], computing returns for trading
strategies based on the concept of pairs trading is non-trivial due to several
factors that complicate the calculation and interpretation of performance
metrics. First of all, there are days in which a trading position is open and
days in which there is no open position. Moreover, a successfully executed
trade delivers a positive cashflow when opened and a non-negative cashflow
when closed, thus preventing the computation of a standard linear return. To
overcome these issues, it is common to rely on the daily mark-to-market profit
& loss indicator.

Let us assume that the strategy is applied over a time period spanning $n$
days. Following [Burgess, 1999], for each day $t=1,\ldots,n$ we compute a
daily return $r_{t}$ in the following way:

  * •

if there is no open position at $t$, then $r_{t}=0$;

  * •

if there is an open position at $t$, then

(3.1)
$r_{t}=S_{t}\frac{S_{t}-S_{t-1}}{V_{t}^{S}}-c\left|S_{t}-S_{t-1}\right|,$

where $V_{t}^{S}=|\lambda^{0}|+\sum_{i=B,S,W}|\lambda^{i}|F_{t}^{i}$
represents the total market exposure of the portfolio components and $c$
represents constant proportional transaction costs, as explained below.

The overall performance of a given strategy is then evaluated by computing

(3.2) $R=\prod_{t=1}^{n}(1+r_{t})-1,$

which can be interpreted as the return over a time period spanning $n$ days.
As in the empirical analysis reported in Section 4 the length of the time
period will vary, we always report the returns in (3.2) in annual terms:
following the usual convention, we multiply $R$ by $250/n$.

When evaluating the performance of statistical arbitrage strategies,
transaction costs play an important role, as discussed in [Do and Faff, 2012].
Indeed, these strategies typically involve frequent trading, several positions
and possibly considerable margins to be posted. Following a common approach in
the pairs trading literature, we take into account these market frictions by
introducing in (3.1) the parameter $c$ representing proportional transaction
costs. This parameter is difficult to estimate. Focusing on the literature on
pairs trading on USD-denominated crude oil futures, we found a minimum value
of $c=20$bps in [Alizadeh and Nomikos, 2008] and a maximum value of $c=60$bps
in [Do and Faff, 2012]. In the present work, we consider trading strategies
that include the Shanghai futures, which is denominated in CNY and therefore
incurs into currency exchange costs. When including the Shanghai futures in
our statistical arbitrage strategies, we assume $c=80$bps. This represents a
significant amount of transaction costs, which will lead to rather
conservative performance measures for our trading strategies.

Unlike pure arbitrage opportunities, statistical arbitrage strategies are not
immune from financial risk. For this reason, we shall also consider the Sharpe
ratio, a classical risk-adjusted performance measure, computed as follows:

(3.3)
$SR=\frac{\frac{1}{n}\sum_{t=1}^{n}r_{t}}{\sqrt{\frac{1}{n}\sum_{t=1}^{n}\left(r_{t}-\frac{1}{n}\sum_{v=1}^{n}r_{v}\right)^{2}}}.$

For the sake of comparability, analogously to the case of the returns in
(3.2), we shall always express the Sharpe Ratio of a given trading strategy in
annual terms. In Section 4.2.1, we will perform an additional analysis by also
computing the value-at-risk of the different statistical arbitrage strategies
described in Section 3.1.

##  4\. Empirical analysis

This section is divided into four sections. In Section 4.1, we perform the
cointegration analysis among the three futures prices over a fixed time
period. In Section 4.2, we evaluate and compare the performance of the
statistical arbitrage strategies introduced in Section 3.1. In Section 4.3, we
verify the robustness of our findings by testing the statistical arbitrage
strategies over different time periods and with respect to different levels of
transaction costs. Finally, in Section 4.4, we replace the Shanghai futures
contract with the Dubai futures and evaluate the performance of the
corresponding statistical arbitrage strategies.

###  4.1. Analysis of statistical arbitrage strategies over a fixed time
window

####  4.1.1. Data description

We consider daily and weekly futures prices at one month maturity for the
Brent, the WTI and the Shanghai futures, converted from CNY to USD. While the
time series for the Brent and the WTI start in the eighties, the Shanghai has
been trading only since March 2018. Therefore, our dataset spans from $t_{0}=$
03/26/2018 to $T=$ 06/30/2023, including 1374 daily observations and 274
weekly ones. The weekly futures prices are displayed in Figure 1.

Figure 1. Weekly futures prices of the three contracts. The black dashed line
corresponds to 07/01/2022, which separates the training sample from the test
sample.

We split the full dataset into a training sample and a test sample. We label
by $t_{B}$ the first date of the test sample, which in our analysis is set at
07/01/2022, represented by a black dashed line in Figure 1. Therefore, the
test sample includes a year of observations, corresponding to 261 (resp. 52)
daily (resp. weekly) observations.111The dataset under analysis is relatively
short and recent and, as such, significantly affected by the market turmoil of
the COVID-19 pandemics. However, given that the Shanghai futures is very
recent, this cannot be avoided.

####  4.1.2. Unit roots, VAR and cointegration analysis

We now study the dynamic relationship among the three futures prices by
performing a cointegration analysis based on weekly data, following the usual
steps of cointegration analysis (see, e.g., [Guidolin and Pedio, 2018, Chapter
4]). The choice of using weekly data is typical in the literature as daily
data tend to be more volatile and noisy due to several factors such as
intraday trading, news releases, and other short-term fluctuations. In the
Supplementary Materials, we report the results of the cointegration analysis
performed on daily data: while our main findings are not affected by the
choice of weekly or daily data, we verify that the performances of our
statistical arbitrage strategies are significantly better when weekly data are
used for the cointegration analysis.

Before assessing the presence of one or more cointegrating relationships among
the three time series of futures prices, we test whether each of them is
integrated of some order. In order to do so, we run standard unit
root/stationarity tests on each of the three series individually. In
particular, we run the ADF and the PP tests ([Said and Dickey, 1984],
[Phillips and Perron, 1988]) for the presence of a unit root. Moreover, we
account for possible structural breaks at an unknown time by running also the
breakpoint Perron unit root test ([Perron, 1997]). The p-values of these tests
are reported in Table 1.

| Brent | Shanghai | WTI
---|---|---|---
ADF | 0.9434 | 0.8507 | 0.9422
PP | 0.8790 | 0.8386 | 0.8711
bP | 0.9548 | 0.9402 | 0.9416
Table 1. P-values of the Augmented Dickey-Fuller (ADF), the Phillips-Perron
(PP) and the breakpoint Perron (bP) tests. The number of lag length for the
ADF and the bP is chosen according to the Schwartz IC.

As can seen from the table, there is a strong evidence that the three series
contain (at least) a unit root. The same tests on the first differences of the
three series deliver p-values all equal to 0.0000, indicating that the series
contain only one unit root. Therefore, we can conclude that the three series
are integrated of order one.

Given the similarities among the three crude oil futures contracts, it is
natural to study whether the stochastic trend found in the three series is
common. Among the two traditional approaches to cointegration, the one of
[Engle and Granger, 1987] and the one of [Johansen, 1988] and [Johansen,
1995], we opt for the latter. Indeed, the Engle and Granger approach would
require to choose a priori a reference futures contract among the three under
consideration and there is no objective way to do so. In order to assess any
cointegrating relationship by means of the Johansen procedure, we first need
to frame them in a vector autoregressive model of order $p$ (VAR($p$)
henceforth). We choose $p$ by minimizing the Schwartz IC, which is known to
deliver the most parsimonious model. On our test sample, we obtain that the
optimal number of lags $p$ is equal to two222This result on the optimal $p$ is
quite stable even when considering different test samples, namely when
changing $t_{0}$ and/or $t_{B}$. In the few remaining cases, the optimal $p$
turns out to be equal either to one or three.. The resulting (reduced-form)
VAR(2) can be written as

$\mathbf{F}_{t}=\mathbf{a}_{0}+\mathbf{A}_{1}\mathbf{F}_{t-1}+\mathbf{A}_{2}\mathbf{F}_{t-2}+\mathbf{u}_{t}$

where $\mathbf{F}_{t}=\left[F_{t}^{B}\quad F_{t}^{S}\quad
F_{t}^{W}\right]^{\top}{}$ are the three variables of interest,
$\mathbf{a}_{0}\in\mathbb{R}^{3}$, $\mathbf{A}_{i}\in\mathbb{R}^{3\times 3}$,
$i=1,2$, are the parameters of the model and $\mathbf{u}_{t}$ is a
tridimensional white noise process with covariance matrix
$\mathbf{\Sigma}\in\mathbb{R}^{3\times 3}$. The VAR(2) model is estimated by
ordinary least squares and the resulting point estimates and standard errors
are

$\left[\begin{array}[]{c}\underset{}{F_{t}^{B}}\\\ \underset{}{F_{t}^{S}}\\\
\underset{}{F_{t}^{W}}\end{array}\right]=\left[\begin{array}[]{c}\underset{(1.095)}{1.546}\\\
\underset{(0.883)}{1.111}\\\
\underset{(1.093)}{1.353}\end{array}\right]+\left[\begin{array}[]{ccc}\underset{(0.249)}{0.614^{**}}&\underset{(0.125)}{0.088}&\underset{(0.234)}{{0.147}}\\\
\underset{(0.201)}{0.510^{**}}&\underset{(0.101)}{0.555^{***}}&\underset{(0.189)}{{-0.100}}\\\
\underset{(0.249)}{-0.221}&\underset{(0.125)}{0.111}&\underset{(0.233)}{{0.968^{***}}}\end{array}\right]\left[\begin{array}[]{c}\underset{}{F_{t-1}^{B}}\\\
\underset{}{F_{t-1}^{S}}\\\ \underset{}{F_{t-1}^{W}}\end{array}\right]$
$\qquad+\left[\begin{array}[]{ccc}\underset{(0.258)}{0.266}&\underset{(0.113)}{-0.157}&\underset{(0.238)}{0.033}\\\
\underset{(0.208)}{-0.221}&\underset{(0.091)}{0.157^{*}}&\underset{(0.192)}{0.078}\\\
\underset{(0.257)}{0.242}&\underset{(0.112)}{-0.214^{*}}&\underset{(0.237)}{0.100}\end{array}\right]\left[\begin{array}[]{c}\underset{}{F_{t-2}^{B}}\\\
\underset{}{F_{t-2}^{S}}\\\
\underset{}{F_{t-2}^{W}}\end{array}\right]+\left[\begin{array}[]{c}\underset{}{u_{t}^{B}}\\\
\underset{}{u_{t}^{S}}\\\ \underset{}{u_{t}^{W}}\end{array}\right]$

where ∗∗∗, ∗∗ and ∗ represent, respectively, a statistical significance at the
1%, 5% and 10% level. Inspecting the roots of the characteristic equation,
this VAR(2) appears to satisfy the stability condition. We can now perform the
cointegration test. We adopt the _trace test_ of [Johansen, 1988] and allow
for a constant in the cointegrating relationship.

r | h | Trace test stat. | Crit. val. | pValue
---|---|---|---|---
0 | 1 | 36.107 | 35.193 | 0.039
1 | 0 | 10.648 | 20.262 | 0.603
2 | 0 | 1.250 | 9.164 | 0.915
SBIC($r=1$): 12.9533
Table 2. Results of [Johansen, 1988] trace test. $h$ is the rejection decision
of the trace test with null hypothesis “there exist less than or equal to $r$
cointegrating relationships”.

Looking at the results reported in Table 2, the test indicates the presence of
only one cointegrating relationship. Moreover, in our dataset we find that
this result is robust also when considering different time windows. We can
therefore estimate the following vector error correction model (VECM
henceforth):

$\Delta\mathbf{F}_{t}=\boldsymbol{\pi}_{0}+\mathbf{\Pi}_{0}\mathbf{F}_{t-1}+\mathbf{\Pi}_{1}\Delta\mathbf{F}_{t-1}+\mathbf{\Pi}_{2}\Delta\mathbf{F}_{t-2}+\boldsymbol{\varepsilon}_{t}$

where $\Delta\mathbf{F}_{t}:=\mathbf{F}_{t}-\mathbf{F}_{t-1}$ is the first
difference operator and $\boldsymbol{\pi}_{0}\in\mathbb{R}^{3}$,
$\mathbf{\Pi}_{i}\in\mathbb{R}^{3\times 3}$, $i=0,1,2$, $c\in\mathbb{R}$ are
the parameters of the model and $\boldsymbol{\varepsilon}_{t}$ is another
tridimensional white noise process with covariance matrix
$\mathbf{\Sigma}_{u}\in\mathbb{R}^{3\times 3}$.

Since we cannot reject that $rk\left[\mathbf{\Pi}_{0}\right]=1$, the matrix
$\mathbf{\Pi}_{0}$ can be decomposed as
$\mathbf{\Pi}_{0}=\boldsymbol{\alpha}\boldsymbol{\beta}^{T}$ with
$\boldsymbol{\beta}\in\mathbb{R}^{3}$. Rewriting $\boldsymbol{\pi}_{0}$ as
$\boldsymbol{\pi}_{0}=\boldsymbol{\alpha}c_{0}$, with $c_{0}\in\mathbb{R}$, we
can derive the expression of the resulting cointegration relationship,
$c_{0}+\boldsymbol{\beta}^{\top}{}\mathbf{F}_{t}$, which is our spread process
$S$, and the adjustment coefficients $\boldsymbol{\alpha}$ that measure the
speed of convergence of each series to the long-run relationship. Therefore,
in our test sample the spread process results to be given by

(4.1) $S_{t}=F_{t}^{B}-0.6982F_{t}^{S}-0.3402F_{t}^{W}+0.4322$

and it is displayed in Figure 2.

Figure 2. Resulting spread process within the test sample, weekly
observations.

As can be seen from this figure, the spread process exhibits a stationary
behavior. Coherently, the ADF and PP tests reject the hull hypothesis of the
presence of a unit root while the KPSS test ([Kwiatkowski et al., 1992])
cannot reject the null hypothesis of stationarity of the series.

Figure 3. ACF and PACF of the spread process over the training sample, weekly
observations.

Moreover, as can be seen from Figure 3, both the ACF and the PACF of $S$
decline to zero. In particular, only the first two lags of the PACF seem to be
statistically different from zero. This indicates some degree of mean-
reversion of $S$ that empirically justifies modelling $S$ by means of the
stochastic model introduced in Section 2.1.

The estimated adjustment coefficients of the VECM are
$\boldsymbol{\alpha}=\left[0.0640\quad 0.4226\quad 0.0684\right]^{\top}$. It
is interesting to observe that the Shanghai futures price realigns to the
long-run equilibrium much faster than the other two futures prices. As a high
speed of adjustment in one or more of the portfolio constituents is a
desirable property when constructing statistical arbitrage strategies, this
provides a preliminary indication of the potential profitability of the
Shanghai futures in the context of crude oil futures portfolios.

####  4.1.3. Filtering the spread

We now estimate the parameters of the OU-HMM for $S$. More precisely: we first
determine the most likely number $N$ of different states of the latent Markov
chain X; we then compute daily filtered estimates of X; finally, we estimate
the model parameters $\boldsymbol{\gamma}$, $\boldsymbol{\alpha}$,
$\boldsymbol{\eta}$ in (2.3) and the transition matrix $\boldsymbol{\Pi}$ of
$X$ by means of the filter-based EM algorithm described in Section 2.2. We
implement the filtering-estimation procedure on daily data. The choice of
daily data is motivatied by two reasons: first, we want to rely on a large
dataset for the estimation of the spread process, in order to capture its
short-term variations; second, we want to possibly implement trading
strategies on a daily basis rather than on a weekly basis.

The daily spread process is determined by using the weights reported in (4.1).
We estimate the whole model for $N\in\\{1,2,3\\}$ and compute the information
criteria and the mean square error of the three resulting estimates at the end
of the training sample. Although most of the times there is full agreement
among the different information criteria on the optimal value of $N$, we
follow the Schwartz IC. According to this criterion, we find that in our
training sample the best model entails $N=2$ states of the Markov chain X.

Figure 4 illustrates the time evolution of the estimates of the model
parameters $\boldsymbol{\gamma}$, $\boldsymbol{\alpha}$, $\boldsymbol{\eta}$
and of the probabilities $\pi_{ii}$, $i=1,2$, of the Markov chain X to remain
in the same state.

Figure 4. Recursively updated estimates of $\gamma_{i}$, $\alpha_{i}$,
$\eta_{i}$ and $\pi_{ii}$, $i=1,2$, for the $N=2$ OU-HMM.

Due to the well-known slow convergence of the EM algorithm, the estimates turn
out to be unstable at the beginning of the training sample. Afterwards, they
stabilize and only major changes in the underlying observations have an impact
on the parameter estimates.

| state 1 | state 2 |  |  | state 1 | state 2
---|---|---|---|---|---|---
$\gamma$ | 0.1632 | 0.0869 |  | $\beta$ | 0.8490 | 0.3233
$\alpha$ | 0.8077 | 0.7313 |  | $a$ | 53.3830 | 78.2361
$\eta$ | 1.4312 | 1.4832 |  | $\xi$ | 25.0833 | 27.2027
$\pi_{11}$ | 0.7138 |  | $\pi_{22}$ | 0.6640
Table 3. Point estimates of the functions of the OU-HMM (both in discrete and
continuous time) at the end of the training sample.

For sake of completeness, the parameter estimates at the end of the training
sample for both the discrete and continuous-time versions of the OU-HMM
presented in Section 2 are reported in Table 3. As we can see, the main
feature that differentiates the two regimes is the value of $\gamma$: the two
regimes display almost the same level of persistence $\alpha$, the same level
of volatility $\eta$ but a different mean level $\gamma$. We point out that
across the other training samples considered in the following section this
difference is even more striking and most of the times we observe
$\gamma_{1}>0$ and $\gamma_{2}<0$.

Figure 5 displays the one-step ahead forecast of $S$ along with the realized
value across the test sample.

Figure 5. One-step ahead forecasts and actual values of $y$ across the test
sample.

As can be seen from this figure, the filtering algorithm produces fairly good
forecasts of the realizations of the spread.

###  4.2. Testing the strategies

We now analyze over the test sample (from $t_{B}=$ 07/01/2022 until $T=$
06/30/2023) the five different statistical arbitrage strategies described in
Section 3.1: the plain Vanilla (PV), the probability interval (ProbI), the
prediction interval (PredI), the realized increment (RI), the predicted
increment (PI), along with the buy-and-hold S&P500 strategy. For each strategy
we compute the opening/closing signals, the annualized return according to
(3.2) and the annualized Sharpe ratio according to (3.3).

Figure 6 illustrates how the signals differ across the five statistical
arbitrage strategies.

Figure 6. Opening/closing signals of the five statistical arbitrage strategies
described in Section 3.1 over the test sample of Section 4.1

We immediately notice that the number of trades over the test sample varies
significantly across the different strategies. Excluding the PV strategy that
entails always an open position as soon as $S_{t}\neq 0$ (which is always the
case), we document a minimum of five trades with the ProbI strategy and a
maximum of twenty-nine trades with the PredI strategy. Moreover, also the
length of the open positions varies considerably: indeed, we generate
positions that are opened on a given day and closed on the following day
together with positions that remain open for several months.

| PV | ProbI | PredI | RI | PI | S&P
---|---|---|---|---|---|---
$R$ | 13.03% | 56.37% | 72.77% | -1.43% | 22.41% | 16.43 %
$SR$ | 0.6082 | 1.0849 | 1.0108 | 0.4440 | 0.6836 | 0.8868
Table 4. Annualized performances with $c=80$bps over the test sample of
Section 4.1 including the Brent, the Shanghai and the WTI futures contracts.

Table 4 summarizes the annualized return $R$ and Sharpe ratio $SR$, as
described in Section 3.2, assuming $c=80$bps. As mentioned in Section 3.2,
this relatively high level of transaction costs yieds a conservative
assessment of the performance of our statistical arbitrage strategies
involving the Shanghai futures. We find that only the PV and RI strategies
deliver a return smaller than the S&P strategy, which entails no transaction
costs. On the contrary, the ProbI, PredI and PI strategies outperforms the S&P
strategy and deliver significant returns. In particular, it turns out that the
forward-looking strategies PredI and PI that are based on the stochastic model
introduced in Section 2 outperform the simpler ProbI and RI strategies that
rely only on historical observations of the spread. We can therefore conclude
that our OU-HMM seems to have promising applications for optimal investment in
international crude oil futures markets.

####  4.2.1. Risk analysis

Statistical arbitrage strategies may involve a considerable amount of risk. In
particular, when opening a position, the investor does not know when the
futures prices will realign again or whether convergence will ever be reached.
As a consequence, the investor might find herself stuck in a costly position
to maintain (due to margin calls) and that does not generate any profit for a
prolonged period of time. In this section, we assess the riskiness of our
statistical arbitrage strategies by means of a Monte Carlo-based analysis.

First of all, we compute the out-of-sample daily Value-at-Risk (VaR
henceforth) of the strategies described in Section 3.1. Moreover, we exploit
the simulations that are required for VaR estimation in order to derive the
empirical distribution of the annualized returns of our strategies. This is
done by implementing the following successive steps:

Figure 7. A simulated out-of-sample trajectory of the spread.

  1. (1)

we estimate the cointegration relationship over the training sample as
described in Section 4.1.2;

  2. (2)

we filter the resulting spread over the training sample following Section
4.1.3;

  3. (3)

at the breaking date $t_{B}$ we retrieve the optimal estimated number of
states $N$ of the Markov chain and the related point estimates of the
parameters of our OU-HMM model;

  4. (4)

for $i=1,...,NSim$:

     * •

we simulate over the test sample (i.e. from $t_{B}$ to $T$) a trajectory of
the latent Markov chain X with the optimal $N$ found before;

     * •

we simulate over the test sample the spread process following (2.3) and the
realizations of X simulated at the previous step using the parameters obtained
at the end of the training sample. Figure 7 shows a simulated trajectory of
the spread. From a qualitative point of view, the out-of-sample trajectory
look similar to the realized one;

     * •

we compute the series of the returns according to (3.1) of the different
statistical arbitrage strategies described in Section 3.1;

     * •

we compute the daily VaR at the confidence levels 1%, 5% and 10%;

     * •

we compute the overall performance of the strategies in terms of annualized
return (3.2) and of Sharpe ratio (3.3);

  5. (5)

we average out the VaRs obtained along all the $NSim$ trajectories;

  6. (6)

we compute the Kernel Density Estimator (KDE henceforth) of the annualized
returns obtained along all the $NSim$ trajectories.

| PV | ProbI | PredI | RI | PI | S&P
---|---|---|---|---|---|---
$VaR_{99\%}$ | -20.21% | -10.86% | -10.16% | -7.93% | -7.62% | -2.80%
$VaR_{95\%}$ | -5.88% | -5.49% | -5.80% | -3.66% | -3.69% | -1.72%
$VaR_{90\%}$ | -4.04% | -3.35% | -3.87% | -1.94% | -2.12% | -1.29%
Table 5. Monte Carlo estimates of the daily VaRs of daily returns over the
test sample of Section 4.1 including the Brent, the Shanghai and the WTI
futures contract ($NSim=1000$; the average t-statistics is equal to almost 6).

Table 5 displays the Monte Carlo estimates of the daily VaRs at the usual
confidence levels. As expected, the attractive performances of statistical
arbitrage opportunities involving the Shanghai futures reported in Table 4 are
mitigated by a considerable level of risk, if compared to the simple S&P
strategy. This implies that an arbitrageur willing to engage in statistical
arbitrage strategies involving the Shanghai futures has to be ready to face
possibly severe drops of the market value of her portfolio which may call for
extra margins to be posted.

Figure 8. Monte Carlo-based KDEs of annual returns over the test sample of
Section 4.1 including the Brent, the Shanghai and the WTI futures contract
($NSim=1000$).

Figure 8 displays the Monte Carlo-based KDEs of the annual returns of the five
strategies considered over the test sample. As we can see, all of the
densities seem to be positively skewed, with a sizeable amount of probability
mass associated to positive annual returns. In particular, this provides
evidence that statistical arbitrage strategies that exploit the OU-HMM model
introduced in Section 2 yield profitable investment opportunities.

###  4.3. Robustness of the results

The empirical results reported and discussed so far depend on the time window
under consideration, as well as on the level of transaction costs. In order to
assess the robustness of our findings, we now systematically repeat the
analysis carried out in the previous sections changing the transaction costs
parameter $c$, the starting date $t_{0}$ of the sample and the breaking date
$t_{B}$ that separates the training sample from the test one.

####  4.3.1. Sensitivity of the performances with respect to $c$

Figure 9. Performances of the strategies with respect to $c$; $t_{0}=$
03/26/2018, $t_{B}=$ 07/01/2022.

Figure 9 displays the performances of the statistical arbitrage strategies as
a function of $c$. Clearly, both the annualized return $R$ and the Sharpe
ratio $SR$ are declining in $c$ (but for the buy-and-hold S&P strategy). The
slope of this decline is proportional to the frequency of opening/closing a
position: the performances of strategies PV, PredI and PI that prescribe the
opening of several trades (see Figure 6) are the most sensible to $c$ while
strategies ProbI and RI are less affected by transaction costs. Anyway, we
find that over the considered time period, statistical arbitrage strategies
that include the Shanghai futures are quite robust with respect to the overall
level of transaction costs. In particular, the benchmark S&P strategy becomes
preferable only when $c$ approaches 1%.

####  4.3.2. Sensitivity of the performances with respect to $t_{0}$

Figure 10. Performances of the strategies with respect to $t_{0}$; $c=80$bps,
$t_{B}=$ 07/01/2022. Missing datapoints are due to the rejection of
cointegration among the futures contract over the related time window.

Figure 10 depicts the performances of our strategies as a function of the
starting date $t_{0}$ of our sample. As the breaking date $t_{B}$ is fixed, if
$t_{0}$ moves forward the sample shrinks and the number of observations in the
training sample decreases. This is expected to worsen the quality of the
filter-based estimates derived according to the algorithm described in Section
2.2, due to the availability of a smaller number of observations. As we see
from Figure 10, reducing the training sample seems to weaken the cointegrating
relationship among the series and, as a consequence, it lowers the performance
of our strategies. In particular, from a risk-adjusted viewpoint, the Sharpe
ratio of the simple S&P strategy seems to outperform almost always the
proposed statistical arbitrage strategies. This indicates that, in order to be
profitable, our statistical arbitrage strategies require a minimum sample size
for the model to be accurately estimated.

####  4.3.3. Sensitivity of the performances with respect to $t_{B}$

Figure 11. Performances of the strategies with respect to $t_{B}$; $c=80$bps,
$t_{0}=$ 03/26/2018.

Figure 11 shows the performances of our strategies as a function of the
breaking date $t_{B}$. As $t_{B}$ moves backward in time, the training sample
shrinks and the test sample becomes longer. In principle, the impact of this
on the performance of statistical arbitrage strategies is twofold: on one end,
reducing the size of the training sample should weaken the estimate of the
cointegration relationship; on the other side, a larger test sample provides
more opportunities for profitable trades. In our dataset, a small decline in
the performance is registered for $t_{B}$ before than March-July 2022.
Bringing $t_{B}$ forward should have ambiguous effects as well: the
cointegration is estimated over a larger training sample but there is less
time in the test sample to take advantage of that. In our dataset we document
a small increase of the performances (especially in terms of Sharpe ratio) for
$t_{B}$ after March-July 2022.

###  4.4. Replacing the Shanghai futures with the Dubai

In order to appreciate the profitability of the Shanghai futures for our
statistical arbitrage strategies, we now replicate the analysis of Section 4.1
replacing the Shanghai crude oil futures by the Dubai futures.

The performances of the statistical arbitrage strategies involving the Dubai
futures are computed using $c=20$bps and are displayed in Table 6.

| PV | ProbI | PredI | RI | PI | S&P
---|---|---|---|---|---|---
$R$ | 6.45% | 7.21% | 7.22% | -0.00% | -0.01% | 16.43 %
$SR$ | 0.6195 | 0.7115 | 0.6852 | 0.0363 | -0.0299 | 0.8868
Table 6. Annualized performances with $c=20$bps over the test sample of
Section 4.1 including the Brent, the Dubai and the WTI futures contract.

These performances are mildly lower than the ones registered in [Alizadeh and
Nomikos, 2008], where only the Brent and the WTI are considered for the pairs
trading strategies and the transaction costs are the same of our example. This
might be consistent with the findings of [Do and Faff, 2010] that document how
net profits from pairs trading strategies tend to decline over time.
Coherently with this fact, the simple buy-and-hold strategy on the S&P
outperforms all other strategies over the test sample.

These findings indicate that, despite the well-established cointegrating
relationship among the Brent, the WTI and the Dubai crude oil futures prices,
the profitability of statistical arbitrage strategies based on cointegration
seems to have been almost exhausted, even in the presence of very mild
transaction costs. On the contrary, there seems to be room for profitable
investment opportunities when including the newly introduced Shanghai futures.

Table 7 displays the Monte Carlo estimates of the daily VaRs of the strategies
involving the Brent, the WTI and the Dubai crude oil futures.

| PV | ProbI | PredI | RI | PI | S&P
---|---|---|---|---|---|---
$VaR_{99\%}$ | -1.14% | -1.05% | -1.13% | -0.77% | -0.87% | -2.80%
$VaR_{95\%}$ | -0.61% | -0.52% | -0.60% | -0.31% | -0.39% | -1.72%
$VaR_{90\%}$ | -0.39% | -0.31% | -0.38% | -0.15% | -0.21% | -1.29%
Table 7. Monte Carlo estimates of the daily VaRs of daily returns over the
test sample of Section 4.1 including the Brent, the Dubai and the WTI futures
contract ($NSim=1000$; the average t-statistics is equal to 7).

Figure 12. Monte Carlo-based KDEs of annual returns over the test sample of
Section 4.1 including the Brent, the Dubai and the WTI futures contract
($NSim=1000$).

Overall, we find that the level of risk faced when trading the Dubai futures
contract instead of the Shanghai one decreases coherently with the lower
performance. Figure 12 shows the KDEs of the annual returns of the strategies
involving the Dubai futures. As we can see, most of the probability mass is
now close or even below zero, thus indicating once more than the profitability
of statistical arbitrage opportunities involving the established crude oil
futures benchmarks has almost vanished. For sake of completeness, Figure 13
shows the performance of the strategies outlined in Section 3.1 as a function
of the transaction costs parameter $c$.

Figure 13. Performances of the strategies with respect to $c$; $t_{0}=$
03/26/2018, $t_{B}=$ 07/01/2022.

##  5\. Conclusions

In this work, we investigated statistical arbitrage strategies that generalize
pairs trading strategies, involving the two established crude oil futures
contracts, the Brent and the WTI, together with the Shanghai crude oil futures
more recently introduced in 2018. We documented that the time series of the
three futures prices are cointegrated. In order to capture the variability of
the resulting cointegration spread, we have proposed a mean-reverting
stochastic process with regime switching modulated by a latent Markov chain.
The model is estimated by relying on a filter-based version of the EM
algorithm, thereby ensuring that the model is dynamically tuned to the market
information. The estimates have been used to build statistical arbitrage
strategies that turned out to be remarkably profitable, even under a
conservative level of transaction costs.

Our findings suggest that the Shanghai crude oil futures, which has already
become the benchmark for the domestic Chinese crude oil market, represents a
valuable investment opportunity for international investors. Moreover, our
results demonstrate that, despite the well established and strong
cointegration among the Brent, the WTI and the Dubai crude oil futures,
profitable statistical arbitrages in these three assets seem to have been
exhausted. This is not the case for the Shanghai futures, which yields
interesting investment performances, even if subject to a considerable amount
of financial risk. The introduction of the Shanghai futures is quite recent
and, therefore, more extensive empirical analysis will be needed in order to
assess the specific features of this asset in the context of global financial
markets.

## References

  * [Alizadeh and Nomikos, 2008] Alizadeh, A. and Nomikos, N. (2008).  Performance of statistical arbitrage in petroleum futures markets.  Journal of Energy Markets, 1(2):3–33.
  * [Avellaneda and Lee, 2010] Avellaneda, M. and Lee, J.-H. (2010).  Statistical arbitrage in the US equities market.  Quantitative Finance, 10(7):761–782.
  * [Baviera and Baldi, 2019] Baviera, R. and Baldi, T. S. (2019).  Stop-loss and leverage in optimal statistical arbitrage with an application to energy market.  Energy Economics, 79:130–143.
  * [Bertram, 2010] Bertram, W. (2010).  Analytic solutions for optimal statistical arbitrage trading.  Physica A, 389(11):2234–2243.
  * [Bock and Mestel, 2009] Bock, M. and Mestel, R. (2009).  A regime-switching relative value arbitrage rule.  In Fleischmann, B., Borgwardt, K.-H., Klein, R., and Tuma, A., editors, Operations Research Proceedings 2008, pages 9–14, Berlin, Heidelberg. Springer.
  * [Bollinger, 2001] Bollinger, J. (2001).  Bollinger on Bollinger bans.  McGraw-Hill, New York (NY).
  * [Bondarenko, 2003] Bondarenko, O. (2003).  Statistical arbitrage and securities prices.  Review of Financial Studies, 16(3):875–919.
  * [Burgess, 1999] Burgess, A. N. (1999).  A Computational Methodology for Modelling the Dynamics of Statistical Arbitrage.  PhD thesis, London Business School.
  * [Caporin et al., 2019] Caporin, M., Fontini, F., and Talebbeydokhti, E. (2019).  Testing persistence of WTI and Brent long-run relationship after the shale oil supply shock.  Energy Economics, 79:21–31.
  * [Cerqueti and Fanelli, 2021] Cerqueti, R. and Fanelli, V. (2021).  Long memory and crude oil’s price predictability.  Annals of Operations Research, 299:895–906.
  * [Cerqueti et al., 2019] Cerqueti, R., Fanelli, V., and Rotundo, G. (2019).  Long run analysis of crude oil portfolios.  Energy Economics, 79:183–205.
  * [Cotter et al., 2022] Cotter, J., Eyiah-Donkor, E., and Potì, V. (2022).  Commodity futures return predictability and intertemporal asset pricing.  Journal of Commodity Markets, page 100289.
  * [Cummins and Bucca, 2012] Cummins, M. and Bucca, A. (2012).  Quantitative spread trading on crude oil and refined products markets.  Quantitative Finance, 12(12):1857–1875.
  * [Do and Faff, 2010] Do, B. and Faff, R. (2010).  Does simple pairs trading still work?  Financial Analysts Journal, 66(4):83–95.
  * [Do and Faff, 2012] Do, B. and Faff, R. (2012).  Are pairs trading profits robust to trading costs?  Journal of Financial Research, 35(2):261–287.
  * [Dunis et al., 2006] Dunis, C., Laws, J., and Evans, B. (2006).  Trading futures spreads: an application of correlation and threshold filters.  Applied Financial Economics, 16(12):903–914.
  * [Elliott et al., 1995] Elliott, R. J., Aggoun, L., and Moore, J. B. (1995).  Hidden Markov Models: Estimation and Control.  Springer, New York.
  * [Elliott and Bradrania, 2018] Elliott, R. J. and Bradrania, R. (2018).  Estimating a regime switching pairs trading model.  Quantitative Finance, 18(5):877–883.
  * [Elliott et al., 2005] Elliott, R. J., Van der Hoek, J., and Malcolm, W. P. (2005).  Pairs trading.  Quantitative Finance, 5(3):271–276.
  * [Endres and Stübinger, 2019] Endres, S. and Stübinger, J. (2019).  A flexible regime switching model with pairs trading application to the s&p 500 high-frequency stock returns.  Quantitative Finance, 19(10):1727–1740.
  * [Engle and Granger, 1987] Engle, R. and Granger, C. (1987).  Co-integration and error correction: representation, estimation and testing.  Econometrica, 55:251–276.
  * [Erlwein et al., 2010] Erlwein, C., Benth, F. E., and Mamon, R. (2010).  HMM filtering and parameter estimation of an electricity spot price model.  Energy Economics, 32(5):1034–1043.
  * [Erlwein and Mamon, 2009] Erlwein, C. and Mamon, R. (2009).  An online estimation scheme for a Hull-White model with HMM-driven parameters.  Statistical Methods and Applications, 18:87–107.
  * [Fontana and Runggaldier, 2010] Fontana, C. and Runggaldier, W. J. (2010).  Credit risk and incomplete information: filtering and EM parameter estimation.  International Journal of Theoretical and Applied Finance, 13(5):683–715.
  * [Galay, 2019] Galay, G. (2019).  Are crude oil markets cointegrated? Testing the co-movement of weekly crude oil spot prices.  Journal of Commodity Markets, 16:100088.
  * [Gatev et al., 2006] Gatev, E., Goetzmann, W. N., and Rouwenhorst, K. G. (2006).  Pairs trading: performance of a relative-value arbitrage rule.  Review of Financial Studies, 19(3):797–827.
  * [Grimm et al., 2020] Grimm, S., Erlwein-Sayer, C., and Mamon, R. (2020).  Discrete-time implementation of continuous-time filters with application to regime-switching dynamics estimation.  Nonlinear Analysis: Hybrid Systems, 35:100814.
  * [Guidolin and Pedio, 2018] Guidolin, M. and Pedio, M. (2018).  Essentials of time series for financial applications.  Academic Press, San Diego (CA).
  * [Hammoudeh et al., 2008] Hammoudeh, S., Ewing, B., and Thompson, M. (2008).  Threshold cointegration analysis of crude oil benchmarks.  The Energy Journal, 29(4):79–95.
  * [Huang and Huang, 2020] Huang, X. and Huang, S. (2020).  Identifying the comovement of price between China‘s and international crude oil futures: a time-frequency perspective.  International Review of Financial Analysis, 72:101562.
  * [Ji et al., 2022] Ji, Q., Zhang, D., and Zhao, Y. (2022).  Intra-day co-movements of crude oil futures: China and the international benchmarks.  Annals of Operations Research, 313(1):77–103.
  * [Johansen, 1988] Johansen, S. (1988).  Statistical analysis of cointegration vectors.  Journal of Economic Dynamics and Control, 12(2):231–254.
  * [Johansen, 1995] Johansen, S. (1995).  Likelihood-based inference in cointegrated vector autoregressive models.  Okford University Press (UK).
  * [Kristoufek and Vosvrda, 2014] Kristoufek, L. and Vosvrda, M. (2014).  Commodity futures and market efficency.  Energy Econmics, 42:50–57.
  * [Kwiatkowski et al., 1992] Kwiatkowski, D., Phillips, P., Schmidt, P., and Shin, Y. (1992).  Testing the null hypothesis of stationarity against the alternative of a unit root: How sure are we that economic time series have a unit root?  Journal of Econometrics, 54(1):159–178.
  * [Lanza et al., 2005] Lanza, A., Manera, M., and Giovannini, M. (2005).  Modeling and forecasting cointegrated relationships among heavy oil and product prices.  Energy Economics, 27:831–848.
  * [Lee and Papanicolaou, 2016] Lee, S. and Papanicolaou, A. (2016).  Pairs trading of two assets with uncertainty in co-integration’s level of mean reversion.  International Journal of Theoretical and Applied Finance, 19(8):1650054.
  * [Niu et al., 2023] Niu, J., Ma, C., and Chang, C. (2023).  The arbitrage strategy in the crude oil futures market of Shanghai international energy exchange.  Economic Change and Restructuring, (56):1201–1223.
  * [Perron, 1997] Perron, P. (1997).  Further evidence on breaking trend functions in macroeconomic variables.  Journal of Econometrics, 80(2):355–385.
  * [Phillips and Perron, 1988] Phillips, P. and Perron, P. (1988).  Testing for a unit root in time series regression.  Biometrika, 75(2):335–346.
  * [Rein et al., 2021] Rein, C., Rüschendorf, L., and Schmidt, T. (2021).  Generalized statistical arbitrage concepts and related gain strategies.  Mathematical Finance, 31(2):563–594.
  * [Said and Dickey, 1984] Said, S. and Dickey, D. (1984).  Testing for unit roots in autoregressive-moving average models of unknown order.  Biometrika, 71(3):599–607.
  * [Sarmento and Horta, 2021] Sarmento, S. and Horta, N. (2021).  A Machine Learning Based Pairs Trading investment strategy.  Springer, Princeton (NJ).
  * [Shao and Hua, 2022] Shao, M. and Hua, Y. (2022).  Price discovery efficiency of China‘s crude oil futures: evidence from the Shanghai crude oil futures market.  Energy Economics, 112:106172.
  * [Tenyakov and Mamon, 2017] Tenyakov, A. and Mamon, R. (2017).  A computing platform for pairs-trading online implementation via a blended Kalman-HMM filtering approach.  Journal of Big Data, 4:46.
  * [Vidyamurthy, 2004] Vidyamurthy, G. (2004).  Pairs Trading: Quantitative Methods and Analysis.  Wiley, Hoboken (NJ).
  * [Wang et al., 2022] Wang, J., Qiu, S., and Yick, H. Y. (2022).  The influence of the Shanghai crude oil futures on the global and domestic oil markets.  Energy, 245:123271.
  * [Yang and Zhou, 2020] Yang, J. and Zhou, Y. (2020).  Return and volatility transmission between China‘s and international crude oil futures markets: a first look.  Journal of Futures Markets, 40(6):860–884.

A HIDDEN MARKOV MODEL FOR STATISTICAL ARBITRAGE
IN INTERNATIONAL CRUDE OIL FUTURES MARKETS
SUPPLEMENTARY MATERIALS

We parallel here the main results of Section 4 using daily data also for the
cointegration analysis. For sake of brevity, we omit here the value of the
several test statistics and point estimates of the analysis carried out in
Subsection 4.1.2, as they are also qualitatively similar, and we just comment
on the results.

We first replicate the analysis of Subsection 4.1. Therefore, we set
$t_{0}=03/26/2018$, $t_{B}=07/01/2022$ and $T=06/30/2023$ and we consider
daily data for estimating the cointegration relationship. The training (resp.
test) sample consists of 1113 (resp. 261) observations.

The ADF, the PP and the bP tests suggest, again, the presence of at least one
unit root in each of the three series of futures prices. Testing the first
differences of the series enables us to conclude that each of them actually
contains precisely one unit root. Looking for the optimal lag length of a VAR
model involving the three series of daily futures prices delivers $p=4$.
Running the Johansen cointegration analysis we cannot reject the presence of
precisely one cointegrating relationship. The resulting spread process here is
given by

(5.1) $S_{t}=F_{t}^{B}-0.5880F_{t}^{S}-0.4619F_{t}^{W}+0.3158,$

which does not differ that much from the one obtained using weekly data. This
spread passes all the unit root/stationarity tests. Although globally smaller
than the ones obtained using weekly data, the adjustment coefficients of the
estimated VECM display the same characteristic observed using weekly data: in
particular, the one for the Shanghai futures is equal to 0.098164 while the
ones for the Brent and the WTI are equal to -0.0117 and 0.0197, respectively.
Therefore, the Shanghai futures proves again to revert much faster to the
long-run equilibrium than the other two futures and this is potentially
profitable when setting up statistical arbitrage strategies.

We now estimate the parameters of the OU-HMM for $S$ in (5.1). Interestingly,
over this precise time window and using daily data, the optimal number of
states of the underlying Markov chain turns out to be equal to three.

Figure 14. Recursively updated estimates of $\gamma_{i}$, $\alpha_{i}$,
$\eta_{i}$, $i=1,\dots,3$ and $\pi_{ij}$, $i,j=1,\dots,3$, for the $N=3$ OU-
HMM.

Figure 14 displays the recursively updated estimates of the parameters of the
OU-HMM. As when the spread was generated by weekly data, the main feature that
differentiates the three states is the constant $\gamma$. While with $N=2$ we
usually get a positive and a negative $\gamma$, with $N=3$ we also have a
third state in which $\gamma$ is almost zero.

Figure 15. One-step ahead forecasts and actual values of $y$ across the test
sample.

As before, the one-step ahead forecasts of the model track the realized values
reasonably well, as can be seen in Figure 15.

We now move to the analysis of the performances of the statistical arbitrage
strategies described in Subsection 3.1 considering the estimates obtained from
daily data only. Table 8 parallels Table 4 and shows that the performances we
obtain now are overall smaller than the ones obtained using weekly data for
the cointegration analysis. In particular, only the ProbI and the PredI
strategies outperform the passive S&P investment. As before, strategies
involving the modelling of the spread (PredI and PI) seem to perform better
than the ones based on its historical analysis (ProbI and RI).

| PV | ProbI | PredI | RI | PI | S&P
---|---|---|---|---|---|---
$R$ | 0.52% | 17.64% | 44.43% | -5.69% | -3.67% | 16.43 %
$SR$ | 0.3848 | 0.5817 | 1.0690 | 0.2700 | -0.1346 | 0.8868
Table 8. Annualized performances with $c=80$bps over the test sample of Subsection 3.1 considering the spread (5.1) obtained from daily data only. | PV | ProbI | PredI | RI | PI | S&P
---|---|---|---|---|---|---
$R$ | -49.81% | -45.76% | -30.94% | -46.92% | -28.02% | 16.43 %
$SR$ | -5.5138 | -4.8878 | -3.0734 | -5.7164 | -3.6193 | 0.8868
Table 9. Annualized performances with $c=20$bps over the test sample of
Subsection 3.1 including the Brent, the Dubai and the WTI futures contract
using daily data only.

The performances of this analysis if the Shanghai futures is replaced by the
Dubai crude oil futures are displayed in Table 9. As we can see, now
performances do not only worsen with respect to the ones in Table 6 but they
get also all negative. This is further evidence that the Shanghai futures
seems way much more profitable if used in a pairs trading setting than the
Dubai futures.

We now replicate also the VaR analysis carried out in Subsection 4.2.1. Table
10 parallels Table 5. As we can see, the VaRs are now slightly larger (which
translates into slightly less risky strategies) than the ones obtained before.
However, the gains in terms of higher VaRs seem to be less than proportional
than the losses in terms of annualized performances.

| PV | ProbI | PredI | RI | PI | S&P
---|---|---|---|---|---|---
$VaR_{99\%}$ | -5.56% | -5.48% | -5.55% | -4.09% | -2.69% | -2.80%
$VaR_{95\%}$ | -2.91% | -2.70% | -2.87% | -1.61% | -1.65% | -1.72%
$VaR_{90\%}$ | -1.82% | -1.50% | -1.69% | -0.70% | -0.78% | -1.29%
Table 10. Monte Carlo estimates of the daily VaRs of daily returns over the
test sample of Subsection 3.1 considering the spread (5.1) obtained from daily
data only ($NSim=1^{\prime}000$; the average t-statistics is equal to 5).

Finally, we replicate the analysis carried out in Subsection 4.3, analyzing
the sensitivities of the performances with respect to $c$, $t_{0}$ and $t_{B}$
when only daily data are used. The results are displayed in Figures 16, 17 and
18. As we can see across different levels of transaction costs parameter $c$,
starting date of the training sample, $t_{0}$, and of the testing sample,
$t_{B}$, the performances are always lower than the ones we obtain when we
stick to weekly data for our cointegration analysis.

Figure 16. Performances of the strategies considering the spread (5.1)
obtained from daily data only with respect to $c$; $t_{0}=$ 03/26/2018,
$t_{B}=$ 07/01/2022.

Figure 17. Performances of the strategies considering the spread (5.1)
obtained from daily data only with respect to $t_{0}$; $c=80$bps, $t_{B}=$
07/01/2022. Missing datapoints are due to the rejection of cointegration among
the futures contract over the related time window.

Figure 18. Performances of the strategies considering the spread (5.1)
obtained from daily data only with respect to $t_{B}$; $c=80$bps, $t_{0}=$
03/26/2018.