Fit Univariate or Multivariate Autoregressive Model

DESCRIPTION:
Fits a model of the form x(t) = a(1)x(t-1) + ... + a(p)x(t-p) + e(t).

USAGE:
ar(x, aic=T, order.max=<<see below>>, method="yule-walker")

REQUIRED ARGUMENTS:
x:
a univariate or multivariate time series, or a vector, or a matrix with columns representing univariate components of a multivariate time series. Missing values are allowed only at the beginning or end of series.

OPTIONAL ARGUMENTS:
aic:
logical flag: if TRUE, use the Akaike information criterion to choose the best order not greater than order.max. If FALSE, order.max will be the order of the fitted model.
order.max:
the maximum order of autoregression to fit to the time series. The default value is 10 * log10(n.used/ncol(x)).
method:
either the string "yule-walker" or "burg" depending on whether Yule-Walker equations or Burg's algorithm is to be used in estimating the autoregression coefficients. Only the first character is necessary.

VALUE:
a list with the following components:
order:
the order of the autoregression fitted. If aic=TRUE, then this is the order less than or equal to order.max which minimizes the AIC, otherwise, it is order.max.
ar:
the autoregressive coefficients. These are in an order by "nser" by "nser" array. where "nser" is the number of univariate components of x. If order is 0, ar will have dimensions 1 by "nser" by "nser" and will be filled with zeros. The first level of the first dimension corresponds to one observation back in time, the second level corresponds to two observations back, etc. The second dimension corresponds to the predicted series, and the third corresponds to the predicting series.
var.pred:
the estimated prediction variance ("nser" by "nser" matrix) of the process with AR coefficients ar.
aic:
vector of the values of Akaike's information criterion for orders 0 through order.max. These have the minimum value subtracted from all of them so the minimum is always zero.
n.used:
the length of the nonmissing part of the input time series.
order.max:
the input or default value of order.max.
partialacf:
vector of the partial autocorrelation coefficients, useful for choosing the proper order of the autoregression. This is computed only if method="yule-walker".
resid:
a time series representing the estimate of e[t] for the model above. The estimates contained in ar are used in the forward direction on the series with mean(s) removed.
method:
the string "yule-walker" or "burg".
series:
the name used for the input x. It includes transformations.

DETAILS:
If method="yule-walker", the autocovariance matrices of the time series are estimated and Whittle's recursion (a multivariate extension of the Levinson-Durbin method) is used to estimate the autoregressive coefficients. The coefficients can also be estimated by using Burg's algorithm, based on the principle of maximum entropy. The output may be used in spec.ar to estimate the spectrum of the process, or acf.plot to produce a plot of the partial autocorrelation function.

The estimation is performed using the sample mean of each univariate series as the estimate of the mean. Remember that the coefficients in ar are for the series with the mean(s) removed.


REFERENCES:
The chapter "Analyzing Time Series" of the S-PLUS Guide to Statistical and Mathematical Analysis.

SEE ALSO:
acf.plot , ar.burg , ar.yw , arima.mle , spec.ar .

EXAMPLES:
a <- ar(log(lynx))
acf.plot(a, conf=T)             # Take a look at the partial correlations
tsplot(a$aic)
# Fit an AR(11) to this time series
lynx.fit <- ar(log(lynx), aic=F, order=11, method="b")

# function to predict using an ar model: # ahead gives the number of predictions to make

pred.ar <- function(series, ar.est, ahead = 1) { order <- ar.est$order series <- as.matrix(series) pred.out <- array(NA, dim = c(order + ahead, ncol(series)), dimnames = list(NULL, dimnames(series)[[2]])) mean.ser <- apply(series, 2, mean) ser.cent <- sweep(series, 2, mean.ser) pred.out[seq(order), ] <- ser.cent[rev(nrow(series) - seq(order) + 1), ] for(i in (order + 1):nrow(pred.out)) { pred.out[i, ] <- apply(aperm(ar.est$ar, c(1, 3, 2)) * as.vector(pred.out[i - seq(order), ]), 3, sum) } sweep(pred.out[ - seq(order), , drop = F], 2, mean.ser, "+") }

sun.ar <- ar(sunspots) tsplot(sun.ar$resid, main="Residuals from AR(27) Fit of Sunspots") # the structure in the residuals means the model is inadequate