I'm trying to model and forecast a time series that is cyclic rather than seasonal (i.e. there are seasonal-like patterns, but not with a fixed period). This should be possible to do using an ARIMA model, as mentioned in Section 8.5 of Forecasting: principles and practice:

The value of $p$ is important if the data show cycles. To obtain cyclic forecasts, it is necessary to have $p\geq 2$ along with some additional conditions on the parameters. For an AR(2) model, cyclic behaviour occurs if $\phi^2_1+4\phi_2<0$.

What are these additional conditions on the parameters in the general ARIMA(p,d,q) case? I haven't been able to find them anywhere.

  • 1
    $\begingroup$ Have you looked into complex roots of the polynomial $\phi(B)$ at all? It seems like this may be what the quotation is referring to. $\endgroup$ – Jason Apr 10 '18 at 5:38

Some graphical intuition

In AR models, cyclic behavior come from complex conjugate roots to the characteristic polynomial. To first give intuition, I've plotted the impulse response functions below to two example AR(2) models.

  1. A persistent process with complex roots.
  2. A persistent process with real roots.

For $j=1\,\ldots, p$, Roots of the characteristic polynomial are $\frac{1}{\lambda_j}$ where $\lambda_1, \ldots, \lambda_p$ are eigenvalues of the $A$ matrix I define below . With a complex conjugate eigenvalues $\lambda = r e^{i \omega t}$ and $\bar{\lambda}= r e^{-i \omega t} $, the $r$ controls the damping (where $r \in [0, 1)$) and $\omega$ controls the frequency of the cosine wave.

Detailed AR(2) example

Let's assume we have the AR(2):

$$ y_t = \phi_1 y_{t-1} + \phi_2 y_{t-2} + \epsilon_t$$

You can write any AR(p) as a VAR(1). In this case, the VAR(1) representation is:

$$ \underbrace{\begin{bmatrix} y_t \\ y_{t-1} \end{bmatrix} }_{X_t} = \underbrace{\begin{bmatrix} \phi_1 & \phi_2 \\ 1 & 0 \end{bmatrix}}_{A} \underbrace{\begin{bmatrix} y_{t-1} \\ y_{t-2} \end{bmatrix}}_{X_{t-1}} + \underbrace{\begin{bmatrix}\epsilon_t \\ 0 \end{bmatrix}}_{U_t}$$ Matrix $A$ governs the dynamics of $X_t$ and hence $y_t$. The characteristic equation of matrix $A$ is: $$ \lambda^2 - \phi_1 \lambda - \phi_2 = 0$$ The eigenvalues of $A$ are: $$ \lambda_1 = \frac{\phi_1 + \sqrt{\phi_1^2 + 4\phi_2}}{2} \quad \quad \lambda_2 = \frac{\phi_1 - \sqrt{\phi_1^2 + 4\phi_2}}{2} $$ The eigenvectors of $A$ are: $$ \mathbf{v}_1 = \begin{bmatrix} \lambda_1 \\ 1 \end{bmatrix} \quad \quad \mathbf{v}_2 = \begin{bmatrix} \lambda_2 \\ 1 \end{bmatrix} $$

Note that $\operatorname{E}[X_{t+k} \mid X_t, X_{t-1}, \ldots] = A^kX_t$. Forming the eigenvalue decomposition and raising $A$ to the $k$th power. $$A^k = \begin{bmatrix} \lambda_1 & \lambda_2 \\ 1 & 1 \end{bmatrix} \begin{bmatrix} \lambda_1^k & 0 \\ 0 & \lambda_2^k \end{bmatrix} \begin{bmatrix} \frac{1}{\lambda_1 - \lambda_2} & \frac{-\lambda_2}{\lambda_1 - \lambda_2}\\ \frac{-1}{\lambda_1 - \lambda_2}& \frac{\lambda_1}{\lambda_1 - \lambda_2} \end{bmatrix} $$

A real eigenvalue $\lambda$ leads to decay as you raise $\lambda^k$. Eigenvalues with non-zero imaginary components leads to cyclic behavior.

Eigenvalues with imaginary component case: $\phi_1^2 + 4\phi_2 < 0$

In the AR(2) context, we have complex eigenvalues if $\phi_1^2 + 4\phi_2 < 0$. Since $A$ is real, they must come in pairs that are complex conjugates of each other.

Following Chapter 2 of Prado and West (2010), let $$ c_t = \frac{\lambda}{\lambda - \bar{\lambda}} y_t - \frac{\lambda \bar{\lambda}}{\lambda - \bar{\lambda}} y_{t-1}$$

You can show the forecast $\operatorname{E}[y_{t+k} \mid y_t, y_{t-1}, \ldots]$ is given by:

\begin{align*} \operatorname{E}[y_{t+k} \mid y_t, y_{t-1}, \ldots] &= c_t \lambda ^k + \bar{c}_t \bar{\lambda}^k \\ &= a_t r^k \cos(\omega k + \theta_t) \end{align*}

Speaking loosely, adding the complex conjugates cancels out their imaginary component leaving you with a single damped cosine wave in the space of real numbers. (Note we must have $0 \leq r < 1$ for stationarity.)

If you want to find $r$, $\omega$, $a_t$, $\theta_t$, start by using Euler's formula that $re^{i \theta} = r \cos \theta + r \sin \theta$, we can write:

$$ \lambda = re^{i\omega} \quad \bar{\lambda} = re^{-i\omega} \quad r = |\lambda| = \sqrt{ -\phi_2}$$ $$\omega = \operatorname{atan2}\left( \operatorname{imag} \lambda, \operatorname{real} \lambda \right) = \operatorname{atan2}\left( \frac{1}{2} \sqrt{-\phi_1^2-4\phi_2}, \frac{1}{2} \phi_1\right) $$

$$ a_t = 2 |c_t| \quad \quad \theta_t = \operatorname{atan2}\left( \operatorname{imag} c_t, \operatorname{real} c_t \right) $$


Note Confusing terminology warning! Relating the characteristic polynomial of A to characteristic polynomial of AR(p)

Another time-series trick is to use the lag operator to write the AR(p) as:

$$ \left(1 - \phi_1 L - \phi_2L^2 - \ldots - \phi_pL^p \right)y_t = \epsilon_t$$

Replace lag operator $L$ with some variable $z$ and people often refer to $1 - \phi_1 z - \ldots - \phi_pz^p$ as the characteristic polynomial of the AR(p) model. As this answer discusses, this is exactly the characteristic polynomial of $A$ where $z = \frac{1}{\lambda}$. The roots $z$ are the reciprocals of the eigenvalues. (Note: for the model to be stationary you want $|\lambda| < 1$, that is inside the unit cirlce, or equivalently $|z|>1$, that is outside the unit circle.)


Prado, Raquel and Mike West, Time Series: Modeling, Computation, and Inference, 2010

  • $\begingroup$ I am surprised I am the only up vote at the moment. Good answer! $\endgroup$ – Taylor Jan 29 '19 at 16:46
  • $\begingroup$ @Taylor It's an old, inactive question. :) $\endgroup$ – Matthew Gunn Jan 29 '19 at 16:47

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.