Econometría II: VAR

Departamento de Economía

Carlos A. Yanes G.

2023-10-24

Paquetes con que se trabaja la sesión

Los paquetes que se van a utilizar en la sesión de hoy son:

Note

Para trabajar en esta ocasión vamos a usar los paquetes de :

library(pacman)
p_load(TSstudio, tidyverse, stats, urca, forecast, ggfortify, ggplot2, tseries, fpp2, vars, BVAR, lmtest)
  • Es importante resaltar los paquetes de Vars y BVAR, ya que con el estimaremos los modelos y criterios. Algo de lo que está explicado acá tiene mucho del libro (Lütkepohl and Krätzig 2004) y adicionalmente de un gran texto como (Hamilton 2020)

Preambulo

Modelo VAR

Definición

Un modelo VAR (Vector Autoregression) es un modelo de ecuaciones simultáneas formado por un sistema de ecuaciones de forma reducida. Es un modelo muy útil cuando existe evidencia de simultaneidad entre un grupo de variables, y que sus relaciones se transmiten a lo largo de un determinado número de períodos.

  • \[B^{-1}BX_t=B^{-1}G_0+B^{-1}G_1X_{t-1}+B^{-1}\varepsilon_t\]
  • \[\begin{aligned} e_{xt}=&-\beta_{21}e_{yt}+ \varepsilon_{xt}\\ e_{xt}+&\beta_{21}e_{yt}= \varepsilon_{xt} \end{aligned}\]

Datos

Investigación de Engle y Granger del mercado laboral

  • \(e\) = log del índice de empleo
  • \(Prod\) = Productividad laboral \(\rightarrow\) Diferencia de los logaritmos de e y el PIB
  • \(U\) = Tasa de desempleo
  • \(rw\) = Salario real - Log del indice de salario real
  • Frecuencia Trimestral
  • Primer trimestre de 1980

Datos a trabajar

Code
data(Canada) # Viene con el paquete vars
Datos=Canada
head(Datos)
               e     prod       rw    U
1980 Q1 929.6105 405.3665 386.1361 7.53
1980 Q2 929.8040 404.6398 388.1358 7.70
1980 Q3 930.3184 403.8149 390.5401 7.47
1980 Q4 931.4277 404.2158 393.9638 7.27
1981 Q1 932.6620 405.0467 396.7647 7.37
1981 Q2 933.5509 404.4167 400.0217 7.13
Code
plot(Datos, nc = 2, xlab = "")

Pruebas

Estacionariedad

Code
#"e"
adf.test(Datos[,1], alternative="stationary")

    Augmented Dickey-Fuller Test

data:  Datos[, 1]
Dickey-Fuller = -2.148, Lag order = 4, p-value = 0.5152
alternative hypothesis: stationary
Code
#"prod"
adf.test(Datos[,2], alternative="stationary")

    Augmented Dickey-Fuller Test

data:  Datos[, 2]
Dickey-Fuller = -2.8725, Lag order = 4, p-value = 0.218
alternative hypothesis: stationary

    Augmented Dickey-Fuller Test

data:  Datos[, 3]
Dickey-Fuller = -2.0558, Lag order = 4, p-value = 0.553
alternative hypothesis: stationary

    Augmented Dickey-Fuller Test

data:  Datos[, 4]
Dickey-Fuller = -2.5988, Lag order = 4, p-value = 0.3303
alternative hypothesis: stationary

Sobre estacionariedad

  • Note que ninguna serie logró ser estacionaria. Debemos entonces aplicar o recurrir a la primera diferencia e intentar corregir y volver a realizar las pruebas concernientes a ella.

  • Si el VAR no es estacionario tendremos problemas en la estimación

Corrección

Primera Diferencia

Code
# Datos diferenciados
DatosD=diff(Datos,differences=1)
#"e"
adf.test(DatosD[,1], alternative="stationary")

    Augmented Dickey-Fuller Test

data:  DatosD[, 1]
Dickey-Fuller = -3.2668, Lag order = 4, p-value = 0.08274
alternative hypothesis: stationary

    Augmented Dickey-Fuller Test

data:  DatosD[, 2]
Dickey-Fuller = -3.2361, Lag order = 4, p-value = 0.08775
alternative hypothesis: stationary

    Augmented Dickey-Fuller Test

data:  DatosD[, 3]
Dickey-Fuller = -3.2352, Lag order = 4, p-value = 0.08789
alternative hypothesis: stationary

    Augmented Dickey-Fuller Test

data:  DatosD[, 4]
Dickey-Fuller = -3.7325, Lag order = 4, p-value = 0.02698
alternative hypothesis: stationary

Selección de orden

Selección de orden de VAR

  • El modelo debe ser elegido de acuerdo a un conjunto de criterios entre ellos los tradicionales o ya conocidos AIC, BIC, se añaden los de Función de pronostico de error FP y el de Hannan and Quin HQ

  • Como siempre, aquel modelo con mejor grupo de criterios es el modelo a elegir. Sin embargo se puede optar por un solo criterio o competencia. P.e solo escoger quien mejor criterio Akaike tenga.

Code
Seleccion=VARselect(DatosD, lag.max = 5, type = "none")
Seleccion$criteria
                  1            2           3            4            5
AIC(n) -5.490494928 -5.701482671 -5.61414426 -5.753612581 -5.575571554
HQ(n)  -5.296970039 -5.314432894 -5.03356960 -4.979513026 -4.607947110
SC(n)  -5.007067476 -4.734627768 -4.16386191 -3.819902774 -3.158434295
FPE(n)  0.004127288  0.003350696  0.00368202  0.003247088  0.003970927
Code
Seleccion$selection
AIC(n)  HQ(n)  SC(n) FPE(n) 
     4      2      1      4 

De forma MCO

Modelo “a mano”

Code
t=2:nrow(DatosD)
ma1=lm(DatosD[t,"e"]~-1+DatosD[(t-1),"e"]+DatosD[(t-1),"prod"]+DatosD[(t-1),"U"]+DatosD[(t-1),"rw"])
ma2=lm(DatosD[t,"prod"]~-1+DatosD[(t-1),"e"]+DatosD[(t-1),"prod"]+DatosD[(t-1),"U"]+DatosD[(t-1),"rw"])
ma3=lm(DatosD[t,"U"]~-1+DatosD[(t-1),"e"]+DatosD[(t-1),"prod"]+DatosD[(t-1),"U"]+DatosD[(t-1),"rw"])
ma4=lm(DatosD[t,"rw"]~-1+DatosD[(t-1),"e"]+DatosD[(t-1),"prod"]+DatosD[(t-1),"U"]+DatosD[(t-1),"rw"])
Code
summary(ma1)

Call:
lm(formula = DatosD[t, "e"] ~ -1 + DatosD[(t - 1), "e"] + DatosD[(t - 
    1), "prod"] + DatosD[(t - 1), "U"] + DatosD[(t - 1), "rw"])

Residuals:
     Min       1Q   Median       3Q      Max 
-0.78224 -0.21426  0.01546  0.30546  1.22327 

Coefficients:
                        Estimate Std. Error t value Pr(>|t|)    
DatosD[(t - 1), "e"]     0.87946    0.11545   7.617 5.16e-11 ***
DatosD[(t - 1), "prod"]  0.20108    0.06149   3.270   0.0016 ** 
DatosD[(t - 1), "U"]     0.26501    0.18679   1.419   0.1599    
DatosD[(t - 1), "rw"]   -0.01211    0.03928  -0.308   0.7586    
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 0.3919 on 78 degrees of freedom
Multiple R-squared:  0.7036,    Adjusted R-squared:  0.6884 
F-statistic: 46.28 on 4 and 78 DF,  p-value: < 2.2e-16
Code
summary(ma2)

Call:
lm(formula = DatosD[t, "prod"] ~ -1 + DatosD[(t - 1), "e"] + 
    DatosD[(t - 1), "prod"] + DatosD[(t - 1), "U"] + DatosD[(t - 
    1), "rw"])

Residuals:
    Min      1Q  Median      3Q     Max 
-1.8131 -0.2911  0.1616  0.4790  1.7160 

Coefficients:
                        Estimate Std. Error t value Pr(>|t|)  
DatosD[(t - 1), "e"]    -0.11538    0.20252  -0.570   0.5705  
DatosD[(t - 1), "prod"]  0.28159    0.10786   2.611   0.0108 *
DatosD[(t - 1), "U"]    -0.57537    0.32765  -1.756   0.0830 .
DatosD[(t - 1), "rw"]    0.06871    0.06891   0.997   0.3218  
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 0.6874 on 78 degrees of freedom
Multiple R-squared:  0.1656,    Adjusted R-squared:  0.1228 
F-statistic:  3.87 on 4 and 78 DF,  p-value: 0.0064
Code
summary(ma3)

Call:
lm(formula = DatosD[t, "U"] ~ -1 + DatosD[(t - 1), "e"] + DatosD[(t - 
    1), "prod"] + DatosD[(t - 1), "U"] + DatosD[(t - 1), "rw"])

Residuals:
     Min       1Q   Median       3Q      Max 
-0.80709 -0.11656  0.00891  0.19765  0.91726 

Coefficients:
                         Estimate Std. Error t value Pr(>|t|)    
DatosD[(t - 1), "e"]    -0.359606   0.094209  -3.817 0.000269 ***
DatosD[(t - 1), "prod"] -0.119164   0.050175  -2.375 0.020008 *  
DatosD[(t - 1), "U"]     0.001336   0.152416   0.009 0.993029    
DatosD[(t - 1), "rw"]    0.106337   0.032054   3.317 0.001382 ** 
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 0.3198 on 78 degrees of freedom
Multiple R-squared:  0.4844,    Adjusted R-squared:  0.4579 
F-statistic: 18.32 on 4 and 78 DF,  p-value: 1.203e-10
Code
summary(ma4)

Call:
lm(formula = DatosD[t, "rw"] ~ -1 + DatosD[(t - 1), "e"] + DatosD[(t - 
    1), "prod"] + DatosD[(t - 1), "U"] + DatosD[(t - 1), "rw"])

Residuals:
    Min      1Q  Median      3Q     Max 
-1.8290 -0.3895 -0.0143  0.8508  3.9638 

Coefficients:
                        Estimate Std. Error t value Pr(>|t|)    
DatosD[(t - 1), "e"]     0.82561    0.29064   2.841  0.00574 ** 
DatosD[(t - 1), "prod"] -0.22389    0.15479  -1.446  0.15208    
DatosD[(t - 1), "U"]     1.07252    0.47022   2.281  0.02528 *  
DatosD[(t - 1), "rw"]    0.53018    0.09889   5.361 8.17e-07 ***
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 0.9865 on 78 degrees of freedom
Multiple R-squared:  0.5519,    Adjusted R-squared:  0.5289 
F-statistic: 24.01 on 4 and 78 DF,  p-value: 5.726e-13

Estimación Directa

Modelo VAR (1)

Code
m1 = VAR(DatosD[,c("e","prod","U","rw")], p = 1, type = "none")
summary(m1)

VAR Estimation Results:
========================= 
Endogenous variables: e, prod, U, rw 
Deterministic variables: none 
Sample size: 82 
Log Likelihood: -228.741 
Roots of the characteristic polynomial:
0.8414 0.7058 0.0892 0.0892
Call:
VAR(y = DatosD[, c("e", "prod", "U", "rw")], p = 1, type = "none")


Estimation results for equation e: 
================================== 
e = e.l1 + prod.l1 + U.l1 + rw.l1 

        Estimate Std. Error t value Pr(>|t|)    
e.l1     0.87946    0.11545   7.617 5.16e-11 ***
prod.l1  0.20108    0.06149   3.270   0.0016 ** 
U.l1     0.26501    0.18679   1.419   0.1599    
rw.l1   -0.01211    0.03928  -0.308   0.7586    
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1


Residual standard error: 0.3919 on 78 degrees of freedom
Multiple R-Squared: 0.7036, Adjusted R-squared: 0.6884 
F-statistic: 46.28 on 4 and 78 DF,  p-value: < 2.2e-16 


Estimation results for equation prod: 
===================================== 
prod = e.l1 + prod.l1 + U.l1 + rw.l1 

        Estimate Std. Error t value Pr(>|t|)  
e.l1    -0.11538    0.20252  -0.570   0.5705  
prod.l1  0.28159    0.10786   2.611   0.0108 *
U.l1    -0.57537    0.32765  -1.756   0.0830 .
rw.l1    0.06871    0.06891   0.997   0.3218  
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1


Residual standard error: 0.6874 on 78 degrees of freedom
Multiple R-Squared: 0.1656, Adjusted R-squared: 0.1228 
F-statistic:  3.87 on 4 and 78 DF,  p-value: 0.0064 


Estimation results for equation U: 
================================== 
U = e.l1 + prod.l1 + U.l1 + rw.l1 

         Estimate Std. Error t value Pr(>|t|)    
e.l1    -0.359606   0.094209  -3.817 0.000269 ***
prod.l1 -0.119164   0.050175  -2.375 0.020008 *  
U.l1     0.001336   0.152416   0.009 0.993029    
rw.l1    0.106337   0.032054   3.317 0.001382 ** 
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1


Residual standard error: 0.3198 on 78 degrees of freedom
Multiple R-Squared: 0.4844, Adjusted R-squared: 0.4579 
F-statistic: 18.32 on 4 and 78 DF,  p-value: 1.203e-10 


Estimation results for equation rw: 
=================================== 
rw = e.l1 + prod.l1 + U.l1 + rw.l1 

        Estimate Std. Error t value Pr(>|t|)    
e.l1     0.82561    0.29064   2.841  0.00574 ** 
prod.l1 -0.22389    0.15479  -1.446  0.15208    
U.l1     1.07252    0.47022   2.281  0.02528 *  
rw.l1    0.53018    0.09889   5.361 8.17e-07 ***
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1


Residual standard error: 0.9865 on 78 degrees of freedom
Multiple R-Squared: 0.5519, Adjusted R-squared: 0.5289 
F-statistic: 24.01 on 4 and 78 DF,  p-value: 5.726e-13 



Covariance matrix of residuals:
            e      prod         U       rw
e     0.15220  0.015207 -0.079134 -0.02804
prod  0.01521  0.465760 -0.008735  0.03991
U    -0.07913 -0.008735  0.100794  0.09348
rw   -0.02804  0.039909  0.093481  0.94001

Correlation matrix of residuals:
            e     prod        U       rw
e     1.00000  0.05712 -0.63891 -0.07413
prod  0.05712  1.00000 -0.04032  0.06031
U    -0.63891 -0.04032  1.00000  0.30370
rw   -0.07413  0.06031  0.30370  1.00000

Modelo VAR(1) testeo

Code
plot(m1, names = "e")

Code
plot(m1, names = "prod")

Code
plot(m1, names = "U")

Code
plot(m1, names = "rw")

Test de fuerza

Code
ser=serial.test(m1, lags.pt = 16, type = "PT.asymptotic")
ser

    Portmanteau Test (asymptotic)

data:  Residuals of VAR object m1
Chi-squared = 212.45, df = 240, p-value = 0.8995
Code
nor=normality.test(m1)
nor$jb.mul
$JB

    JB-Test (multivariate)

data:  Residuals of VAR object m1
Chi-squared = 6.8343, df = 8, p-value = 0.5546


$Skewness

    Skewness only (multivariate)

data:  Residuals of VAR object m1
Chi-squared = 4.7921, df = 4, p-value = 0.3093


$Kurtosis

    Kurtosis only (multivariate)

data:  Residuals of VAR object m1
Chi-squared = 2.0422, df = 4, p-value = 0.728
Code
stability(m1)
$e

Empirical Fluctuation Process: OLS-based CUSUM test 

Call: efp(formula = formula, data = data, type = type, h = h, dynamic = dynamic, 
    rescale = rescale)


$prod

Empirical Fluctuation Process: OLS-based CUSUM test 

Call: efp(formula = formula, data = data, type = type, h = h, dynamic = dynamic, 
    rescale = rescale)


$U

Empirical Fluctuation Process: OLS-based CUSUM test 

Call: efp(formula = formula, data = data, type = type, h = h, dynamic = dynamic, 
    rescale = rescale)


$rw

Empirical Fluctuation Process: OLS-based CUSUM test 

Call: efp(formula = formula, data = data, type = type, h = h, dynamic = dynamic, 
    rescale = rescale)
Code
plot(stability(m1), nc = 2)

Code
arch1 = arch.test(m1, lags.multi = 6)
arch1$arch.mul

    ARCH (multivariate)

data:  Residuals of VAR object m1
Chi-squared = 637.69, df = 600, p-value = 0.1389
Code
plot(arch1, names = "rw")

Raices de Polinomio

Code
roots(m1)
[1] 0.84141809 0.70575263 0.08920161 0.08920161

Pronostico

Pronostico de variables

Impulso respuesta

Función Impulso-Respuesta

Code
irf.emp <- vars::irf(m1, impulse = "prod", response = "U",n.ahead = 20, boot = TRUE)
plot(irf.emp) 
  • Por un choque positivo de un 1% en la productividad el desempleo se reduce hasta en un 8% en el primer periodo

Descomposición de Varianza

Descomposición

Code
des.varmod <- vars::fevd(m1, n.ahead = 10)
plot(des.varmod, nc=2)

Gracias por su atención!!

Referencias

Hamilton, James D. 2020. Time Series Analysis. Princeton university press.
Lütkepohl, Helmut, and Markus Krätzig. 2004. Applied Time Series Econometrics. Cambridge university press.