ARIMA model: definition

• Autoregressive integrated moving-average models (ARIMA): Let \(d \in \mathbb{N}\), then series \(\{X_t\}\) is an ARIMA\((p,d,q)\) process if \[ Y_t = (1 - B)^d X_t \] is a causal ARMA\((p,q)\) process.

• Difference equation (DE) for an ARIMA\((p,d,q)\) process \[ \phi^*(B) X_t = \phi(B) (1-B)^d X_t = \theta(B) Z_t, \quad \{Z_t\} \sim \textrm{WN}(0, \sigma^2) \]
  • \(\phi(z)\): polynomial of degree \(p\), and \(\phi(z) \neq 0\) for \(|z| \leq 1\)
  • \(\theta(z)\): polynomial of degree \(q\)
  • \(\phi^*(z) = \phi(z) (1-z)^d\): has a zero of order \(d\) at \(z=1\)

• An ARIMA process with \(d > 0\) is NOT stationary!

ARIMA mean is not dertermined by the DE

• \(\{X_t\}\) is an ARIMA\((p,d,q)\) process. We can add an arbitrary polynomial trend of degree \(d-1\) \[ W_t = X_t + A_0 + A_1 t + \cdots + A_{d-1} t^{d-1} \] with \(A_0, \ldots, A_{d-1}\) being any random variables, and \(\{W_t\}\) still satisfies the same ARIMA\((p,d,q)\) difference equation

• In other words, the ARIMA DE determines the second-order properties of \(\{(1-B)^d X_t\}\) but not those of \(\{X_t\}\)

  • For parameter estimation: \(\boldsymbol\phi\), \(\boldsymbol\theta\), and \(\sigma^2\) are estimated based on \(\{(1-B)^d X_t\}\) rather than \(\{X_t\}\)
  • For forecast, we need additional assumptions

Fit data using ARIMA processes

• Whether to fit a finite time series using
  • non-stationary models (such as ARIMA), or
  • directly using stationary models (such as ARMA)?
• If the fitted stationary ARMA model’s \(\phi(\cdot)\) have zeros very close to unit circles, then fitting an ARIMA model is better
  • Parameter estimation is stable
  • The differenced series may only need a low-order ARMA
• Limitation of ARIMA: only permits data to be nonstationary in a very special way
  • Non-stationary: can have zeros anywhere on the unit circle \(|z|=1\)
  • ARIMA model: only has a zero of multiplicity \(d\) at the point \(z=1\)

Transformation and Identification Techniques

Unit Root Test

library(tseries);
## Note: the lag k here is the AR order p
adf.test(dowj, k = 2);

## 
##  Augmented Dickey-Fuller Test
## 
## data:  dowj
## Dickey-Fuller = -1.3788, Lag order = 2, p-value = 0.8295
## alternative hypothesis: stationary

Forecast ARIMA models

Seasonal ARIMA (SARIMA) Model: definition

• Suppose \(d, D\) are non-negative integers. \(\{X_t\}\) is a seasonal ARIMA\((p, d, q)\) \(\times\) \((P, D, Q)_s\) process with period \(s\) if the differenced series \[ Y_t = (1-B)^d (1-B^s)^D X_t \] is a causal ARMA process defined by \[ \phi(B) \Phi(B^s) Y_t = \theta(B) \Theta(B^s) Z_t, \quad \{Z_t\} \sim \textrm{WN}(0, \sigma^2) \] where \[\begin{align*} \phi(z) & = 1 - \phi_1 z - \cdots - \phi_p z^p, \quad \Phi(z) = 1 - \Phi_1 z - \cdots - \Phi_P z^P \\ \theta(z) & = 1 + \theta z + \cdots + \theta_q z^q, \quad \Theta(z) = 1 + \Theta z + \cdots + \Theta_Q z^Q \end{align*}\]

• \(\{Y_t\}\) is causal if and only if neither \(\phi(z)\) or \(\Phi(z)\) has zeros inside the unit circle

• Usually, \(s=12\) for monthly data

Special case: seasonal ARMA (SARMA)

• Between-year model: for monthly data, each one of the 12 time series is generated by the same ARMA\((P, Q)\) model \[ \Phi(B^{12}) Y_t = \Theta(B^{12}) U_t, \quad \{U_{j+12t}, t \in \mathbb{Z}\} \sim \textrm{WN}(0, \sigma^2_U) \]

• SARMA\((P, Q)\) with period \(s\): in the above between-year model, the period 12 can be changed to any positive integer \(s\)
  • If \(\{U_{t}, t \in \mathbb{Z}\}\sim \textrm{WN}(0, \sigma^2_U)\), then the ACVF \(\gamma(h) = 0\) unless \(h\) divides \(s\) evenly. But this may not be ideal for real life application! E.g., this Feb is correlated with last Feb, but not this Jan.
• General SARMA\((p, q)\) \(\times\) \((P, Q)\) with period \(s\): incorporate dependency between the 12 series by letting \(\{U_t\}\) be ARMA: \[ \Phi(B^{s}) Y_t = \Theta(B^{s}) U_t, \quad \phi(B) U_t = \theta(B)Z_t, \quad \{Z_t\} \sim \textrm{WN}(0, \sigma^2) \]
  • Equivalent DE for the general SARMA: \[ \phi(B) \Phi(B^s) Y_t = \theta(B) \Theta(B^s) Z_t \]

Fit a SARIMA Model

• Period \(s\) is known

1. Find \(d\) and \(D\) to make the difference series \(\{Y_t\}\) to look stationary

2. Examine the sample ACF and sample PACF of \(\{Y_t\}\) at lags being multiples of \(s\), to find orders \(P, Q\)

3. Use \(\hat{\rho}(1), \ldots, \hat{\rho}(s-1)\) to find orders \(p, q\)

4. Use AICC to decide among competing order choices

5. Given orders \((p, d, q, P, D, Q)\), estimate MLE of parameters \((\phi, \theta, \Phi, \Theta, \sigma^2)\)

Regression with ARMA errors: OLS estimation

• Linear model with ARMA errors \(\mathbf{W} = (W_1, \ldots, W_n)'\): \[ Y_t = \mathbf{x}_t' \boldsymbol\beta + W_t, \quad t = 1, \ldots, n, \quad \{W_t\} \sim \text{causal ARMA}(p,q) \]
  • Note: each row is indexed by a different time \(t\)!
  • Error covariance \(\boldsymbol\Gamma_n = E(\mathbf{W} \mathbf{W}')\)
• Ordinary least squares (OLS) estimator \[ \hat{\boldsymbol\beta}_{\text{OLS}} = (\mathbf{X}'\mathbf{X})^{-1}\mathbf{X}'\mathbf{Y} \]
  • Estimated by minimizing \((\mathbf{Y} - \mathbf{X}\boldsymbol\beta)' (\mathbf{Y} - \mathbf{X}\boldsymbol\beta)\)
  • OLS is unbiased, even when errors are dependent!

Regression with ARMA errors: GLS estimation

• Generalized least squares (GLS) estimator \[ \hat{\boldsymbol\beta}_{\text{GLS}} = (\mathbf{X}'\boldsymbol\Gamma_n^{-1}\mathbf{X})^{-1}\mathbf{X}'\boldsymbol\Gamma_n^{-1}\mathbf{Y} \]

  • Estimated by minimizing the weighted sum of squares \[ (\mathbf{Y} - \mathbf{X}\boldsymbol\beta)' \boldsymbol\Gamma_n^{-1} (\mathbf{Y} - \mathbf{X}\boldsymbol\beta) \]

  • Covariance \[ \textrm{Cov}\left( \hat{\boldsymbol\beta}_{\text{GLS}} \right) = (\mathbf{X}'\boldsymbol\Gamma_n^{-1}\mathbf{X})^{-1} \]

  • GLS is the best linear unbiased estimator, i.e., for any vector \(\mathbf{c}\) and any unbiased estimator \(\hat{\boldsymbol\beta}\), we have \[ \textrm{Var}(\mathbf{c}' \hat{\boldsymbol\beta}_{\text{GLS}}) \leq \textrm{Var}(\mathbf{c}' \hat{\boldsymbol\beta}) \]

When \(\{W_t\}\) is an AR\((p)\) process

• We can apply \(\phi(B)\) to each side of the regression equation and get uncorrelated, zero-mean, constant-variance errors \[ \phi(B) \mathbf{Y} = \phi(B) \mathbf{X} \boldsymbol\beta + \phi(B) \mathbf{W} = \phi(B) \mathbf{X} \boldsymbol\beta + \mathbf{Z} \]

• Under the transformed target variable \[ Y_t^* = \phi(B) Y_t, \quad t = p+1, \ldots, n \] and the transformed design matrix \[ X_{t, j}^* = \phi(B) X_{t, j}, \quad t = p+1, \ldots, n, \quad j = 1, \ldots, k \] the OLS estimator is the best linear unbiased estimator

• Note: after the transformation, the regression sample size reduces to \(n-p\)

Regression with ARMA errors: MLE

• MLE of \(\boldsymbol\beta, \boldsymbol\phi, \boldsymbol\theta, \sigma^2\) can be estimated by maximizing the Gaussian likelihood with error covariance \(\boldsymbol\Gamma_n\)

• An iterative scheme

  1. Compute \(\hat{\boldsymbol\beta}_{\text{OLS}}\) and regression residuals \(\mathbf{Y} - \mathbf{X}\hat{\boldsymbol\beta}_{\text{OLS}}\)

  2. Based on the estimated residuals, compute MLE of the ARMA\((p, q)\) parameters

  3. Based on the fitted ARMA model, compute \(\hat{\boldsymbol\beta}_{\text{GLS}}\)

  4. Compute regression residuals \(\mathbf{Y} - \mathbf{X}\hat{\boldsymbol\beta}_{\text{GLS}}\), and return to Step 2 until estimators stabilize

• Asymptotic properties of MLE: If \(\{W_t\}\) is a causal and invertible ARMA, then
  • MLEs are asymptotically normal
  • Estimated regression coefficients are asymptotically independent of estimated ARMA parameters

References

• Brockwell, Peter J. and Davis, Richard A. (2016), Introduction to Time Series and Forecasting, Third Edition. New York: Springer

• Weigt, George (2018) ITSM-R Reference Manual http://www.eigenmath.org/itsmr-refman.pdf

• R package tseries https://cran.r-project.org/web/packages/tseries/index.html