Canadian Lynx Series

A first time series model of the Canadian lynx data was fitted by P.A.P. Moran [a13] in 1953. He observed that the cycle is very asymmetrical with a sharp and large peak and a relatively smooth and small trough. The log transformation gives a series which appears to vary symmetrically about the mean. As the actual population of lynx is not exactly proportional to the number caught, a better representation would perhaps be obtained by incorporating an additional "error of observation" in the model, thereby resulting in a more complicated model. The log transformation substantially reduces the effect of ignoring this error of observation; therefore, after Moran, nearly all the time series analysis of the lynx data in the literature have used the log-transformed data. Let

$$ X _ {t} = $$

$$ = { \mathop{\rm log} } _ {10 } ( \textrm{ number recorded as trapped in year 1820 } + t ) $$

( $ t = 1 \dots 114 $). Because of the apparently slow damping in the amplitude of the sample correlogram, Moran discarded the idea of a sinusoidal-shape model and proposed an $ { \mathop{\rm AR} } ( 2 ) $- model.

In 1977, M.J. Campbell and A.M. Walker [a2] believed that an appropriate model of the lynx data should be, in some sense, "between" a pure harmonic model and a pure auto-regression. Subsequently this led them to combining a harmonic component with an $ { \mathop{\rm AR} } ( 2 ) $- model. Two models with frequencies $ 9.5 $, $ { \mathop{\rm CW} } ( 9.5 ) $, and $ 9.63 $, $ { \mathop{\rm CW} } ( 9.63 ) $, were recommended. At about the same time, an $ { \mathop{\rm AR} } ( 11 ) $- model based on the Akaike information criterion was fitted [a20]. In the discussion of the above two papers, D.R. Cox [a4] suggested a polynomial model.

In 1979, R.J. Bhansali [a1] used a mixed spectrum analysis to analyze the lynx data.

Using the Canadian lynx data as a case study, H. Tong and K.S. Lim [a24] fitted a class of non-linear models called the self-exciting threshold auto-regressive model (SETAR model) to the log-transformed lynx data. They demonstrated that this model has interesting features in non-linear oscillations, such as jump resonance, amplitude-frequency dependency, limit cycles, subharmonics, and higher harmonics. Later, in 1981, it was discovered that the self-exciting threshold auto-regressive model also generates chaos [a11]. In the discussion of [a24], T. Subba Rao and M.M. Gabr [a16] proposed a subset bilinear model, $ { \mathop{\rm SBL} } ( 11 ) $, to the first one hundred log-transformed lynx data and an $ { \mathop{\rm SBL} } ( 9 ) $ to the first one hundred original lynx data.

In 1981, they used the maximum-likelihood estimation coupled with the Akaike information criterion to fit a subset bilinear model $ { \mathop{\rm SBL} } ( 12 ) $ to the first one hundred log-transformed lynx data [a5]. Their model was able to produce small values of the noise variance and the mean-squared errors of the one-step-ahead predictions, but it failed to detect the inherited behaviour of the data.

V. Haggan and T. Ozaki [a8] fitted an exponential auto-regressive model, $ { \mathop{\rm EXPAR} } ( 11 ) $, another class of non-linear models, to the mean-deleted log-transformed lynx data in 1981. Ozaki [a15] felt that the almost symmetric series generated by this model was unsatisfactory. Subsequently, in 1982, he fitted two more exponential auto-regressive models to the full set of log-transformed lynx data with mean deleted. One of them, $ { \mathop{\rm EXPAR} } ( 2 ) $, could reproduce the asymmetric limit cycle structure of the lynx data; the other, $ { \mathop{\rm EXPAR} } ( 9 ) $, with smaller variance of fitted residuals, was believed to be more appropriate for forecasting.

D.F. Nicholls and D.G. Quinn [a14], in 1982, fitted another new class of time series models, called a random coefficient auto-regressive model (RCA model) to the first one hundred log-transformed lynx data, using a maximum-likelihood method or the conditional least-squares method.

In 1984, Haggan, S.M. Heravi and M.B. Priestley [a7] fitted a state-dependent model (SDM) of $ { \mathop{\rm AR} } ( 2 ) $ to the log-transformed lynx data.

Using the revised computer program in [a21], a $ { \mathop{\rm SETAR} } ( 2;5,2 ) $ was fitted to the first one hundred log-transformed lynx data, and a $ { \mathop{\rm SETAR} } ( 2;7,2 ) $ was fitted to the full set of log-transformed data as follows:

$$ X _ {t} = 0.546 + 1.032X _ {t - 1 } - 0.173X _ {t - 2 } + $$

$$ + 0.171X _ {t - 3 } - 0.431X _ {t - 4 } + 0.332X _ {t - 5 } - $$

$$ - 0.284X _ {t - 6 } + 0.210X _ {t - 7 } + \epsilon _ {t} ^ {( 1 ) } $$

if $ X _ {t - 2 } \leq 3.116 $, and

$$ X _ {t} = 2.632 + 1.492X _ {t - 1 } - 1.324X _ {t - 2 } + \epsilon _ {t} ^ {( 2 ) } $$

if $ X _ {t - 2 } > 3.116 $. Here, $ { \mathop{\rm var} } ( \epsilon _ {t} ^ {( 1 ) } ) = 0.0258 $, $ { \mathop{\rm var} } ( \epsilon _ {t} ^ {( 2 ) } ) = 0.0505 $( pooled variance equals $ 0.0360 $).

This model was able to describe the biological features of the Canadian lynx data such as:

1) its cyclical behaviour of about 9–10 years per cycle;

2) the rise periods exceed the descent periods in the cycles;

3) the delay parameter of $ 2 $ in $ X _ {t - 2 } \leq 3.116 $ is associated with the biological cycle that a Canadian lynx is fully grown in the autumn of its second year and births of kittens (1–4 per litter) take place about 63 days after breeding in March–April;

4) the threshold estimate, $ 3.116 $, lies in the vicinity of the anti-mode of the histogram of the lynx data, which implies that there is insufficient information in the data to model more precisely the functional form of the dynamics over the state space near the sample mean.

A comparative study of some of the above models was carried out in [a12]. The $ { \mathop{\rm SETAR} } $ models were ranked to be the best among the models considered. For an extensive discussion of the lynx data, see [a22].

B.Y. Thanoon (1988) [a18] fitted the two subset $ { \mathop{\rm SETAR} } $( $ { \mathop{\rm SSETAR} } $) models which produced limit cycles with two subcycles giving an average period of $ 9.5 $ years. He commented that the $ { \mathop{\rm SSETAR} } $ detected the inherited behaviour of the data better than the full $ { \mathop{\rm SETAR} } $ model in terms of the auto-covariance function.

In 1989, R.S. Tsay [a25] fitted a two-thresholds $ { \mathop{\rm SETAR} } ( 3;1,7,2 ) $ when proposing a new procedure for testing and building $ { \mathop{\rm TAR} } $ models. Around the same time, techniques from multivariate analysis were applied, [a23]; namely, the principal coordinate analysis and dendograms to twelve time series models reported in the literature.

In 1991, G.H. Yu and Y.C. Lin [a29] suggested a subset auto-regressive model, $ { \mathop{\rm SAR} } ( 1,3,9,12 ) $, to the log-transformed lynx data when proposing a method for selecting a best $ { \mathop{\rm SAR} } $ model automatically.

Applying the cross-validatory approach, in 1992, B. Cheng and Tong [a3] found the embedding dimension of the lynx data to be $ 3 $.

Using the lynx data as an example, in 1993, J. Geweke and N. Terui [a6] proposed a Bayesian approach for deriving the exact posterior distributions of the delay and threshold parameters.

In 1994, T. Teräsvirta [a17] fitted a logistic smooth transition auto-regressive model ( $ { \mathop{\rm LSTAR} } ( 11 ) $), which had a limit cycle of 77 years with eight subcycles of lengths 9 and 10 years. The lynx data has been used for the non-parametric identification of non-linear time series in selecting significant lags [a19]. The lynx data for 1821–1924 has been used [a27] to estimate $ f _ {m} ( \cdot ) $ and $ \lambda _ {m} ( \cdot ) $, the $ m $- step Lyapunov-like index, where $ \lambda _ {m} ( x ) = { {df _ {m} } / {dx } } $, and the last ten data to check the predicted values. The data was also used [a28] for subset selection in non-parametric stochastic regression. The subset of lags $ 1 $, $ 3 $, and $ 6 $ are selected from the original lynx data.

In 1995, C. Kooperberg, C.J. Stone and Y.K. Truong [a9] used their automatic procedure to estimate the mixed spectral distribution of the log-transformed lynx data and found the lynx cycle to be $ 9.5 $ years.

In 1996, D. Lai [a10] used a BDS statistic to test the residuals from the models of Moran, $ { \mathop{\rm SETAR} } ( 2;2,2 ) $, $ { \mathop{\rm SETAR} } ( 2;7,2 ) $, $ { \mathop{\rm SETAR} } ( 3;1,7,2 ) $, $ { \mathop{\rm EXPAR} } ( 11 ) $, $ { \mathop{\rm EXPAR} } ( 2 ) $, and Cox's polynomial model. He concluded that Tong's $ { \mathop{\rm SETAR} } ( 2;7,2 ) $ was found to be the best. C. Wong and R. Kohn [a26] used a Bayesian approach for estimating non-parametrically an additive auto-regressive model for the lynx data.

References[edit]