MAX_ARMA Function
Exact maximum likelihood estimation of the parameters in a univariate ARMA (autoregressive, moving average) time series model.
Usage
result = MAX_ARMA (w, p, q)
Input Parameters
w—Array of length number of observations containing the time series.
p—Non-negative number of autoregressive parameters.
q—Non-negative number of moving average parameters.
Returned Value
result—Array of length 1 + p + q with the estimated constant, AR and MA parameters.
Input Keywords
Double—If present and nonzero, double precision is used.
N_predict—Maximum lead time for forecasts. N_predict must be greater than zero. Forecast and N_predict must be used together.
Confidence—Value in the exclusive interval (0, 100) used to specify the confidence level of the forecasts. Typical choices for Confidence are 90.0, 95.0, and 99.0. Default: Confidence = 95.0
Backward_origin—Maximum backward origin. Backward_origin must be greater than or equal to zero and less than or equal to N_ELEMENTS(w) – MAX(p,q).
Forecasts at origins N_ELEMENTS(w)Backward_origin through N_ELEMENTS(w) are generated. Default: Backward_origin = 0
Init_ar—If specified, Init_ar is an array of length p containing preliminary estimates of the autoregressive parameters. If p = 0, then the corresponding arguments are ignored.
Init_ma—If specified, Init_ma is an array of length q containing preliminary estimates of the moving average parameters. If q = 0, then the corresponding arguments are ignored.
Maxit—Maximum number of estimation iterations. Default: Maxit = 300
Output Keywords
Avar—Estimate of innovation variance.
Log_likelihood—Value of –2 * (ln(likelihood)) for the fitted model.
Forecast—An array of length N_predict by (Backward_origin + 3) containing the forecasts up to N_predict steps ahead and the information necessary to obtain confidence intervals. Forecast and N_predict must be used together.
W_mean—(Input/Output) Contains the estimate on input and an update of the mean on output. Default: Time series w is centered about its sample mean.
Discussion
MAX_ARMA is derived from the maximum likelihood estimation algorithm described by Akaike, Kitagawa, Arahata and Tada (1979), and the XSARMA routine published in the TIMSAC-78 Library.
Using the notation developed in the Time Domain Methodology at the beginning of this chapter, the stationary time series Wt with mean μ can be represented by the nonseasonal autoregressive moving average (ARMA) model by the following relationship:
φ  (B)(Wtμ ) = θ (B)at
where:
B is the backward shift operator defined by:
and:
MAX_ARMA estimates the AR coefficients and the MA coefficients using maximum likelihood estimation.
MAX_ARMA checks the initial estimates for both the autoregressive and moving average coefficients to ensure that they represent a stationary and invertible series respectively. If:
are the initial estimates for a stationary series then all (complex) roots of the following polynomial will fall outside the unit circle:
If:
are initial estimates for an invertible series then all (complex) roots of the polynomial:
will fall outside the unit circle.
Initial values for the AR and MA coefficients can be supplied by keywords Init_ar and Init_ma. Otherwise, estimates are computed internally by the method of moments. MAX_ARMA computes the roots of the associated polynomials. If the AR estimates represent a non-stationary series, MAX_ARMA issues a warning message and replaces Init_ar with initial values that are stationary. If the MA estimates represent a non-invertible series, MAX_ARMA issues a terminal error, and new initial values have to be sought.
MAX_ARMA also validates the final estimates of the AR coefficients to ensure that they too represent a stationary series. This is done to guard against the possibility that the internal log-likelihood optimizer converged to a non-stationary solution. If non-stationary estimates are encountered, MAX_ARMA issues a fatal error message.
For model selection, the ARMA model with the minimum value for AIC might be preferred:
MAX_ARMA can also handle white noise processes, i.e. ARMA(0,0) Processes.
Example 1
Consider the Wolfer Sunspot data (Anderson 1971, p. 660) consisting of the number of sunspots observed each year from 1770 through 1869. In this example, MAX_ARMA is used to fit the following ARMA(2,1) model:
with the sample mean of the time series {wt}.
For these data, MAX_ARMA calculated the following model:
Defining the overall constant , we obtain the following equivalent representations:
and:
 
n_obs = 100  ; Total number of observations
p = 2        ; Number of autoregressive parameters
q = 1        ; Number of moving average parameters
w = FLTARR(n_obs)  ; The array containing the time series
z = STATDATA(2)    ; Get Wolfer sunspot data
FOR i=0L, n_obs-1 DO w(i) = z(21+i,1)
 
parameters = MAX_ARMA(w, p, q, $
                      Maxit=12000, $
                      Avar=avar, $
                      Log_likelihood=log_likeli)
 
PRINT, parameters(1), parameters(2), Format= $
  "('AR estimates are ',  F8.4, ' and ', F8.4, '.')"
PRINT, parameters(3), $
  Format="('MA estimate is ', F8.4, '.')"
PRINT, parameters(0), $
  Format="('Constant estimate is ', F8.4, '.')"
PRINT, log_likeli, $
  Format="('-2*ln(Maximum Log Likelihood) = ', F8.4, '.')"
PRINT, avar, Format="('White noise variance = ', F8.4, '.')"
Output
AR estimates are   1.2273 and  -0.5626.
MA estimate is  -0.3808.
Constant estimate is  15.7508.
-2*ln(Maximum Log Likelihood) = 539.5843.
White noise variance = 214.5020.
Example 2
This is the same as example 1, but now initial values for the AR and MA parameters are explicitly given.
n_obs = 100  ; Total number of observations
p = 2        ; Number of autoregressive parameters
q = 1        ; Number of moving average parameters
w = FLTARR(n_obs)  ; The array containing the time series
z = STATDATA(2)    ; Get Wolfer sunspot data
FOR i=0L, n_obs-1 DO w(i) = z(21+i,1)
 
; Initial values for AR
init_ar = [1.244e0, -0.575e0]
; Initial values for MA
init_ma = [-0.1241e0]
 
parameters = MAX_ARMA(w, p, q, $
                      Maxit=12000, $
                      Avar=avar, $
                      Log_likelihood=log_likeli, $
                      Init_ar=init_ar, $
                      Init_ma=init_ma)
 
PRINT, parameters(1), parameters(2), Format= $
  "('AR estimates are ',  F8.4, ' and ', F8.4, '.')"
PRINT, parameters(3), $
  Format="('MA estimate is ', F8.4, '.')"
PRINT, parameters(0), $
  Format="('Constant estimate is ', F8.4, '.')"
PRINT, log_likeli, $
  Format="('-2*ln(Maximum Log Likelihood) = ', F8.4, '.')"
PRINT, avar, Format="('White noise variance = ', F8.4, '.')"
Output
AR estimates are   1.2273 and  -0.5623.
MA estimate is  -0.3804.
Constant estimate is  15.7373.
-2*ln(Maximum Log Likelihood) = 539.5843.
White noise variance = 214.5052.