VAR, VARMA and VMA Models 8-1
Revised Chapter 8 in Specifying and Diagnostically Testing Econometric Models (Edition 3)
© by Houston H. Stokes 6 January 2013. All rights reserved. Preliminary Draft
Chapter 8
Time Series Analysis Part II: VAR, VARMA and VMA Models
8.0 Introduction
8.1 The Relationship Between VAR, VARMA and Structural Models
8.2 The Structural VAR Formulation
Table 8.1 Matlab Program to Invert zero order bivariate and trivariate VAR(2) Matrix
Table 8.2 Beta Inverse for bivariate and trivariate VAR Structural Models
8.3 Identification of VAR and VARMA Models
8.4 Testing Series for Nonlinearity With Hinich Tests
Table 8.3 Number of False Positives as a Function of M
8.5 Examples
Table 8.5 Estimated Coefficients for Gas Furnace Data
Table 8.7 Estimating VAR Models for the Gas Furnace Data
Table 8.8 Forecast Validation Macro
Table 8.9 Out-of-Sample Forecast Diagnostics
Table 8.10 Out-of-Sample Forecasts for GASIN
Table 8.11 Out-of-Sample Forecasts for GASOUT
8.6 Testing for Coefficient and Variance Changes
Table 8.12 Two Period Stock Watson Test of Money Interest Model – 100 Replications
Table 8.13 Stock-Watson Model estimated over a range of dates – For RES (1979) Data.
Table 8.14 Code to perform Plots using Extract File
Figure 8.1 Moving Stock-Watson values For Money Interest Rate Series
Figure 8.2 Moving Stock-Watson Differences for Money Interest Series Part 1
Figure 8.3 Moving Stock-Watson Differences for Money Interest Series Part 2
Table 8.15 Moving Stock Watson Model with Critical Values
Table 8.16 Plotting of Frankle Data
Figure 8.5 Moving Plots of Frankle Data for Differences 1
Figure 8.6 Moving Plots of Frankle Data for Differences 2
8.7 Variance Decomposition of a VAR Model
Table 8.17 VARRESDC and ST_RES
Table 8.18 Test Program to Illustrate a Variance Decomposition
8.8 Alternate IRF definitions
Table 8.19 IRF a subroutine to calculate the IRF of a VAR model
Table 8.20 Illustrating the calculation of the IRF of a VAR model
8.9 Testing a VAR Model for Stability
Table 8.21 VARSTAB a subroutine to test a VAR model for stability
Table 8.22 Code to illustrate VARSTAB using B34S and Stata
8.10 Conclusion
Time Series Analysis Part II: VAR, VARMA and VMA Models
8.0 Introduction
The B34S commands btiden and btest control the identification and estimation of VAR and VARMA models. The basic code for this section was obtained from the Department of Statistics, University of Wisconsin, and is documented in Tiao et al (1979)[1]and Tiao and Box (1981). After first discussing how VAR and VARMA models relate to structural models and the concept of Granger (1969) causality, VAR and VARMA model identification is reviewed. These models are used to test for Granger (1969) causality. The Hinich (1982) nonlinearity test is shown to be an added diagnostic tool to determine whether the model is correctly specified. The VAR and VARMA models are illustrated, using the gas furnace data, which is also used in chapter 12 to illustrate the frequency approach to the VAR model estimation proposed by Geweke (1982a, 1982b, 1982c, 1984).
8.1 The Relationship Between VAR, VARMA and Structural Models
Quenouille (1957) argued that a dynamic system of simultaneous equations could be written and estimated in the form
(8.1-1)
where is a row vector of the t th observation on k series , and G(B) and D(B) are k by k polynomial matrices in which each element, , is itself a polynomial vector in the lag operator B. We assume that has been suitably differenced to achieve stationarity. Assuming k=3, then (8.1-1) can be written in expanded form as
(8.1-2)
If for
(8.1-3)
then equation (8.1-2) reduces to three ARIMA models of the form of equation (7.1-9), where
(8.1-4)
and
(8.1-5)
and where are the AR and MA polynomials of the vector of . In this case, the three series, , are independent. We can say that series is exogenous with respect to series and and that is exogenous with respect to if for j>i
(8.1-6)
In such a situation, a transfer function of the form of equation (7.2-1) can be estimated. If, on the other hand for j > i
(8.1-7)
and/or
(8.1-8)
then there is feedback in the system and a transfer function model of the form of equation (7.2-1) is not the appropriate way to proceed. Note that using the above formulation there are no zero order terms for and for and zero order terms = 1 for . Such a formulation has been called the standard VAR formulation by Enders (2004, 264). Section 8.2 outlines how this formulation relates to a structural formulation.
The above discussion has highlighted the fact that the VARMA model of the form of equation (8.1-1) is a very general functional form of which the transfer function and the ARIMA model are increasingly more special cases. If we assume that and or that D(B) is a matrix of degree 0 in B, then equation (8.1-1) reduces to the VAR form of the model
(8.1-9)
where . In the more general case, a VARMA model, such as equation (8.1-1), can be written as a VAR model, provided that D(B) is invertible. Here
.(8.1-10)
In a like manner, equation (8.1-1) can be written in VMA form as
,(8.1-11)
where
or that G(B) is a matrix of degree zero in B. In general,
(8.1-12)
if G(B) is invertible.[2] It is to be stressed that provided invertibility conditions are satisfied, equations (8.1-1), (8.1-9) and (8.1-11) are alternative forms. Usually, equation (8.1-1), the VARMA form, is the most parsimonious representation. Equation (8.1-9), the VAR form, is usually estimated first as a way to identify the order of the VARMA model. Sims (1980) advocated estimating the model in the form of equation (8.1-9) and calculating R(B) in equation (8.1-11) as . The pattern in the elements of the polynomials in R(B), element by element, would trace the effects of "shocks" on the variables in . For example, the term R21(B) measures the effect of an unexplained shock "innovation" in the first series on the second series.
The concept of Granger (1969) causality is related to the econometric concept of exogeneity. A series is said to Granger cause a series if, and only if, a model that predicts as a function of only its past has a greater sum of squares of the error term than a model that predicts as a function of its own past and the past of . Thus, in equation (8.1-2) if and is exogenous, we can say that Granger (1969) causes. Using the Tiao-Box (1981) approach, the assumptions needed to estimate equation (8.1-9) include the following:
- What variables to place in .
- The orders of the differencing in to make the series stationary.
- The maximum degree of any element in the matrix Q(B).
These assumptions are substantially less that those needed to estimate a simultaneous equations model of the form of equation (4.1-1) but more stringent than Engle - Granger (1987), who suggested cointegrated models where some or all the variables in the vector could be nonstationary as long as the k vectors of residuals in are stationary. These cointegration models are discussed further in 12.4.
Zellner and Palm (1974) recommend rewriting equation (8.1-1) as
, (8.1-13)
where G*(B) and │G(B)│ are the adjoint and determinant of G(B). Equation (8.1-13) quickly reduces to
,(8.1-14)
which expresses the VARMA model in the form of a restricted, seemingly unrelated, autoregressive model with correlated moving average errors that are not restricted to be only in period t as was the case with the SUR model. This can be seen if we rewrite the ith equation of (8.1-14) as
,(8.1-15)
VAR, VARMA and VMA Models 8-1
where is the i th row of . Following Zellner's and Palm's (1974) classic paper showing the relationship between VARMA models and structural equations models, we next can partition into the k1 endogenous and k2 exogenous variables and rewrite equation (8.1-1) as
, (8.1-16)
where are nowki by kj matrices.The assumptions of structural equation modeling assume that the k1 variables in vector are endogenous and that the variables in the vector are exogenous. This restriction implies that
.(8.1-17)
The assumption insures that shocks from the endogenous variables do not impact the exogenous variables while the assumption insures that in the structural equation (8.1-18) below only the exogenous variables enter the model, not their shocks. Further detail on this is contained in section 8.2 below. Further insight into this can be obtained by assuming the VAR form of the Model in (8.1-9) and setting since x is exogenous by assumption. If (8.1-12) is used to form R(B), then as can be demonstrated by the following code.[3]
/; Illustrates that if VAR is lower triangular,
/; VMA will be lower triangular
b34sexec matrix;
* problem if var lower T => vma will be lower T;
a=array(:1.,.0, .0,1.,-.3,.4,.0,-.1,.2 .3 .0 .1);
ia=index(2,2,3);
nterms=20;
call echooff;
call polymdisp(:display a ia);
call polyminv(a,ia,ainv,iainv,nterms);
call names(all);
call print(%p,%det);
call polymdisp(:display ainv iainv);
call polymdisp(:display %adj %iadj);
call polymmult(a ia ainv iainv test itest);
call polymdisp(:display test itest);
call polymdisp(:extract ainv iainv vec1 index(1,1,0));
call polymdisp(:extract ainv iainv vec2 index(2,1,0));
call polymdisp(:extract ainv iainv vec3 index(1,2,0));
call polymdisp(:extract ainv iainv vec4 index(2,2,0));
call names(all);
call tabulate(vec1,vec2,vec3,vec4);
b34srun;
that inverts
And forms for elements (1,1), (2.1), (1,2) and (2,2) of R(B):
Order R(1,1) R(2,1) R(1,2) R(2,2)
0 1.000 0.000 0.000 1.000
1 0.3000 -0.4000 0.000 0.1000
2 -0.1100 -0.4600 0.000 -0.9000E-01
3 -0.9300E-01 -0.5200E-01 0.000 -0.1900E-01
4 -0.5900E-02 0.1110 0.000 0.7100E-02
5 0.1683E-01 0.4656E-01 0.000 0.2610E-02
6 0.6229E-02 -0.1141E-01 0.000 -0.4490E-03
7 -0.1497E-02 -0.1334E-01 0.000 -0.3059E-03
8 -0.1695E-02 -0.1463E-02 0.000 0.1431E-04
9 -0.2090E-03 0.2315E-02 0.000 0.3202E-04
10 0.2763E-03 0.9699E-03 0.000 0.1771E-05
11 0.1247E-03 -0.1823E-03 0.000 -0.3025E-05
12 -0.1785E-04 -0.2480E-03 0.000 -0.4796E-06
13 -0.3029E-04 -0.3684E-04 0.000 0.2545E-06
14 -0.5518E-05 0.3859E-04 0.000 0.7341E-07
15 0.4403E-05 0.1884E-04 0.000 -0.1811E-07
16 0.2425E-05 -0.2081E-05 0.000 -0.9153E-08
17 -0.1533E-06 -0.4383E-05 0.000 0.8960E-09
18 -0.5309E-06 -0.8963E-06 0.000 0.1005E-08
19 -0.1286E-06 0.6070E-06 0.000 0.1089E-10
A major contribution of VARMA model building is that maintained structural equations assumptions are now subject to formal testing. Sims (1980) makes this point forcefully. From equations (8.1-16) and (8.1-17) we obtain the implied structural equations
(8.1-18)
and
.(8.1-19)
Assuming that is invertible, equation (8.1-18) can be transformed to the final form set of k1 equations, where all endogenous variables are functions of the exogenous variables
(8.1-20)
or the reduced form set of k1 equations, where the endogenous variables are a function of the exogenous variables and the lagged endogenous variables
, (8.1-21)
where
(8.1-22)
and
.(8.1-23)
VAR, VARMA and VMA Models 8-1
The above discussion shows that the structural equations assumptions place testable restrictions on the form of and. The revealed structure of these polynomial matrices will impose restrictions on the reduced form model given in equation (8.1-21) and final form model given in equation (8.1-20). Rather than seeing time series and structural equations modeling as distinct approaches, Zellner and Palm (1974) showed how these econometric techniques were related. Econometric practice has not been the same since. The steps to identify a VARMA model will be discussed in the next section.
Equation (8.1-1) writes the series in only in terms of their lags and the lags of the other series. Contemporaneous effects are captured in the error process since, in general, the covariance matrix is not diagonal.[4] Granger and Newbold (1977, 223) discuss the possible transformations. Equation (8.1-1) assumes the off diagonal elements of order 0 in G(B) and D(B) are zero but that the covariance of the error terms was not necessarily zero. The covariance between the error process will be zero if there is no instantaneous causality in the system. Granger and Newbold (1977, 223) argue that this form of the model, which they call model A, can be transformed into a model in which the off diagonal elements in G(B) are not zero, yet the covariance matrix, , is diagonal. This transformation requires calculation of a k by k P matrix defined such that
,(8.1-24)
where off diagonal elements in P measure the contemporaneous relationship between the variables in the vector. The btiden command, which is used to estimate the VAR model, and the btest command, which is used to estimate the VMA and VARMA models, automatically print P, P-1 and the diagonal elements of as well as. These values can be used to transform the estimated model to what Granger and Newbold (1977, 224) call model C:
,(8.1-25)
since .[5]
In a manner similar to equation (7.1-4) the VAR, VARMA and VMA models can be broken down into two factors or matrices on the left and two factors or matrices on the right. In terms of equation (8.1-1), we write and where and are the regular VAR and VMA matrices of order k by k and and are the seasonal VAR and seasonal VMA matrices of order k by k. Using this setup equation (8.1-1) becomes
.(8.1-26)
VAR, VARMA and VMA Models 8-1
Following logic similar to (7.1-23 and (7.1-24), if there are seasonal elements of the model, estimating (8.1-26) in place of (8.1-1) will result in fewer parameters. The disadvantage of (8.1-26) is that it is inherently more complex and some of the cross terms may not be needed.
8.2 The Structural VAR Formulation
The VAR formulations discussed in the prior section do not contain contemporaneous terms in and in equation (8.1-1).[6] This formulation is ideally suited for Granger (1969) causality tests.If the assumption of no contemporaneous terms is relaxed, we have the structural form of the VAR, which was discussed in Enders (2004, 264-266)(2010, 297-349) in the context of a 2 variable model. After first discussing this example, the results will be generalized using the symbolic capability in Matlab. Generalizing Enders first order structural VAR model to order k and defining produces
(8.2-1)
(8.2-2)
Combining (8.2-1) and (8.2-2) and moving contemporaneous terms to the left
(8.2-3)
or in more compact notation
(8.2-4)
Provided exists, the standard VAR form, (8.1-9) becomes
(8.2-5)
In terms of (8.1-9)
(8.2-6)
and . The importance of this derivation is that it shows that while by assumption the elements of structural form of the model shown in equations (8.2-1) and (8.2-2) , are not contemporaneously related, this is not true for the elements of the standard form or VAR model shown in equation (8.1-9). To understand what is involved, the Matlab program in Table 8.1 calculates for two and three variable models and for a 2 variable model, the answers of which are shown in Table 8.2. The error vector for the standard VAR form is calculated in Table 8.2 and reduces to
(8.2-7)
(8.2-8)
Table 8.1 Matlab Program to Invert zero order bivariate and trivariate VAR(2) Matrix
disp('Inversion of structural bivariate zero order VAR matrix')
syms b12b21g111g121g211g221e1e2
%%
B=[1, b12;b21,1]
inv(B)
disp(' ')
disp('Covariance of standard form VAR Error term')
e_struc=[e1;e2]
E=inv(B)*e_struc
disp(' ')
disp('Effect of Lag 1 variables on left hand side')
gamma1=[g111, g121; g211,g221]
A1=inv(B)*gamma1
%% Inversion of structural trivariate matrix
disp(' ')
disp('Inversion of structural trivariate zero order VAR Matrix')
syms b12b13b21b23b31b32
B=[ 1, b12, b13;
b21, 1, b23;
b31, b32, 1]
inv(B)
Table 8.2 Beta Inverse for bivariate and trivariate VAR Structural Models
Inversion of structural bivariate zero order VAR matrix
B =
[ 1, b12]
[ b21, 1]
ans =
[ -1/(-1+b12*b21), b12/(-1+b12*b21)]
[ b21/(-1+b12*b21), -1/(-1+b12*b21)]
Covariance of standard form VAR Error term
e_struc =
e1
e2
E =
-1/(-1+b12*b21)*e1+b12/(-1+b12*b21)*e2
b21/(-1+b12*b21)*e1-1/(-1+b12*b21)*e2
Effect of Lag 1 variables on left hand side
gamma1 =
[ g111, g121]
[ g211, g221]
A1 =
[ -1/(-1+b12*b21)*g111+b12/(-1+b12*b21)*g211,
-1/(-1+b12*b21)*g121+b12/(-1+b12*b21)*g221]
[ b21/(-1+b12*b21)*g111-1/(-1+b12*b21)*g211,
b21/(-1+b12*b21)*g121-1/(-1+b12*b21)*g221]
Inversion of structural trivariate zero order VAR Matrix
B =
[ 1, b12, b13]
[ b21, 1, b23]
[ b31, b32, 1]
ans =
[ (-1+b23*b32)/(-1+b23*b32+b12*b21-b21*b13*b32-b31*b12*b23+b31*b13),
(b12-b13*b32)/(-1+b23*b32+b12*b21-b21*b13*b32-b31*b12*b23+b31*b13),
-(b12*b23-b13)/(-1+b23*b32+b12*b21-b21*b13*b32-b31*b12*b23+b31*b13)]
[ (b21-b23*b31)/(-1+b23*b32+b12*b21-b21*b13*b32-b31*b12*b23+b31*b13),
(-1+b31*b13)/(-1+b23*b32+b12*b21-b21*b13*b32-b31*b12*b23+b31*b13),
-(-b23+b13*b21)/(-1+b23*b32+b12*b21-b21*b13*b32-b31*b12*b23+b31*b13)]
[ -(b21*b32-b31)/(-1+b23*b32+b12*b21-b21*b13*b32-b31*b12*b23+b31*b13),
-(-b32+b12*b31)/(-1+b23*b32+b12*b21-b21*b13*b32-b31*b12*b23+b31*b13),
(-1+b12*b21)/(-1+b23*b32+b12*b21-b21*b13*b32-b31*b12*b23+b31*b13)]
From (8.2-7) and (8.2.8) we can show that the variances of the errors of the standard VAR form are time-independent but are potentially correlated. This can be seen if we note that
(8.2-9)
(8.2-10)
and
(8.2-11)
(8.2-12)
Rewrite (8.2-5) as
(8.2-12)
Solving for from the second equation of (8.2-12)
(8.2-13)
and substituting in the first gives
(8.2-14)
which when explicitly solved for produces a constrained ARIMA model with error terms from the two standard VAR equations.
(8.2-15)
For the corresponding equations are
(8.2-16)
(8.2-17)
and
(8.2-18)
The above derivation[7] shows how a constrained ARIMA model with its error decomposed is based on the underlying VAR parameters. Stability of the general structural VAR model (8.2-5) implies that
(8.2-19)
converges. Note that from (8.1-1). The -1.0 times the elements of make up the non zero order elements of Q
The importance of the above discussion is that it highlights, using a simple example, some important aspects of the relationships between the VAR model in both its standard form, its structural form and the transfer function and generalized ARIMA model.
8.3 Identification of VAR and VARMA Models
The first step of VARMA model building is to make all series stationary so that the calculated ACF and CCF can be interpreted. Assume two series, GASIN and GASOUT, from the gas furnace data discussed in chapter 7. The btiden command
b34sexec btiden nstder=3.0$
title=('calculation of acf and ccf on gas data')$
seriesn var=gasin title=('input gas from b-j book')$
seriesn var=gasout title=('output gas from b-j book')$
iden lagrho=36 isacf isccf ivalue$ b34seend$
will list and plot 36 autocorrelations and cross correlations. If the ivalue parameter is left off, these are only plotted in terms of + or - for significant and . for not significant terms. The nstder parameter on the btidensentence allows specification of a value different from the usual 2 SE convention. The ACF and CCF plots are useful to test for stationarity. "Line up" plots can optionally be made for the raw, transformed or differenced data to assist in spotting lead lag relationships visually. Hinich (1982) tests for Gaussianity and nonlinearity, which are discussed in section 8.3, are also possible options.[8]Once the series are deemed stationary, the ESTVAR sentence is used to estimate VAR models of the form of equation (8.1-9) for increasingly longer lags. The setup listed below will run a VAR model of up to order 6 on the gas furnace data.
b34sexec btiden $
Title=('calculation of var model on gas data')$
seriesn var=gasin title=('input gas from b-j book')$
seriesn var=gasout title=('output gas from b-j book')$
estvar p=6 output=normal lagrho=36$ b34seend$
If the keyword ilarf is set on the estvar sentence, only the last VAR fit is printed. If output=brief is set, the estimated coefficients are not given and only a summary table is printed.
Inspection of the cross correlations and autocorrelations of the residuals determines whether the order of the VAR model, k, is sufficient to summarize all the information in the series. Assume that P=12, yet after the 6th-order VAR is estimated, there are no spikes in the cross correlations of the residuals. The researcher next should set P=6 and rerun the command. The answers for the 6th- order VAR model will not be exactly the same as were found the first time. This surprising result is due to the btiden command deleting observations based on the maximum lag of the final VAR model estimated. Using this data set for all VAR estimations allows comparisons of the effects of adding another lag. Hence, once the maximum lag is determined, it is a good plan to rerun the model with that lag specified as the maximum. The ilarf option is especially useful on this second run to save paper.
VAR, VARMA and VMA Models 8-1
The next step in estimating a VARMA model is to delete the nonsignificant VAR coefficients and rerun the model in constrained VAR form with the btest command. Inspection of the cross correlations of the residuals indicates how to specify the MA part of the model. As a general rule, the same principles used in the transfer function identification should be used.
Before we proceed, it is important to discuss some of the diagnostic aids available in the estvar sentence of the btiden command. As each order of the VAR model is estimated, output includes the residual covariance matrix S(j), the residual correlation matrix RS(j), the eigenvalues and eigenvectors of S(j) and the determinant and reciprocal of S(j). Zero eigenvalues in S(J) indicate that, rather than having k independent series, we actually have a process driven by less than k innovation series. Tiao and Box (1981) recommend calculating, which is distributed as chi-square with degrees of freedom, where and