STAT ______________ 758, Fall 2008, Statistical Lab 3 _________ SIMULATING ARMA MODELS
Goals: Learn how to use to simulate AR, MA, and ARMA processes and perform essential ARMA analysis. Learn how to estimate ARMA and its essential characteristics.
Assignments: 1. Simulate causal and invertible ARMA(2,2) model of your choice with a fixed Normal iid noise (simulate and store the noise), of length 50. 2. Estimate the four model parameters, acf, pacf, acvf, compare them with the respective theoretical values/functions. Compare the estimated residuals with the true ones, find the correlation between the two sets. Discuss the quality of the estimation. 3. Repeat the analysis of 1-2 increasing the length of the time series. Find such a length N0 that the four parameters are estimated with less than 5% relative error. 4. Estimate the model parameters assuming that the true order of the model is ARMA(1,1). Compare the quality of estimation, using the variance of the estimated residuals and comparing the estimated and true residuals. Discuss the results. Report: All assignments require a printed report. Your goal is to support your conclusions using theoretical and/or practical arguments using (do not forget about nice graphs.) Due date: Thursday, October 23, 2008 in class. commands used in these handouts are collected in R-files Lab3.R, which is available on the course Web site http://unr.edu/homepage/zal/STAT758.htm
Page 1
STAT 758 Fall 2008 Lab 3
1. Introduction
ARMA models play an important role in time series analysis, mainly because of their flexibility in fitting a given autocorrelation function. It can be shown that one can find an ARMA model that will exactly match an arbitrary number of the first values of any acf. In particular, if the acf of a time series has only a finite number of non-zero values, it can be fit exactly with a MA(q) process, and if the pacf of a time series has only a finite number of non-zero values, it can be fit exactly by an AR(p) process. The first step in studying the -techniques of ARMA modeling is to learn how to simulate ARMA time series, how to find acf and pacf of a given ARMA model, and how to find an MA representation of a given ARMA model. These topics will be considered in Sections 2-7. The next step is to estimate ARMA parameters from the observed time series; this will be discussed in Sections 8-9.
2. Simulating MA(q)
In , there are two ways of modeling MA(q) processes. The first method is based on the command X<-filter(Z,filter,sides=1) It uses the following model representation:
X t = 1Z t + 2 Zt -1 + ... + q+1Z t -q .
Inputs: Z is the time series of innovations; it should be a white noise in order to obtain a proper MA model. Caution: does not check whether Z is a white noise and will work with any time series! filter is a vector with MA coefficients. For example, filter=c(1,2) corresponds to a MA(1) process
X t = Z t + 2Z t -1 .
sides = 1 for MA models. Output is an object of the time series class. No check is made whether an MA model is invertible. The second method of MA simulation will be considered in Section 4 below.
3. Simulating AR(p)
In , there are two ways of modeling AR(p) processes. The first method is based on the command X<-filter(Z,filter,method=`r')
Page 2
STAT 758 Fall 2008 Lab 3
It uses the following model representation:
X t = Z t + 1 X t -1 + ... + p X p .
Inputs: Z is the time series of innovations; it should be a white noise in order to obtain a proper AR model. Caution: does not check whether Z is a white noise and will work with any time series! filter is a vector with regression coefficients. filter=c(0.9) corresponds to an AR(1) process For example,
X t = Z t + 0.9 X t -1 .
method = `r' for AR models. Output is an object of the time series class. No check is made whether an AR model is causal, the output may diverge if it is not. The second method of AR(p) simulation will be considered in Section 4 below.
4. Simulating ARMA(p,q)
The most general way of simulating ARMA models is implemented in the command X<-arima.sim(n, model, ...) It uses the following model representation:
X t = 1 X t -1 + ... + p X p + Z t + 1Z t -1 + ... + q Z t -q
Basic inputs: n is the length of the output time series model is a list with model parameters. Notice that here we explicitly assume a unit coefficient by Zt. For example, X<-arima.sim(n=100,list(ar=c(0.9),ma=c(0.2))) produces an ARMA(1,1) model
X t = 0.9 X t -1 + Z t + 0.2 Z t -1 .
Output is an object of the time series class. No check is made whether an MA part of the model is invertible. The AR part is checked for stationarity and causality, an error message will appear if an attempt is made to simulate a noncausal model. The command arima.sim allows several other inputs that can be very useful in applied modeling.
Page 3
STAT 758 Fall 2008 Lab 3
Additional sd inputs: is the standard deviation of the normal white noise (sd=1 by default) rand.gen is the name of a function used to generate innovations (rand.gen = rnorm by default) innov is a time series of innovations (useful to model several models with the same white noise) arima.sim can be used to generate MA models (by omitting ar part), AR models (by omitting ma part), or white noises (by submitting an empty list of parameters).
5. Theoretical ACF and PACF computation
allows to compute the theoretical values of the autocorrelation function for stationary ARMA processes (do not confuse this with estimation of the acf from a given time series). This is done with the command (mind the capital letters) Acf<-ARMAacf(ar=...,ma=...,lag.max=r) Example: Acf<-ARMAacf(ar=c(.2),ma=c(.3,.2),lag.max=4) computes the first five (for lags 0,1,2,3,4) values of the acf of the process
X t = 0.2 X t -1 + Z t + 0.3Z t -1 + 0.2 Z t - 2 ,
which is (0)=1. 0, (1)=0.498, (2)=0.248, (3)= 0.050, and (4)= 0.010. You also can compute the theoretical partial autocorrelation function using the command Acf<-ARMAacf(ar=...,ma=...,lag.max=r,pacf=TRUE)
6. MA() representation of ARMA model
A moving average representation of an ARMA model,
X t = i Z t - k ,
k =0
is calculated by the command ARMAtoMA(ar=..., ma=..., lag.max=r) Notice that the calculated coefficients do not include the first 1, for example: P<-ARMAtoMA(ar=0.9, lag.max=2)
Page 4
STAT 758 Fall 2008 Lab 3
results in P=(0.9, 0.81). Also, there is no check for causality, so the resulting series may diverge.
7. Variance computations
Notably, the command ARMAtoMA() can be used to compute the variance of a particular ARMA model. For that we compute a long enough series P of MA coefficients (say, r = 100), and then sum their squares plus 1: sum(P^2)+1 Example: P<-ARMAtoMA(ar=0.9, lag.max=100) V<-sum(P^2)+1 V [1] 5.263158 We know that the theoretical variance for AR(1) model with coefficient is given by
V =
1 1 = = 5.263158 . 2 1- 1 - 0.9 2
Thus, the numerical approach gives a very accurate approximation.
8. ARMA estimation
Estimation of the ARMA parameters is done with help the command arima(X, order = c(p, d, q)). We will discuss here a simple use of his command, and will proceed with more advanced applications in the next Labs. The triplet (p,d,q) specifies the order of the model. We will only vary the AR order p and the MA order q, keeping the order of differencing d = 0. Also, we will only apply here the estimation procedure to the true ARMA time series with known values of p and q.
9. ACF, ACVF, and PACF estimation
Estimation of ACF, ACVF, and PACF for any time series (not necessarily ARMA) is done using the command > acf(x,lag.max = r,type = c("correlation", "covariance", "partial"),plot = TRUE)
Page 5
STAT 758 Fall 2008 Lab 3
Appendix 1: Lab3.R
#====================...