A study on the nonlinear relationship

between WTI crude oil spot and futures prices

Eun-Young Kim*

		< Abstract >
	.
Keywords：Threshold Vector Error Correction Models, Nonlinearity, Structural Breaks, Regimes, Crude Oil Spot and Futures Prices

金融工學硏究

* Ph.D. in Economics, Researcher of Statistical Research Institute in Pusan National University

E-mail：keyivy@empal.com

Ⅰ. Introduction

Crude oil spot and futures prices and their volatilities have been under the great attention from the industries and academia. Unfortunately the studies on these data have been limited to conventional methodologies for analysis such as cointegration analysis, VAR(Vector AutoRegressive) models, VECM(Vector Error Correction Models) and so on. So I'd like to widen the established area just one step forward through the application of multiple-regime bivariate TVECM(Threshold Vector Error Correction Models). Threshold models has been devised by Tong(1978) and TVECM has been created by adding threshold variables to the ordinary VECM. This transformation allows me to use basis as a threshold variable and investigate investors' arbitrage behaviors under different market situations. Of course, there is a constraint for this, because a threshold variable in TVECM is an error correction term, so any two variables in this term should be cointegrated. This constraint can be however easily overcome in the fields of economics and finance, because there are a lot of pairs which are cointegrated. That is, threshold models including TVECM are said to have strong applicability.

To guarantee quite different market situations, I used rather long-period raw data for dividing it into several subsamples. For this, I used the ‘dynamic programming algorithm’ presented by Bai and Perron(2003). On the 5 subsample periods, I also performed many preliminary tests before the application of multiple-regime bivariate TVECM. To add some interest on this study, I set more than 2 regimes, indeed 3 regimes(lower, middle, and upper regimes). Hence I characterized the markets of 5 subsamples as a backwardation market or a contango marekt and applied 3-regime bivariate TVECM on every subsample for getting some hints about investors' arbitrage behaviors under relatively different market conditions.

This paper is organized as follows: I explained some important methodologies for analysis in this study in section 2. I provided the structure and nature of the data in section 3. I reported all the test results including 3-regime bivariate TVECM in section 4, and section 5 concludes this study.

Ⅱ. Research Methodologies

1. Nonlinearity Test(BDS Test)

Brock, Dechert, and Scheinkman(1987) suggested a new test statistic for checking whether the I.I.D. assumption of any time series can be satisfied or not. This test statistic was shown to be used for detecting serial dependence of any time series in the residual analysis(Kim et al., 2003). Two different hypotheses used in hypothesis testing are as follows:

Null hypothesis([eqn: H _{0}]): This time series is I.I.D..

Alternative hypothesis([eqn: H _{1}]): This time series is not I.I.D..

BDS test statistic can be calculated like this. There is a time series with [eqn: N] observations, which should be the first difference of the natural logarithms of raw data in a time series.

[eqn: LEFT { x _{i} RIGHT } = LEFT [ x _{1} ,`x _{2} ,` CDOTS ,`x _{N} RIGHT ]]

You can select [eqn: m], embedding dimension and embed the time series into [eqn: m]-dimensional vectors by taking each [eqn: m] successive points in this time series. This procedure turns the series of scalars into a series of vectors with overlapping elements.

[eqn: x _{1} ^{m} = LEFT ( x _{1} ,`x _{2} ,` CDOTS ,`x _{m} RIGHT )# x _{2} ^{m} = LEFT ( x _{2} ,`x _{3} ,` CDOTS ,`x _{m+1} RIGHT )# ```````````````````````````````` VDOTS # x _{N-m} ^{m} = LEFT ( x _{N-m} ,`x _{N-m+1} ,` CDOTS ,`x _{N} RIGHT )]

You can compute the correlation integral, which measures the spatial correlation among the points, by adding the number of pairs of points [eqn: LEFT ( i,`j RIGHT )], where [eqn: 1 LEQ i LEQ N] and [eqn: 1 LEQ j LEQ N], in the [eqn: m]-dimensional space which are “close” in the sense that the points are within a radius or tolerance [eqn: epsilon] of each other.

[eqn: C _{epsilon ,m} = {1} over {N _{m} LEFT ( N _{m} -1 RIGHT )} sum _{i != j} ^{} I _{i,j; epsilon }]                          (1)

where [eqn: I _{i,j; epsilon } =1] if [eqn: LEFT ∥ x _{i} ^{m} -x _{j} ^{m} RIGHT ∥ LEQ epsilon]

[eqn: I _{i,j; epsilon } =0] if [eqn: LEFT ∥ x _{i} ^{m} -x _{j} ^{m} RIGHT ∥ > epsilon]

Brock, Dechert, and Scheinkman(1987) proved that if this time series is I.I.D.

[eqn: C _{epsilon ,m} APPROX LEFT [ C _{epsilon ,1} RIGHT ] ^{m}]                              (2)

If [eqn: {N} over {m}] is greater than 200, the values of [eqn: {epsilon } over {sigma }] are from 0.5 to 2(Lin, 1997) and [eqn: m] can have its values between 2 and 5(Brock, Dechert, and Scheinkman, 1988) and the quantity [eqn: LEFT [ C _{epsilon ,m} - LEFT ( C _{epsilon ,1} RIGHT ) ^{m} RIGHT ]] has an asymptotic normal distribution with mean zero and variance [eqn: V _{epsilon m}] defined as follows:

[eqn: V _{epsilon ,m} =4 LEFT [ K ^{m} +2 sum _{j=1} ^{m-1} K ^{m-j} C _{epsilon } ^{2j} + LEFT ( m-1 RIGHT ) ^{2} -m ^{2} KC _{epsilon } ^{2m-2} RIGHT ]]          (3)

where [eqn: K=K _{epsilon } = {6} over {N _{m} LEFT ( N _{m} -1 RIGHT ) LEFT ( N _{m} -2 RIGHT )} sum _{i<j<N} ^{} h _{i,j,N; epsilon }]

[eqn: h _{i,j,N; epsilon } = {LEFT [ I _{i,j; epsilon } I _{j,N; epsilon } +I _{i,N; epsilon } I _{N,j; epsilon } +I _{j,i; epsilon } I _{i,N; epsilon } RIGHT ]} over {3}]

Consequently, the BDS test statistic can be formulated like this:                                    [eqn: BDS _{epsilon ,m} = {sqrt {N} LEFT [ C _{epsilon ,m} - LEFT ( C _{epsilon ,1} RIGHT ) ^{m} RIGHT ]} over {sqrt {V _{epsilon ,m}}}]                       (4)

It has a limiting standard normal distribution under the null hypothesis of I.I.D. as  [eqn: N -> INF] and obtains its critical values using the standard normal distribution. Due to the type of alternative hypothesis BDS test becomes a two-tailed test and the rejection regions are placed at the two tails of probability density curve. Therefore if you can reject the null hypothesis under a given significance level [eqn: alpha], the time series is not I.I.D., so this data should be fitted with nonlinear models.

2. The Detection of Structural Change Points

1) Tests on structural change

The most important classes of tests on structural change are the tests from the generalized fluctuation(Kuan and Hornik, 1995) and the tests based on [eqn: F] statistics(Hansen, 1992; Andrews, 1993; Andrews and Ploberger, 1994). The former includes the CUSUM and MOSUM tests and the fluctuation test, while the latter includes the Chow and the sup[eqn: F] test.

① Generalized fluctuation tests

These tests fit a model to the given data and derive an empirical process, that captures the fluctuation either in residuals or in estimates. The limiting processes are known for these empirical processes, so the boundaries can be calculated and the probability to cross them under the null hypothesis(no structural change) is [eqn: alpha]. If the empirical process path crosses theses boundaries, the fluctuation is improbably large and consequently the null hypothesis should be rejected at significance level [eqn: alpha]. Hence we can infer the presence of structural change at this time point.

② [eqn: F] statistics-based tests

The tests based on [eqn: F] statistics can be used alternatively to investigate whether the null hypothesis of “no structural change” holds. While the generalized fluctuation tests are fit for various patterns of structural changes, the [eqn: F] statistics-based tests first suggested by Chow(1960) and followed by Andrews(1993) and Andrews and Ploberger(1994) are only for testing against a single shift alternative.

Chow(1960) test was designed to perform this function with the prior knowledge of change point itself. It has been regarded as the major drawback of “Chow test”, but Andrews and Ploberger contributed to overcome the weak point of Chow test and develop this type of tests by allowing for no prior knowledge of change point.

Bai and Perron(1998, 2003) went several more steps forward by extending the cases of structural changes to no break versus [eqn: l] breaks and to [eqn: l] versus [eqn: l+1] breaks. [eqn: l] can be an arbitrary but fixed number.

2) The methodology of Bai and Perron(2003)

Bai and Perron(2003) dealt with the problem of structural break point estimation and presented an efficient algorithm called ‘the dynamic programming algorithm’ to obtain global minimizers of the sum of squared residuals. This algorithm uses at most least-squares operations of order [eqn: O LEFT ( T ^{2} RIGHT )] for any number of structural changes [eqn: m], while a standard grid procedure does least-squares operations of order [eqn: O LEFT ( T ^{m} RIGHT )]. With this algorithm they successfully performed the sequential detection of structural change points by increasing the number of breaks one at a time.

They considered the following multiple linear regression with [eqn: m] structural breaks ([eqn: m+1] subsamples):

[eqn: y _{t} =x _{t} ` prime ` beta +z _{t} ` prime delta _{j} +u _{t`} `````` LEFT ( t=T _{j-1} +1,` CDOTS ,`T _{j} RIGHT )]                 (5)

for [eqn: j=1,` CDOTS ,`m+1]. [eqn: y _{t}] is the observed response variable at time [eqn: t], [eqn: x _{t} LEFT ( p TIMES 1 RIGHT )] and [eqn: z _{t} LEFT ( q TIMES 1 RIGHT )] are vectors of covariates, [eqn: beta] and [eqn: delta _{j} LEFT ( j=1,` CDOTS ,`m+1 RIGHT )] are corresponding vectors of coefficients and [eqn: u _{t}] is the disturbance at time [eqn: t]. The break points[eqn: LEFT ( T _{1} ,` CDOTS ,`T _{m} RIGHT )] are treated as unknown([eqn: T _{0} =0] and [eqn: T _{m+1} =T]). They tried to estimate the unknown regression coefficients together with the break points when [eqn: T] observations on [eqn: LEFT ( y _{t} ,`x _{t} ,`z _{t} RIGHT )] are in their hands.

Before using dynamic programming algorithm, the construction of the triangular matrix of sums of squared residuals for all possible segments should be made. The vertical axis of this matrix indicates the initial date of a segment and the horizontal axis the ending date, so each element means an estimated sum of squared residuals corresponding to the associated segment. The global sum of squared residuals for any [eqn: m]-partition[eqn: LEFT ( T _{1} ,` CDOTS ,`T _{m} RIGHT )] and for any value of [eqn: m] should be a particular combination of these [eqn: T LEFT ( T+1 RIGHT ) /2] sums of squared residuals. The estimates of break dates [eqn: LEFT ( hat{T} ` _{1} ,` CDOTS ,` hat{T} ` _{m} RIGHT )] correspond to this linear combination with a minimal value.   They wrote that in practice, less than [eqn: T LEFT ( T+1 RIGHT ) /2] segments were possible in this paper. For this procedure, some minimum distance [eqn: h] between each breaks may be imposed which can be quite helpful for reducing the number of segments along with other factors.

If the sum of squared residuals of the relevant segments have been calculated and stored, a dynamic programming algorithm can be applied. This method exams optimal one-break partitions sequentially. Let [eqn: SSR LEFT ( LEFT { T _{r,n} RIGHT } RIGHT )] be the sum of squared residuals associated with the optimal partition having [eqn: r] breaks using the first [eqn: n] observations. The optimal partition solves the following recursive problem:

[eqn: SSR LEFT ( LEFT { T _{r,n} RIGHT } RIGHT ) = pile{min#mh LEQ j LEQ T-h} LEFT [ SSR LEFT ( LEFT { T _{m-1,`j} RIGHT } RIGHT ) +SSR LEFT ( j+1`,`T RIGHT ) RIGHT ]]     (6)

This procedure begins by evaluating the optimal one-break partition for all subsamples that allow a possible break ranging from observations [eqn: h] to [eqn: T-mh]. The entire process is composed of the following 4 steps. The first step is to store a set of [eqn: T- LEFT ( m+1 RIGHT ) h+1] optimal one-break partitions with their associated sum of squared residuals. The second step is to make a search for optimal partitions with two breaks. For each of the possible ending dates, the procedure looks at which one-break partition can be input to achieve a minimal sum of squared residuals. The result is a set of [eqn: T- LEFT ( m+1 RIGHT ) h+1] optimal two breaks partitions. This method continues sequentially until a set of [eqn: T- LEFT ( m+1 RIGHT ) h+1] optimal [eqn: LEFT ( m-1 RIGHT )] breaks partitions are obtained. The third step is to check which of these optimal [eqn: LEFT ( m-1 RIGHT )] breaks partitions yields an overall minimal sum of squared residuals when combined with an additional segment. The fourth step is to create a single optimal [eqn: m] breaks partition.

They also constructed confidence intervals for the parameters [eqn: beta] and [eqn: delta], and the break dates. They proposed their own test statistic for multiple breaks as follows:

[eqn: F _{T} LEFT ( lambda _{1} ,` CDOTS ,` lambda _{k} ;`q RIGHT ) = {1} over {T} LEFT ( {T- LEFT ( k+1 RIGHT ) q-p} over {kq} RIGHT ) hat{delta } ` prime R` prime LEFT ( R hat{V} ` LEFT ( hat{delta } RIGHT ) R` prime RIGHT ) ^{-1} R hat{delta }]         (7)

This is a sup[eqn: F] type test of no structural break[eqn: LEFT ( m=0 RIGHT )] versus [eqn: m=k] breaks. [eqn: LEFT ( T _{1} ,` CDOTS ,`T _{k} RIGHT )] is a partition such that [eqn: T _{i} = LEFT [ T lambda _{i} RIGHT ] LEFT ( i=1,` CDOTS ,`k RIGHT )], where [eqn: lambda _{i} = {T _{i}} over {T}] can be called ‘break fractions’. [eqn: R] is the conventional matrix such that [eqn: LEFT ( R delta RIGHT ) prime = LEFT ( delta _{1} prime - delta _{2} prime ,` CDOTS ,` delta _{k} prime - delta _{k+1} prime ` RIGHT )]. [eqn: hat{V} ` LEFT ( hat{delta } RIGHT )] is the estimate of the variance covariance matrix of [eqn: hat{delta }] which is robust to serial correlation and heteroscedasticity. According to Andrews(1993) and others, sup[eqn: F _{T} LEFT ( k;q RIGHT ) =F _{T} LEFT ( hat{lambda } ` _{1} ,` CDOTS ,` hat{lambda } ` _{k} ;`q RIGHT )] holds. Here [eqn: hat{lambda } ` _{1} ,` CDOTS ,` hat{lambda } ` _{k}] minimize the global sum of squared residuals which has the same meaning of maximizing the [eqn: F]-test assuming spherical errors. In this papar, I'll use ‘dynamic programming algorithm’ presented by Bai and Perron(2003) to detect the multiple structural break points.

3. Threshold Vector Error Correction Model(TVECM)

To understand the structure and nature of bivariate threshold vector error correction models, it could be better to see first a much simpler Threshold model as an example such as SETAR(Self-Exciting Threshold AutoRegressive) models.

If [eqn: LEFT { x _{t} RIGHT }] follows a [eqn: k]-regime SETAR model with threshold variable [eqn: x _{t-d}],

[eqn: x _{t} = phi _{0} ^{LEFT ( i RIGHT )} + phi _{1} ^{LEFT ( i RIGHT )} x _{t-1} + CDOTS + phi _{p} ^{LEFT ( i RIGHT )} x _{t-p} + epsilon _{t} ^{LEFT ( i RIGHT )}]  ([eqn: i=1,` CDOTS ,`k])            (8)

where [eqn: x _{t-d}] is a threshold variable([eqn: tau _{i-1} LEQ x _{t-d} < tau _{i}]) and [eqn: tau _{i}]'s are integers such as [eqn: - INF = tau _{0} < tau _{1} < CDOTS < tau _{k-1} < tau _{k} = INF]. [eqn: k] and [eqn: d] are all positive integers and [eqn: d] is called ‘a delay parameter’. [eqn: LEFT ( i RIGHT )] means regimes and [eqn: p] is a lag of AR(AutoRegressive) models. [eqn: LEFT { epsilon _{t} ^{LEFT ( i RIGHT )} RIGHT }] is assumed to be I.I.D. with mean 0 and variance [eqn: sigma _{i} ^{2}] and independent  for different [eqn: LEFT ( i RIGHT )] regimes.

Let's take an example of SETAR models with [eqn: k=2], [eqn: d=1], [eqn: tau =0.5], and [eqn: p=2].

[eqn: x _{t} = phi _{0} ^{LEFT ( 1 RIGHT )} + phi _{1} ^{LEFT ( 1 RIGHT )} x _{t-1} + phi _{2} ^{LEFT ( 1 RIGHT )} x _{t-2} +a _{t} ^{LEFT ( 1 RIGHT )}]   [eqn: LEFT ( x _{t-1} <0.5 RIGHT )]              (9)

[eqn: x _{t} = phi _{0} ^{LEFT ( 2 RIGHT )} + phi _{1} ^{LEFT ( 2 RIGHT )} x _{t-1} + phi _{2} ^{LEFT ( 2 RIGHT )} x _{t-2} +a _{t} ^{LEFT ( 2 RIGHT )}]   [eqn: LEFT ( x _{t-1} GEQ 0.5 RIGHT )]              (10)

In this SETAR model, there are two regimes in the form of AR models which can be explained with linear models, but the entire model is explicitly nonlinear. This is one of attractions of threshold models.

The threshold vector error correction models I'll use in this paper is the vector error correction models with a threshold variable as an error correction term. This means that if you have cointegrated variables which are composed of error correction term, you can investigate the short-run dynamics of two variables of interest based on the different regimes(situations) caused by a threshold variable.

The general bivariate threshold vector error correction models can be specified as follows:

[eqn: TRIANGLE f _{t} = alpha _{i} ^{f} LEFT ( b _{t-d} RIGHT ) + sum _{j=1} ^{p} beta _{i,t-j} ^{f} TRIANGLE f _{t-j} + sum _{j=1} ^{q} gamma _{i,t-j} ^{f} TRIANGLE s _{t-j} + epsilon _{t} ^{f}]                  (11)

[eqn: TRIANGLE s _{t} = alpha _{i} ^{s} LEFT ( b _{t-d} RIGHT ) + sum _{j=1} ^{p} beta _{i,t-j} ^{s} TRIANGLE f _{t-j} + sum _{j=1} ^{q} gamma _{i,t-j} ^{s} TRIANGLE s _{t-j} + epsilon _{t} ^{s}]   ([eqn: i=1,` CDOTS ,`k])   (12)

where [eqn: b _{t-d} (=f _{t-d} -s _{t-d}], [eqn: tau _{i-1} LEQ b _{t-d} < tau _{i}]) is a threshold variable. In this paper, [eqn: f _{t}] and [eqn: s _{t}] are natural logarithmic values of [eqn: F _{t}](the closing price of the nearby WTI futures contract at time [eqn: t]) and [eqn: S _{t}](the spot price at time [eqn: t]), respectively. [eqn: tau _{i}]'s are integers such as [eqn: - INF = tau _{0} < tau _{1} < CDOTS < tau _{k-1} < tau _{k} = INF]. [eqn: k] and [eqn: d] are all positive integers and [eqn: d] is called ‘a delay parameter’. [eqn: LEFT ( i RIGHT )] means regimes, [eqn: DELTA] is a difference operator, and [eqn: p] and [eqn: q] are  lags of AR(AutoRegressive) models. [eqn: LEFT { epsilon _{t} ^{LEFT ( i RIGHT )} RIGHT }] is assumed to be I.I.D. with mean 0 and variance [eqn: sigma _{i} ^{2}] and independent  for different [eqn: LEFT ( i RIGHT )] regimes.

Let's take an example of bivariate threshold vector error correction models with [eqn: k=3], [eqn: d=1], [eqn: p=q=2], [eqn: tau _{1} < tau _{2}], and three regimes(regime 1: lower, regime 2: middle, regime 3: upper). There can be 6 equations as follows:

① lower regime[eqn: LEFT ( b _{t-1} < tau _{1} RIGHT )]

[eqn: TRIANGLE f _{t} = alpha _{1} ^{f} LEFT ( b _{t-1} RIGHT ) + beta _{1,t-1} ^{f} TRIANGLE f _{t-1} + beta _{1,t-2} ^{f} TRIANGLE f _{t-2} + gamma _{1,t-1} ^{f} TRIANGLE s _{t-1} + gamma _{1,t-2} ^{f} TRIANGLE s _{t-2} + epsilon _{t} ^{f}]

[eqn: TRIANGLE s _{t} = alpha _{1} ^{s} LEFT ( b _{t-1} RIGHT ) + beta _{1,t-1} ^{s} TRIANGLE f _{t-1} + beta _{1,t-2} ^{s} TRIANGLE f _{t-2} + gamma _{1,t-1} ^{s} TRIANGLE s _{t-1} + gamma _{1,t-2} ^{s} TRIANGLE s _{t-2} + epsilon _{t} ^{s}]

② middle regime[eqn: LEFT ( tau _{1} LEQ b _{t`-1} < tau _{2} RIGHT )]

[eqn: TRIANGLE f _{t} = alpha _{2} ^{f} LEFT ( b _{t-1} RIGHT ) + beta _{2,t-1} ^{f} TRIANGLE f _{t-1} + beta _{2,t-2} ^{f} TRIANGLE f _{t-2} + gamma _{2,t-1} ^{f} TRIANGLE s _{t-1} + gamma _{2,t-2} ^{f} TRIANGLE s _{t-2} + epsilon _{t} ^{f}]

[eqn: TRIANGLE s _{t} = alpha _{2} ^{s} LEFT ( b _{t-1} RIGHT ) + beta _{2,t-1} ^{s} TRIANGLE f _{t-1} + beta _{2,t-2} ^{s} TRIANGLE f _{t-2} + gamma _{2,t-1} ^{s} TRIANGLE s _{t-1} + gamma _{2,t-2} ^{s} TRIANGLE s _{t-2} + epsilon _{t} ^{s}]

③ upper regime[eqn: LEFT ( tau _{2} LEQ b _{t-1} RIGHT )]

[eqn: TRIANGLE f _{t} = alpha _{3} ^{f} LEFT ( b _{t-1} RIGHT ) + beta _{3,t-1} ^{f} TRIANGLE f _{t-1} + beta _{3,t-2} ^{f} TRIANGLE f _{t-2} + gamma _{3,t-1} ^{f} TRIANGLE s _{t-1} + gamma _{3,t-2} ^{f} TRIANGLE s _{t-2} + epsilon _{t} ^{f}]

[eqn: TRIANGLE s _{t} = alpha _{3} ^{s} LEFT ( b _{t-1} RIGHT ) + beta _{3,t-1} ^{s} TRIANGLE f _{t-1} + beta _{3,t-2} ^{s} TRIANGLE f _{t-2} + gamma _{3,t-1} ^{s} TRIANGLE s _{t-1} + gamma _{3,t-2} ^{s} TRIANGLE s _{t-2} + epsilon _{t} ^{s}]

Ⅲ. Data

The raw data for this research are WTI crude oil closing prices of the nearby futures contracts traded in NYMEX and corresponding spot prices. The sample period is from January 2, 1986 to December 30, 2010 and data type is daily time series. I considered the trading days when both futures and spot markets opened and the final data set has 6,244 observations in total. I downloaded these data from the homepage of the Energy Information Administration and edited them for the purpose of this analysis.

<Figure 1> shows the plot of WTI crude oil closing prices of the nearby futures contracts and basis([eqn: F _{t} -S _{t}]), and <Figure 2> shows the same plot in the case of log-transformed data. The log-transformed data will be used in the following analysis.

<Figure 1> WTI crude oil closing prices of the nearby futures contracts and basis(raw data)

<Figure 2> WTI crude oil closing prices of the nearby futures contracts and basis(log-transformed data)

Ⅳ. Empirical Analysis

1. Stationarity test

The unit root tests were performed on the log-transformed data([eqn: f _{t} =ln LEFT ( F _{t} RIGHT )] and [eqn: s _{t} =ln LEFT ( S _{t} RIGHT )]) and their first difference data for checking stationarity of these time series. The results of Augmented Dickey Fuller test and Phillips Perron test are reported in <Table 1>.

<Table 1> Unit root test results

	[eqn: f _{t}]	[eqn: TRIANGLE f _{t}]	[eqn: s _{t}]	[eqn: TRIANGLE s _{t}]
ADF	-2.9416	-17.2893[eqn: {} ^{***}]	-3.0074	-17.1793[eqn: {} ^{***}]
PP	-17.6701	-5978.482[eqn: {} ^{***}]	-18.3086	-5618.572[eqn: {} ^{***}]
Notes. [eqn: {} ^{***}] indicates that the null hypothesis can be rejected at the 1% significance level.

The null hypotheses of the above unit root tests are ‘this time series has a unit root’. The test results mean that log-transformed time series is not stationary but their first difference time series become stationary. This codition guarantees that the two log-transformed time series can be tested for cointegration as well.

2. Nonlinearity test

BDS test as a kind of nonlinearity tests was performed on the entire sample period. The data don't have any linear dependence because they have been log-transformed and first differenced. The null and alternative hypotheses are as follows:

[eqn: H _{0}]: This time series is I.I.D..

[eqn: H _{1}]: This time series is not I.I.D..

<Table 2> shows the nonlinearity test results.

<Table 2> Nonlinearity test results

• in the case of [eqn: TRIANGLE f _{t}]

embedding dimension[eqn: LEFT ( m RIGHT )]	[eqn: epsilon / sigma]	BDS test statistic
2	0.0128	11.7992[eqn: {} ^{***}]
	0.0256	13.7647[eqn: {} ^{***}]
	0.0384	16.6785[eqn: {} ^{***}]
	0.0511	19.2069[eqn: {} ^{***}]
3	0.0128	16.4283[eqn: {} ^{***}]
	0.0256	17.9193[eqn: {} ^{***}]
	0.0384	21.0021[eqn: {} ^{***}]
	0.0511	23.5035[eqn: {} ^{***}]
Notes. [eqn: {} ^{***}] indicates that the null hypothesis can be rejected at the 1% significance level.

• in the case of [eqn: TRIANGLE s _{t}]

embedding dimension[eqn: LEFT ( m RIGHT )]	[eqn: epsilon / sigma]	BDS test statistic
2	0.0131	12.8271[eqn: {} ^{***}]
	0.0263	15.4343[eqn: {} ^{***}]
	0.0394	17.8306[eqn: {} ^{***}]
	0.0526	18.6221[eqn: {} ^{***}]
3	0.0131	17.0067[eqn: {} ^{***}]
	0.0263	19.7777[eqn: {} ^{***}]
	0.0394	22.6229[eqn: {} ^{***}]
	0.0526	23.4999[eqn: {} ^{***}]
Notes. [eqn: {} ^{***}] indicates that the null hypothesis can be rejected at the 1% significance level.

As a result of BDS test, the null hypothesis can be rejected even at the 1% significance level. It can be therefore valid to fit nonlinear models to these data.

3. The detection of structural break points

According to dynamic programming algorithm presented by Bai and Perron(2003), the structural break point [eqn: m=4] has been chosen by the criterion of the minimum RSS(Residual Sum of Squares). Another criterion for selecting the appropriate case, the minimum BIC has been added. The values of RSS and BIC are reported in <Table 3> and the division of the entire sample period into 5 subsamples according to these two criteria is summarized in <Table 4>.

<Table 3> The values of RSS and BIC according to [eqn: m]

[eqn: m]	0	1	2	3	4	5
RSS	2116.16	438.93	275.90	265.28	247.36	267.26
BIC	10980.30	1177.49	-1703.79	-1931.36	-2350.49	-1849.86

<Table 4> The subsamples divided by test results

	period 1	period 2	period 3	period 4	period 5
period	1986/01/02~ 1989/09/25	1989/09/26~ 1993/06/17	1993/06/18~ 1999/09/07	1999/09/08~ 2004/08/05	2004/08/06~ 2010/11/30
# of obs.	936	936	1562	1223	1587
ratio(%)	15.00%	15.00%	25.02%	19.59%	25.42%

4. Descriptive statistics for subsamples

<Table 5> reports descriptive statistics for 5 subsamples.

<Table 5> The descriptive statistics for 5 subsamples

• in the case of [eqn: f _{t}]

	period 1	period 2	period 3	period 4	period 5
mean	2.8346	3.6070	2.8888	3.3604	4.2244
maximum	3.2798	3.6993	3.2817	3.7935	4.9787
minimum	2.3437	2.7279	2.3721	2.8593	3.5225
standard deviation	0.1600	0.1583	0.1753	0.1686	0.2721
skewness	-0.4254	1.8138	-0.4061	-0.2649	0.3147
kurtosis	2.9804	6.6926	3.1008	3.2896	3.1500
Jarque-Bera	28.2426 (0.0000)	1044.9990 (0.0000)	43.5874 (0.0000)	18.5805 (0.0000)	27.6792 (0.0000)

• in the case of [eqn: s _{t}]

	period 1	period 2	period 3	period 4	period 5
mean	2.8362	3.0671	2.8900	3.3620	4.2226
maximum	3.2783	3.7153	3.2790	3.7930	4.9789
minimum	2.3273	2.7363	2.3814	2.8622	3.4105
standard deviation	0.1605	0.1611	0.1757	0.1695	0.2738
skewness	-0.4275	1.7448	-0.3985	-0.2818	0.2939
kurtosis	2.9865	6.5079	3.0985	3.2709	3.1834
Jarque-Bera	28.5231 (0.0000)	954.8354 (0.0000)	41.9819 (0.0000)	19.9211 (0.0000)	25.0723 (0.0000)

Stationarity test(unit root test), nonlinearity test(BDS test), cointegration test, and normality test(Jarque-Bera test) were performed on two log-transformed time series for 5 subsamples. All the results were not reported on this paper for the save of space. Test results said that ① log-transformed futres and spot time series are not stationary but they turned into stationary ones after first difference, ② They proved to be fit for nonlinear modelling in every 5 subsample, ③ They are cointegrated in each subsamples, ④ according to the results of Jarque-Bera test, they don't follow a normal distribution in every subsample, as you can check this result with skewness and kurtosis as well.

5. The application of bivariate TVECM

The entire sample period has been divided into 5 subsamples via the detection of strcutural break points[eqn: LEFT ( m=4 RIGHT )]. This allows me to investigate whether there is a significant difference between investor's behaviors under the different market situations using bivariate TVECM, so the first thing we should do is to define the character of the market. For this, ‘basis’ can be used as a criterion. The market can be classified as ‘backwardation(negative basis)’ if spot prices are above futures ones in a specific period, and the market can be classified as ‘contango(positive basis)’ if futures prices are above spot ones. <Table 6> summarizes this characterization of markets in 5 subsamples.

<Table 6> The characterization of subsamples

		period 1	period 2	period 3	period 4	period 5
period		1986/01/02 1989/09/25	1989/09/26~ 1993/06/17	1993/06/18~ 1999/09/07	1999/09/08~ 2004/08/05	2004/08/06~ 2010/11/30
# of obs.		936	936	1562	1223	1587
basis	(+)	381(40.71)	358(38.25)	593(37.96)	535(43.74)	881(55.51)
	0	61( 6.52)	74( 7.91)	151( 9.67)	49( 4.01)	93( 5.86)
	(-)	494(52.78)	504(53.85)	818(52.37)	639(52.25)	613(38.63)
characteristic of the market		backwardation	backwardation	backwardation	backwardation	contango

The estimation results of 3-regime bivariate TVECM has been reported in <Table 7>.

<Table 7> Bivariate TVECM estimation results

model	[eqn: TRIANGLE f _{t} = alpha _{i} ^{f} ecm _{t-1} + beta _{i} ^{f} LEFT ( L RIGHT ) TRIANGLE f _{t-1} + gamma _{i} ^{f} LEFT ( L RIGHT ) TRIANGLE s _{t-1} + epsilon _{t} ^{f}] [eqn: TRIANGLE s _{t} = alpha _{i} ^{s} ecm _{t-1} + beta _{i} ^{s} LEFT ( L RIGHT ) TRIANGLE f _{t-1} + gamma _{i} ^{s} LEFT ( L RIGHT ) TRIANGLE s _{t-1} + epsilon _{t} ^{s}]

period 1		lower regime	middle regime	upper regime
[eqn: TRIANGLE f _{t}]	[eqn: alpha]	-0.6545[eqn: {} ^{***}]	5.9471	-0.1529
	[eqn: beta _{t-1}]	0.4998[eqn: {} ^{*}]	0.1972	0.1408
	[eqn: beta _{t-2}]	1.2851[eqn: {} ^{***}]	0.7915	0.0318
	[eqn: gamma _{t-1}]	-0.6979[eqn: {} ^{**}]	-0.5782	-0.0069
	[eqn: gamma _{t-2}]	-1.5076[eqn: {} ^{***}]	-0.1582	-0.0364
[eqn: TRIANGLE s _{t}]	[eqn: alpha]	-0.2353	5.2431	0.3615[eqn: {} ^{**}]
	[eqn: beta _{t-1}]	-0.0662	0.1297	0.0803
	[eqn: beta _{t-2}]	1.4276[eqn: {} ^{***}]	0.6399	-0.0689
	[eqn: gamma _{t-1}]	-0.1236	-0.5077	0.0293
	[eqn: gamma _{t-2}]	-1.5641[eqn: {} ^{***}]	0.0127	0.0735

period 2		lower regime	middle regime	upper regime
[eqn: TRIANGLE f _{t}]	[eqn: alpha]	0.0559	4.9412	-0.0615
	[eqn: beta _{t-1}]	-0.7215	0.1202	0.0805
	[eqn: beta _{t-2}]	0.4017	0.7907	-0.0220
	[eqn: beta _{t-3}]	-0.9686[eqn: {} ^{**}]	-0.3959	-0.1468
	[eqn: gamma _{t-1}]	0.4469	-0.6200	0.0578
	[eqn: gamma _{t-2}]	-0.5912	-0.1590	0.0112
	[eqn: gamma _{t-3}]	0.6064[eqn: {} ^{*}]	0.0587	0.1737
[eqn: TRIANGLE s _{t}]	[eqn: alpha]	0.3216	5.1327	0.3813[eqn: {} ^{**}]
	[eqn: beta _{t-1}]	-1.1178[eqn: {} ^{**}]	0.0382	0.0712
	[eqn: beta _{t-2}]	0.7237	0.6777	-0.0809
	[eqn: beta _{t-3}]	-0.7836[eqn: {} ^{*}]	-0.3645	-0.0107
	[eqn: gamma _{t-1}]	0.8449[eqn: {} ^{*}]	-0.5484	0.0393
	[eqn: gamma _{t-2}]	-0.8265[eqn: {} ^{*}]	-0.0379	0.0838
	[eqn: gamma _{t-3}]	0.4080	-0.0069	0.0415

period 3		lower regime	middle regime	upper regime
[eqn: TRIANGLE f _{t}]	[eqn: alpha]	-0.0109	0.4406	0.4841
	[eqn: beta _{t-1}]	0.0806	0.1063	-0.2973
	[eqn: beta _{t-2}]	-0.0244	0.0752	-0.2087
	[eqn: beta _{t-3}]	0.1460	0.1197	0.1535
	[eqn: gamma _{t-1}]	-0.1488	-0.0638	0.3387
	[eqn: gamma _{t-2}]	0.0162	-0.0815	0.0662
	[eqn: gamma _{t-3}]	-0.2072[eqn: {} ^{*}]	-0.1003	0.1408
[eqn: TRIANGLE s _{t}]	[eqn: alpha]	0.6960[eqn: {} ^{***}]	0.9986	1.2548[eqn: {} ^{***}]
	[eqn: beta _{t-1}]	-0.1579	0.0341	-0.3487
	[eqn: beta _{t-2}]	-0.3224[eqn: {} ^{*}]	0.1022	-0.5176
	[eqn: beta _{t-3}]	0.2438[eqn: {} ^{*}]	0.1298	0.2288
	[eqn: gamma _{t-1}]	0.0709	-0.0071	0.3642
	[eqn: gamma _{t-2}]	0.2502	-0.0878	0.3234
	[eqn: gamma _{t-3}]	-0.2509[eqn: {} ^{*}]	-0.1088	0.1158

period 4		lower regime	middle regime	upper regime
[eqn: TRIANGLE f _{t}]	[eqn: alpha]	0.1278	0.8938	0.6457[eqn: {} ^{*}]
	[eqn: beta _{t-1}]	-0.0528	0.0882	-0.4663[eqn: {} ^{*}]
	[eqn: beta _{t-2}]	-0.1437	0.0114	-0.3437
	[eqn: gamma _{t-1}]	-0.0062	-0.0431	0.5224
	[eqn: gamma _{t-2}]	0.1327	-0.0347	0.2638
[eqn: TRIANGLE s _{t}]	[eqn: alpha]	0.8936[eqn: {} ^{***}]	1.3790	1.3965[eqn: {} ^{***}]
	[eqn: beta _{t-1}]	-0.3422[eqn: {} ^{**}]	0.0113	-0.5145[eqn: {} ^{*}]
	[eqn: beta _{t-2}]	-0.4872[eqn: {} ^{***}]	0.0311	-0.6162[eqn: {} ^{*}]
	[eqn: gamma _{t-1}]	0.2506[eqn: {} ^{*}]	0.0214	0.5462[eqn: {} ^{*}]
	[eqn: gamma _{t-2}]	0.4142[eqn: {} ^{***}]	-0.0327	0.4811

period 5		lower regime	middle regime	upper regime
[eqn: TRIANGLE f _{t}]	[eqn: alpha]	-0.0109	0.4406	0.4841
	[eqn: beta _{t-1}]	0.0806	0.1063	-0.2973
	[eqn: beta _{t-2}]	-0.0244	0.0752	-0.2087
	[eqn: beta _{t-3}]	0.1460	0.1197	0.1535
	[eqn: gamma _{t-1}]	-0.1488	-0.0638	0.3387
	[eqn: gamma _{t-2}]	0.0162	-0.0815	0.0662
	[eqn: gamma _{t-3}]	-0.2072[eqn: {} ^{*}]	-0.1003	0.1408
[eqn: TRIANGLE s _{t}]	[eqn: alpha]	0.6960[eqn: {} ^{***}]	0.9986	1.2548[eqn: {} ^{***}]
	[eqn: beta _{t-1}]	-0.1579	0.0341	-0.3487
	[eqn: beta _{t-2}]	-0.3224[eqn: {} ^{*}]	0.1022	-0.5176
	[eqn: beta _{t-3}]	0.2438[eqn: {} ^{*}]	0.1298	0.2288
	[eqn: gamma _{t-1}]	0.0709	-0.0071	0.3642
	[eqn: gamma _{t-2}]	0.2502	-0.0878	0.3234
	[eqn: gamma _{t-3}]	-0.2509[eqn: {} ^{*}]	-0.1088	0.1158

Ⅴ. Conclusion

To analyze investor's arbitrage behaviors under different market conditions, I divided the entire sample period(January 2, 1986 ~ December 30, 2010) into 5 subsamples, which represent relatively different characters of the market(4 backwardations and 1 contango). For the division of the sample, ‘dynamic programming algorithm’ presented by Bai and Perron(2003) has been applied. Other preliminary tests including stationarity test, nonlinearity test, normaliaty test, and cointegration test have been used at the same time. The main methodology for this study was 3-regime bivariate TVECM with the basis as a threshold variable. 3 regimes are lower, middle, and upper ones. In a lower regime, basis is below the lower bound of threshold band. On the contrary, in an upper regime, basis is above the upper bound of threshold band. A middle regime is therefore between the two. So we can expect that there can be an interaction of WTI crude oil futures and spot prices in both lower and upper regimes, but not in a middle regime. The statistical significance of coefficient [eqn: alpha]can reflect this expectation. Coefficient [eqn: alpha]'s in both lower and upper regimes are proven to be statistically significant, but any coefficient [eqn: alpha] in any middle regime don't have statistical significance at all. This means that there aren't any arbitrage trading in middle regimes.

Coefficient [eqn: alpha] can be an indicator for the type of investor's arbitrage behavior. If a investor is in a contango market[eqn: LEFT ( f _{t} >s _{t} RIGHT )], he'll sell his futures contracts and buy spots. This moves spot prices upward and futures prices downward. That is, there can be a positive relationship between basis and spot prices, and we can expect [eqn: alpha ^{s}] to be statistically significant and positive. Let's look at the other side. If a investor is in a backwardation market[eqn: LEFT ( s _{t} >f _{t} RIGHT )], he'll sell his spots and buy futures contracts. This moves futures prices upward and spot prices downward. That is, there can be a negative relationship between basis and futures prices, and we can expect [eqn: alpha ^{f}] to be statistically significant and negative. These expectations can be checked through the estimation results of 3-regime bivariate TVECM. The above analysis results are relatively rich compared with those of ordinary VECM.

This study can be extended into prediction area afterwards. I plan to compare in-sample prediction of linear models with that of this nonlinear models for proving that nonlinear models have a higher predictability via methodology proposed by Diebold and Mariano(1995).

References

Andrews D.W.K. (1993), “Tests for parameter instability and structural change with unknown change point,” Econometrica 61, pp.821-856.

Andrews D.W.K. and Ploberger W. (1994), “Optimal tests when a nuisance parameter is present only under the alternative,” Econometrica 62, pp.1383-1414.

Bai J. and Perron P. (1998), “Estimating and testing linear models with multiple structural changes,” Econometrica 66, pp.47-78.

Brock W.A., Dechert W.D., and Scheinkman J.A. (1987), “A test for independence based on the correlation dimension,” Working paper, Department of Economics, University of Wisconsin, Madison.

Diebold F.X. and Mariano R.S. (1995), “Comparing predictive Accuracy,” Journal of Business and Economic Statistics 13(3), pp.253-263.

Hansen B.E. (1992), “Tests for parameter instability in regressions with I(1) processes,” Journal of Business and Economic Statistics 10, pp.321-335.

Kim H.S., Kang D.S., and Kim J.H. (2003), “The BDS statistic and residual test,” Stochastic Environmental Research and Risk Assessment No. 17, pp.104-115.

Kuan C.M. and Hornik K. (1995), “The generalized fluctuation test: A unifying view,”  Econometric Reviews 14, pp.135-161.

Lin, K. (1997), “The ABC's of BDS,” Journal of Computational Intelligence in Finance, 97(July/August), pp.23-26.

Ploberger W., Kramer W. (1992), “The CUSUM test with OLS residuals,” Econometrica 60(2), pp.271-285.

Tong, H. (1978) “On a threshold model,” In Pattern recognition and signal    processing, (ed. C. H. Chen). Sijthoff and Noordhoff, Amsterdam.