Vector Autoregressive Modelling A brief introduction to vars package - - PowerPoint PPT Presentation

▶
vector autoregressive modelling
SMART_READER_LITE
LIVE PREVIEW

Vector Autoregressive Modelling A brief introduction to vars package - - PowerPoint PPT Presentation

May 06, 2023 339 likes •686 views

Vector Autoregressive Modelling A brief introduction to vars package Jilber Urbina R User Group Barcelona May 2012. Jilber Urbina (RUG BCN) vars Package May 2012. 1 / 33 Outline Introduction to VAR models 1 Estimation 2 Restricted VARs

slide-1
SLIDE 1

Vector Autoregressive Modelling

A brief introduction to vars package Jilber Urbina

R User Group Barcelona

May 2012.

Jilber Urbina (RUG BCN) vars Package May 2012. 1 / 33

slide-2
SLIDE 2

Outline

1

Introduction to VAR models

2

Estimation

3

Restricted VARs

4

Diagnostic testing

5

Causality Analysis

6

Forecasting

7

Impulse Response Function

8

Forecast Error Variance Decomposition This presentation is entirely based on «Using the vars package» downloadable from http://rss.acs.unt.edu/Rdoc/library/vars/doc/vars.pdf

Jilber Urbina (RUG BCN) vars Package May 2012. 2 / 33

slide-3
SLIDE 3

Introduction

VAR stands for Vector Autoregressive which is a natural extension to univariate Autoregressive models (AR). In its basic form, a VAR consists of a set of K endogenous variables yt = (y1t, . . . , ykt, . . . , yKt) for k = 1, . . . , K. The VAR(p)-process is then de…ned as: yt = A1yt1 + + Apytp + ut =

p

∑

i=1

Aiyti + ut where Ai are (K K) coe¢cients matrices for i = 1, . . . , p and ut is a Kdimensional white noise process with time invariant positive de…nite covariance matrix E(utu0

t) = Σu.

A VAR process is said to be stable if and only if all its characteristic roots lie inside the complex unit circle, such roots are obtained by: det (IK A1z . . . Apzp) 6= 0 for jzj 1. where z are the eigenvalues.

Jilber Urbina (RUG BCN) vars Package May 2012. 3 / 33

slide-4
SLIDE 4

Data

Before estimating the VAR model, we are going to say some words about the data set we are going to work with. Data set: Canada consists of canadian macroeconomic time series published by the OECD. The sample range is from the 1stQ 1980 until 4thQ 2000. Variables: e: Employment prod: Productivity rw: Real Wage U: Unemployment Rate.

Jilber Urbina (RUG BCN) vars Package May 2012. 4 / 33

slide-5
SLIDE 5

Data

Canada is available from vars package In order to access the data we have to execute the following lines: install.packages(’vars’) library(vars) data(Canada) ?Canada

Jilber Urbina (RUG BCN) vars Package May 2012. 5 / 33

slide-6
SLIDE 6

Estimation

Empirically speaking the …rst step when dealing with VARs is to determine the lag length (p) to be included in the estimation process. This is a decision making problem where its solution is based on some information

  • criteria. An optimal lag-order can be determined according to an

information criterion or the …nal prediction error of a VAR(p) with the function VARselect(). > args(VARselect) function (y, lag.max = 10, type = c("const", "trend", "both","none"))

Jilber Urbina (RUG BCN) vars Package May 2012. 6 / 33

slide-7
SLIDE 7

Lag length selection.

> VARselect(Canada, lag.max = 5, type = "const") $selection AIC(n) HQ(n) SC(n) FPE(n) 3 2 2 3 $criteria

1 2 3 4 5 AIC(n)

  • 5.919117819
  • 6.45220283
  • 6.499021907
  • 6.247207996
  • 6.027766024

HQ(n)

  • 5.726859935
  • 6.06768706
  • 5.922248254
  • 5.478176460
  • 5.066476603

SC(n)

  • 5.439229647
  • 5.49242648
  • 5.059357389
  • 4.327655306
  • 3.628325161

FPE(n) 0.002976003 0.00175206 0.001685528 0.002201523 0.002811116

We’ll work with 2 lags.

Jilber Urbina (RUG BCN) vars Package May 2012. 7 / 33

slide-8
SLIDE 8

Estimation

Now we’re ready for estimating the 4-variables second order VAR, VAR(2). The VAR() function allows us to estimate the model, its arguments are: > args(VAR) function (y, p = 1, type = c("const", "trend", "both", "none")) The VAR(2) is estimated with the function VAR() and as deterministic regressors a constant is included. > var.2c <- VAR(Canada, p = 2, type = "const")

Jilber Urbina (RUG BCN) vars Package May 2012. 8 / 33

slide-9
SLIDE 9

Summary

The function summary returns a list object of class varest and the list elements are explained in detail in the function’s help …le. Let us now focus on two methods, namely summary and plot. The summary method simply applies the summary.lm method to the lm objects contained in varresult. > summary(var.2c)

VAR Estimation Results: ========================= Endogenous variables: e, prod, rw, U Deterministic variables: const Sample size: 82 Log Likelihood: -175.819 Roots of the characteristic polynomial: 0.995 0.9081 0.9081 0.7381 0.7381 0.1856 0.1429 0.1429 Call: VAR(y = Canada, p = 2, type = "const")

Jilber Urbina (RUG BCN) vars Package May 2012. 9 / 33

slide-10
SLIDE 10

Estimation: Summary

Estimation results for equation prod: =====================================

prod = e.l1 + prod.l1 + rw.l1 + U.l1 + e.l2 + prod.l2 + rw.l2 + U.l2 + const Estimate Std. Error t value Pr(>jtj) e.l1

  • 0.17277

0.26977

  • 0.640

0.52390 prod.l1 1.15043 0.10995 10.464 3.57e-16 *** rw.l1 0.05130 0.09934 0.516 0.60710 U.l1

  • 0.47850

0.36470

  • 1.312

0.19362 e.l2 0.38526 0.25688 1.343 0.18346 prod.l2

  • 0.17241

0.11881

  • 1.451

0.15104 rw.l2

  • 0.11885

0.09985

  • 1.190

0.23778 U.l2 1.01592 0.37285 2.725 0.00808 *** const

  • 166.77552

100.43388

  • 1.661

0.10109

Residual standard error: 0.6525 on 73 degrees of freedom Multiple R-Squared: 0.9787, Adjusted R-squared: 0.9764 F-statistic: 419.3 on 8 and 73 DF, p-value: < 2.2e-16

Jilber Urbina (RUG BCN) vars Package May 2012. 10 / 33

slide-11
SLIDE 11

Plot

plot(var.2c)

Jilber Urbina (RUG BCN) vars Package May 2012. 11 / 33

slide-12
SLIDE 12

roots(): VAR stability

As it was said before, a VAR model is stable (stationary) if and only if the roots of its characteristic polynomial are all equal or less than 1. > summary(var.2c)$roots > roots(var.2c) [1] 0.9950338 0.9081062 0.9081062 0.7380565 0.7380565 0.1856381 0.1428889 [8] 0.1428889 > all(roots(var.2c) <=1) [1] TRUE Although, the …rst eigenvalue is pretty close to unity, for the sake of simplicity, we assume a stable VAR(2)-process with a constant as deterministic regressor. Therefore, the estimated VAR model is stable, hence it is stationary.

Jilber Urbina (RUG BCN) vars Package May 2012. 12 / 33

slide-13
SLIDE 13

restrict(): Restricted VARs

It is obvious that not all regressors enter signi…cantly. With the function restrict() the user has the option to re-estimate the VAR either by signi…cance (argument method = "ser") or by imposing zero restrictions manually (argument method = "manual"). In the former case, each equation is re-estimated separately as long as there are t-values that are in absolut value below the threshhold value set by the function’s argument

  • thresh. In the latter case, a restriction matrix has to be provided that

consists of 0/1 values, thereby selecting the coe¢cients to be retained in the model. The function’s arguments are therefore: > args(restrict) function (x, method = c("ser", "manual"), thresh = 2, resmat = NULL)

Jilber Urbina (RUG BCN) vars Package May 2012. 13 / 33

slide-14
SLIDE 14

restrict(): Restricted VARs

Let’s re-estimate the model leaving out all non-signi…cant coe¢cients at 5% level.

Jilber Urbina (RUG BCN) vars Package May 2012. 14 / 33

slide-15
SLIDE 15

B(): Coe¢cients of the restricted VAR

Jilber Urbina (RUG BCN) vars Package May 2012. 15 / 33

slide-16
SLIDE 16

Normality

Jarque-Bera test for univariate case: JB = n 6

  • s2 + 1

4(k 3)2

  • χ2

2

where s is the skewness, k stands for kurtosis, and n is the sample size. Jarque-Bera test for multivariate case: JBmv = s2

3 + s2 4 χ2 2K

where s2

3 and s2 4 are computed according to:

s2

3

= Tb0

3b3

6 s2

4

= T 24(b4 3K )0(b4 3K )

with b3 and b4 are the third and fourth non-central moment vectors of the standardized residuals ^

us

t = ˜

P (ˆ ut

_

^ ut) where ˜ P is the lower triangular

choleski matrix.

Jilber Urbina (RUG BCN) vars Package May 2012. 16 / 33

slide-17
SLIDE 17

normality(): Jarque-Bera test for normality. Univariate

In package ‘vars’ the functions for diagnostic testing are normality.test(), serial() and stability(). > args(normality.test) function (x, multivariate.only = TRUE) > normality.test(var.2c, multivariate.only = TRUE) $e JB-Test (univariate) data: Residual of e equation Chi-squared = 0.1535, df = 2, p-value = 0.9261 $prod JB-Test (univariate) data: Residual of prod equation Chi-squared = 4.2651, df = 2, p-value = 0.1185

Jilber Urbina (RUG BCN) vars Package May 2012. 17 / 33

slide-18
SLIDE 18

normality(): Jarque-Bera test for normality. Univariate

$rw JB-Test (univariate) data: Residual of rw equation Chi-squared = 0.3348, df = 2, p-value = 0.8459 $U JB-Test (univariate) data: Residual of U equation Chi-squared = 0.5664, df = 2, p-value = 0.7534

Jilber Urbina (RUG BCN) vars Package May 2012. 18 / 33

slide-19
SLIDE 19

normality(): Jarque-Bera test for normality. Multivariate

> normality.test(var.2c, multivariate.only = TRUE) $JB JB-Test (multivariate) data: Residuals of VAR object var.2c Chi-squared = 5.094, df = 8, p-value = 0.7475

Jilber Urbina (RUG BCN) vars Package May 2012. 19 / 33

slide-20
SLIDE 20

Serial autocorrelation

To test for multivariate serial autocorrelation I will use the Breusch-Godfrey LM-statistic which is based upon the following auxiliary regressions: ^ ut = A1yt1 + . . . + Apytp + B1^ ut1 + . . . + Bh^ uth + εt where the null hypothesis is: H0 : B1 = . . . = Bh = 0 H1 : 9Bi 6= 0 for i = 1, 2, . . . , h. The test statistic is de…ned as: LMh = T(K tr( ˜ Σ1

R ˜

Σe)) χ2

(hK 2)

where ˜ ΣR and ˜ Σe assign the residual covariance matrix of the restricted and unrestricted model, respectively.

Jilber Urbina (RUG BCN) vars Package May 2012. 20 / 33

slide-21
SLIDE 21

serial.test()

> args(serial.test) function (x, lags.pt = 16, lags.bg = 5, type = c("PT.asymptotic", "PT.adjusted", "BG", "ES")) > var2c.serial <- serial.test(var.2c, lags.pt = 16, lags.bg = 5, type=’BG’) Breusch-Godfrey LM test data: Residuals of VAR object var.2c Chi-squared = 92.6282, df = 80, p-value = 0.1581 We cannot reject the null hypothesis: Residuals of the system are not autocorrelated.

Jilber Urbina (RUG BCN) vars Package May 2012. 21 / 33

slide-22
SLIDE 22

plot()

Jilber Urbina (RUG BCN) vars Package May 2012. 22 / 33

slide-23
SLIDE 23

stability()

The stability of the regression relationships in a VAR(p) can be assessed with the function stability(). An empirical ‡uctuation process is estimated for each regression by passing the function’s arguments to the efp()-function contained in the package strucchange. > args(stability) function (x, type = c("OLS-CUSUM", "Rec-CUSUM", "Rec-MOSUM", "OLS-MOSUM", "RE", "ME", "Score-CUSUM", "Score-MOSUM", "fluctuation"), h = 0.15, dynamic = FALSE, rescale = TRUE) var2c.stab <- stability(var.2c, type="OLS-CUSUM") plot(var2c.stab, plot.type = "single")

Jilber Urbina (RUG BCN) vars Package May 2012. 23 / 33

slide-24
SLIDE 24

plot()

Jilber Urbina (RUG BCN) vars Package May 2012. 24 / 33

slide-25
SLIDE 25

Causality Analysis.

Very often (almots always in economics) researchers are interested in the detection of causalities between variables. The most common one is the Granger-Causality test (Granger 1969). We can conduct causality test using causality() function. > args(causality) function (x, cause = NULL, vcov. = NULL, boot = FALSE, boot.runs = 100) The function causality() is now applied for investigating if the real wage and productivity is causal to employment and unemployment. > causality(var.2c, cause = c("rw", "prod")) $Granger Granger causality H0: prod rw do not Granger-cause e U data: VAR object var.2c F-Test = 3.4529, df1 = 8, df2 = 292, p-value = 0.0008086 $

Jilber Urbina (RUG BCN) vars Package May 2012. 25 / 33

slide-26
SLIDE 26

Forecasting

A predict-method for objects with class attribute varest is available. > var.f10 <- predict(var.2c, n.ahead=10, ci=0.95) > var.f10 $ fcst $ e fcst lower upper CI [1,] 962.6557 961.9446 963.3668 0.7111044 [2,] 963.6538 962.3422 964.9654 1.3116050 [3,] 964.6932 962.8261 966.5603 1.8670903 [4,] 965.6882 963.3092 968.0671 2.3789396 [5,] 966.5814 963.7240 969.4387 2.8573301 [6,] 967.3460 964.0344 970.6576 3.3116112 [7,] 967.9769 964.2302 971.7236 3.7467269 [8,] 968.4827 964.3193 972.6461 4.1633974 [9,] 968.8798 964.3200 973.4396 4.5597736 [10,] 969.1877 964.2546 974.1208 4.9330961

Jilber Urbina (RUG BCN) vars Package May 2012. 26 / 33

slide-27
SLIDE 27

plotting the forecast using plot()

plot(var.f10, plot.type = "single")

Jilber Urbina (RUG BCN) vars Package May 2012. 27 / 33

slide-28
SLIDE 28

plotting the forecast using fanchart()

fanchart(var.f10, plot.type = "single")

Jilber Urbina (RUG BCN) vars Package May 2012. 28 / 33

slide-29
SLIDE 29

Impulse -Response Analysis.

The impulse response analysis is based upon the Wold moving average representation of a VAR(p)-process. It is used to investigate the dynamic interactions between the endogenous variables. > args(irf) function (x, impulse = NULL, response = NULL, n.ahead = 10,

  • rtho = TRUE, cumulative = FALSE, boot = TRUE, ci =

0.95,runs = 100, seed = NULL, ...) > Canada2 <- Canada[, c(3, 1, 4, 2)] > colnames(Canada2) [1] "rw" "e" "U" "prod" > var.2c.alt <- VAR(Canada2, p = 2, type = "const") > irf.rw.eU <- irf(var.2c.alt, impulse = "rw", response = c("e",’U’))

Jilber Urbina (RUG BCN) vars Package May 2012. 29 / 33

slide-30
SLIDE 30

Plot Impulse -Responses

> plot(irf.rw.eU)

Jilber Urbina (RUG BCN) vars Package May 2012. 30 / 33

slide-31
SLIDE 31

Forecast Error Variance Decomposition (FEVD)

The FEVD allows the user to analyse the contribution of variable j to the h-step forecast error variance of variable k. > args(fevd) function (x, n.ahead = 10, ...) NULL > var2c.fevd <- fevd(var.2c, n.ahead = 5) > var2c.fevd$e e prod rw U [1,] 1.0000000 0.00000000 0.000000000 0.000000000 [2,] 0.9633815 0.02563062 0.004448081 0.006539797 [3,] 0.8961692 0.06797131 0.013226872 0.022632567 [4,] 0.8057174 0.11757589 0.025689192 0.051017495 [5,] 0.7019003 0.16952744 0.040094324 0.088477959

Jilber Urbina (RUG BCN) vars Package May 2012. 31 / 33

slide-32
SLIDE 32

plot (FEVD)

Jilber Urbina (RUG BCN) vars Package May 2012. 32 / 33

slide-33
SLIDE 33

It’s over, thank you for attending this presentation.

Vector Autoregressive Modelling

A brief introduction to vars package

Jilber Urbina RUG-BCN

Jilber Urbina (RUG BCN) vars Package May 2012. 33 / 33