Systematic Oil Time Spread Trading : Relative Analysis Between Oil Time Spread & Inventory

Systematic Oil Time Spread Trading : Relative Analysis Between Oil Time Spread & Inventory

Systematic Oil Time Spread Trading : We provide empirical results to an oil time spread model described in a paper by I. Bouchouev. To understand the model and research for the oil data, we analyze the model input variable and the model calibration procedure. We also analyze and backtest the potential trading signals for oil time spread based on inventory data. 

1. Introduction  

In the commodity market, oil time spread, which is defined by the price difference between the front month and second month futures contracts, is an extremely important object that represents the dynamics of short-term oil prices. It helps define the oil futures price term structure, contango and backwardation. By the standard storage theory, oil time spread is linked with the oil storage by convenience yield and storage cost, and it is widely known that the spread is negatively correlated with the inventory, but its exact formulation has been difficult to establish explicitly. In the paper “A Stylized Model of the Oil Price Squeeze” [1], I. Bouchouev proposes a novel approach to model the oil time spread under the assumption that the oil inventory follows a mean reversion stochastic process and the boundary conditions are modelled by point impulse functions, while the price dynamics are determined by the standard diffusion equation for derivatives pricing. In this paper we will present some empirical results for this model. We will provide data analysis for model input variable selection, study for model calibration, and analyze potential trading signals for oil time spread.

2. Model Overview  

In this section we will present a summary for the time spread model described in the paper by I. Bouchouev [1].

x is oil inventories, it follows the mean reversion stochastic process

dx=k(\overline{x}-x)dt+\sigma dz\text{, }x_{\min}\le x\le x_{\max}

 f(t,T;x) is the spread function, defined by


where S is the spot price, F is the future price. In practice, S could be viewed as the front month futures price and F is the second month futures price, which means f(t,T;x) represents the time spread. 

f(t,T;x) should solve the following partial differential equation 

\frac{\partial f}{\partial t}+\frac{\sigma^2}{2}\frac{\partial^2f}{\partial x^2}+k(\overline{x}-x)\frac{\partial f}{\partial x}=0
x_{\min}\le x\le x_{\max}\text{ ,    }0<t<T

with boundary conditions

f(T,T;x)=\delta (x-x_{min})-\delta(x-x_{max})=
+\infty \text{,   }&x=x_{min}\\
-\infty  \text{,   }&x=x_{max}\\
0  \text{,   } & x_{min}< x< x_{max}

Point impulse function (Dirac’s delta function) is defined as

+\infty \text{, }&x=x_0\\
0    \text{, }& x\neq x_0

and with additional normalization property that


Below is the closed form solution for f(t,T;x).

If k=0, we have

f(t,T;x)=\frac{1}{\sqrt{2\pi(T-t)}\sigma} \Big(exp(-\frac{(x-x_{min})^2}{2\sigma^2(T-t)})-exp(-\frac{(x-x_{max})^2}{2\sigma^2(T-t)})\Big)

If k≠0, we define y=(x-\overline{x})e^{-k(T-t)}\, \tau=\int_t^Te^{-2k(T-s)}ds, if k is time invariant,   could also be written as \tau=\frac{1}{2k}(1-e^{-2k(T-t)}), then we have

f(t,T;x)=\frac{1}{\sqrt{2\pi\tau}\sigma} \Big(exp(-\frac{(y+\bar{x}-x_{min})^2}{2\sigma^2\tau})-exp(-\frac{(y+\bar{x}-x_{max}))^2}{2\sigma^2\tau})\Big)

  Notice that the unit of variable x will affect the model parameters. For example, if the unit of x is in thousands barrels, then x, xmin and xmax should all have the unit in thousands barrels. For simplicity, we provide examples of the behavior of the model curve with different combinations of parameters under the assumption that x is in percentage units. In practice, since the 0% boundary might never be achieved due to the minimum maintenance level, while the 100% capacity boundary might be approximately reached when the market is in “super contango”, we assume xmin=15% as an approximate minimum maintenance level and xmax=100% as the maximum capacity. 

Figure 1: T-t=1/12, =0.5, x=0.5

Figure 2: T-t=1/12, =0.5, k=1

Figure 3: T-t=1/12, x=0.5, k=1

Figure 4: =0.5, x=0.5, k=1

3. Empirical Results  

3.1 The Data  

We use the data from the EIA official website. The data includes daily and weekly WTI crude oil prices, Cushing weekly and monthly inventory, PADD2 weekly and monthly inventory, working storage capacity by PADD districts and production by PADD districts. However, the storage capacity data is only available for years after 2011, which is a restriction for model calibration. In order to obtain more data points for model calibration, we also model the Cushing inventory capacity for years between 2008 and 2011 by applying a sigmoid function modelled with market capacity. The market capacity is a slowly increasing function with time and it is updated semi-annually, while the sigmoid function satisfies such feature on its right tail.

Systematic Oil Time Spread Trading : Relative Analysis Between Oil Time Spread & Inventory

We provide an example graph for the fitted capacity curve. Such capacity modelling is only applied for calibration, if there is market capacity data available, we use the market capacity, otherwise we use the fitted capacity to obtain more data points. The capacity modelling ensures us to cover two periods with “super contango” in the calibration. In order to exclude the capacity modelling influence in data analysis for other sections such as input variable selection and trading signals analysis, we restrict the data period between 2011/03/01 and 2020/07/21, when we have the capacity data. For model calibration, the data period will be restricted between 2008/01/02 and 2020/07/21 and we use model capacity for periods before 2011. We will use the weekly data for input variable selection and model calibration since they are less noisy. For the trading signals analysis section, we use the daily crude oil prices data.

Figure 5: Capacity Modelling

3.2 Input Variable Selection  

Input variable selection is important to model calibration. In the oil market, there is many inventory data available, such as U.S. overall inventory data, PADD districts inventory data and Cushing inventory data, but not all of them are suitable for model input. In this section we present a general method for input variable selection, which is based on cointegration analysis between candidate inventory variables and oil time spread.

  Before we proceed to cointegration analysis, we need to pre-process the inventory data. There are some common techniques for standardizing oil inventory variables such as dividing the raw inventory data by the working capacity, subtracting the long term inventory average, or adjusting by the historical minimum and maximum inventory with a certain length of period lookback. In the paper “What Does the Convenience Yield Curve Tell Us about the Crude Oil Market?” by R. Alquist, G. Bauer and A. Diez de los Rios [2], there is another method which could de-trend the raw inventory by its production. We reproduce the algorithm described in [2] and compare it with other common techniques.

Systematic Oil Time Spread Trading : Relative Analysis Between Oil Time Spread & Inventory

  After pre-processing the inventory variable, we conduct cointegration analysis. Assuming we have It and yt as the inventory variable and spread variable. The cointegration result is the residual, denoted as et, obtained by applying linear regression between It and yt, where e_t=y_t-kI_t-c. Here k and c are the linear regression coefficients with the highest R2. We also apply the augmented Dickey–Fuller (ADF) test to each cointegration result et to check for the stationarity. If et pass the ADF test, we check the cointegration correlation, defined as r=R2, which is the square root of corresponding R2. Inventory variable with higher r means that it might carry much more information, and thus becoming a more proper input variable for the model. We also need to test the feasibility of the input variable. We provide a statistics summary table for the candidate inventory variables.

  In the table, r is the cointegration correlation, ADF is the ADF test statistic, p is the p-value, W means that the statistic is obtained by the weekly data while M is for the monthly data. We also present explanations to each candidate variable. Cushing I is the Cushing raw inventory data, Cushing De-Trend is the cushing inventory data de-trended with its production by the algorithm described in [2], PADD2 I is the raw PADD2 inventory data, PADD2 De-Trend is the PADD2 inventory data de-trended with its production by the algorithm described in [2], Cushing I/C is the Cushing inventory data divided by its working capacity, Cushing I – 3Y Avg is the Cushing inventory data subtracted with its 3 years moving average, Cushing I – 2Y Avg is the Cushing inventory data subtracted with its 2 years moving average, Cushing MinMax 2Y is the Cushing inventory data adjusted with its 2 years lookback historical minimum and maximum, Cushing MinMax 3Y is the Cushing inventory data adjusted with its 3 years lookback historical minimum and maximum. Such minimum and maximum with m years lookback adjustment is defined by 

\text{m Years MinMax Adjusted Inventory}=\frac{\text{Inventory - m Years Min}}{\text{m Years Max - m Years Min}}

  The critical value for weekly data is -3.44 at 1% significance level and -2.87 at 5% significance level. The critical value for monthly data is -3.49 at 1% significance level and -2.89 at 5% significance level.

Inventory Variable r(W) r(M) ADF(W) ADF(M) p(W) p(M) 
Cushing I 0.580.60 -5.20 -7.31 8.7e-06 1.3e-10 
Cushing De-Trend 0.610.60 -5.43 -7.60 2.9e-06 2.4e-11 
PADD2 I 0.48 0.49 -4.76 -6.92 6.4e-05 1.1e-09 
PADD2 De-Trend 0.58 0.58 -5.22 -7.61 7.8e-06 2.2e-11 
Cushing I/C 0.61 0.61 -5.34 -7.4 4.5e-06 7.5e-11 
Cushing I-3Y Avg 0.62 0.62 -5.53 -7.54 1.8e-06 3.3e-11 
Cushing I-2Y Avg 0.65 0.65 -5.76 -7.87 5.6e-07 5.1e-12 
Cushing MinMax 2Y 0.62 0.60 -5.58 -7.86 1.4e-06 5.4e-12 
Cushing MinMax 3Y 0.60 0.59 -5.46 -7.44 2.5e-06 5.9e-11 

Table 1: Cointegration Correlation Statistics Summary Table

  After analyzing the statistics summary table and input feasibility for each candidate variable, we decide to use weekly Cushing I/C as the final model input variable. The weekly Cushing I/C data maintains the raw data information after adjustment and its cointegration correlation is close to the highest cointegration correlation in the table. Furthermore, inventory over capacity could be viewed as the objective inventory level in percentage unit and thus helping reduce the model parameters complexity. We could easily define the parameters xmin and xmax between 0% and 100%.  

3.3 Model Calibration  

We use the WTI crude oil prices data for the time spread variable, which is the price difference between the front month and second month futures contracts. Based on the import variable selection result, we use Cushing I/C as the input variable for the model. Both variables are within weekly frequency. We apply the Huber loss function as the loss function since it could help reduce the effect of extreme outliers while maintaining their influence to the calibration. The Huber loss function is defined by 

L_{\delta} (y,\hat{y})=
\frac{1}{2} (y-\hat{y})^2  \text{ ,}& |y-\hat{y}|\leq \delta\\
\delta |y-\hat{y}|-\frac{1}{2} \delta^2   \text{,}& otherwise\\

Here y is the observed market time spread while y is the model spread defined by \hat{y}=f(t,T;x)=f(t,T;x | \sigma, k, \bar{x}, x_{min}, x_{max}). For simplicity, we use the same assumption for xmin and xmax as above, where xmin=15% and xmax=100%. Since the dynamics for parameter k is difficult to measure, we also keep it as a constant k=k0, i.e. k0=1. The time variable is also fixed to T-t=1/12 since the time spread lag is one month. Thus, given a certain penalty parameter , the calibration becomes solving the problem below 

\min_{\bar{x}, \sigma} \sum_{i=1}^{n} L_{\delta} (y_i,\hat{y}_i |  k=k_0, x_{min}=15\%, x_{max}=100\%)

We provide a sample calibration graph and a summary table for parameters estimation with different choices of under two cases, k0=1 and k0=1.5.

Figure  6: =0.5,T-t=1/12,=0.42, x=0.63, k=1

Table 2: Parameter Estimation with different choices of and k0=1

Table 3: Parameter Estimation with different choices of and k0=1.5

  From the calibration graph we could observe that the model spread curve actually captures the macro dynamics of the market time spread. It maintains the negative relationship between time spread and inventory and it has a very nice non-linearity. The parameters estimation for x and is relatively stable with proper choice of . Intuitively, x is the long-term mean inventory level while represents the annualized implied volatility for inventory. Although the variance of the market spread is very large and many data points are away from the model spread curve, the calibration result is pretty good for constant parameters setting and one factor modelling. The market data is also affected by lots of other noises which increases the difficulty to do a precise matching. Within such observation, we could also obtain some intuition for potential trading signals for the time spread using inventory data.

3.4 Trading Signals Analysis  

The skewness of the model curve indicates that the mean reversion feature of inventory is stronger at both ends than in the middle, which could potentially benefit us if we trade time spread, when the current inventory level is within some specific regions, based on long-term inventory mean reversion signal, such as longing when inventory excesses its long-term average and shorting when inventory is below its long-term average. We could deduce that there might exist some flip points, such that, when inventory passes the flip boundaries and enters into the specific regions, the long-term inventory mean reversion signal might be strong enough to beat other noises and generate profits. When the inventory does not pass the flip points, such mean reversion feature might be distorted by other noises and we could introduce a short-term inventory momentum signal to ensure that we have a systematic trading strategy. 

Figure  7: Intuition

3.4.1 Notations  

t: The discrete time points for the data.

It: Inventory at t.

Ct: Capacity at t.

St=It/Ct: Inventory level at t.

MA(t,k): The inventory moving average level with k weeks lookback at t.

Mt=St-MA(t,k1): Inventory subtracted by its long-term moving average with k1 weeks lookback at t. 

M’t=St-MA(t,k2): Inventory subtracted by its short-term moving average with k2 weeks lookback at t.

MR: The long-term inventory mean reversion signal.

STMT: The short-term inventory momentum signal.

3.4.2 Flip Boundaries  

We introduce an extra parameter , and we could apply MR signal when |Mt|≥ϵ otherwise we could apply STMT signal. Using as the flip boundary is one of the common techniques but we need to take care about its robustness.

3.4.3 Signals  

Long-term Mean Reversion Signal (MR): If Mt0, then signal=+1, else signal=-1. 

Short-term Momentum Signal (STMT): If M’tM’t-1, then signal=-1, else signal=+1. 

Combining signals and flip boundaries together, we have the following strategy: If |Mt|≥ϵ, we apply MR signal, then if Mt0, we have signal=+1, else signal=-1. If |Mt|≤ϵ, we apply STMT signal, then if M’tM’t-1, we have signal=-1, else signal=+1. Notice that the strategy guarantees that we apply exactly one of the signal when we rebalance our position, either MR or STMT, based on the filp points. 

3.4.4 Trading Settings  

The inventory reports are released each Wednesday, we enter the market on Wednesday upon observing the new inventory data, and we hold the same position until the next inventory reports releasing date. That means we rebalance our position each week, but the profits and losses are calculated on a daily basis. For the rolling, we roll our position 5 days before the front month contract expiration date. We do not consider transaction cost since there are ways to execute the time spread with no transaction cost in the industry. Our benchmark is the consistent short signal, which means holding short position throughout the whole backtesting period. The short benchmark works because the market is under structured contango for recent years. The short benchmark Sharpe ratio is 0.50. 

3.4.5 Backtesting Results  

We present a summary Sharpe ratio table and a return over maximum drawdown (RoD) table with =0.1 for different combinations of long-term moving average lookback k1 and short-term moving average lookback k2. Notice that k1 affects the MR signal and k2 affects the STMT signal.

  From the Sharpe ratio table, we could see the overall performance beats the short benchmark. The maximum Sharpe ratio is 0.94 and is achieved with k1=28 and k2=7. However, when k1 increases to 32, the Sharpe ratio drops from 0.94 to 0.65, and the RoD drops from 0.43 to 0.30, which is not stable. On the other hand, the Sharpe ratio around k1=48 and k2=8 looks relatively stable and it achieves a local maximum of RoD in its threshold, which is 0.40. We are interested in knowing its robustness with a slight change in . Intuitively, cannot be too large or too small, as it will reduce our strategy to apply one single signal only, either STMT or MR , both cases might lead to an extreme and unstable Sharpe ratio or RoD by accident. We restrict in a range of [0.05,0.15]. We provide another sample graph for the Sharpe ratio with in the range and fixed k1=48 and k2=8, and we could see that the Sharpe ratio around =0.1 is relatively stable.   

Table 4: Sharpe Ratio Summary Table =0.1

Table 5: RoD Summary Table =0.1

Figure 8: Sharpe Ratio with k1=48, k2=8

Systematic Oil Time Spread Trading : Relative Analysis Between Oil Time Spread & Inventory

4. Conclusions  

The calibration result of the time spread model gives us a high level overview of the possibility to model the oil time spread using stochastic mean reversion model for oil inventory. However, such modelling under constant parameters setting could only capture the macro dynamics of the market time spread. Part of the difficulty is that the market data includes many unknown noises and some data are missing which also adds restrictions to the calibration. More research could be conducted to study the model parameters k and to improve the model calibration.

  The backtesting results for our strategy are very impressive since we obtain some steady Sharpe ratio between 0.7 and 0.9 within a threshold of proper choices of as flip boundaries, and the strategy beats the short benchmark. The potential risks such as drawdown measured by RoD in the backesting is not as good as the Sharpe ratio. But as a one factor strategy, such backtesting performance is reasonable and fair enough to convince us that the performance is not a result from data mining.

Systematic Oil Time Spread Trading : Relative Analysis Between Oil Time Spread & Inventory

Written by Yifan Wei and Xingjia Gu
Advisor: Ilia Bouchouev


[1]  I. Bouchouev, 2020. “A Stylized Model of the Oil Price Squeeze.” Working Paper, Pentathlon Investments.

[2]  R. Alquist, G. Bauer and A. Diez de los Rios, 2014. “What Does the Convenience Yield Curve Tell Us about the Crude Oil Market?” Staff Working Paper 14-42, Bank of Canada.

Systematic Oil Time Spread Trading : Relative Analysis Between Oil Time Spread & Inventory

Courant Institute of Mathematical Sciences, New York University,

New York, NY 10012-1185 USA

British Petroleum’s Head Quant Cetin Karakus Talks With Rebellion Research CEO Alexander Fleiss on Machine Learning​ & Big Data​ in the Energy​ Market

The Automation of the Oil and Gas Industry