Econometría II: SVAR

Departamento de Economía

Carlos A. Yanes G.

2023-11-07

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, readr, lmtest)
  • Es importante resaltar el paquete vars, ya que con ellos trabajaremos. Algo de lo que se muestra en está el libro de (Enders 2012) y desde luego el de (Hamilton 2020)

Preambulo

Investigación en Macroeconomía

  • En macroeconomía el estudio de variables tales como PIB y el Desempleo se traducen en indicadores de como va la economía de un país.
  • Anteriormente vimos los modelos VAR de forma reducida, a continuación la idea es mirar una parte concerniente a lo propiamente teórico y con esto estimar los choques contemporáneos.
  • Los modelos VAR estructurales o SVAR tienen una apreciación mas completa.
  • Blanchard y Quah plantearon un ejercicio para evaluar la dinámica del PIB y el desempleo a la luz del modelo teórico de oferta y demanda agregada. En este se establece que los choques de oferta generan efectos permanentes sobre el PIB, mientras que los choques de demanda, asociados al desempleo, solo afectan al producto de manera transitoria.

Datos

Investigación de Blanchard y Quah (1989)

  • \(PIB\) = Producto interno Bruto
  • \(NACIONAL\) = Tasa de desempleo
  • Frecuencia Anual
  • Enero de 1975

Datos a trabajar

Code 
svardata <- read_csv("svardata.csv")
head(svardata)
# A tibble: 6 × 3
  Date       NACIONAL    PIB
  <date>        <dbl>  <dbl>
1 1975-01-01    10.5  179472
2 1976-01-01    10.4  188118
3 1977-01-01     9.4  195921
4 1978-01-01     8.2  212502
5 1979-01-01     8.9  223941
6 1980-01-01     9.10 233119
Code 
pib=ts(svardata$PIB,start = 1975) # Serie de tiempo 1
td=ts(svardata$NACIONAL,start = 1975) # Serie de tiempo 2
Code 
par(mfrow=c(1,2))
dlpib=diff(log(pib)) # Diferencia logaritmica del PIB
plot(dlpib) # Gráfico
plot(td)

Adiciones

  • Tenemos un desempleo con dato (atípico), se sugiere corregir ya que se debe a una crisis como tal en 1999.
  • Solo se hace una regresión y se extrae el error o residuo.
  • Utilizamos esa nueva variable y tenemos un modelo VAR original reducido.
  • Eliminamos un dato de la tasa de desempleo (ya que tenemos la diferencia de la otra serie)
  • Seleccionamos el mejor modelo VAR
Code 
# Dummy de crisis: 1 para los años entre 1997-2003, 0 en otro caso
dum_crisis=ifelse(time(td)>=1997&time(td)<=2003,1,0)
td_nueva=resid(lm(td~dum_crisis)) # regresion de la TD en funcion de dummy
plot(td_nueva)

Code 
Y1 = na.omit(cbind(dlpib,td_nueva=td_nueva[-1]))
# Dummy 1999: 1 para año 1999, 0 en otro caso
dum_99=ifelse(time(td)==1999,1,0)[-1]
VARselect(y = Y1,exogen = dum_99,type = "const")
$selection
AIC(n)  HQ(n)  SC(n) FPE(n) 
     9      1      1      1 

$criteria
                 1             2             3             4             5
AIC(n) -7.32502122 -7.1179909957 -6.9697499784 -6.9522779543 -7.0103185784
HQ(n)  -7.20355872 -6.9357972518 -6.7268249866 -6.6486217145 -6.6459310906
SC(n)  -6.95858724 -6.5683400322 -6.2368820270 -6.0361930151 -5.9110166513
FPE(n)  0.00066058  0.0008177027  0.0009604595  0.0009987325  0.0009746776
                  6             7            8             9           10
AIC(n) -6.845572099 -7.1700928353 -7.101975566 -7.3592862646 -7.310976927
HQ(n)  -6.420453363 -6.6842428516 -6.555394335 -6.7519737850 -6.642933200
SC(n)  -5.563053184 -5.7043569325 -5.453022676 -5.5271163861 -5.295590061
FPE(n)  0.001207706  0.0009369616  0.001107064  0.0009813198  0.001245549
Code 
# Estimaci?n VAR en forma reducida
var1=VAR(y = Y1,exogen = cbind(dum99=dum_99),type = "const")
summary(var1)

VAR Estimation Results:
========================= 
Endogenous variables: dlpib, td_nueva 
Deterministic variables: const 
Sample size: 41 
Log Likelihood: 39.348 
Roots of the characteristic polynomial:
0.5226 0.168
Call:
VAR(y = Y1, type = "const", exogen = cbind(dum99 = dum_99))


Estimation results for equation dlpib: 
====================================== 
dlpib = dlpib.l1 + td_nueva.l1 + const + dum99 

             Estimate Std. Error t value Pr(>|t|)    
dlpib.l1     0.215900   0.139358   1.549 0.129836    
td_nueva.l1  0.002181   0.001620   1.346 0.186361    
const        0.030040   0.005972   5.030 1.28e-05 ***
dum99       -0.074340   0.018304  -4.061 0.000243 ***
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1


Residual standard error: 0.01753 on 37 degrees of freedom
Multiple R-Squared: 0.3954, Adjusted R-squared: 0.3463 
F-statistic: 8.064 on 3 and 37 DF,  p-value: 0.0002931 


Estimation results for equation td_nueva: 
========================================= 
td_nueva = dlpib.l1 + td_nueva.l1 + const + dum99 

            Estimate Std. Error t value Pr(>|t|)    
dlpib.l1      6.7386    11.3465   0.594 0.556200    
td_nueva.l1   0.4746     0.1319   3.599 0.000930 ***
const        -0.3997     0.4862  -0.822 0.416314    
dum99         5.4914     1.4903   3.685 0.000729 ***
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1


Residual standard error: 1.427 on 37 degrees of freedom
Multiple R-Squared: 0.4184, Adjusted R-squared: 0.3712 
F-statistic: 8.872 on 3 and 37 DF,  p-value: 0.0001466 



Covariance matrix of residuals:
              dlpib  td_nueva
dlpib     0.0003074 -0.002972
td_nueva -0.0029717  2.037629

Correlation matrix of residuals:
           dlpib td_nueva
dlpib     1.0000  -0.1187
td_nueva -0.1187   1.0000

Testeo de fuerza

Code 
Portmanteau.var = function(modelo,rezago.max){
  Estad=c()
  p.valor=c()
  Estad.ajus=c()
  p.valor.=c()
  for (i in 1:rezago.max) {
    Estad[i]=serial.test(x = modelo,lags.pt = i,type = "PT.asymptotic")$serial$statistic
    p.valor[i]=serial.test(x = modelo,lags.pt = i,type = "PT.asymptotic")$serial$p.value
    Estad.ajus[i]=serial.test(x = modelo,lags.pt = i,type = "PT.adjusted")$serial$statistic
    p.valor.[i]=serial.test(x = modelo,lags.pt = i,type = "PT.adjusted")$serial$p.value
  }
  resultado=cbind(Lag=1:rezago.max,Estad,p.valor,Estad.ajus,p.valor.)
  print("Prueba Portmanteau")
  return(resultado)
}

# LM Breush-Pagan

LMBreush.var = function(modelo,rezago.max){
  Estadistico=c()
  p.valor=c()
  for (i in 1:rezago.max) {
    Estadistico[i]=serial.test(x = modelo,lags.bg = i,type = "BG")$serial$statistic
    p.valor[i]=serial.test(x = modelo,lags.bg = i,type = "BG")$serial$p.value
  }
  resultado=cbind(Lag=1:rezago.max,Estadistico,p.valor)
  print("Prueba LM Breush Godfrey")
  return(resultado)
}  # crear función
Code 
Portmanteau.var(modelo = var1,rezago.max = 10)
[1] "Prueba Portmanteau"
      Lag      Estad   p.valor Estad.ajus  p.valor.
 [1,]   1  0.4211294 0.0000000  0.4316577 0.0000000
 [2,]   2  2.8371994 0.5854286  2.9716287 0.5625844
 [3,]   3  9.9673936 0.2673219 10.6647329 0.2214289
 [4,]   4 14.4722973 0.2715698 15.6566532 0.2074728
 [5,]   5 18.4865382 0.2961856 20.2284275 0.2101081
 [6,]   6 22.2570616 0.3267263 24.6453265 0.2153383
 [7,]   7 31.3745662 0.1432737 35.6399644 0.0594590
 [8,]   8 31.8450668 0.2808616 36.2245257 0.1369691
 [9,]   9 34.5207621 0.3482289 39.6527602 0.1656761
[10,]  10 34.9669221 0.5175780 40.2428429 0.2878904
Code 
LMBreush.var(modelo = var1,rezago.max = 10)
[1] "Prueba LM Breush Godfrey"
      Lag Estadistico   p.valor
 [1,]   1    1.116245 0.8916863
 [2,]   2    7.748019 0.4584635
 [3,]   3   16.632373 0.1639572
 [4,]   4   21.460207 0.1614913
 [5,]   5   22.893208 0.2940687
 [6,]   6   27.071559 0.3011286
 [7,]   7   34.854043 0.1741560
 [8,]   8   36.077934 0.2836010
 [9,]   9   39.633837 0.3111235
[10,]  10   43.472752 0.3257545
  • Siempre la clave es NO rechazar Ho. \[H_0: \varepsilon_t \; \sim (0, \sigma^2)\]
  • Eso lo tenemos con valores P-values mayores a 0.10 para ambas pruebas.
  • Recuerde que para la Portmanteau test se lee a partir del siguiente rezago del orden del modelo.

SVAR

Modelo SVAR

  • El modelo Estructural debe tener su matriz de restricción.
  • Esa matriz hay que crearla.
  • No debe olvidar el algebra de la diapositiva de VAR
Code 
var1$varresult %>%
    stargazer(type = "text", no.space = TRUE,
              column.labels = colnames(Y1),
              dep.var.labels.include = FALSE)

==========================================================
                                  Dependent variable:     
                              ----------------------------
                                  dlpib           td      
                                   (1)            (2)     
----------------------------------------------------------
dlpib.l1                          0.216          6.739    
                                 (0.139)       (11.347)   
td_nueva.l1                       0.002        0.475***   
                                 (0.002)        (0.132)   
const                            0.030***       -0.400    
                                 (0.006)        (0.486)   
dum99                           -0.074***      5.491***   
                                 (0.018)        (1.490)   
----------------------------------------------------------
Observations                        41            41      
R2                                0.395          0.418    
Adjusted R2                       0.346          0.371    
Residual Std. Error (df = 37)     0.018          1.427    
F Statistic (df = 3; 37)         8.064***      8.872***   
==========================================================
Note:                          *p<0.1; **p<0.05; ***p<0.01
Code 
B0 <- diag(2)
diag(B0) <- NA
B0[2, 1] <- NA
B0[2, 2] <- NA
B0
     [,1] [,2]
[1,]   NA    0
[2,]   NA   NA
Code 
mod_svar <- SVAR(var1, estmethod = "direct", Amat = B0, hessian = TRUE, method = "BFGS")
summary(mod_svar)

SVAR Estimation Results:
======================== 

Call:
SVAR(x = var1, estmethod = "direct", Amat = B0, hessian = TRUE, 
    method = "BFGS")

Type: A-model 
Sample size: 41 
Log Likelihood: 35.138 
Method: direct 
Number of iterations: 64 
Convergence code: 0 

Estimated A matrix:
          dlpib td_nueva
dlpib    57.037   0.0000
td_nueva -6.814   0.7055

Estimated standard errors for A matrix:
         dlpib td_nueva
dlpib    6.299  0.00000
td_nueva 8.940  0.07791

Estimated B matrix:
         dlpib td_nueva
dlpib        1        0
td_nueva     0        1

Covariance matrix of reduced form residuals (*100):
           dlpib td_nueva
dlpib    0.03074   0.2969
td_nueva 0.29690 203.7601

Impulso-Respuesta

IMR

Code 
par(mfrow=c(4,4), cex=.5, mar = c(4,4,2,1))
mod_svar %>% irf(n.ahead = 40) %>% plot(plot.type = "single")

Code 
par(mfrow = c(2,2), cex = .5, mar = c(4,4,2,1))
mod_svar %>% irf(n.ahead = 40, response = "td_nueva") %>% plot(plot.type = "single")

  • Mire que el desempleo sobre ella misma reduce su nivel hasta llegar a cero en los primeros periodos. Se prevee que si el choque de demanda aumenta en un 1% se reduzca en los siguientes periodos. El desempleo no es persistente para el caso colombiano.

  • Con respecto al PIB note que al principio un incremento de un 1% el desempleo para el periodo (2) llega a su máximo casi un 20% pero luego decae totalmente (tenemos un frente de largo plazo)

Varianza

Varianza

Code 
#### FEVD ####
par(mar = c(4,5,2,1))
mod_fevd <- mod_svar %>% fevd(n.ahead = 40) 
mod_fevd %>% plot(addbars = 3)

  • Mire como poco o casi nada tiene que ver por lo pronto una con la otra. Todas dependen de si mismas.
  • En los siguientes periodos se nota mas la necesidad de la otra. Se debe revizar otras formas o restricciones.

Gracias por su atención!!

Referencias

Enders, Walter. 2012. “Applied Econometric Time Series.” Privredna Kretanja i Ekonomska Politika 132: 93.
Hamilton, James D. 2020. Time Series Analysis. Princeton university press.

Universidad del Norte