Capítulo 5 Modelos Multivariados: VAR, Cointegración, ARDL y Filtro de Kalman
En este capítulo removeremos el supuesto de que el análisis es univariado, ya que introduciremos la posibilidad de que los procesos generadores de datos compartan información entre dos o más series. Como primera aproximación desarrollaremos el concepto de Causalidad de Granger. Mediante esta metodología discutiremos cuándo dos series se causan estadísticamente. Posteriormente, introduciremos una técnica más sofisticada conocida como la metodología de Vectores Autorregresivos (VAR), la cual es una generalización de los procesos Autorregresivos (AR) que analizamos en los primeros capítulos. Finalmente, introduciremos la técnica de cointegración y de rezagos distribuidos (ARDL) para los casos en que las series que analicemos sean procesos no estacionarios.
A partir de este punto, asumiremos que las series empleadas son estacionarias en sus primeras diferencias y solo nos preocuparemos por su estacionariedad en los casos particulares de Cointegración y los modelos ARDL.
5.1 Causalidad de Granger
Hasta ahora hemos supuesto que una serie puede ser explicada únicamente con la información contenida en ella misma. No obstante, en adelante trataremos de analizar el caso en el que buscamos determinar relaciones entre variables y cómo el comportamiento de una serie influye en las demás. Una de las relaciones más importantes entre variables es la de causalidad. En este caso analizaremos el procedimiento de Granger (1969), conocido como causalidad de Granger. En adelante asumiremos que las series involucradas son débilmente estacionarias.
La exposición de esta sección –las definiciones de causalidad y de causalidad instantánea en términos de varianzas de predicción, la clasificación de las relaciones causales posibles y la formulación de la prueba– sigue a Kirchgässner, Wolters y Hassler (2013, cap. 3).
Sean \(X\) y \(Y\) dos series débilmente estacionarias. Definamos a \(I_t\) un conjunto de toda la información disponible hasta el momento \(t\). Asimismo, digamos que \(\overline{X}_t\) y \(\overline{Y}_t\) son los conjuntos de toda la información disponible (actual y pasada) de \(X\) y \(Y\), respectivamente. Es decir: \[\begin{eqnarray*} \overline{X}_t & := & \{ X_t, X_{t-1}, X_{t-2}, \ldots \} \\ \overline{Y}_t & := & \{ Y_t, Y_{t-1}, Y_{t-2}, \ldots \} \\ I_t & := & \overline{X}_t \cup \overline{Y}_t \end{eqnarray*}\]
Adicionalmente, definamos \(\sigma^2(*)\) como la varianza del término de error estimado de una regresión dada. Dicho lo anterior, digamos que:
Existe Causalidad de Granger o \(X\) causa a \(Y\) si y solo si, una regresión lineal da como resultado que:
\[\begin{equation} \sigma^2 (Y_{t+1} | I_t) < \sigma^2 (Y_{t+1} | I_t \setminus \overline{X}_t) \end{equation}\]
Es decir, que la varianza del error de pronóstico de \(Y_{t+1}\) que se obtiene utilizando toda la información disponible es MENOR que la que se obtiene al excluir del conjunto de información la historia completa de \(X\). En palabras: el pasado de \(X\) aporta información útil para predecir a \(Y\) que no está contenida en el propio pasado de \(Y\).
Existe Causalidad de Granger Instantánea o \(X\) causa de forma instantánea a \(Y\) si y solo si, una regresión lineal da como resultado:
\[\begin{equation} \sigma^2 (Y_{t+1} | \{ I_t, X_{t+1} \}) < \sigma^2 (Y_{t+1} | I_t) \end{equation}\]
Nótese que aquí el conjunto de información se amplía con el valor contemporáneo de \(X\): la causalidad instantánea no dice que \(X\) anticipe a \(Y\), sino que ambas se mueven juntas dentro del mismo periodo una vez descontado su pasado. Es una noción simétrica –si \(X\) causa instantáneamente a \(Y\), entonces \(Y\) causa instantáneamente a \(X\)–, de modo que, a diferencia de la causalidad en el primer sentido, no permite hablar de una dirección. Esa simetría tiene una consecuencia práctica que conviene tener presente: la causalidad instantánea suele reflejar que ambas series responden a un mismo factor común no incluido en el análisis, o simplemente que la frecuencia de los datos es demasiado baja para separar el orden de los acontecimientos.
Ambas definiciones se aplican, por supuesto, intercambiando los papeles de \(X\) y \(Y\). Como la causalidad en el primer sentido sí tiene dirección, conviene organizar los resultados posibles preguntando por separado si cada serie causa a la otra. Eso da cuatro configuraciones, que resumimos en el Cuadro 5.1 con la notación habitual en la literatura.
| ¿\(X\) causa a \(Y\)? | ¿\(Y\) causa a \(X\)? | Configuración | Notación |
|---|---|---|---|
| No | No | Las series no se anticipan entre sí | \((X, Y)\) |
| Sí | No | \(X\) precede a \(Y\) | \((X \rightarrow Y)\) |
| No | Sí | \(Y\) precede a \(X\) | \((X \leftarrow Y)\) |
| Sí | Sí | Retroalimentación | \((X \leftrightarrow Y)\) |
A estas cuatro configuraciones hay que superponer la causalidad instantánea, que puede estar presente o ausente en cada una de ellas –por ejemplo, dos series pueden no anticiparse en absoluto y aun así moverse juntas dentro del mismo periodo, caso que se denota \((X - Y)\)–. Cruzando ambas dimensiones se obtiene la clasificación completa de ocho casos excluyentes que presentan Kirchgässner, Wolters y Hassler (2013, 98). En el trabajo aplicado, sin embargo, lo que suele reportarse es la prueba de precedencia temporal del Cuadro 5.1, que es la que desarrollamos a continuación.
Por lo anterior, representaremos mediante un \(AR(p)\) con variables exógenas lo siguiente: \[\begin{equation} A(L) \begin{bmatrix} Y_t \\ X_t \end{bmatrix} = \begin{bmatrix} a_{11}(L) & a_{12}(L) \\ a_{21}(L) & a_{22}(L) \end{bmatrix} \begin{bmatrix} Y_t \\ X_t \end{bmatrix} = \begin{bmatrix} V_t \\ U_t \end{bmatrix} \tag{5.1} \end{equation}\]
O en su versión \(MA(q)\) con variables exógenas: \[\begin{equation} \begin{bmatrix} Y_t \\ X_t \end{bmatrix} = B(L) \begin{bmatrix} V_t \\ U_t \end{bmatrix} = \begin{bmatrix} b_{11}(L) & b_{12}(L) \\ b_{21}(L) & b_{22}(L) \end{bmatrix} \begin{bmatrix} V_t \\ U_t \end{bmatrix} \end{equation}\]
Para determinar el test de causalidad utilizaremos una especificación similar a la de la ecuación (5.1). Para probar si \(X\) causa a \(Y\), consideraremos la siguiente regresión: \[\begin{equation} Y_t = \alpha_0 + \sum^{k_1}_{k = 1} a^k_{11} Y_{t-k} + \sum^{k_2}_{k = k_0} a^k_{12} X_{t-k} + U_{1,t} \end{equation}\]
Donde \(k_0 = 1\) y, en general, se asume que \(k_1 = k_2\). Asimismo, el valor de estas constantes se puede determinar con el criterio de Akaike (o cualquier otro criterio de información). No obstante, algunos autores sugieren que una buena práctica es considerar valores de \(k_1\) y \(k_2\) que recorran los valores 4, 8, 12 y 16.
Dicho lo anterior, el test de causalidad de Granger se establece con una prueba F, en la cual se prueba la siguiente hipótesis nula: \[\begin{equation} H_0: a^1_{12} = a^2_{12} = \ldots = a^{k2}_{12} = 0 \end{equation}\]
library(ggplot2)
library(dplyr)
library(stats)
library(MASS)
library(strucchange)
library(zoo)
library(sandwich)
library(urca)
library(lmtest)
library(vars)
#
load("BD/Datos_Ad.RData")
#
INPC <- ts(Datos_Ad$INPC_Ad,
start = c(2000, 1),
freq = 12)
DLINPC <- diff(log( ts(Datos_Ad$INPC_Ad, start = c(2000, 1), freq = 12) ))
TC <- ts(Datos_Ad$TC_Ad,
start = c(2000, 1),
freq = 12)
DLTC <- diff(log( ts(Datos_Ad$TC_Ad, start = c(2000, 1), freq = 12) ))
CETE28 <- ts(Datos_Ad$CETE28_Ad,
start = c(2000, 1),
freq = 12)
DLCETE28 <- diff(log( ts(Datos_Ad$CETE28_Ad, start = c(2000, 1), freq = 12) ))
## Fecha final de la muestra, calculada a partir de los datos para que el
## texto no dependa de la vintage de la base.
fin_num <- end(INPC)
meses_es <- c("enero", "febrero", "marzo", "abril", "mayo", "junio",
"julio", "agosto", "septiembre", "octubre", "noviembre",
"diciembre")
fin_txt <- paste(meses_es[fin_num[2]], "de", fin_num[1])Ejemplo. Consideremos como variables analizadas al Índice Nacional de Precios al Consumidor (\(INPC_t\)), al Tipo de Cambio (\(TDC_t\)) y al rendimiento anual de los Cetes a 28 días (\(CETE28_t\)), todas desestacionalizadas para el periodo de enero de 2000 a mayo de 2026. Dado que la metodología de Granger supone que las series son estacionarias, utilizaremos las diferencias logarítmicas de cada una de las tres series (es decir, utilizaremos una transformación del tipo \(ln(X_t) - ln(X_{t-1})\)). La Figura 5.1 muestra las series en su transformación de diferencias logarítmicas.
par(mfrow=c(3, 1))
plot(DLINPC, xlab = "Tiempo",
main = "Diferencias Logarítmicas del INPC",
col = "darkgreen")
plot(DLTC, xlab = "Tiempo",
main = "Diferencias Logarítmicas del Tipo de Cambio",
col = "darkblue")
plot(DLCETE28, xlab = "Tiempo",
main = "Diferencias Logarítmicas de los Cetes a 28 días",
col = "darkred")
Figura 5.1: Series en diferencias logarítmicas dadas por las siguientes expresiones: \(DLINPC_t = ln(INPC_t) - ln(INPC_{t-1})\), \(DLTC_t = ln(TC_t) - ln(TC_{t-1})\) y \(DLCETE28_t = ln(CETE28_t) - ln(CETE28_{t-1})\).
Por simplicidad, en el Cuadro 5.2 se muestra el resultado de aplicar el test de Granger a diferentes especificaciones, con rezagos 4, 8, 12 y 16, sólo para la serie de Tipo de Cambio en diferencias logarítmicas. En cada una de las pruebas se compara el modelo considerado como regresor a la variable que es candidata de causar, respecto del modelo sin considerar a dicha variable.
rezagos_g <- c(4, 8, 12, 16)
granger_rows <- lapply(rezagos_g, function(p) {
gt <- grangertest(DLTC ~ DLINPC, order = p)
pv <- gt$`Pr(>F)`[2]
data.frame(Rezagos = p,
Fstat = round(gt$F[2], 4),
Prob = round(pv, 5),
Signif = ifelse(pv < 0.001, "***",
ifelse(pv < 0.01, "**",
ifelse(pv < 0.05, "*", ""))))
})
# Valores p de la prueba del cuadro y de dos contrastes de referencia que se
# comentan en el texto (se calculan aquí para que la discusión no dependa de
# números escritos a mano).
pv_g <- vapply(granger_rows, function(d) d$Prob, numeric(1))
pv_tc_inpc <- grangertest(DLINPC ~ DLTC, order = 8)$`Pr(>F)`[2]
pv_inpc_ce <- grangertest(DLCETE28 ~ DLINPC, order = 4)$`Pr(>F)`[2]
knitr::kable(
do.call(rbind, granger_rows),
col.names = c("Rezagos", "Estadística $F$", "Probabilidad ($>F$)", "Signif."),
caption = "Prueba de si $DLINPC_t$ Granger causa a $DLTC_t$.",
align = c("c", "c", "c", "c"),
booktabs = TRUE,
escape = FALSE
)| Rezagos | Estadística \(F\) | Probabilidad (\(>F\)) | Signif. |
|---|---|---|---|
| 4 | 0.9404 | 0.44080 | |
| 8 | 0.6918 | 0.69877 | |
| 12 | 0.6752 | 0.77498 | |
| 16 | 0.6344 | 0.85484 |
De acuerdo con el Cuadro 5.2, con esta muestra no se rechaza la hipótesis nula en ninguna de las cuatro especificaciones: los valores \(p\) van de 0.441 a 0.855, muy por encima de cualquier nivel de significancia convencional. Es decir, la inflación no aporta información para predecir la tasa de depreciación cambiaria más allá de lo que ya aporta la propia historia del tipo de cambio.
Este resultado merece tres comentarios, porque ilustra bien los límites
de la prueba. Primero, un no rechazo no demuestra la ausencia de
relación: sólo indica que, con esta muestra y esta especificación, no hay
evidencia suficiente en contra de la hipótesis nula. Segundo, la prueba
es estrictamente bivariada, de modo que ignora la posibilidad de que la
relación entre inflación y tipo de cambio esté mediada por terceras
variables –como la tasa de interés–, que es justamente la limitación
que motiva el enfoque multivariado de la siguiente sección. Tercero, la
conclusión depende del par de variables que se examine: al aplicar la
misma prueba con cuatro rezagos, la inflación sí resulta causal en el
sentido de Granger sobre la tasa de los Cetes (valor \(p\) de
0.0151), y en el sentido inverso el tipo de cambio
sobre la inflación se acerca a la significancia con ocho rezagos (valor
\(p\) de 0.066). Por eso la práctica recomendable es
reportar todas las combinaciones y no sólo la que confirma la hipótesis
de interés; el resto de los pares se puede replicar con el código de
este capítulo cambiando las variables de grangertest().
5.2 Definición y representación del Sistema o Modelo de Vectores Autorregresivos (VAR(p))
En esta sección ampliaremos la discusión planteada en el apartado anterior. En el sentido de que en la sección pasada nuestra discusión se limitó al análisis de causalidad entre dos variables a la vez, que si bien es posible extenderlo a más variables, es un procedimiento limitado a casos particulares por las siguientes razones.
El procedimiento de causalidad de Granger supone que es posible identificar un sistema de ecuaciones que debe conformarse una vez que se ha identificado el sentido de la causalidad. Así, el proceso anterior necesita del conocimiento previo de las relaciones que existen entre las variables.
Adicionalmente, no resuelve el problema más general que está relacionado con cómo identificar la causalidad cuando se tienen múltiples variables con múltiples sentidos de causalidad. En esta sección analizaremos una mejor aproximación al problema de cómo identificar la causalidad múltiple. Por lo tanto, como mecanismo para solucionar el problema planteado, analizaremos el caso de un Sistema o Modelo de Vectores Autorregresivos conocido como VAR.
El primer supuesto del que partiremos es que existe algún grado de endogeneidad entre las variables consideradas en el análisis. Adicionalmente, el segundo supuesto que estableceremos es que requerimos que las variables que tengamos consideradas sean estacionarias.
Por lo anterior, diremos que un modelo de Vectores Autorregresivos (VAR) es un procedimiento fundado en el supuesto de que las variables consideradas son estacionarias. Así, hasta este momento hemos pasado de modelos univariados a modelos multivariados, pero no hemos podido dejar de asumir que las series son estacionarias.
En lo subsecuente asumiremos que las series empleadas son estacionarias y sólo lo demostraremos cuando, en su caso, sea necesario. Esto no significa que el lector deba asumir estacionariedad. Por el contrario, siempre debe probar que las series son estacionarias antes de iniciar la implementación de cualquier técnica de series de tiempo.
El desarrollo del modelo \(VAR(p)\) que sigue –su representación, las condiciones de estabilidad, las matrices de autocovarianzas y los criterios de información para elegir el orden– sigue a Lütkepohl (2005, caps. 2 y 4).
Ahora bien, iniciaremos con el establecimiento de la representación del proceso. Digamos que tenemos un proceso estocástico \(\mathbf{X}_t\) estacionario vectorial de dimensión \(k\): \[\begin{equation*} \mathbf{X}_t = \begin{bmatrix} X_{1t} \\ X_{2t} \\ \vdots \\ X_{kt} \end{bmatrix} \end{equation*}\]
Para cualquier \(i = 1, 2, \ldots, p\): \[\begin{equation*} \mathbf{X}_{t-i} = \begin{bmatrix} X_{1t-i} \\ X_{2t-i} \\ \vdots \\ X_{kt-i} \end{bmatrix} \end{equation*}\]
Donde cada \(X_{kt}\) en \(\mathbf{X}_t\) es una serie de tiempo por sí misma. De esta forma, la expresión reducida del modelo o el proceso \(VAR(p)\) estará dado por: \[\begin{equation} \mathbf{X}_t = \boldsymbol{\delta} + A_1 \mathbf{X}_{t-1} + A_2 \mathbf{X}_{t-2} + \ldots + A_p \mathbf{X}_{t-p} + \mathbf{U}_{t} \tag{5.2} \end{equation}\]
Donde cada uno de las \(A_i\), \(i = 1, 2, \ldots, p\), son matrices cuadradas de dimensión \(k\) y \(\mathbf{U}_t\) representa un vector de dimensión \(k \times 1\) con los residuales en el momento del tiempo \(t\) que son, por individual, un proceso puramente aleatorio. También se incorpora un vector de términos constantes denominado como \(\mathbf{\delta}\), el cual es de dimensión \(k \times 1\).
Así, la ecuación (5.2) supone la siguiente estructura del vector \(\boldsymbol{\delta}\): \[\begin{equation*} \boldsymbol{\delta} = \begin{bmatrix} \delta_{1} \\ \delta_{2} \\ \vdots \\ \delta_{k} \end{bmatrix} \end{equation*}\]
También, la ecuación (5.2) supone que cada matriz \(A_i\), \(i = 1, 2, \ldots, p\) está definida de la siguiente forma: \[\begin{equation*} \mathbf{A}_i = \begin{bmatrix} a^{(i)}_{11} & a^{(i)}_{12} & \ldots & a^{(i)}_{1k} \\ a^{(i)}_{21} & a^{(i)}_{22} & \ldots & a^{(i)}_{2k} \\ \vdots & \vdots & \ddots & \vdots \\ a^{(i)}_{k1} & a^{(i)}_{k2} & \ldots & a^{(i)}_{kk} \end{bmatrix} \end{equation*}\]
Donde \(i = 1, 2, \ldots, p\).
Retomando la ecuación (5.2) y considerando que podemos ocupar el operador rezago \(L^j\) de forma análoga al caso del modelo \(AR(p)\), pero aplicado a un vector, tenemos las siguientes ecuaciones: \[\begin{eqnarray} \mathbf{X}_t - A_1 \mathbf{X}_{t-1} - A_2 \mathbf{X}_{t-2} - \ldots - A_p \mathbf{X}_{t-p} & = & \boldsymbol{\delta} + \mathbf{U}_{t} \nonumber \\ \mathbf{X}_t - A_1 L \mathbf{X}_{t} - A_2 L^2 \mathbf{X}_{t} - \ldots - A_p L^p \mathbf{X}_{t} & = & \boldsymbol{\delta} + \mathbf{U}_{t} \nonumber \\ (I_k - \mathbf{A_1} L - \mathbf{A_2} L^2 - \ldots - \mathbf{A_p} L^p) \mathbf{X}_t & = & \boldsymbol{\delta} + \mathbf{U}_{t} \nonumber \\ \mathbf{A}(L) \mathbf{X}_t & = & \boldsymbol{\delta} + \mathbf{U}_{t} \tag{5.3} \end{eqnarray}\]
Adicionalmente, requeriremos que dado que \(\mathbf{U}_t\) es un proceso puramente aleatorio, este debe cumplir con las siguientes condiciones:
El valor esperado del término de error es cero: \[\begin{equation} \mathbb{E}[\mathbf{U}_t] = 0 \end{equation}\]
Existe una matriz de varianzas y covarianzas entre los términos de error contemporáneos dada por: \[\begin{eqnarray} \mathbb{E}[\mathbf{U}_t \mathbf{U}_t'] & = & \mathbb{E} \left[ \begin{bmatrix} U^{(t)}_{1} \\ U^{(t)}_{2} \\ \vdots \\ U^{(t)}_{k} \end{bmatrix} \begin{bmatrix} U^{(t)}_{1} & U^{(t)}_{2} & \ldots & U^{(t)}_{k} \end{bmatrix} \right] \nonumber \\ & = & \mathbb{E} \begin{bmatrix} U^{(t)}_{1} U^{(t)}_{1} & U^{(t)}_{1} U^{(t)}_{2} & \ldots & U^{(t)}_{1} U^{(t)}_{k} \\ U^{(t)}_{2} U^{(t)}_{1} & U^{(t)}_{2} U^{(t)}_{2} & \ldots & U^{(t)}_{2} U^{(t)}_{k} \\ \vdots & \vdots & \ldots & \vdots \\ U^{(t)}_{k} U^{(t)}_{1} & U^{(t)}_{k} U^{(t)}_{2} & \ldots & U^{(t)}_{k} U^{(t)}_{k} \end{bmatrix} \nonumber \\ & = & \begin{bmatrix} \sigma^2_1 & \sigma_{12} & \ldots & \sigma_{1k} \\ \sigma_{21} & \sigma^2_2 & \ldots & \sigma_{2k} \\ \vdots & \vdots & \ldots & \vdots \\ \sigma_{k1} & \sigma_{k2} & \ldots & \sigma^2_k \end{bmatrix} \nonumber \\ & = & \mathbf{\Sigma}_{UU} \tag{5.4} \end{eqnarray}\]
La matriz de varianzas y covarianzas no contemporáneas es nula. Es decir, que para todo \(t \neq s\): \[\begin{eqnarray} \mathbb{E} [\mathbf{U}_t \mathbf{U}_s'] & = & \mathbb{E} \left[ \begin{bmatrix} U^{(t)}_{1} \\ U^{(t)}_{2} \\ \vdots \\ U^{(t)}_{k} \end{bmatrix} \begin{bmatrix} U^{(s)}_{1} & U^{(s)}_{2} & \ldots & U^{(s)}_{k} \end{bmatrix} \right] \nonumber \\ & = & \mathbb{E} \begin{bmatrix} U^{(t)}_{1} U^{(s)}_{1} & U^{(t)}_{1} U^{(s)}_{2} & \ldots & U^{(t)}_{1} U^{(s)}_{k} \\ U^{(t)}_{2} U^{(s)}_{1} & U^{(t)}_{2} U^{(s)}_{2} & \ldots & U^{(t)}_{2} U^{(s)}_{k} \\ \vdots & \vdots & \ldots & \vdots \\ U^{(t)}_{k} U^{(s)}_{1} & U^{(t)}_{k} U^{(s)}_{2} & \ldots & U^{(t)}_{k} U^{(s)}_{k} \end{bmatrix} \nonumber \\ & = & \mathbf{0} \tag{5.5} \end{eqnarray}\]
Las ecuaciones (5.4) y (5.5) significan que los residuales \(\mathbf{U}_t\) pueden estar correlacionados entre ellos solo en el caso de que la información sea contemporánea, pero no tienen información en común entre residuales de otros periodos.
Al igual que en el caso del modelo o especificación \(AR(p)\), en la especificación del modelo \(VAR(p)\) existen condiciones de estabilidad. El proceso es estable si el determinante del polinomio matricial asociado a \(\mathbf{A}(L)\) en la ecuación (5.3) no se anula sobre el círculo unitario ni en su interior: \[\begin{equation} Det[I_k - A_1 z - A_2 z^2 - \ldots - A_p z^p] \neq 0 \quad \text{para todo} \quad |z| \leq 1 \end{equation}\]
Es decir, todas las raíces del polinomio determinante deben ubicarse fuera del círculo unitario o, equivalentemente, todos los eigenvalores de la matriz compañera del sistema deben tener módulo menor que 1.
La ecuación (5.3) puede ser reescrita en una forma similar a un proceso de MA. Al respecto, de forma similar a la siguiente ecuación podemos construir un modelo \(VARMA(p,q)\), el cual no estudiamos en este libro.
Retomando el primer planteamiento, podemos escribir: \[\begin{eqnarray} \mathbf{X}_t & = & \mathbf{A}^{-1}(L) \boldsymbol{\delta} + \mathbf{A}^{-1}(L) \mathbf{U}_t \nonumber \\ & = & \boldsymbol{\mu} + \boldsymbol{\beta}(L) \mathbf{U}_t \tag{5.6} \end{eqnarray}\]
Donde \(\boldsymbol{\mu}\) es un vector de \(k \times 1\) constantes y \(\boldsymbol{\beta}(L)\) es una matriz que depende de \(L\).
Por el lado de las matrices que representan la autocovarianza, éstas resultan de resolver lo siguiente: \[\begin{equation} \Gamma_X(\tau) = E[(\mathbf{X}_t - \mu)(\mathbf{X}_{t-\tau} - \mu)'] \end{equation}\]
Ahora, sin pérdida de generalidad digamos que la especificación VAR(p) en la ecuación (5.2) no tiene constante, por lo que \(\boldsymbol{\delta} = \mathbf{0}\), lo que implica que \(\boldsymbol{\mu} = \mathbf{0}\). De esta forma las matrices de autocovarianza resultan de: \[\begin{eqnarray*} \Gamma_{\mathbf{X}}(\tau) & = & E[(\mathbf{X}_t)(\mathbf{X}_{t-\tau})'] \\ & = & \mathbf{A_1} E[(\mathbf{X}_{t-1})(\mathbf{X}_{t-\tau})'] + \mathbf{A_2} E[(\mathbf{X}_{t-2})(\mathbf{X}_{t-\tau})'] \\ & & + \ldots + \mathbf{A_p} E[(\mathbf{X}_{t-p})(\mathbf{X}_{t-\tau})'] + E[(\mathbf{U}_t(\mathbf{X}_{t-\tau})'] \end{eqnarray*}\]
Finalmente, al igual que en el caso \(AR(p)\), requerimos de una métrica que nos permita determinar el número de rezagos óptimo \(p\) en el \(VAR(p)\). Así, establecemos criterios de información similares a los del \(AR(p)\) dados por:
1.Final Prediction Error (FPE): \[\begin{equation} FPE(p) = \left[ \frac{T + kp + 1}{T - kp - 1} \right]^k |\mathbf{\Sigma}_{\hat{U}\hat{U}}(p)| \end{equation}\]
Akaike Criterion (AIC): \[\begin{equation} AIC(p) = ln|\mathbf{\Sigma}_{\hat{U}\hat{U}}(p)| + (k + p k^2) \frac{2}{T} \end{equation}\]
Hannan - Quinn Criterion (HQ): \[\begin{equation} HQ(p) = ln|\mathbf{\Sigma}_{\hat{U}\hat{U}}(p)| + (k + p k^2) \frac{2ln(ln(T))}{T} \end{equation}\]
Schwartz Criterion (SC): \[\begin{equation} SC(p) = ln|\mathbf{\Sigma}_{\hat{U}\hat{U}}(p)| + (k + p k^2) \frac{ln(T)}{T} \end{equation}\]
Donde la matriz de varianzas y covarianzas contemporáneas estará dada por: \[\begin{equation*} \mathbf{\Sigma}_{\hat{U}\hat{U}}(p) = \mathbb{E} \left[ \begin{bmatrix} U^{(t)}_{1} \\ U^{(t)}_{2} \\ \vdots \\ U^{(t)}_{k} \end{bmatrix} \begin{bmatrix} U^{(t)}_{1} & U^{(t)}_{2} & \ldots & U^{(t)}_{k} \end{bmatrix} \right] \end{equation*}\]
Ejemplo. Ahora veamos un ejemplo de estimación de \(VAR(p)\). Para el ejemplo utilizaremos las series de INPC, Tipo de Cambio, rendimiento de los Cetes a 28 días, el IGAE y el Índice de Producción Industrial de los Estados Unidos, todas desestacionalizadas y para el período de enero de 2000 a mayo de 2026. Dado que el supuesto de estacionariedad sigue presente en nuestro análisis, emplearemos cada una de las series en su versión de diferencias logarítmicas. Las Figuras 5.2 y 5.3 muestran las series referidas.
library(ggplot2)
library(dplyr)
library(stats)
library(MASS)
library(strucchange)
library(zoo)
library(sandwich)
library(urca)
library(lmtest)
library(vars)
#
load("BD/Datos_Ad.RData")
#
DLINPC <- diff(log( ts(Datos_Ad$INPC_Ad, start = c(2000, 1), freq = 12) ))
DLTC <- diff(log( ts(Datos_Ad$TC_Ad, start = c(2000, 1), freq = 12) ))
DLCETE28 <- diff(log( ts(Datos_Ad$CETE28_Ad, start = c(2000, 1), freq = 12) ))
DLIGAE <- diff(log( ts(Datos_Ad$IGAE_Ad, start = c(2000, 1), freq = 12) ))
DLIPI <- diff(log( ts(Datos_Ad$IPI_Ad, start = c(2000, 1), freq = 12) ))
Datos <- data.frame(cbind(DLINPC, DLTC, DLCETE28, DLIGAE, DLIPI))
Datos <- ts(Datos,
start = c(2000, 2), freq = 12)plot(Datos, plot.type = "s",
col = c("darkgreen", "darkblue", "darkred", "black", "purple"),
main = "Series en Diferencias logarítmicas",
xlab = "Tiempo", ylab = "Variacion")
legend("bottomright", c("INPC", "TC", "CETES28", "IGAE", "IPI"),
cex = 0.6, lty = 1:1,
col = c("darkgreen", "darkblue", "darkred", "black", "purple"))
Figura 5.2: Series en diferencias logarítmicas (Forma 1)
plot(Datos, plot.type = "m",
col = "darkgreen",
main = "Series en Diferencias logarítmicas", xlab = "Tiempo")
Figura 5.3: Series en diferencias logarítmicas (Forma 2)
Dicho lo anterior, a continuación mostraremos la tabla que resume el valor de los distintos criterios de información para una especificación de un \(VAR(p)\) con constante. Nótese que es posible especificar un \(VAR(p)\) con tendencia, siempre que exista evidencia de que algunas de las series sean estacionarias alrededor de una tendencia. Caso que no aplica hasta este momento, ya que nuestro análisis de estacionariedad es claro respecto a la media constante (más adelante aportaremos la evidencia de esto), lo cual elimina la posibilidad de incluir una tendencia.
En el Cuadro 5.3 reportamos el número de rezagos propuesto a partir de cada criterio de información y en el Cuadro 5.4 reportamos los resultados de aplicar una prueba de criterios de información para diferentes valores de rezagos. Del cual se concluye que el número óptimo de rezagos es 2 (según el criterio AIC y el FPE) y 1 (según el criterio HQ y el SC). Cuando los criterios no coinciden –lo que es frecuente, porque el AIC penaliza menos la inclusión de parámetros y por ello tiende a seleccionar órdenes mayores– la práctica habitual en el análisis de modelos VAR es privilegiar el AIC o el FPE: subparametrizar la dinámica deja autocorrelación en los residuales y distorsiona las funciones de impulso-respuesta, mientras que un rezago adicional sólo cuesta eficiencia. La decisión debe validarse después con las pruebas de diagnóstico de los residuales.
vs <- vars::VARselect(Datos, lag.max = 12, type = "const")
knitr::kable(
as.data.frame(t(vs$selection)),
col.names = c("AIC", "HQ", "SC", "FPE"),
caption = paste("Número de rezagos según cada criterio de información para",
"modelos VAR($p$) con constante: series $DLINPC_t$, $DLTC_t$,",
"$DLCETE28_t$, $DLIGAE_t$ y $DLIPI_t$."),
align = rep("c", 4),
booktabs = TRUE
)| AIC | HQ | SC | FPE |
|---|---|---|---|
| 2 | 1 | 1 | 2 |
crit_df <- as.data.frame(t(vs$criteria))
crit_df <- cbind(Rezagos = seq_len(nrow(crit_df)), crit_df)
knitr::kable(
crit_df,
col.names = c("Rezagos", "AIC", "HQ", "SC", "FPE"),
caption = paste("Criterios de información para modelos VAR($p$) con constante:",
"$DLINPC_t$, $DLTC_t$, $DLCETE28_t$, $DLIGAE_t$ y $DLIPI_t$."),
align = c("c", "r", "r", "r", "r"),
digits = 4,
booktabs = TRUE
)| Rezagos | AIC | HQ | SC | FPE |
|---|---|---|---|---|
| 1 | -44.1233 | -43.9766 | -43.7565 | 0 |
| 2 | -44.2162 | -43.9472 | -43.5437 | 0 |
| 3 | -44.1813 | -43.7900 | -43.2032 | 0 |
| 4 | -44.1290 | -43.6154 | -42.8451 | 0 |
| 5 | -44.0397 | -43.4038 | -42.4502 | 0 |
| 6 | -43.9921 | -43.2339 | -42.0969 | 0 |
| 7 | -44.0111 | -43.1307 | -41.8102 | 0 |
| 8 | -43.9984 | -42.9957 | -41.4919 | 0 |
| 9 | -43.9569 | -42.8320 | -41.1447 | 0 |
| 10 | -43.9477 | -42.7004 | -40.8298 | 0 |
| 11 | -43.8614 | -42.4919 | -40.4378 | 0 |
| 12 | -43.8141 | -42.3223 | -40.0849 | 0 |
De esta forma, justificamos la estimación de un \(VAR(2)\). Los resultados del mismo se reportan en los siguientes cuadros, en los que se muestra el resultado de una de las ecuaciones. Los resultados restantes se encuentran en el código de R mostrado más abajo. Primero mostraremos los resultados de las raíces del polinomio característico en el Cuadro 5.5, seguido de un cuadro para la ecuación del IGAE en el Cuadro 5.6 (por simplicidad se omiten las otras cuatro ecuaciones del VAR(2)), y del Cuadro 5.7 con la matriz \(\mathbf{\Sigma}_{\hat{U}\hat{U}}\) estimada del VAR.
# Se llama a la función con el prefijo del paquete (vars::) porque otros
# paquetes que se cargan más adelante en el libro --por ejemplo MTS, en el
# capítulo de volatilidad-- también exportan una función llamada VAR() con
# argumentos distintos, y el caché de knitr puede adjuntarlos antes de
# ejecutar este bloque.
VAR_p <- vars::VAR(Datos, p = 2, type = "const")
# La salida completa de summary(VAR_p) reporta las cinco ecuaciones del
# sistema; en el texto se muestran únicamente los cuadros que se discuten.
# Para revisarla, ejecute la siguiente línea en su sesión de R:
# summary(VAR_p)rts <- roots(VAR_p)
rts_mat <- matrix(round(rts, 4), nrow = 2, byrow = TRUE)
knitr::kable(
rts_mat,
col.names = rep("", ncol(rts_mat)),
caption = "Módulos de las raíces del polinomio característico del VAR(2).",
booktabs = TRUE
)| 0.5208 | 0.5208 | 0.5114 | 0.4713 | 0.4169 |
| 0.3541 | 0.3541 | 0.1639 | 0.0946 | 0.0946 |
coef_IGAE <- summary(VAR_p)$varresult$DLIGAE$coefficients
coef_df <- as.data.frame(coef_IGAE)
lbl_map <- c(
DLINPC.l1 = "$DLINPC_{t-1}$", DLTC.l1 = "$DLTC_{t-1}$",
DLCETE28.l1 = "$DLCETE28_{t-1}$", DLIGAE.l1 = "$DLIGAE_{t-1}$",
DLIPI.l1 = "$DLIPI_{t-1}$", DLINPC.l2 = "$DLINPC_{t-2}$",
DLTC.l2 = "$DLTC_{t-2}$", DLCETE28.l2 = "$DLCETE28_{t-2}$",
DLIGAE.l2 = "$DLIGAE_{t-2}$", DLIPI.l2 = "$DLIPI_{t-2}$",
const = "$\\delta_4$"
)
rownames(coef_df) <- lbl_map[rownames(coef_df)]
coef_df$Signif <- ifelse(coef_df[, 4] < 0.001, "***",
ifelse(coef_df[, 4] < 0.01, "**",
ifelse(coef_df[, 4] < 0.05, "*",
ifelse(coef_df[, 4] < 0.1, ".", ""))))
knitr::kable(
coef_df,
col.names = c("Coef.", "Error Est.", "Estad. $t$", "Prob. $(>|t|)$", "Signif."),
caption = "Ecuación $DLIGAE_t$ estimada en el VAR(2).",
align = c("r", "r", "r", "r", "c"),
digits = 4,
booktabs = TRUE,
escape = FALSE
)| Coef. | Error Est. | Estad. \(t\) | Prob. \((>|t|)\) | Signif. | |
|---|---|---|---|---|---|
| \(DLINPC_{t-1}\) | 0.7285 | 0.3534 | 2.0614 | 0.0401 |
|
| \(DLTC_{t-1}\) | -0.1788 | 0.0296 | -6.0303 | 0.0000 | *** |
| \(DLCETE28_{t-1}\) | 0.0086 | 0.0148 | 0.5823 | 0.5608 | |
| \(DLIGAE_{t-1}\) | -0.0144 | 0.0876 | -0.1640 | 0.8698 | |
| \(DLIPI_{t-1}\) | 0.4252 | 0.1127 | 3.7741 | 0.0002 | *** |
| \(DLINPC_{t-2}\) | -0.7362 | 0.3560 | -2.0678 | 0.0395 |
|
| \(DLTC_{t-2}\) | 0.0434 | 0.0308 | 1.4075 | 0.1603 | |
| \(DLCETE28_{t-2}\) | -0.0050 | 0.0147 | -0.3390 | 0.7349 | |
| \(DLIGAE_{t-2}\) | -0.1606 | 0.0834 | -1.9268 | 0.0549 | . |
| \(DLIPI_{t-2}\) | -0.1025 | 0.1143 | -0.8966 | 0.3707 | |
| \(\delta_4\) | 0.0016 | 0.0017 | 0.9377 | 0.3491 |
sig_mat <- summary(VAR_p)$covres
vnames <- c("$DLINPC$", "$DLTC$", "$DLCETE28$", "$DLIGAE$", "$DLIPI$")
colnames(sig_mat) <- vnames
rownames(sig_mat) <- vnames
knitr::kable(
sig_mat,
caption = "Matriz $\\mathbf{\\Sigma}_{\\hat{U}\\hat{U}}$ estimada del VAR(2).",
digits = 4,
booktabs = TRUE,
escape = FALSE,
format.args = list(scientific = TRUE, digits = 3)
)| \(DLINPC\) | \(DLTC\) | \(DLCETE28\) | \(DLIGAE\) | \(DLIPI\) | |
|---|---|---|---|---|---|
| \(DLINPC\) | 0e+00 | 0e+00 | 0.0e+00 | 0e+00 | 0e+00 |
| \(DLTC\) | 0e+00 | 7e-04 | 2.0e-04 | 0e+00 | 0e+00 |
| \(DLCETE28\) | 0e+00 | 2e-04 | 2.4e-03 | 0e+00 | 1e-04 |
| \(DLIGAE\) | 0e+00 | 0e+00 | 0.0e+00 | 2e-04 | 1e-04 |
| \(DLIPI\) | 0e+00 | 0e+00 | 1.0e-04 | 1e-04 | 1e-04 |
Finalmente, en el Cuadro 5.8 reportamos las pruebas de diagnóstico del VAR(2): correlación serial de los residuales (Breusch- Godfrey multivariada), normalidad (Jarque-Bera multivariada) y efectos ARCH. Los tres contrastes rechazan su hipótesis nula. La correlación serial residual indica que dos rezagos no bastan para agotar la dinámica del sistema –lo que sugiere estimar un orden mayor al que eligieron el AIC y el FPE, o revisar la inclusión de variables relevantes–; el rechazo de la normalidad se explica en buena medida por los valores atípicos de 2008-2009 y de 2020, que podrían tratarse con variables dicotómicas; y el rechazo de la homocedasticidad apunta a la presencia de heterocedasticidad condicional, es decir, al tipo de estructura que modelaremos explícitamente en el Capítulo 6. Vale la pena insistir en la lectura de conjunto: los diagnósticos no invalidan el ejercicio, pero sí advierten que los errores estándar y las bandas de confianza de las funciones de impulso-respuesta deben tomarse con reserva.
diag_s2 <- serial.test(VAR_p, lags.bg = 2, type = "BG")
diag_s4 <- serial.test(VAR_p, lags.bg = 4, type = "BG")
diag_s6 <- serial.test(VAR_p, lags.bg = 6, type = "BG")
diag_n <- normality.test(VAR_p)
diag_a <- arch.test(VAR_p, lags.multi = 6)
diagn_fila <- function(nombre, stat, pval) {
data.frame(
Prueba = nombre,
Estadistico = round(stat, 3),
pvalor = round(pval, 4),
Conclusion = ifelse(pval < 0.05, "Se rechaza", "No se rechaza")
)
}
knitr::kable(
rbind(
diagn_fila("Corr. serial — $\\chi^2(2)$",
as.numeric(diag_s2$serial$statistic), diag_s2$serial$p.value),
diagn_fila("Corr. serial — $\\chi^2(4)$",
as.numeric(diag_s4$serial$statistic), diag_s4$serial$p.value),
diagn_fila("Corr. serial — $\\chi^2(6)$",
as.numeric(diag_s6$serial$statistic), diag_s6$serial$p.value),
diagn_fila("Normalidad JB — $\\chi^2$",
as.numeric(diag_n$jb.mul$JB$statistic), diag_n$jb.mul$JB$p.value),
diagn_fila("ARCH multivariado — $\\chi^2$",
as.numeric(diag_a$arch.mul$statistic), diag_a$arch.mul$p.value)
),
col.names = c("Prueba", "Estadístico", "p-valor", "Conclusión ($H_0$)"),
caption = "Pruebas de diagnóstico sobre los residuales del VAR(2).",
align = c("l", "r", "r", "l"),
booktabs = TRUE,
escape = FALSE,
row.names = FALSE
)| Prueba | Estadístico | p-valor | Conclusión (\(H_0\)) |
|---|---|---|---|
| Corr. serial — \(\chi^2(2)\) | 71.226 | 0.0259 | Se rechaza |
| Corr. serial — \(\chi^2(4)\) | 130.790 | 0.0210 | Se rechaza |
| Corr. serial — \(\chi^2(6)\) | 226.172 | 0.0001 | Se rechaza |
| Normalidad JB — \(\chi^2\) | 58055.139 | 0.0000 | Se rechaza |
| ARCH multivariado — \(\chi^2\) | 2195.659 | 0.0000 | Se rechaza |
5.3 Análisis de Impulso-Respuesta
Una de las grandes ventajas que aporta el análisis de los modelos VAR es el análisis de Impulso-Respuesta. Dicho análisis busca cuantificar el efecto que tiene en \(\mathbf{X}_t\) una innovación o cambio en los residuales de cualquiera de las variables en un momento definido. Partamos de la ecuación (5.6) de forma que tenemos: \[\begin{eqnarray} \mathbf{X}_t & = & \mathbf{A}^{-1}(L) \delta + \mathbf{A}^{-1}(L) \mathbf{U}_t \nonumber \\ & = & \mu + \mathbf{B}(L) \mathbf{U}_t \nonumber \\ & = & \mu + \Psi_0 \mathbf{U}_t + \Psi_1 \mathbf{U}_{t-1} + \Psi_2 \mathbf{U}_{t-2} + \Psi_3 \mathbf{U}_{t-3} + \ldots \end{eqnarray}\]
Donde \(\Psi_0 = I\) y cada una de las \(\Psi_i = \mathbf{B}_i\), \(i = 1, 2, \ldots\). De esta forma se verifica el efecto que tiene en \(\mathbf{X}_t\) cada una de las innovaciones pasadas. Por lo que el análisis de Impulso-Respuesta cuantifica el efecto de cada una de esas matrices en las que hemos descompuesto a \(\mathbf{B}(L)\).
Ejemplo. Retomando el modelo \(VAR(2)\) anteriormente estimado, en las siguientes figuras reportamos las gráficas de Impulso-Respuesta de la serie \(DLTC_t\) ante cambios en los residuales del resto de las series y de la propia serie.
IR_DLINPC <- irf(VAR_p, n.ahead = 12, boot = TRUE,
ci = 0.95, response = "DLINPC")
# Las respuestas puntuales y sus bandas se pueden inspeccionar con
# print(IR_DLINPC) o graficar con plot(IR_DLINPC).
Figura 5.4: Impulso - Respuesta en \(DLTC_t\)
Figura 5.5: Impulso - Respuesta en \(DLTC_t\)
Figura 5.6: Impulso - Respuesta en \(DLTC_t\)
Figura 5.7: Impulso - Respuesta en \(DLTC_t\)
Figura 5.8: Impulso - Respuesta en \(DLTC_t\)
Los resultados muestran que la respuesta de \(DLTC_t\) ante impulsos en los términos de error fue estadísticamente significativa sólo para algunos de los casos y en periodos cortos de tiempo. El resto de los resultados de Impulso-Respuesta se puede replicar con el código disponible en el repositorio de GitHub del proyecto.
5.4 Identificación de los Choques Estructurales en los Modelos VAR
En su trabajo seminal “Macroeconomics and Reality”, Sims (1980) cuestionó los modelos macroeconométricos utilizados en ese momento. Antes de la aparición de los modelos de referencia modernos, como los Vectores Autorregresivos (VAR) y los modelos estocásticos dinámicos de equilibrio general (DSGE), se estimaban modelos de ecuaciones simultáneas derivados de la escuela keynesiana, especificados con \(n\) cantidad de ecuaciones. Por lo tanto, para entender el surgimiento de los modelos VAR en la macroeconometría, es necesario partir del contexto en el que se escribió el trabajo seminal de Sims.
A finales de los años 70 e inicios de los 80, comenzó una crisis en torno a los modelos estructurales tradicionales. Los modelos macroeconómicos de gran escala –como los de la Cowles Commission– utilizaban restricciones teóricas para identificar relaciones estructurales que, en muchos casos, eran impuestas de manera arbitraria. Estos modelos estaban altamente parametrizados y basados en supuestos teóricos fuertes –como las expectativas adaptativas–, lo que los hacía sensibles a errores específicos. De esta situación surge en parte la Crítica de Lucas, en la cual Robert Lucas (1976) mostró que este tipo de modelos fallaba al capturar cambios en el comportamiento cuando las políticas cambiaban, debido a que no trataban adecuadamente las expectativas racionales.
Ante este panorama, Sims propuso los modelos VAR, lo que eventualmente le valió el Premio Nobel de Economía. El objetivo de estos modelos era ofrecer una alternativa más flexible y empíricamente realista, sin imponer supuestos estructurales fuertes a priori. De esta manera, los VAR permiten modelar un conjunto de variables macroeconómicas como funciones de sus propios rezagos, sin imponer una estructura teórica rígida. Además, todos los choques y relaciones se tratan inicialmente de forma simétrica y empírica. Así, los modelos VAR en su forma reducida son modelos completamente ateóricos.
Sin embargo, aunque los modelos VAR permiten capturar la dinámica conjunta entre variables sin imponer restricciones teóricas fuertes a priori, tienen una limitación que no es menor: los choques estimados a partir de la forma reducida -forma que hemos trabajado hasta ahora- están correlacionados entre sí. Es decir, los errores \(u_t\) de la forma reducida reflejan combinaciones lineales de múltiples choques estructurales y esto nos impide identificar la interacción contemporánea de las variables. Los choques \(u_t\) no se pueden interpretar directamente porque no son ortogonales. En este sentido, para encontrar choques estructurales, debemos encontrar choques \(w_t\) que sí sean ortogonales y tengan una interpretación económica.
Comencemos con la forma estructural del VAR para comprender sus parámetros: \[\begin{equation} B_0 y_t = B_1 y_{t-1} + \ldots + B_p y_{t-p} + w_t \end{equation}\]
donde \(w_t\) es un término de error con media cero y no correlacionado en el tiempo, también conocido como innovación estructural o shock estructural. Además, se asume que el término de error es incondicionalmente homocedástico a menos de que se indique lo contrario. La matriz \(B_0\) es no singular y establece la interacción contemporánea entre las variables del modelo. De este modo, el modelo puede escribirse de manera compacta como: \[\begin{equation} B(L) y_t = w_t \end{equation}\]
donde \(B(L) \equiv B_0 - B_1 L - B_2 L^2 - \ldots - B_p L^p\) es el polinomio autorregresivo en rezagos. La matriz de covarianza del término de error estructural es normalizada tal que: \[\begin{equation} \mathbb{E}(w_t w_t') \equiv \Sigma_w = I_K. \end{equation}\]
Con esto, sabemos que existen tantos shocks estructurales como variables en el modelo. Asimismo, los shocks estructurales, por definición, son mutuamente no correlacionados lo que implica que \(\Sigma_w\) es diagonal. Por último, también sabemos que si normalizamos la varianza de todos los shocks estructurales a uno, no implica una pérdida de generalidad siempre y cuando los elementos diagonales \(B_0\) permanezcan sin restricciones.
Sin embargo, para que el modelo \(B_0 y_t = B_1 y_{t-1} + \ldots + B_p y_{t-p} + w_t\) sea considerado como un VAR estructural no es suficiente con que los elementos de \(w_t\) no se encuentren correlacionados, sino que los choques también deben ser económicamente interpretables. Para esto, derivamos la forma reducida de este modelo VAR estructural de tal modo que \(y_t\) es una función de los rezagos de \(y_t\) únicamente. Si multiplicamos ambos lados de la representación estructural del VAR por \(B_0^{-1}\) \[\begin{equation} B_0^{-1} B_0 y_t = B_0^{-1} B_1 y_{t-1} + ... + B_0^{-1} B_p y_{t-p} + B_0^{-1} w_t \end{equation}\]
de tal modo que puede ser representado como \[\begin{equation} y_t = A_1 y_{t-1} + ... + A_p y_{t-p} + u_t \end{equation}\]
donde \(A_i = B_0^{-1} B_i\), \(i = 1,...,p\), y \(u_t = B_0^{-1} w_t\). Así, las innovaciones en forma reducida de \(u_t\) son un promedio ponderado de los shocks estructurales \(w_t\). De manera compacta, el modelo puede ser expresado como: \[\begin{equation} A(L) y_t = u_t \end{equation}\]
donde \(A(L) = I_K - A_1L - A_2 L^2 - \ldots - A_p L^p\) es el polinomio autorregresivo en rezagos. Los métodos de estimación estándar permiten obtener estimaciones de los parámetros de la forma reducida \(A_i\), para \(i=1,...p\), de las innovaciones reducidas \(u_t\) y de su matriz de covarianza \(\mathbb{E}(u_t u_t') \equiv \Sigma_u\).
Sin embargo, el debate que se ha generado desde entonces es, una vez estimada la forma reducida, de qué manera se puede recuperar la representación estructural del modelo VAR. Es decir, queremos recuperar la matriz \(B_0\) -o bien, su inversa- que contiene las relaciones estructurales contemporáneas de las variables, recordando que \(u_t = B_0^{-1}w_t\).
El método para recuperar dicha matriz ha generado gran debate en el análisis macroeconométrico, dando pie incluso a tesis doctorales al respecto. Con esto en mente, ahora veremos algunos de los métodos clásicos para recuperar las relaciones estructurales de las variables. Nos concentraremos en los siguientes:
1. Identificación Recursiva
2. Restricción de largo plazo
3. Imposición de restricciones de signo
Nota: Los paquetes disponibles en R son un poco limitados en cuanto a los métodos de identificación de las innovaciones estructurales relacionados con la teoría económica. En los trabajos actuales de macroeconometría suele ocuparse, mayoritariamente, el método de Identificación Recursiva. Para los métodos de Restricción de largo plazo y de Imposición de restricciones de signo, el repositorio del proyecto incluye un script de Matlab complementario.
5.4.1 VAR Recursivo - Identificación Recursiva
Una de las formas más populares para recuperar las innovaciones estructurales es mediante la identificación recursiva, un caso particular de las restricciones de corto plazo que se implementa utilizando la descomposición de Cholesky.
La idea detrás de este método es simple e intuitiva: queremos recuperar las innovaciones estructurales \(w_t\) a partir de las innovaciones de forma reducida \(u_t\). Para lograr esto, buscamos ortogonalizar los errores de forma reducida, lo que en este contexto significa transformarlos en un conjunto de choques no correlacionados contemporáneamente. Es decir, deseamos que \(\mathbb{E}[w_t w_t'] = I\).
Para lograr esto, definimos una matriz \(P\) de tamaño \(K \times K\), triangular inferior y con diagonal principal positiva, de tal modo que:
\[\begin{equation} \Sigma_u = \mathbb{E}[u_t u_t'] = P P' \end{equation}\]
La matriz \(P\) es conocida como la descomposición de Cholesky inferior de \(\Sigma_u\). Por ejemplo, en un VAR con tres variables, \(P\) tendría la forma:
\[ P = \begin{bmatrix} p_{11} & 0 & 0 \\ p_{21} & p_{22} & 0 \\ p_{31} & p_{32} & p_{33} \end{bmatrix} \]
Ahora, si recordamos que en el modelo estructural los errores están relacionados por \(u_t = B_0^{-1}w_t\) y que
\[\begin{equation} \Sigma_u = B_0^{-1} \mathbb{E}[w_t w_t'] B_0^{-1'} = B_0^{-1} B_0^{-1'} \end{equation}\]
entonces una solución válida al problema de identificación estructural es asumir que \(B_0^{-1} = P\). Dado que \(P\) es triangular inferior, contiene exactamente \(\frac{K(K-1)}{2}\) ceros impuestos, lo que cumple la condición de orden para identificar todos los elementos libres de la matriz \(B_0^{-1}\). Nótese que la inversa de una matriz triangular inferior es también triangular inferior, de modo que si \(B_0^{-1}\) es triangular inferior, entonces \(B_0\) lo es igualmente.
Sin embargo, es importante destacar que esta ortogonalización de los errores solo es válida si la estructura recursiva impuesta por la matriz \(P\) se puede justificar con fundamentos económicos. El método impone que la primera variable del sistema responde solo a su propio shock contemporáneo, la segunda puede responder al shock de la primera, la tercera a los dos primeros, y así sucesivamente. Por eso se dice que el modelo estructural resultante es recursivo: se impone una estructura jerárquica de causalidad contemporánea, en lugar de inferirla directamente de los datos.
En este sentido, la descomposición de Cholesky no descubre las relaciones estructurales sino que las impone. Por lo tanto, el orden de las variables en el VAR sí importa, ya que determina la interpretación de los shocks estructurales. Por esta razón, se recomienda justificar el orden mediante teoría económica. Por ejemplo, en aplicaciones del Banco de México, el orden se suele establecer de acuerdo con la hipótesis de una economía pequeña y abierta, ordenando las variables de la más exógena a la más endógena.
Veamos un ejemplo de la identificación recursiva mediante la
descomposición de Cholesky. Utilizaremos la base USA que acompaña al
paquete svars (Lange et al. 2023) –el mismo sistema trimestral de tres
variables para Estados Unidos que se emplea habitualmente en la
literatura de VAR estructurales–: la brecha del producto (x), la inflación anualizada del
deflactor del PIB (pi) y la tasa de fondos federales (i). El orden en
que se declaran las variables es precisamente el que impone la estructura
recursiva: los choques de demanda afectan contemporáneamente a la
inflación y a la tasa de política, la inflación afecta a la tasa de
política, y la tasa de política no afecta contemporáneamente a las otras
dos –el supuesto habitual de que la política monetaria opera con
rezago–.
# install.packages("svars")
library(svars)
library(ggplot2)
library(ggfortify)
# x = Porcentaje de desviación logarítmica del PIB real con respecto a la
# estimación del producto potencial
# pi = Crecimiento anualizado trimestre a trimestre del deflactor del PIB
# i = Tasa de interés de los fondos federales
var.reducido <- vars::VAR(USA, lag.max = 10, ic = "AIC")
var.estructural <- id.chol(var.reducido)
# Matriz de impacto estructural estimada B_0^{-1} (triangular inferior)
knitr::kable(
var.estructural$B,
col.names = c("Choque en $x$", "Choque en $\\pi$", "Choque en $i$"),
caption = paste("Matriz de impacto contemporáneo $B_0^{-1}$ estimada",
"por identificación recursiva (descomposición de",
"Cholesky) para el sistema de tres variables de",
"Estados Unidos."),
digits = 4,
booktabs = TRUE,
escape = FALSE
)| Choque en \(x\) | Choque en \(\pi\) | Choque en \(i\) | |
|---|---|---|---|
| x | 0.6834 | 0.0000 | 0.0000 |
| pi | -0.0364 | 1.0727 | 0.0000 |
| i | 0.2245 | 0.1818 | 0.7672 |
Los ceros por arriba de la diagonal del Cuadro 5.9 no son un resultado de la estimación, sino las restricciones que impone el método. Con la matriz de impacto ya identificada, las funciones de impulso-respuesta estructurales se pueden calcular junto con bandas de confianza obtenidas por bootstrap –utilizamos el wild bootstrap, robusto a heterocedasticidad–, como se muestra en la Figura 5.9.
set.seed(1234)
boot.svar <- wild.boot(var.estructural, n.ahead = 24, nboot = 500, nc = 1)
plot(boot.svar)
Figura 5.9: Funciones de impulso-respuesta estructurales del VAR recursivo con bandas de confianza al 95% obtenidas por wild bootstrap (500 réplicas)
La lectura de la Figura 5.9 es la habitual en la literatura macroeconométrica: un choque contractivo de política monetaria –un aumento no anticipado de la tasa de fondos federales– reduce la brecha del producto y, con mayor rezago, la inflación, mientras que un choque de demanda positivo eleva simultáneamente el producto, la inflación y la respuesta de política.
5.4.2 VAR Estructural - Restricción de largo plazo
Otra manera de identificar las innovaciones estructurales de un modelo VAR, es mediante el método propuesto por Blanchard y Quah (1989) en “The Dynamic Effects of Aggregate Demand and Supply Disturbances”.
La intuición de este método es simple: Blanchard y Quah querían medir el efecto de choques de oferta agregada \(w_t^{AS}\) y demanda agregada \(w_t^{AD}\) sobre el Producto y el Desempleo, imponiendo restricciones en la respuesta acumulada, es decir, en el efecto acumulado permanente (horizonte \(h \rightarrow \infty\)). La lógica económica detrás de este método es, básicamente, que el choque de demanda agregada no tiene efectos de largo plazo sobre el nivel del PIB real.
Ahora bien, sean:
- \(ur_t\): tasa de desempleo en EE.UU.
- \(gdp_t\): logaritmo del PIB real de EE.UU.
Y definimos el vector \(z_t\)
\[\begin{equation} z_t = \begin{bmatrix} \Delta gdp_t \\ ur_t \end{bmatrix} \sim I(0) \end{equation}\]
Aunque \(gdp_t \sim I(1)\), su primera diferencia es estacionaria (\(I(0)\)), lo cual justifica que \(z_t \sim I(0)\). Por otro lado, podemos suponer que este vector es generado por un VAR reducido:
\[\begin{equation} A(L) z_t = u_t \end{equation}\]
donde \(A(L)=I_2-A_1L-\ldots-A_pL^p\), y \(u_t \sim (0, \Sigma_u)\) es ruido blanco. Si lo representamos como su forma estructural, tenemos que:
\[\begin{equation} B(L) z_t = w_t \end{equation}\]
donde \(B(L) = B_0-B_1L-\ldots-B_pL^p = B_0A(L)\), y si recordamos a qué es igual \(w_t\), tenemos que:
\[\begin{equation} w_t = B_0 u_t \rightarrow u_t = B_0^{-1} w_t \rightarrow \Sigma_u = B_0^{-1} (B_0^{-1})' \end{equation}\]
Ahora, pasemos a la representación MA estructural para entender dónde se impone la restricción. Desde el modelo MA, tenemos que:
\[\begin{equation} z_t = B(L)^{-1} w_t = \Theta (L) w_t \end{equation}\]
Esto implica que, dado que \(z_t \sim I(0)\), el efecto de un choque estructural sobre \(z_t\) se desvanece con el tiempo y, por tanto, tanto \(\Delta gdp_t\) como \(ur_t\) eventualmente vuelven a su nivel original tras un choque. Sin embargo, el nivel del PIB real -y no su crecimiento- no necesariamente vuelve a su nivel inicial, por lo que la suma acumulada de las respuestas al choque nos da su efecto permanente. Ahora, pasemos a la matriz de efectos acumulados de largo plazo, definida como:
\[\begin{equation} \Theta (1) = \sum_{i=0}^{\infty} \Theta_i = B(1)^{-1} \end{equation}\]
Dado que queremos que el PIB real vuelva a su tendencia tras un choque de demanda -recordemos que no tiene efectos de largo plazo- implica imponer un cero en la posición \(\theta_{12}\) de la matriz \(\Theta(1)\):
\[\begin{equation} \Theta (1) = \begin{bmatrix} \theta_{11} (1) & 0 \\ \theta_{21} (1) & \theta_{22} (1) \end{bmatrix} \end{equation}\]
el hecho de que \(\theta_{12}(1) = 0\) implica que el choque de demanda no afecta permanentemente al PIB real y, por otro lado, que \(\theta_{11}(1)\) permanezca libre, implica que los choques de oferta sí pueden tener efectos de largo plazo.
Para la identificación, recordemos que:
\[\begin{equation} \Theta (1) = B(1)^{-1} = A(1)^{-1} B_0^{-1} \rightarrow B_0^{-1} = A(1) \Theta(1) \end{equation}\]
De este modo, una vez que conocemos \(A(1)\) y \(\Theta(1)\) podemos recuperar \(B_0^{-1}\), recordando que la restricción en \(\Theta(1)\) es equivalente a imponer una restricción sobre \(B_0\).
Para estimar \(\Theta(1)\), sabemos que:
\[\begin{equation} \Sigma_u = B_0^{-1} (B_0^{-1})' = \Theta(1) \Theta(1)' \end{equation}\]
y, entonces, si usamos \(A(1)\):
\[\begin{equation} A(1)^{-1} \Sigma_u A(1)^{-1'} = \Theta(1) \Theta(1)' \end{equation}\]
lo cual implica que podemos calcular el lado izquierdo usando solo parámetros estimados del VAR reducido y luego identificar \(\Theta(1)\) imponiendo que tenga forma triangular inferior, y aplicar la descomposición de Cholesky para tener finalmente que:
\[\begin{equation} B_0^{-1} = A(1) \Theta (1) \end{equation}\]
5.4.3 VAR Estructural - Imposición de restricciones de signo
Una alternativa a la identificación de las innovaciones estructurales además de la Identificación recursiva es la Imposición de restricciones de signo. Sigue un poco la misma lógica: queremos imponer restricciones de signo en vez de proponer una matriz triangular inferior en la matriz \(B_0^{-1}\). Este método ofrece una alternativa menos restrictiva que imponer ceros exactos. Para explicar este método, tomemos el siguiente modelo ejemplo.
Consideremos un modelo bivariado de un mercado de bienes con:
- Un choque de demanda \(w_t^{demanda}\)
- Un choque de oferta \(w_t^{oferta}\)
en donde nuestras variables observables son:
- Precio \(p_t\)
- Cantidad \(q_t\)
Y los errores en forma reducida están dados por:
\[\begin{equation} u_t = \begin{bmatrix} u_t^q \\ u_t^p \end{bmatrix} = B_0^{-1} w_t \quad \text{donde} \quad w_t = \begin{bmatrix} w_t^{oferta} \\ w_t^{demanda} \end{bmatrix} \end{equation}\]
Es sencillo: la interpretación económica es que los efectos de los choques \(u_t\) dependen de la pendiente de las curvas de demanda y oferta.
En un método tradicional -restricción de exclusión-, se asume que la oferta es vertical en el corto plazo, lo cual implica que los choques de demanda no afectan la cantidad contemporáneamente:
\[\begin{equation} \begin{bmatrix} u_t^q \\ u_t^p \end{bmatrix} = \begin{bmatrix} \ast & 0 \\ \ast & \ast \end{bmatrix} \begin{bmatrix} w_t^{oferta} \\ w_t^{demanda} \end{bmatrix} \end{equation}\]
En donde los asteriscos representan coeficientes libres y el cero es una restricción de exclusión.
En cambio, si tomamos la teoría económica básica en el enfoque alternativo de Imposición de restricción de signo esperaríamos que, ante un choque de oferta positivo -la curva de oferta se desplaza a la derecha-, la cantidad aumente y el precio caiga, mientras que, ante un choque de demanda positiva -la curva de demanda se desplaza a la derecha-, la cantidad aumenta y el precio aumenta. De este modo, tenemos que:
\[\begin{equation} \begin{bmatrix} u_t^q \\ u_t^p \end{bmatrix} = \underbrace{ \begin{bmatrix} + & + \\ - & + \end{bmatrix} }_{B_0^{-1}} \begin{bmatrix} w_t^{\text{oferta}} \\ w_t^{\text{demanda}} \end{bmatrix} \end{equation}\]
Donde \(+\) indica signo positivo estricto y \(-\) negativo estricto.
La diferencia clave con el método de Identificación recursiva mediante la descomposición de Cholesky es que, en estos últimos modelos, los parámetros están puntualmente identificados mientras que con restricciones de signo los parámetros no se identifican exactamente, sino que quedan dentro de un conjunto compatible con las restricciones. De este modo, no hay una única solución, pues hay muchas matrices \(B_0^{-1}\) que cumplen con los signos impuestos. También, se pueden imponer restricciones de signo de la siguiente forma:
\[\begin{equation} \begin{bmatrix} + & + \\ 0 & + \end{bmatrix} \quad \text{o} \quad \begin{bmatrix} + & + \\ - & 0 \end{bmatrix} \end{equation}\]
Estas formas también pueden ser admisibles pero algunas podrían llegar a ser problemáticas si llevan a choques no identificables entre sí.
Se dice que los modelos identificados por signos son más generales que los recursivos pero no son modelos “anidados”, o sea, uno no puede validar o rechazar el otro con datos. Así, las restricciones de signo permiten más flexibilidad, pero imponen otras limitaciones como perder identificación puntual.
Ahora, de manera más puntual, para imponer las restricciones de signo estáticas consideremos un modelo VAR en su forma estructural:
\[\begin{equation} B_0 y_t = B_1 y_{t-1} + ... + B_p y_{t-p} + w_t \end{equation}\]
y normalizamos la matriz de varianzas-covarianzas del término de error estructural \(w_t\) tal que:
\[\begin{equation} \mathbb{E}(w_t w_t') \equiv \Sigma_w = I_K. \end{equation}\]
Ahora, sea \(u_t = P \eta_t\), donde \(u_t\) es la innovación del VAR en forma reducida y \(P\) es la descomposición de Cholesky de \(\Sigma_u\). Por construcción, los choques \(\eta_t\) son mutuamente no correlacionados y tienen varianza unitaria. Por supuesto, no hay razón para que estos choques correspondan a choques estructurales económicamente interpretables, como los choques de oferta y demanda en el modelo bivariado mencionado anteriormente. Sin embargo, podemos buscar soluciones candidatas \(w_t^*\) para los choques estructurales desconocidos \(w_t\) construyendo un gran número de combinaciones de los choques \(\eta_t\) de la forma:
\[\begin{equation} w_t^* = Q' \eta_t, \end{equation}\]
donde \(Q'\) es una matriz ortogonal cuadrada tal que \(Q'Q = QQ' = I_K\) y \(u_t = PQ \eta_t = PQw_t^*\). Por lo tanto, cada solución candidata \(w_t^*\) consiste en choques no correlacionados con varianza unitaria. Que una de estas soluciones candidatas \(w_t^*\) sea una solución admisible para el choque estructural desconocido \(w_t\), dado el vector de parámetros en forma reducida, depende de si la matriz de impacto estructural implicado \(PQ\) satisface las restricciones de signo mantenidas sobre \(B_0^{-1}\). De este modo, solo conservamos las soluciones que satisfagan las restricciones de signo y se descartan las demás. El hecho de repetir este procedimiento nos permite caracterizar el conjunto de todos los modelos estructurales que son consistentes con las restricciones de signo mantenidas y los parámetros en forma reducida. Conocer \(PQ\) permite construir todos los coeficientes estructurales de respuesta al impulso que se derivan de las estimaciones de los parámetros en forma reducida.
Los dos enfoques comunes para construir las matrices ortogonales \(Q\) están basados en:
Matrices de rotación de Givens.
Transformación de Householder.
La paquetería del código de MatLab, VAR Toolbox, ocupa el enfoque de Matrices de rotación de Givens.
5.5 Cointegración
Hasta ahora hemos usado el supuesto de que las series son estacionarias para el conjunto de técnicas \(ARMA(p,q)\) y \(VAR(p)\). No obstante, dado que relajamos el supuesto de estacionariedad (incluyendo la estacionariedad en varianza) y que establecimos una serie de pruebas para determinar cuándo una serie es estadísticamente estacionaria, ahora podemos plantear una técnica llamada Cointegración. Para esta técnica consideraremos sólo series que son \(I(1)\) y reconoceremos que se originó con los trabajos de Engle y Granger (1987), Stock (1987) y Johansen (1988).
5.5.1 Definición y propiedades del proceso de cointegración
Cointegración puede ser caracterizada o definida en palabras sencillas como que dos o más variables tienen una relación común estable en el largo plazo. Es decir, estas no suelen tomar caminos o trayectorias diferentes, excepto por períodos de tiempo transitorios y eventuales. A continuación, utilizaremos la definición de Engle y Granger (1987) de cointegración.
Sea \(\mathbf{Y}\) un vector de k-series de tiempo, decimos que los elementos en \(\mathbf{Y}\) están cointegrados en un orden (d, c), es decir, \(\mathbf{Y} \sim CI(d, c)\), si todos los elementos de \(\mathbf{Y}\) son series integradas de orden d, I(d), y si existe al menos una combinación lineal no trivial \(\mathbf{Z}\) de esas variables que es de orden I(d - c), donde \(d \geq c > 0\), si y sólo si: \[\begin{equation} \boldsymbol{\beta}_i' \mathbf{Y}_t = \mathbf{Z}_{it} \sim I(d-c) \end{equation}\]
Donde \(i = 1, 2, \ldots, r\) y \(r < k\).
A los diferentes vectores \(\boldsymbol{\beta}_i\) se les denomina como vectores de cointegración. El rango de la matriz de vectores de cointegración \(r\) es el número de vectores de cointegración linealmente independientes. En general diremos que los vectores de la matriz de cointegración \(\boldsymbol{\beta}\) tendrán la forma de: \[\begin{equation} \boldsymbol{\beta}' \mathbf{Y}_t = \mathbf{Z}_t \end{equation}\]
Antes de continuar hagamos algunas observaciones. Si todas las variables de \(\mathbf{Y}\) son I(1) y \(0 \leq r < k\), diremos que las series no cointegran si \(r = 0\). Si esto pasa, entonces, como demostraremos más adelante, la mejor opción será estimar un modelo VAR(p) en diferencias. Adicionalmente, asumiremos que \(c = d = 1\), por lo que la relación de cointegración, en su caso, generará combinaciones lineales \(\mathbf{Z}\) estacionarias.
5.5.2 Cointegración para modelos de más de una ecuación o para modelos basados en Vectores Autorregresivos
Sean \(Y_1, Y_2, \ldots, Y_k\) las series que forman \(\mathbf{Y}\), todas ellas I(1); entonces los siguientes casos son posibles:
Si \(r = 1\) entonces existe una única relación de cointegración y el caso puede analizarse con el enfoque de Engle y Granger.
Si \(r > 1\) entonces se trata de un caso de cointegración múltiple, que requiere el enfoque de Johansen.
Por lo anterior, en este libro analizaremos el caso de Cointegración de Johansen, que cubre ambas situaciones. Ahora plantearemos la forma de estimar el proceso de cointegración. El primer paso para ello es determinar un modelo VAR(p) con las k-series no estacionarias (series en niveles)–en este punto se vuelve fundamental caracterizar las series a través de pruebas de raíces unitarias–. Elegimos el valor de \(p\) mediante el uso de los criterios de información. De esta forma tendremos una especificación similar a: \[\begin{equation} \mathbf{Y}_t = \sum_{j=1}^p \mathbf{A}_j \mathbf{Y}_{t-j} + \mathbf{D}_t + \mathbf{U}_t \tag{5.7} \end{equation}\]
Donde \(\mathbf{U}_t\) es un término de error k-dimensional puramente aleatorio; \(\mathbf{D}_t\) contiene los componentes determinísticos de constante y tendencia, y \(\mathbf{A}_i\), \(i = 1, 2, \ldots, p\), son matrices de \(k \times k\) coeficientes. Notemos que el VAR(p) involucrado en este caso, a diferencia del VAR anteriormente estudiado, puede incluir un término de tendencia. Esto en razón de que hemos relajado el concepto de estacionariedad.
Si reescribimos la ecuación (5.7) en su forma de Vector Corrector de Errores (VEC, por sus siglas en inglés), restamos \(\mathbf{Y}_{t-1}\) en ambos lados y usamos la identidad \(\mathbf{Y}_{t-j} = \mathbf{Y}_{t-1} - \sum_{i=1}^{j-1} \Delta \mathbf{Y}_{t-i}\) válida para \(j \geq 2\): \[\begin{eqnarray} \Delta \mathbf{Y}_t & = & \sum_{j=1}^p \mathbf{A}_j \mathbf{Y}_{t-j} + \mathbf{D}_t - \mathbf{Y}_{t-1} + \mathbf{U}_t \nonumber \\ & = & (\mathbf{A}_1 - \mathbf{I}) \mathbf{Y}_{t-1} + \mathbf{A}_2 \mathbf{Y}_{t-2} + \ldots + \mathbf{A}_p \mathbf{Y}_{t-p} + \mathbf{D}_t + \mathbf{U}_t \nonumber \\ & = & \left( \sum_{j=1}^{p} \mathbf{A}_j - \mathbf{I} \right) \mathbf{Y}_{t-1} - \sum_{j=2}^{p} \mathbf{A}_j \sum_{i=1}^{j-1} \Delta \mathbf{Y}_{t-i} + \mathbf{D}_t + \mathbf{U}_t \nonumber \\ & = & \left( \sum_{j=1}^{p} \mathbf{A}_j - \mathbf{I} \right) \mathbf{Y}_{t-1} + \sum_{j=1}^{p-1} \mathbf{A}^*_j \Delta \mathbf{Y}_{t-j} + \mathbf{D}_t + \mathbf{U}_t \nonumber \\ & = & - \left( \mathbf{I} - \sum_{j=1}^{p} \mathbf{A}_j \right) \mathbf{Y}_{t-1} + \sum_{j=1}^{p-1} \mathbf{A}^*_j \Delta \mathbf{Y}_{t-j} + \mathbf{D}_t + \mathbf{U}_t \nonumber \\ \Delta \mathbf{Y}_t & = & - \Pi \mathbf{Y}_{t-1} + \sum_{j=1}^{p-1} \mathbf{A}^*_j \Delta \mathbf{Y}_{t-j} + \mathbf{D}_t + \mathbf{U}_t \tag{5.8} \end{eqnarray}\]
Donde \(\mathbf{A}_j^* = - \sum_{i=j+1}^p \mathbf{A}_i\), para \(j = 1, 2, \ldots, p-1\), y la matriz \(\Pi\) captura todas las relaciones de largo plazo entre las variables. Para que exista cointegración, \(\Pi\) debe tener rango reducido \(r\), con \(0 \leq r < k\). Por lo tanto, tenemos que dicha matriz en la ecuación (5.8) se puede factorizar como: \[\begin{equation} \Pi_{(k \times k)} = \Gamma_{(k \times r)} \boldsymbol{\beta}_{(r \times k)}' \tag{5.9} \end{equation}\]
Donde \(\boldsymbol{\beta}_{(r \times k)}' \mathbf{Y}_{t-1}\) son \(r\) combinaciones linealmente independientes que son estacionarias.
Una advertencia sobre la convención de signos: aquí definimos \(\Pi = \mathbf{I} - \sum_j \mathbf{A}_j\), de modo que el término de corrección de error aparece con signo negativo en la ecuación (5.8). Buena parte de la literatura –y la salida de los paquetes de cómputo– define en cambio \(\Pi^{*} = - \Pi = \boldsymbol{\alpha} \boldsymbol{\beta}'\), con lo que el VEC se escribe \(\Delta \mathbf{Y}_t = \boldsymbol{\alpha} \boldsymbol{\beta}' \mathbf{Y}_{t-1} + \ldots\) y los coeficientes de ajuste de \(\boldsymbol{\alpha}\) deben ser negativos para que las desviaciones del equilibrio de largo plazo se corrijan. Con nuestra convención, \(\Gamma = - \boldsymbol{\alpha}\) y los signos se invierten. Al interpretar resultados conviene verificar cuál de las dos convenciones utiliza el programa empleado.
Dada la ecuación (5.8) podemos establecer la aproximación de Johansen (1988) que se realiza mediante una estimación por Máxima Verosimilitud de la ecuación: \[\begin{equation} \Delta \mathbf{Y}_t + \Gamma \boldsymbol{\beta}' \mathbf{Y}_{t-1} = \sum_{j=1}^{p-1} \mathbf{A}^*_j \Delta \mathbf{Y}_{t-j} + \mathbf{D}_t + \mathbf{U}_t \end{equation}\]
Donde una vez estimado el sistema: \[\begin{equation} \boldsymbol{\beta} = [v_1, v_2, \ldots, v_r] \end{equation}\]
Cada \(v_i\), \(i = 1, 2, \ldots, r\) es un vector propio que está asociado con los \(r\) valores propios positivos, mismos que están asociados con la prueba de hipótesis de cointegración. Dicha hipótesis está basada en dos estadísticas con las que se determina el rango \(r\) de \(\Pi\):
Prueba de Traza: \(H_0 :\) Existen al menos \(r\) valores propios positivos o Existen al menos \(r\) relaciones de largo plazo estacionarias.
Prueba del valor propio máximo o \(\lambda_{max}\): \(H_0 :\) Existen \(r\) valores propios positivos o Existen \(r\) relaciones de largo plazo estacionarias.
Ejemplo. Para ejemplificar el procedimiento de cointegración utilizaremos las series de INPC, Tipo de Cambio, rendimiento de los Cetes a 28 días, IGAE e Índice de Producción Industrial de Estados Unidos, para el mismo periodo de enero de 2000 a mayo de 2026. Quizá el marco teórico de la relación entre las variables no sea del todo correcto, pero dejando de lado ese problema, estimaremos si las 5 series cointegran.
Por principio, probaremos que todas las series son I(1), lo cual es cierto (ver Script para mayores detalles). En las Figuras 5.10, 5.11 y 5.12 se muestran las series en niveles y en diferencias, con lo cual ilustramos como es viable que las series sean I(1).
library(ggplot2)
library(dplyr)
library(stats)
library(MASS)
library(strucchange)
library(zoo)
library(sandwich)
library(urca)
library(lmtest)
library(vars)
#
load("BD/Datos_Ad.RData")
#
### Conversión a series de tiempo
Datos <- ts(Datos_Ad[7: 11],
start = c(2000, 1),
freq = 12)
LDatos <- log(Datos)
DLDatos <- diff(log(Datos, base = exp(1)),
lag = 1,
differences = 1)plot(LDatos,
plot.type = "m", nc = 2,
col = c("darkgreen", "darkblue", "darkred", "orange", "purple"),
#main = "Series en Logaritmos",
xlab = "Tiempo")
Figura 5.10: Series en niveles (logaritmos) para la prueba de Cointegración
plot(DLDatos,
plot.type = "m", nc = 2,
col = c("darkgreen", "darkblue", "darkred", "orange", "purple"),
#main = "Series en Diferencias Logaritmicas",
xlab = "Tiempo")
Figura 5.11: Series en Diferencias Logarítmicas para la prueba de Cointegración
plot(cbind(LDatos, DLDatos),
plot.type = "m", nc = 2,
col = c("darkgreen", "darkblue", "darkred", "orange", "purple"),
#main = "Comparación de Series en Diferencias",
xlab = "Tiempo")
Figura 5.12: Comparación de Series en Diferencias para la prueba de Cointegración
Posteriormente, determinamos cuál es el orden adecuado de un VAR(p) en niveles. En el Cuadro 5.10 mostramos los resultados de los criterios de información para determinar el número de rezagos óptimos, el cual resultó en \(p = 3\) para los criterios AIC y FPE, y \(p = 2\) para los criterios HQ y SC. Por lo tanto, decidiremos utilizar un VAR(3) con tendencia y constante. Note que es posible elegir otros modelos de VAR que incluyan: solo tendencia, solo constante o ninguno de estos elementos.
vs_niv <- vars::VARselect(LDatos, lag.max = 10, type = "both")
crit_niv <- as.data.frame(t(vs_niv$criteria))
tab_niv <- data.frame(
Rezagos = seq_len(nrow(crit_niv)),
AIC = round(crit_niv[, 1], 4),
HQ = round(crit_niv[, 2], 4),
SC = round(crit_niv[, 3], 4),
FPE = formatC(crit_niv[, 4], format = "e", digits = 3))
knitr::kable(
tab_niv,
col.names = c("Rezagos", "AIC", "HQ", "SC", "FPE"),
caption = paste("Criterios de información para diferentes",
"especificaciones de modelos VAR(p) con término",
"constante y tendencia de las series $LINPC_t$,",
"$LTC_t$, $LCETE28_t$, $LIGAE_t$ y $LIPI_t$."),
align = c("c", "r", "r", "r", "r"),
booktabs = TRUE,
row.names = FALSE
)| Rezagos | AIC | HQ | SC | FPE |
|---|---|---|---|---|
| 1 | -43.9249 | -43.7550 | -43.5000 | 8.389e-20 |
| 2 | -44.2918 | -44.0005 | -43.5634 | 5.813e-20 |
| 3 | -44.3099 | -43.8972 | -43.2780 | 5.711e-20 |
| 4 | -44.2537 | -43.7197 | -42.9183 | 6.045e-20 |
| 5 | -44.2256 | -43.5702 | -42.5867 | 6.224e-20 |
| 6 | -44.1297 | -43.3530 | -42.1873 | 6.861e-20 |
| 7 | -44.1073 | -43.2092 | -41.8614 | 7.031e-20 |
| 8 | -44.1614 | -43.1420 | -41.6121 | 6.679e-20 |
| 9 | -44.1168 | -42.9760 | -41.2640 | 7.008e-20 |
| 10 | -44.0863 | -42.8241 | -40.9300 | 7.257e-20 |
El Cuadro 5.10 corresponde a la especificación con constante y tendencia. Conviene verificar la sensibilidad de la elección del orden a los componentes determinísticos incluidos; el siguiente cuadro reporta el número de rezagos que selecciona cada criterio bajo las cuatro especificaciones posibles.
### Selección de rezagos VAR(p) según los componentes determinísticos
tipos <- c(both = "Constante y tendencia", trend = "Sólo tendencia",
const = "Sólo constante", none = "Ninguno")
sel_tipos <- t(sapply(names(tipos), function(tp)
vars::VARselect(LDatos, lag.max = 10, type = tp)$selection))
rownames(sel_tipos) <- tipos
knitr::kable(
sel_tipos,
col.names = c("AIC", "HQ", "SC", "FPE"),
caption = paste("Número de rezagos seleccionado por cada criterio de",
"información según los componentes determinísticos",
"incluidos en el VAR(p) en niveles."),
align = rep("c", 4),
booktabs = TRUE
)| AIC | HQ | SC | FPE | |
|---|---|---|---|---|
| Constante y tendencia | 3 | 2 | 2 | 3 |
| Sólo tendencia | 3 | 2 | 2 | 3 |
| Sólo constante | 3 | 2 | 2 | 3 |
| Ninguno | 3 | 2 | 2 | 3 |
El mismo número de rezagos los utilizaremos para probar la Cointegración, ya sea por una estadística de la Traza o por una del máximo valor propio. Dado que los resultados se sostienen, sólo mostraremos uno de los casos en que las series cointegran y únicamente para el caso de la prueba de la traza (el otro caso está disponible en el código de R disponible abajo). En el Cuadro 5.12 reportamos los resultados del Test de Cointegración para un modelo con 3 rezagos y término constante en la relación de largo plazo.
library(urca)
CA_1 <- ca.jo(LDatos, type = "trace", ecdet = "const", K = 3,
spec = "longrun")
# Primer vector de cointegración, normalizado respecto de su primera
# entrada (LINPC). Se calcula aquí para que el texto, la ecuación de largo
# plazo y la gráfica de residuales usen siempre los valores estimados.
beta_hat <- CA_1@V[, 1]
beta_hat <- beta_hat / beta_hat[1]
# Versión con 4 decimales para el texto y las ecuaciones
beta_txt <- formatC(beta_hat, format = "f", digits = 4)
beta_abs <- formatC(abs(beta_hat), format = "f", digits = 4)
pct_cv <- if (knitr::is_latex_output()) "\\%" else "%"
cv_cols <- paste0(c("10", "5", "1"), pct_cv)
traza_df <- data.frame(
Hip = c("$r \\leq 4$", "$r \\leq 3$", "$r \\leq 2$",
"$r \\leq 1$", "$r = 0$"),
Stat = round(as.numeric(CA_1@teststat), 2),
round(CA_1@cval, 2))
knitr::kable(
traza_df,
col.names = c("Hipótesis nula", "Estadística", cv_cols),
caption = paste("Prueba de la traza para cointegración considerando",
"un VAR(3) con término constante de las series",
"$LINPC_t$, $LTC_t$, $LCETE28_t$, $LIGAE_t$ y $LIPI_t$."),
align = c("c", "c", "c", "c", "c"),
booktabs = TRUE,
escape = FALSE,
row.names = FALSE
)| Hipótesis nula | Estadística | 10% | 5% | 1% |
|---|---|---|---|---|
| \(r \leq 4\) | 4.57 | 7.52 | 9.24 | 12.97 |
| \(r \leq 3\) | 13.73 | 17.85 | 19.96 | 24.60 |
| \(r \leq 2\) | 26.02 | 32.00 | 34.91 | 41.07 |
| \(r \leq 1\) | 47.71 | 49.65 | 53.12 | 60.16 |
| \(r = 0\) | 138.43 | 71.86 | 76.07 | 84.45 |
Los resultados del Cuadro 5.12 indican que la hipótesis nula \(r = 0\) se rechaza al \(5\%\), mientras que la hipótesis \(r \leq 1\) no se rechaza. La secuencia de la prueba de la traza se detiene en ese punto, por lo que concluimos que existe evidencia estadística de exactamente un vector de cointegración (\(r = 1\)). Dicho vector, normalizado respecto de su primera entrada, es: \[\begin{equation} \boldsymbol{\beta} = \left[ \begin{matrix} 1.0000 \\ 0.8379 \\ 1.1526 \\ -3.9087 \\ -6.7655 \\ 44.0974 \\ \end{matrix} \right] \end{equation}\]
Donde las entradas corresponden, en ese orden, a \(LINPC_t\), \(LTC_t\), \(LCETE28_t\), \(LIGAE_t\), \(LIPI_t\) y al término constante. Como \(\boldsymbol{\beta}' \mathbf{Y}_t\) debe ser estacionario, al despejar la primera variable obtenemos la relación de largo plazo estimada: \[\begin{eqnarray*} LINPC_t & = & - 0.8379 LTC_t - 1.1526 LCETE28_t \\ & & + 3.9087 LIGAE_t + 6.7655 LIPI_t \\ & & - 44.0974 \end{eqnarray*}\]
Los signos de los coeficientes de largo plazo deben interpretarse con cautela: como se advierte más abajo, la evidencia de cointegración de este sistema es frágil y la especificación del modelo tiene una motivación principalmente ilustrativa.
### Estimación del VAR(p)
VAR_1 <- vars::VAR(LDatos, p = 3, type = "both")
#summary(VAR_1)
#plot(VAR_1, names = "INPC_Ad")
#plot(VAR_1, names = "TC_Ad")
#plot(VAR_1, names = "CETE28_Ad")
#plot(VAR_1, names = "IGAE_Ad")
#plot(VAR_1, names = "IPI_Ad")
# Cointegration Test:
#ca.jo = function (x, type = c("eigen", "trace"), ecdet = c("none", "const",
#"trend"), K = 2, spec = c("longrun", "transitory"), season = NULL,
#dumvar = NULL)
#summary(ca.jo(LDatos, type = "trace", ecdet = "trend", K = 3, spec = "longrun"))
#summary(ca.jo(LDatos, type = "trace", ecdet = "const", K = 3, spec = "longrun"))
#summary(ca.jo(LDatos, type = "trace", ecdet = "none", K = 3, spec = "longrun"))
# La salida completa de la prueba de Johansen (valores propios, vectores
# de cointegración y coeficientes de ajuste) se obtiene con:
# summary(CA_1)Considerando lo anterior, podemos determinar \(\hat{U}_t\) para esta ecuación de cointegración. En la Figura 5.13 mostramos los residuales estimados. Derivado de la inspección visual, parecería que estos no son estacionarios –cuando, de existir cointegración, deberían serlo–. De esta forma, una prueba deseable es aplicar todas las pruebas de raíces unitarias a esta serie para verificar que es I(0); al hacerlo (el código está disponible en el repositorio del proyecto) se encuentra que es posible que no sea estacionaria, por lo que la evidencia de cointegración debe tomarse con cautela.
# Combinación lineal beta'Y_t (los residuales de la relación de largo
# plazo). Se usa el vector estimado, no valores fijos en el texto.
U <- LDatos %*% beta_hat[1:5] + beta_hat[6]
U <- ts(as.numeric(U), start = start(LDatos), frequency = frequency(LDatos))
plot(U,
main = "Residuales de la Ecuación de Cointegración",
ylab = expression(hat(U)[t]), xlab = "Tiempo",
type = "l",
col = "darkred")
abline(h = mean(U), lty = 2, col = "gray40")
Figura 5.13: Residuales estimados de la ecuación de cointegración, \(\hat{U}_t = \boldsymbol{\hat{\beta}}' \mathbf{Y}_t\)
5.6 Modelos ARDL
5.6.1 Teoría
Una vez que hemos analizado diversas técnicas de series de tiempo, el problema consiste en decidir cuál corresponde a un conjunto de series dado. La Figura 5.14 resume la regla de decisión que se desprende de lo visto hasta aquí: todo depende del orden de integración de las series y, cuando son \(I(1)\), de si cointegran o no. Un esquema equivalente puede consultarse en Shrestha y Bhatta (2018).
Figura 5.14: Selección del método de estimación según el orden de integración de las series y la existencia de cointegración. MCO: mínimos cuadrados ordinarios; VAR: vectores autorregresivos; VEC: vector de corrección de error; ARDL: rezagos distribuidos autorregresivos
Conviene leer el esquema como una guía, no como una receta: los resultados de las pruebas de raíz unitaria rara vez son unánimes –como comprobamos en el Capítulo 4– y, cuando hay duda sobre el orden de integración, el enfoque ARDL tiene la ventaja de seguir siendo válido con órdenes mixtos.
En este caso incorporaremos los modelos de rezagos distribuidos autorregresivos (ARDL, autoregressive distributed lag, por sus siglas en inglés). El procedimiento de Johansen no puede aplicarse directamente cuando las variables incluidas tienen órdenes de integración mixtos –es decir, cuando unas son \(I(0)\) y otras \(I(1)\)–, y es precisamente en ese escenario donde el enfoque ARDL resulta útil, ya que se estima por mínimos cuadrados ordinarios (MCO) sobre la especificación en niveles y diferencias.
Este tipo de modelos toma suficientes rezagos para capturar el mecanismo generador de datos. También es posible llegar a una especificación del mecanismo corrector de errores a partir de una transformación lineal del ARDL.
Consideremos la siguiente ecuación: \[\begin{equation} Y_t = \alpha + \delta X_t + \gamma Z_t + U_t \tag{5.10} \end{equation}\]
Dada la ecuación (5.10) podemos establecer su forma de mecanismo corrector de errores en forma ARDL dada por: \[\begin{eqnarray*} \Delta Y_t & = & \alpha + \sum_{i = 1}^p \beta_i \Delta Y_{t-i} + \sum_{i = 1}^p \delta_i \Delta X_{t-i} + \sum_{i = 1}^p \gamma_i \Delta Z_{t-i} \\ & & + \lambda_1 Y_{t-1} + \lambda_2 X_{t-1} + \lambda_3 Z_{t-1} + U_t \end{eqnarray*}\]
Donde los coeficientes \(\beta_i\), \(\delta_i\), \(\gamma_i\) representan la dinámica de corto plazo y las \(\lambda\)’s la dinámica de largo plazo.
La hipótesis nula de la prueba de límites (bounds test) de Pesaran, Shin y Smith (2001) es conjunta: \(H_0 : \lambda_1 = \lambda_2 = \lambda_3 = 0\), es decir, que ninguno de los niveles rezagados entra en la ecuación y, por lo tanto, que no existe una relación de largo plazo entre las variables. La prueba se contrasta con una estadística \(F\) cuya distribución no es estándar: Pesaran, Shin y Smith (2001) tabulan dos conjuntos de valores críticos –una cota inferior, que supone que todas las variables son \(I(0)\), y una cota superior, que supone que todas son \(I(1)\)–. Si la estadística excede la cota superior, se concluye que existe relación de largo plazo; si queda por debajo de la cota inferior, se concluye que no existe; y si cae entre ambas cotas, el resultado es inconcluso sin información adicional sobre el orden de integración de las series.
En la práctica estimamos una especificación con rezagos distribuidos: \[\begin{equation} Y_t = \alpha + \sum_{i = 1}^p \beta_i Y_{t-i} + \sum_{i = 1}^p \delta_i X_{t-i} + \sum_{i = 1}^p \gamma_i Z_{t-i} + U_t \end{equation}\]
Además de verificar si las series involucradas son estacionarias y decidir el número de rezagos \(p\) mediante criterios de información.
5.6.2 Ejemplo
5.6.2.1 Descripción del problema
Supongamos que queremos modelar el logaritmo del dinero real (M2) como una función del logaritmo del ingreso real (LRY), de la tasa de los bonos (IBO) y de la tasa de los depósitos bancarios (IDE).
El problema es que la aplicación de una regresión de MCO en datos no estacionarios daría lugar a una regresión espuria.
Los parámetros estimados serían consistentes solo si las series estuvieran cointegradas.
5.6.2.2 Importamos Datos desde un dataset de R:
Utilizaremos la base denmark, incluida en el paquete ARDL de R, que
reproduce los datos trimestrales de demanda de dinero de Dinamarca
empleados originalmente por Johansen y Juselius
(1990). Se trata de un dataframe con 55 renglones y 5
variables, del primer trimestre de 1974 al tercero de 1987:
LRM: logaritmo del dinero real (M2)
LRY: logaritmo del ingreso real
LPY: logaritmo del deflactor de precios
IBO: tasa de los bonos
IDE: tasa de los depósitos bancarios
## [1] "LRM" "LRY" "LPY" "IBO" "IDE"5.6.2.3 Procedimiento
1. Selección del orden de rezagos. La función auto_ardl() estima
todas las combinaciones de rezagos hasta un orden máximo y devuelve las
que minimizan el criterio de información (por defecto, el AIC).
models <- auto_ardl(LRM ~ LRY + IBO + IDE, data = denmark, max_order = 5)
# Las 5 mejores combinaciones de rezagos según el criterio AIC
knitr::kable(
head(models$top_orders, 5),
caption = paste("Cinco especificaciones ARDL con el menor criterio de",
"información (AIC). Las columnas indican el número de",
"rezagos de cada variable."),
row.names = FALSE,
digits = 3,
booktabs = TRUE
)| LRM | LRY | IBO | IDE | AIC |
|---|---|---|---|---|
| 3 | 1 | 3 | 2 | -251.026 |
| 3 | 1 | 3 | 3 | -250.114 |
| 2 | 2 | 0 | 0 | -249.627 |
| 3 | 2 | 3 | 2 | -249.109 |
| 3 | 2 | 3 | 3 | -248.186 |
El orden seleccionado es un
\(ARDL(3, 1, 3, 2)\): tres rezagos de
LRM, uno de LRY, tres de IBO y dos de IDE. El Cuadro
5.14 reporta los coeficientes estimados de esa
especificación.
knitr::kable(
coef(summary(BestMod)),
col.names = c("Coef.", "Error Est.", "Estad. $t$", "Prob. $(>|t|)$"),
caption = "Estimación del modelo $ARDL$ seleccionado para $LRM_t$.",
digits = 4,
booktabs = TRUE,
escape = FALSE
)| Coef. | Error Est. | Estad. \(t\) | Prob. \((>|t|)\) | |
|---|---|---|---|---|
| (Intercept) | 2.6202 | 0.5678 | 4.6149 | 0.0000 |
| L(LRM, 1) | 0.3192 | 0.1367 | 2.3358 | 0.0247 |
| L(LRM, 2) | 0.5326 | 0.1324 | 4.0239 | 0.0003 |
| L(LRM, 3) | -0.2687 | 0.1021 | -2.6305 | 0.0121 |
| LRY | 0.6728 | 0.1312 | 5.1295 | 0.0000 |
| L(LRY, 1) | -0.2574 | 0.1472 | -1.7491 | 0.0881 |
| IBO | -1.0785 | 0.3217 | -3.3525 | 0.0018 |
| L(IBO, 1) | -0.1062 | 0.5858 | -0.1813 | 0.8571 |
| L(IBO, 2) | 0.2877 | 0.5691 | 0.5055 | 0.6161 |
| L(IBO, 3) | -0.9947 | 0.3925 | -2.5341 | 0.0154 |
| IDE | 0.1255 | 0.5545 | 0.2263 | 0.8222 |
| L(IDE, 1) | -0.3280 | 0.7213 | -0.4547 | 0.6518 |
| L(IDE, 2) | 1.4079 | 0.5520 | 2.5503 | 0.0148 |
2. Modelo de corrección de error no restringido (UECM). La misma información puede reordenarse en la forma de corrección de error, que separa la dinámica de corto plazo (los términos en diferencias) de la relación de largo plazo (los niveles rezagados).
UECM_BestMod <- uecm(BestMod)
knitr::kable(
coef(summary(UECM_BestMod)),
col.names = c("Coef.", "Error Est.", "Estad. $t$", "Prob. $(>|t|)$"),
caption = "Modelo de corrección de error no restringido (UECM) del $ARDL$ estimado.",
digits = 4,
booktabs = TRUE,
escape = FALSE
)| Coef. | Error Est. | Estad. \(t\) | Prob. \((>|t|)\) | |
|---|---|---|---|---|
| (Intercept) | 2.6202 | 0.5678 | 4.6149 | 0.0000 |
| L(LRM, 1) | -0.4169 | 0.0917 | -4.5479 | 0.0001 |
| L(LRY, 1) | 0.4154 | 0.1176 | 3.5317 | 0.0011 |
| L(IBO, 1) | -1.8917 | 0.3911 | -4.8368 | 0.0000 |
| L(IDE, 1) | 1.2053 | 0.4469 | 2.6971 | 0.0103 |
| d(L(LRM, 1)) | -0.2639 | 0.1019 | -2.5898 | 0.0134 |
| d(L(LRM, 2)) | 0.2687 | 0.1021 | 2.6305 | 0.0121 |
| d(LRY) | 0.6728 | 0.1312 | 5.1295 | 0.0000 |
| d(IBO) | -1.0785 | 0.3217 | -3.3525 | 0.0018 |
| d(L(IBO, 1)) | 0.7070 | 0.4687 | 1.5083 | 0.1395 |
| d(L(IBO, 2)) | 0.9947 | 0.3925 | 2.5341 | 0.0154 |
| d(IDE) | 0.1255 | 0.5545 | 0.2263 | 0.8222 |
| d(L(IDE, 1)) | -1.4079 | 0.5520 | -2.5503 | 0.0148 |
3. Modelo de corrección de error restringido (RECM). Al imponer la
restricción de largo plazo se obtiene el término de corrección de error
(ect), cuyo coeficiente mide la velocidad de ajuste hacia el
equilibrio. Nota: se permite que la constante forme parte de la relación
de corto plazo (caso 2), en lugar de la de largo plazo (caso 3).
RECM_BestMod <- recm(UECM_BestMod, case = 2)
knitr::kable(
coef(summary(RECM_BestMod)),
col.names = c("Coef.", "Error Est.", "Estad. $t$", "Prob. $(>|t|)$"),
caption = "Modelo de corrección de error restringido (RECM) del $ARDL$ estimado.",
digits = 4,
booktabs = TRUE,
escape = FALSE
)| Coef. | Error Est. | Estad. \(t\) | Prob. \((>|t|)\) | |
|---|---|---|---|---|
| d(L(LRM, 1)) | -0.2639 | 0.0901 | -2.9301 | 0.0054 |
| d(L(LRM, 2)) | 0.2687 | 0.0913 | 2.9435 | 0.0052 |
| d(LRY) | 0.6728 | 0.1159 | 5.8047 | 0.0000 |
| d(IBO) | -1.0785 | 0.3002 | -3.5921 | 0.0008 |
| d(L(IBO, 1)) | 0.7070 | 0.4436 | 1.5938 | 0.1183 |
| d(L(IBO, 2)) | 0.9947 | 0.3649 | 2.7258 | 0.0092 |
| d(IDE) | 0.1255 | 0.4829 | 0.2598 | 0.7962 |
| d(L(IDE, 1)) | -1.4079 | 0.4887 | -2.8810 | 0.0062 |
| ect | -0.4169 | 0.0785 | -5.3111 | 0.0000 |
El coeficiente del término de corrección de error es negativo y
altamente significativo (ect \(=\)
-0.4169), lo que indica que
aproximadamente el 42% de la desviación respecto de la relación de
largo plazo se corrige cada trimestre. Su signo negativo es la condición
que garantiza que el sistema regrese al equilibrio: si en un trimestre la
demanda de dinero está por encima de su nivel de largo plazo, en los
trimestres siguientes tenderá a caer.
4. Prueba de límites (bounds test). Antes de interpretar la relación de largo plazo debemos verificar que exista.
##
## Bounds F-test (Wald) for no cointegration
##
## data: d(LRM) ~ L(LRM, 1) + L(LRY, 1) + L(IBO, 1) + L(IDE, 1) + d(L(LRM, 1)) + d(L(LRM, 2)) + d(LRY) + d(IBO) + d(L(IBO, 1)) + d(L(IBO, 2)) + d(IDE) + d(L(IDE, 1))
## F = 5.1168, p-value = 0.004418
## alternative hypothesis: Possible cointegration
## null values:
## k T
## 3 1000La estadística \(F\) de la prueba es 5.1168 con un valor \(p\) de 0.0044, por lo que se rechaza la hipótesis nula de ausencia de relación de nivel: la estadística supera la cota superior de los valores críticos de Pesaran, Shin y Smith (2001) para \(k = 3\) regresores. Existe, por lo tanto, evidencia de una relación de largo plazo entre la demanda de dinero, el ingreso real y las dos tasas de interés, y tiene sentido calcular los multiplicadores de largo plazo.
5. Multiplicadores de largo plazo. Se obtienen como el cociente entre la suma de los coeficientes de cada regresor y uno menos la suma de los coeficientes autorregresivos.
mult <- multipliers(BestMod)
knitr::kable(
mult,
col.names = c("Término", "Estimado", "Error Est.", "Estad. $t$", "Prob. $(>|t|)$"),
caption = "Multiplicadores de largo plazo del modelo $ARDL$ estimado.",
digits = 4,
booktabs = TRUE,
escape = FALSE,
row.names = FALSE
)| Término | Estimado | Error Est. | Estad. \(t\) | Prob. \((>|t|)\) |
|---|---|---|---|---|
| (Intercept) | 6.2857 | 0.7719 | 8.1429 | 0.000 |
| LRY | 0.9965 | 0.1239 | 8.0405 | 0.000 |
| IBO | -4.5381 | 0.5203 | -8.7222 | 0.000 |
| IDE | 2.8915 | 0.9951 | 2.9058 | 0.006 |
Los tres multiplicadores son significativos y sus signos son los que predice la teoría de la demanda de dinero. La elasticidad ingreso de largo plazo es prácticamente unitaria (0.996), resultado consistente con la teoría cuantitativa: en el largo plazo, la demanda real de dinero crece en la misma proporción que el ingreso real. El coeficiente de la tasa de los bonos es negativo (-4.538), como corresponde al costo de oportunidad de mantener dinero, mientras que el de la tasa de los depósitos bancarios es positivo (2.892), pues una mayor remuneración de los depósitos incrementa la demanda del agregado M2, que los incluye. Nótese que, al estar las tasas expresadas en fracciones y las demás variables en logaritmos, estos dos coeficientes son semielasticidades y no elasticidades.
5.6.2.4 Gráfica de la ecuación de cointegración
La Figura 5.15 compara la serie observada de \(LRM_t\) con la relación de largo plazo estimada: la distancia entre ambas líneas es precisamente el término de corrección de error cuyo coeficiente interpretamos arriba.
Datos <- cbind.zoo(LRM = denmark[,"LRM"], Result)
Datos <- xts(Datos)
plot(Datos, legend.loc = "right")
Figura 5.15: Demanda real de dinero observada (\(LRM_t\)) y relación de cointegración estimada por el modelo ARDL
5.7 Filtro Kalman
En secciones pasadas, entendimos la naturaleza del filtrado univariado como, por ejemplo, el Filtro Hodrick-Prescott. Este ha sido utilizado ampliamente en la literatura macroeconométrica, sobre todo para la estimación del Producto Potencial y, de este modo, la Brecha del Producto. Sin embargo, este método presenta dos limitaciones importantes: 1) es univariado, es decir, la tendencia de cada serie se estima únicamente con la información de esa serie, sin aprovechar la información de otras variables relacionadas y, 2) el problema de fin de muestra discutido en el Capítulo 3, que la variante de St-Amant y van Norden (1997) atenúa a costa de introducir dos parámetros adicionales (\(u_{ss}\) y \(\lambda_{ss}\)) cuya elección es discrecional.
De este modo, el Filtro Kalman puede ser una alternativa viable y metodológicamente robusta, sobre todo si queremos trabajar con variables no observables como el Producto Potencial o, incluso, la Tasa de interés neutral, variable clave en la decisión de Política Monetaria de los Bancos Centrales.
En este sentido, el Filtro Kalman nació a partir del artículo seminal de Rudolf E. Kálmán (1960) A New Approach to Linear Filtering and Prediction Problems. La idea principal de Kálmán era describir los procesos estocásticos con un modelo de “estado” y resolver, de alguna manera, el problema de estimación con un ciclo predicción-corrección ejecutado en tiempo real. Así, el Filtro Kalman es un algoritmo de estimación recursivo a partir de mediciones ruidosas. De este modo, también podemos entender al Filtro Kalman como una especie de filtro multivariado. Sus aplicaciones no se han limitado únicamente a la Economía, pues ha sido utilizado para trabajos de la NASA como el pronóstico de la posición y velocidad del Apollo 11 en 1969.
5.7.1 Modelos Estado-Espacio y Teoría Económica
Para hacer uso del Filtro Kalman, necesitamos de un modelo estructural que vincule a nuestras variables (observables y no observables), representado en la Forma Estado-Espacio. En general, casi todos los modelos de series de tiempo pueden ser representados como Modelos Estado-Espacio. El modelo estructural puede basarse en la teoría económica como, por ejemplo, un modelo DSGE. Para efectos prácticos, tomaremos como ejemplo el Modelo de Proyección Trimestral (QPM) del IMF (2020). No nos enfocaremos en el por qué y cómo del Modelo QPM, pues el objetivo principal es entender el Filtro Kalman y su aplicación.
5.7.1.1 Modelo de Proyección Trimestral (QPM) del IMF (2020)
Supongamos una economía pequeña y abierta, en donde:
Demanda Agregada (Curva IS) \[\begin{equation} \hat{y}_t = b_1 \hat{y}_{t-1} - b_2 mci_t + b_3 \hat{y}_t^* + \varepsilon_t^y \end{equation}\] \[\begin{equation} mci_t = b_4 \hat{r}_t + (1-b_4)(-\hat{z}_t) \end{equation}\]
donde:
\(\hat{y}_t =\) Brecha del Producto Doméstico,
\(\hat{y}_{t-1} =\) Brecha del Producto Doméstico rezagada,
\(mci_t =\) Índice de Condiciones Monetarias Real,
\(\hat{y}_t^* =\) Brecha del Producto Foráneo,
\(\varepsilon_t^y =\) Choque de Demanda,
\(\hat{r}_t =\) Brecha de la tasa de interés real y
\(\hat{z}_t =\) Brecha del tipo de cambio real.
Inflación (Curva de Phillips) \[\begin{equation} \pi_t = a_1 \pi_{t-1} + (1-a_1) E\{\pi_{t+1}\} + a_2 rmc_t + \varepsilon_t^\pi \end{equation}\] \[\begin{equation} rmc_t = a_3 \hat{y}_t + (1-a_3) \hat{z}_t \end{equation}\]
donde:
\(\pi_t =\) Inflación,
\(\pi_{t-1} =\) Inflación rezagada,
\(E\{\pi_{t+1}\} =\) Expectativa de inflación,
\(rmc_t =\) Costo Marginal Real,
\(\epsilon_t^\pi =\) Choque de Oferta (cost-push shock),
\(\hat{y}_t =\) Brecha del Producto Doméstica y
\(\hat{z}_t =\) Brecha del tipo de cambio real.
Paridad Descubierta de Tasas de Interés (UIP) \[\begin{equation} S_t = (1 - e_1) \mathbb{E}_t \{ S_{t+1} \} + e_1 \left[ S_{t-1} + \frac{2(\pi_t - \pi_t^{*} + \Delta \bar{z}_t)}{4} + \frac{i_t^{*} - i_t + prem_t}{4} \right] + \varepsilon_s \end{equation}\]
donde:
\(\mathbb{E}_t \{ S_{t+1}\} =\) Expectativas sobre el tipo de cambio nominal futuro,
\(S_{t-1} =\) Tipo de cambio nominal rezagado,
\(\pi_t =\) Inflación,
\(\pi_t^* =\) Inflación foránea,
\(\Delta \bar{z}_t =\) Depreciación del tipo de cambio tendencial,
\(i_t^* =\) Tasa de interés real foránea,
\(i_t =\) Tasa de interés real doméstica,
\(prem_t =\) Premio Riesgo-País y
\(\varepsilon_t =\) Choque del Tipo de cambio.
Función de Reacción de Política Monetaria (Regla de Taylor) \[\begin{equation} i_t = g_1 i_{t-1} + (1 - g_1) \left\{ i_t^n + g_2 \left( \mathbb{E}_t [\pi_{t+N}^4] - \pi_{t+N}^T \right) + g_3 \hat{y}_t \right\} + \varepsilon_t \end{equation}\]
donde:
\(i_t =\) Tasa de interés (o política),
\(i_{t-1} =\) Tasa de interés rezagada,
\(i_t^n =\) Tasa de interés natural o neutral,
\(\mathbb{E}_t [\pi_{t+N}^4] =\) Inflación esperada (YoY) en el periodo \(_{t+N}\),
\(\pi_{t+N}^T =\) Meta de inflación,
\(\hat{y}_t =\) Brecha del Producto Doméstica y
\(\varepsilon_t =\) Choque de política monetaria no sistémico.
Bloque externo
En este caso, todas las variables siguen procesos autorregresivos simples
\[ \begin{array}{l} \text{\textbf{Brecha del producto foráneo}} \\ \hat{y}^{*}_{t} = \rho_{y}\,\hat{y}^{*}_{t-1} + \varepsilon^{\hat{y}^{*}}_{t} \\[1.5ex] \text{\textbf{Tasa natural (tendencia) de interés real foránea}} \\ \tilde{r}^{*}_{t} = \rho_{\tilde r}\,\tilde{r}^{*}_{t-1} + (1 - \rho_{\tilde r})\,\tilde{r}^{*\,\mathrm{SS}} + \varepsilon^{\tilde{r}^{*}}_{t} \\[1.5ex] \text{\textbf{Tasa de interés nominal foránea}} \\ i^{*}_{t} = \rho_{i}\,i^{*}_{t-1} + (1 - \rho_{i})(r^{*}_{t} + \pi^{*}_{t}) + \varepsilon^{i^{*}}_{t} \\[1.5ex] \text{\textbf{Brecha de la tasa de interés real foránea}} \\ \hat{r}^{*}_{t} = r^{*}_{t} - \tilde{r}^{*}_{t} \\[1.5ex] \text{\textbf{Tasa de interés real foránea}} \\ r^{*}_{t} = i^{*}_{t} - \pi^{*}_{t} \\[1.5ex] \text{\textbf{Inflación foránea}} \\ \pi^{*}_{t} = \rho_{\pi}\,\pi^{*}_{t-1} + (1 - \rho_{\pi})\,\pi^{*\,\mathrm{SS}} + \varepsilon^{\pi^{*}}_{t} \end{array} \]
5.7.1.2 Modelos Estado-Espacio
La idea general detrás de los Modelos Estado-Espacio es que una (múltiples) variable observable \(Y_1,...Y_T\) depende de un posible estado no observable \(X_t\) que sigue un proceso estocástico. La relación entre \(Y_t\) y \(X_t\) es descrita por la Ecuación de Observación
\[\begin{equation} Y_t = C X_t + v_t \tag{5.11} \end{equation}\]
donde \(C\) es una matriz de parámetros conocidos y \(v_t\) es un vector de errores observados (si los hay), los cuales son independientes e idénticamente distribuidos con media cero y matriz de covarianza \(R\), es decir, \(v_t \sim iid(0, R)\). De este modo, la ecuación (5.11) describe la relación estática entre las variables observables y no observables.
La ecuación que describe la dinámica de las variables de estado, también conocida como Ecuación de Transición, es definida como
\[\begin{equation} X_t = A X_{t-1} + \varepsilon_t, \tag{5.12} \end{equation}\]
en donde \(A\) es una matriz de parámetros conocidos, \(\varepsilon_t\) es un vector de errores, los cuales son independientes e idénticamente distribuidos con media cero y matriz de covarianza \(Q\), es decir, \(\varepsilon_t \sim iid(0,Q)\). Adicionalmente, se asume que \(v_t\) y \(\varepsilon_t\) son independientes entre sí, es decir, \(Cov(v_t,\varepsilon_t) = 0\).
La Figura 5.16 resume la estructura: los estados evolucionan por su cuenta según la matriz \(A\), cada uno recibe su propio choque \(\varepsilon_t\), y lo único que el investigador observa es la proyección \(C X_t\) contaminada por el error de medición \(v_t\).
Figura 5.16: Estructura de un modelo estado-espacio. La fila superior son los estados no observables, que evolucionan según la ecuación de transición; la inferior, las variables observables, ligadas a los estados por la ecuación de observación
Afortunadamente, el Modelo QPM está construido de tal manera que sus ecuaciones de transición y de observación vinculan las variables observables con las no observables. El modelo completo –con sus 34 variables de transición, 12 choques y las ecuaciones de comportamiento e identidades– se declara en un archivo de texto que el FMI distribuye con su curso; a modo de ilustración, reproducimos únicamente el bloque que corresponde a las cuatro ecuaciones de comportamiento que presentamos arriba:10
%% === Demanda agregada (curva IS) ===
L_GDP_GAP = b1*L_GDP_GAP{-1} - b2*MCI + b3*L_GDP_RW_GAP + SHK_L_GDP_GAP;
MCI = b4*RR_GAP + (1-b4)*(- L_Z_GAP);
%% === Inflación (curva de Phillips) ===
DLA_CPI = a1*DLA_CPI{-1} + (1-a1)*DLA_CPI{+1} + a2*RMC + SHK_DLA_CPI;
RMC = a3*L_GDP_GAP + (1-a3)*L_Z_GAP;
%% === Regla de política monetaria (tipo Taylor, prospectiva) ===
RS = g1*RS{-1} + (1-g1)*(RSNEUTRAL + g2*(D4L_CPI{+4} - D4L_CPI_TAR{+4})
+ g3*L_GDP_GAP) + SHK_RS;
%% === Paridad descubierta de tasas de interés (UIP) ===
L_S = (1-e1)*L_S{+1} + e1*(L_S{-1} + 2/4*(D4L_CPI_TAR - ss_DLA_CPI_RW
+ DLA_Z_BAR)) + (- RS + RS_RW + PREM)/4 + SHK_L_S;
%% Convenciones de nomenclatura del modelo:
%% _GAP desviación cíclica respecto de la tendencia
%% _BAR tendencia (equilibrio) ss_ valor de estado estacionario
%% DLA_ cambio trimestral D4L_ cambio anual
%% _RW variable externa SHK_ residual de la ecuaciónNótese la correspondencia directa entre este bloque y las ecuaciones que
escribimos en la sección anterior: L_GDP_GAP es la brecha del producto
\(\hat{y}_t\), MCI es el índice de condiciones monetarias \(mci_t\),
DLA_CPI es la inflación \(\pi_t\) y RS es la tasa de política \(i_t\). La
notación {-1} y {+1} indica rezagos y adelantos, y es justamente la
presencia de adelantos –las expectativas– lo que impide escribir el
modelo directamente en la forma estado-espacio y obliga al paso que
describimos a continuación.
Sin embargo, la Ecuación de Transición de los estados presentes dependen de los estados pasados, mientras que el Modelo QPM es prospectivo, es decir, las variables dependen de sus expectativas así que la ecuación toma la forma de
\[\begin{equation} F \mathbb{E}_t[X_{t+1}] + G X_t + H X_{t-1} + M \eta_t = 0 \tag{5.13} \end{equation}\]
donde \(F\), \(G\), \(H\), \(M\) son matrices de parámetros estructurales conocidos, \(\mathbb{E}_t\) es el operador de expectativas y \(\eta_t\) es un vector de choques estructurales. Para pasar de este modelo a la Forma Estado-Espacio, tenemos que resolver para Expectativas Racionales. Para esto, proponemos una hipótesis de solución
\[\begin{equation} X_t = A X_{t-1} + P \eta_t \quad \text{con} \quad \varepsilon_t \equiv P \eta_t \tag{5.14} \end{equation}\]
Ahora, suponemos que \(\eta_t\) sigue un proceso autorregresivo, es decir, tiene su propia dinámica:
\[\begin{equation} \eta_{t+1} = S \eta_t + \gamma_{t+1} \quad \text{donde} \quad \gamma_{t+1} \sim iid(0,W) \end{equation}\]
Así, sabemos que
\[\begin{equation} \mathbb{E}_t[X_{t+1}] = A X_t + P \mathbb{E}_t[\eta_{t+1}] = A (A X_{t-1} + P \eta_t) + P S \eta_t \end{equation}\]
De este modo, si sustituimos la ecuación (5.14) en la (5.13), tenemos que las matrices deben de satisfacer la identidad
\[\begin{equation} F [A (A X_{t-1} + P \eta_t) + P (S \eta_t)] + G [A X_{t-1} + P \eta_t] + H X_{t-1} + M \eta_t = 0 \tag{5.15} \end{equation}\]
para todos los valores \(X_{t-1}\) y \(\eta_t\), de tal modo que \(A\) y \(P\) son las soluciones de ecuaciones
\[\begin{equation} F A^2 + G A + H = 0, \quad \quad F A P + F P S + G P + M = 0 \tag{5.16} \end{equation}\]
De esta forma, obtenemos nuestro modelo en la Forma Estado-Espacio.
5.7.2 Recursión del Filtro Kalman
Ya que tenemos nuestro modelo en la Forma Estado-Espacio como se muestra en las ecuaciones (5.11) y (5.12), podemos hacer uso del Filtro Kalman.
Para comenzar, definamos \(Y_{t/t-1}\) y \(X_{t/t-1}\) como los valores predichos de \(Y_t\) y \(X_t\), respectivamente, condicionado a la información del periodo \(_{t-1}\). Además, definamos \(X_{t/t}\) como la estimación de \(X_t\) condicionada a la información del periodo \(_t\).
Como se mencionó en un inicio, el Filtro Kalman utiliza un algoritmo recursivo para cada periodo \(1,2,...,T\). Cada recursión está compuesta de dos pasos. Supongamos que tenemos los datos observados y las estimaciones del estado para el periodo \(_{t-1}\), es decir, conocemos \(Y_{t-1}\) y \(X_{t-1/t-1}\). De este modo, para el periodo \(_t\):
- Paso 1. Predicción del estado y observación de las variables en el periodo \(_t\) dada la información de \(_{t-1}\). Como no tenemos nueva información, el promedio de la predicción se basa en el Modelo Estado-Espacio (5.11) - (5.12):
\[\begin{equation} X_{t/t-1} = A X_{t-1/t-1} \tag{5.17} \end{equation}\]
\[\begin{equation} Y_{t/t-1} = C X_{t/t-1} \tag{5.18} \end{equation}\]
- Paso 2. Actualización de la estimación del estado dada la información en el periodo \(_t\). Una vez que conocemos la nueva información del periodo \(_t\), existirá una discrepancia entre la información actual (\(Y_t\)) y la predicción del modelo (\(Y_{t/t-1}\)). El Filtro Kalman utiliza dicha discrepancia para actualizar la estimación de los estados no observables utilizando la matriz \(K_t\), que se llama Ganancia de Kalman:
\[\begin{equation} X_{t/t} = X_{t/t-1} + K_t (Y_t - Y_{t/t-1}) \tag{5.19} \end{equation}\]
5.7.3 Ganancia de Kalman
Supongamos que los errores siguen una distribución normal, es decir, \(v_t \sim N(0,R), \varepsilon_t \sim N(0,Q)\).
Para la primera estimación (Paso 1) de la Recursión de Kalman, \(P_{t/t-1}\) denota la varianza del error de predicción de las variables de estado (\(X_t\)) dada la información del periodo \(_{t-1}\), o sea, \(P_{t/t-1} = Var_{t-1}[(X_t - X_{t/t-1})]\). De manera similar, \(F_{t/t-1}\) es la varianza del error de predicción de las variables observables (\(Y_t\)) dada la información del periodo \(_{t-1}\), o bien, \(F_{t/t-1} = Var_{t-1}[(Y_t - Y_{t/t-1})]\).
Para la actualización (Paso 2) de la Recursión de Kalman, \(P_{t/t}\) denota la varianza del error de las variables de estado (\(X_t\)) dada la información del periodo \(_t\), es decir, \(P_{t/t} = Var_t[(X_t - X_{t/t})]\).
Así, la distribución multivariada del vector \((Y_t,X_t)\) toma la forma
\[\begin{equation} \begin{pmatrix} Y_t \\ X_t \end{pmatrix} \Bigg| Y_t,...,Y_{t-1} \sim N \Bigg(\begin{pmatrix} Y_{t/t-1} \\ X_{t/t-1} \end{pmatrix}, \begin{pmatrix} F_{t/t-1} & C P_{t/t-1} \\ P_{t/t-1} C' & P_{t/t-1} \end{pmatrix}\Bigg). \tag{5.20} \end{equation}\]
De (5.20), la distribución condicional de \(X_t\) dado \(Y_t\) toma la forma de
\[\begin{equation} X_t \mid Y_t, Y_1, ... , Y_{t-1} \sim N \begin{aligned} ( & X_{t|t-1} + P_{t|t-1} C' F_{t|t-1} {}^{-1} [Y_t - Y_{t|t-1}], \\ & P_{t|t-1} - P_{t|t-1} C' F_{t|t-1} {}^{-1} C P_{t|t-1} ) \end{aligned} \tag{5.21} \end{equation}\]
Y de la ecuación (5.21) obtenemos directamente que
\[\begin{equation} X_{t|t} = X_{t|t-1} + P_{t|t-1} C' F_{t|t-1} {}^{-1} \left[ Y_t - Y_{t|t-1} \right] = X_{t|t-1} + K_t \left[ Y_t - Y_{t|t-1} \right], \tag{5.22} \end{equation}\]
\[\begin{equation} P_{t|t} = P_{t|t-1} - P_{t|t-1} C' F_{t|t-1} {}^{-1} C P_{t|t-1} = (I - K_t C) P_{t|t-1}, \tag{5.23} \end{equation}\]
en donde la Ganancia de Kalman es
\[\begin{equation} K_t = P_{t/t-1} C' F_{t/t-1} {}^{-1}. \tag{5.24} \end{equation}\]
Nótese que \(F_{t/t-1}\) sigue:
\[\begin{equation} F_{t|t-1} = Var_{t-1} \left[( Y_t - Y_{t|t-1} \right)] = Var_{t-1} \left[( C X_t + v_t - C X_{t|t-1} \right)] = C P_{t|t-1} C' + R \end{equation}\]
Por el contrario, \(P_{t/t-1}\) sigue:
\[\begin{align} P_{t|t-1} &= Var_{t-1} \left[( X_t - X_{t|t-1} \right)] \\ &= Var_{t-1} \left[( A X_{t-1} + \varepsilon_t - A X_{t-1|t-1} \right)] \\ &= A P_{t-1|t-1} A' + Q \\ &= A (I - K_{t-1} C) P_{t-1|t-2} A' + Q \end{align}\]
Y si sustituimos \(F_{t/t-1}\) y \(P_{t/t-1}\) en las ecuaciones (5.22) - (5.24), obtenemos las siguientes ecuaciones para el algoritmo del Filtro Kalman
\[\begin{equation} K_t = P_{t/t-1} C' (C P_{t/t-1} C' + R)^{-1} \tag{5.25} \end{equation}\]
\[\begin{equation} X_{t/t} = X_{t/t-1} + K_t [Y_t - Y_{t/t-1}] \tag{5.26} \end{equation}\]
\[\begin{equation} P_{t+1|t} = A (I - K_t C) P_{t/t-1} A' + Q \tag{5.27} \end{equation}\]
De este modo, la Ganancia de Kalman en (5.25) depende de manera positiva de la varianza del error de predicción del estado (\(P_{t/t-1}\)) y de manera negativa de la varianza de la ecuación de observación \(R\). Esto es bastante intuitivo:
Si existen grandes errores en las predicciones del estado \(X_{t/t-1}\) en el Paso 1 de la Recursión de Kalman, es decir, \(P_{t/t-1}\) es grande, se le da más peso a la nueva información observada \(Y_t\), lo cual significa que la Ganancia de Kalman es grande. En particular, esto es cierto para las ecuaciones de transición que no se cumplen de manera estricta, es decir, \(Q\) es grande;
Si la información es ruidosa, o sea, \(R\) es grande, se le da menos peso a la información nueva por lo que \(K_t\) es pequeño.
5.7.4 Algoritmo del Filtro Kalman y del Suavizador de Kalman
Ya que encontramos las expresiones más relevantes, podemos combinar todas para entender cómo funciona el Filtro Kalman.
Sabemos que el filtro es recursivo. De este modo, para poder correrlo necesitamos condiciones iniciales para todos los estados al principio de la muestra (\(X_{0|0}\)), al igual que para la varianza del error de predicción del estado (\(P_{0|0}\)). Para muestras pequeñas, las condiciones iniciales pueden afectar significativamente los resultados del filtro. Si no tenemos información previa, una opción es asumir que las condiciones iniciales de los estados y de la varianza del error de predicción no están lejos de sus valores en estado estacionario. En particular, para modelos estacionarios, se puede establecer que las condiciones iniciales son iguales a los valores en estado estacionario.
Ahora, ya que hemos asignado las condiciones iniciales, el Filtro Kalman corre el algoritmo recursivo para cada periodo de tiempo \(1,2,...,T\). Cada recursión está compuesta de dos pasos:
- Paso 1. Predicción utilizando el Modelo Estado-Espacio:
\[\begin{equation} X_{t/t-1} = A X_{t-1/t-1}, \tag{5.28} \end{equation}\]
\[\begin{equation} Y_{t/t-1} = C X_{t/t-1}. \tag{5.29} \end{equation}\]
- Paso 2. Actualización utilizando la información:
\[\begin{equation} K_t = P_{t/t-1} C' (C P_{t/t-1} C' + R)^{-1}, \tag{5.30} \end{equation}\]
\[\begin{equation} X_{t/t} = X_{t/t-1} + K_t [Y_t - Y_{t/t-1}], \tag{5.31} \end{equation}\]
\[\begin{equation} P_{t+1/t} = A (I - K_t C) P_{t/t-1} A' + Q. \tag{5.32} \end{equation}\]
El Filtro Kalman nos proporciona el estimador lineal óptimo de mínimos cuadrados del estado, dada la información hasta el periodo \(_t\). Si los errores están distribuidos normalmente, es óptimo entre todos los estimadores de \(X_t\) condicionados a la información disponible hasta ese periodo.
Para cada periodo \(_t\), el Filtro Kalman utiliza únicamente la información disponible hasta ese periodo para estimar \(X_t\). Dichas estimaciones se encuentran en \(X_{t/t-1}\) (del Paso 1) y en \(X_{t/t}\) (del Paso 2). Sin embargo, para estimar \(X_t\), comúnmente es deseable utilizar toda la información disponible, hasta el último periodo \(T\). Es decir, podría ser deseable obtener estimaciones \(X_{t/T}\), que se conocen como estimaciones suavizadas de los estados. Precisamente eso es lo que hace el Suavizador de Kalman (Kalman smoother).
Para correr el Suavizador de Kalman, necesitamos primero correr el Filtro Kalman y obtener todas las estimaciones de los estados: \(X_{1/1}, X_{2/2},...,X_{T/T}\).
Después, el algoritmo del Suavizador de Kalman funciona recursivamente hacia atrás, desde \(T-1, T-2,...,1\). Para cada periodo \(_{t+1}\), observamos la diferencia entre las estimaciones suavizadas y las estimaciones filtradas de los estados: \(X_{t+1/T}-X_{t+1/t}\). El Suavizador de Kalman utiliza la diferencia para actualizar las estimaciones filtradas de los estados en el periodo \(_t\), con el fin de obtener las estimaciones suavizadas:
\[\begin{equation} X_{t/T} = X_{t/t} + J_t [X_{t+1/T} - X_{t+1/t}]. \tag{5.33} \end{equation}\]
La matriz \(J_t\) pondera cuánto se corrige la estimación filtrada con la información posterior y está dada por
\[\begin{equation} J_t = P_{t/t} A' P_{t+1/t}^{-1}, \tag{5.34} \end{equation}\]
es decir, la corrección es mayor cuanto más incierta es la estimación filtrada del estado en \(_t\) (\(P_{t/t}\) grande) y menor cuanto más incierta es la predicción del estado en \(_{t+1}\) (\(P_{t+1/t}\) grande). Nótese la analogía con la Ganancia de Kalman de la ecuación (5.30): ambas matrices reparten el peso entre la información propia y la nueva evidencia según sus varianzas relativas.
5.7.5 Ejemplo: la brecha del producto de la economía mexicana
Cerremos el capítulo aplicando la recursión que acabamos de derivar. El problema es el que motivó toda la sección: separar una serie observada en una tendencia y un ciclo, ninguno de los cuales se observa directamente. Tomemos el modelo de tendencia local lineal, que escribe el nivel de la serie como una tendencia \(\mu_t\) con pendiente \(\beta_t\) más un componente transitorio:
\[\begin{eqnarray} y_t & = & \mu_t + \varepsilon_t \nonumber \\ \mu_t & = & \mu_{t-1} + \beta_{t-1} \nonumber \\ \beta_t & = & \beta_{t-1} + \zeta_t \tag{5.35} \end{eqnarray}\]
Con el vector de estado \(X_t = (\mu_t, \beta_t)'\), este modelo tiene exactamente la forma de las ecuaciones (5.12) y (5.11):
\[ A = \begin{pmatrix} 1 & 1 \\ 0 & 1 \end{pmatrix}, \quad C = \begin{pmatrix} 1 & 0 \end{pmatrix}, \quad Q = \begin{pmatrix} 0 & 0 \\ 0 & \sigma^2_\zeta \end{pmatrix}, \quad R = \sigma^2_\varepsilon \]
El código siguiente implementa la recursión tal como la escribimos en las ecuaciones (5.28) a (5.32), seguida del suavizador de la ecuación (5.33). Son unas cuantas líneas, y vale la pena leerlas junto a las fórmulas: cada renglón corresponde a un paso del algoritmo.
filtro_kalman <- function(y, A, C, Q, R, X0, P0) {
n <- length(y); k <- nrow(A)
Xf <- matrix(NA, k, n); Pf <- array(NA, c(k, k, n)) # filtrado
Xp <- matrix(NA, k, n); Pp <- array(NA, c(k, k, n)) # predicción
X <- X0; P <- P0
for (t in seq_len(n)) {
# Paso 1: predicción
X <- A %*% X
P <- A %*% P %*% t(A) + Q
Xp[, t] <- X; Pp[, , t] <- P
# Paso 2: actualización con la nueva observación
F <- C %*% P %*% t(C) + R # varianza del error de predicción
K <- P %*% t(C) %*% solve(F) # ganancia de Kalman
X <- X + K %*% (y[t] - C %*% X) # corrección del estado
P <- P - K %*% C %*% P
Xf[, t] <- X; Pf[, , t] <- P
}
# Suavizador: recursión hacia atrás
Xs <- Xf
for (t in (n - 1):1) {
J <- Pf[, , t] %*% t(A) %*% solve(Pp[, , t + 1])
Xs[, t] <- Xf[, t] + J %*% (Xs[, t + 1] - Xp[, t + 1])
}
list(filtrado = Xf, suavizado = Xs)
}Aplicamos el filtro al logaritmo del IGAE desestacionalizado (escalado por 100, de modo que las desviaciones se leen en por ciento). El único parámetro por fijar es la razón señal-ruido \(q = \sigma^2_\zeta / \sigma^2_\varepsilon\), que determina qué tan flexible es la tendencia.
load("BD/Datos_Ad.RData")
y_igae <- 100*log(ts(Datos_Ad$IGAE_Ad, start = c(2000, 1), frequency = 12))
lambda <- 14400 # mismo valor que usamos con el filtro HP mensual
A <- matrix(c(1, 0, 1, 1), 2, 2)
C <- matrix(c(1, 0), 1, 2)
Q <- matrix(c(0, 0, 0, 1/lambda), 2, 2)
R <- matrix(1, 1, 1)
kf <- filtro_kalman(y_igae, A, C, Q, R,
X0 = matrix(c(y_igae[1], 0), 2, 1),
P0 = diag(1e6, 2))
tendencia <- ts(kf$suavizado[1, ], start = start(y_igae), frequency = 12)
brecha <- y_igae - tendenciaEl resultado conecta este capítulo con el Capítulo
3 de una manera que conviene subrayar: el filtro
de Hodrick-Prescott es el suavizador de Kalman de este modelo, con la
correspondencia \(\lambda = 1/q\). No es una analogía, es una identidad;
podemos comprobarla numéricamente comparando la tendencia que acabamos de
suavizar con la que produce hpfilter():
## [1] 0.00000005402188La diferencia entre ambas es del orden de \(10^{-8}\), es decir, error de redondeo. Visto así, el filtro HP deja de ser una receta y se vuelve un caso particular de un marco mucho más general: cambiando \(A\), \(C\), \(Q\) y \(R\) se pueden estimar tendencias con otras dinámicas, incorporar varias series observables para identificar un mismo estado no observable, o agregar ecuaciones de comportamiento como las del QPM.
La Figura 5.17 muestra la descomposición resultante.
par(mfrow = c(2, 1), mar = c(4, 4, 3, 1))
plot(y_igae, col = "gray40", lwd = 1, xlab = "Tiempo",
ylab = "100 x log(IGAE)",
main = "IGAE desestacionalizado y tendencia estimada")
lines(tendencia, col = "darkblue", lwd = 2)
legend("topleft", c("Observado", "Tendencia (suavizador de Kalman)"),
col = c("gray40", "darkblue"), lwd = c(1, 2), cex = 0.8, bty = "n")
plot(brecha, col = "darkred", lwd = 1.2, xlab = "Tiempo",
ylab = "Desviación (%)",
main = "Brecha del producto estimada")
abline(h = 0, lty = 2, col = "gray40")
Figura 5.17: Descomposición del IGAE desestacionalizado en tendencia y brecha mediante el suavizador de Kalman aplicado al modelo de tendencia local lineal de la ecuación (5.35)
La lectura económica es inmediata. La brecha recoge las dos recesiones del periodo: la de 2008-2009 y, sobre todo, el desplome de la pandemia, con un mínimo de -22.5% en mayo de 2020. Conviene cerrar con la advertencia que ya conocemos del Capítulo 3 y que aquí reaparece: como el suavizador usa toda la muestra, las estimaciones del final se revisan cuando llegan datos nuevos. Es el mismo problema de fin de muestra del filtro HP, y es la razón por la que los bancos centrales estiman la brecha con modelos multivariados –como el QPM– en lugar de con un filtro univariado.
Nota: Para este ejemplo, se ocupa el QPM básico presentado por el IMF (2020) en su curso “Monetary Policy Analysis and Forecasting” disponible en la plataforma edX. El código de Matlab fue realizado por el IMF.
5.7.6 Representación Estado-Espacio de modelos ARMA
Aunque motivamos el Filtro Kalman con un modelo macroeconómico estructural, la Forma Estado-Espacio es mucho más general: prácticamente todos los modelos estudiados en este libro pueden escribirse en ella. Como ilustración, consideremos el proceso \(AR(2)\) del Capítulo 3:
\[\begin{equation} Y_t = a_1 Y_{t-1} + a_2 Y_{t-2} + U_t \end{equation}\]
Si definimos el vector de estado \(X_t = (Y_t, Y_{t-1})'\), el proceso se escribe como
\[\begin{equation} X_t = \begin{pmatrix} a_1 & a_2 \\ 1 & 0 \end{pmatrix} X_{t-1} + \begin{pmatrix} U_t \\ 0 \end{pmatrix}, \quad \quad Y_t = \begin{pmatrix} 1 & 0 \end{pmatrix} X_t, \tag{5.36} \end{equation}\]
que es exactamente la estructura de las ecuaciones (5.12) y (5.11), con \(v_t = 0\) (no hay error de medición). Los componentes de medias móviles se incorporan ampliando el vector de estado con los choques rezagados. Esta representación es la base de, por ejemplo, la función arima() de R, que estima modelos \(ARMA\) por máxima verosimilitud exacta precisamente vía el Filtro Kalman; además, permite manejar observaciones faltantes de forma natural, pues en los periodos sin dato simplemente se omite el paso de actualización.
5.7.7 Estimación por Máxima Verosimilitud
Hasta ahora asumimos que las matrices del sistema (\(A\), \(C\), \(Q\), \(R\)) son conocidas. En la práctica, éstas dependen de un vector de parámetros \(\theta\) que debe estimarse. El propio Filtro Kalman proporciona la herramienta para hacerlo mediante la llamada descomposición del error de predicción. Definamos la innovación del periodo \(_t\) como
\[\begin{equation} u_t = Y_t - Y_{t/t-1}, \end{equation}\]
es decir, la parte de \(Y_t\) que el modelo no pudo anticipar. Bajo el supuesto de normalidad de los errores, \(u_t \sim N(0, F_{t/t-1})\) y las innovaciones no están correlacionadas entre sí, por lo que la función log-verosimilitud de la muestra completa se puede escribir como la suma
\[\begin{equation} ln L(\theta) = -\frac{1}{2} \sum_{t=1}^{T} \left[ n \, ln(2 \pi) + ln |F_{t/t-1}| + u_t' F_{t/t-1}^{-1} u_t \right], \tag{5.37} \end{equation}\]
donde \(n\) es el número de variables observables. Así, para cada valor candidato de \(\theta\) se corre el Filtro Kalman, se acumulan las innovaciones \(u_t\) y sus varianzas \(F_{t/t-1}\), y se evalúa la ecuación (5.37); un algoritmo numérico de optimización busca entonces el \(\theta\) que maximiza \(ln L(\theta)\). Éste es el procedimiento que utilizan los paquetes de cómputo que estiman modelos estado-espacio (en R, por ejemplo, StructTS(), arima() o el paquete KFAS; en Matlab, el toolbox utilizado por el IMF para el QPM). Cuando algún estado es no estacionario (por ejemplo, una tendencia estocástica como el Producto Potencial), las condiciones iniciales se tratan con una inicialización difusa: se asigna a \(P_{0/0}\) una varianza muy grande para reflejar la ignorancia inicial sobre el nivel del estado.
5.8 Resumen del capítulo
En este capítulo extendimos el análisis univariado a sistemas multivariados, introduciendo herramientas para capturar la dinámica conjunta de varias series de tiempo.
Causalidad de Granger. Dos series \(X\) y \(Y\) exhiben causalidad de Granger de \(X\) a \(Y\) cuando la información pasada de \(X\) mejora la predicción de \(Y\) más allá de lo que aporta la historia de \(Y\) sola. Es una noción estadística de precedencia temporal, no de causalidad económica en sentido estricto. El test se implementa comparando las varianzas del error de dos regresiones auxiliares mediante una prueba \(F\).
Modelos VAR(p). El VAR(p) generaliza el AR(p) a sistemas de \(k\) variables, donde cada variable depende de sus propios rezagos y de los rezagos de todas las demás. La condición de estacionariedad requiere que todas las raíces del polinomio matricial característico caigan fuera del círculo unitario. La selección del orden \(p\) se realiza mediante los criterios AIC, SC, HQ y FPE. Las funciones de impulso-respuesta (IRF) cuantifican la reacción dinámica del sistema ante un choque unitario en una de las innovaciones.
VAR Estructural (SVAR). La identificación de los choques estructurales requiere imponer restricciones adicionales. Los tres enfoques principales son: (i) descomposición de Cholesky (identificación recursiva), que impone una estructura triangular inferior en la matriz de impacto contemporáneo; (ii) restricciones de largo plazo de Blanchard y Quah (1989), que restringen el efecto acumulado permanente de ciertos choques; y (iii) restricciones de signo, que imponen solo la dirección del impacto.
Cointegración y VEC. Cuando series no estacionarias \(I(1)\) comparten una tendencia estocástica común, se dice que cointegran. El procedimiento de Johansen determina el número de vectores de cointegración (rango de cointegración) mediante las estadísticas de la traza y del máximo valor propio. El modelo VEC (Vector de Corrección de Error) incorpora tanto la dinámica de corto plazo como el ajuste hacia el equilibrio de largo plazo.
Modelos ARDL. El enfoque de rezagos distribuidos autorregresivos (ARDL) de Pesaran, Shin y Smith (2001) permite estimar relaciones de largo plazo incluso cuando las variables son una mezcla de \(I(0)\) e \(I(1)\). La prueba de límites (bounds test) contrasta si existe una relación de nivel entre las variables.
Filtro de Kalman y Modelos Estado-Espacio. El Filtro de Kalman es un algoritmo recursivo de dos pasos —predicción y actualización— que estima variables de estado no observables a partir de observaciones ruidosas. Sus aplicaciones incluyen la estimación de brechas del producto, tendencias estocásticas y parámetros variables en el tiempo. El Suavizador de Kalman extiende el filtro al usar toda la información muestral para refinar las estimaciones, y la descomposición del error de predicción permite estimar los parámetros del sistema por máxima verosimilitud; los modelos ARMA de los capítulos anteriores son un caso particular de esta representación.
5.9 Ejercicios
Los siguientes ejercicios combinan preguntas teóricas con aplicaciones en R. Los marcados como (Teórico) se resuelven analíticamente; los marcados como (Computacional) utilizan la base BD/Datos_Ad.RData y las funciones de los paquetes vars, urca y ARDL empleadas en este capítulo.
(Teórico) Considere el VAR(1) bivariado \(\mathbf{X}_t = \boldsymbol{\delta} + A_1 \mathbf{X}_{t-1} + \mathbf{U}_t\), con \(A_1 = \begin{bmatrix} 0.5 & 0.3 \\ 0.2 & 0.4 \end{bmatrix}\) y \(\boldsymbol{\delta} = \begin{bmatrix} 1 \\ 1 \end{bmatrix}\).
- Calcule los valores propios de \(A_1\) y determine si el sistema es estable (estacionario).
- Calcule la media incondicional del proceso, \(\boldsymbol{\mu} = (I - A_1)^{-1} \boldsymbol{\delta}\).
- Muestre cómo un VAR(2) bivariado puede reescribirse como un VAR(1) de dimensión 4 (forma companion) y enuncie la condición de estabilidad en términos de la matriz acompañante.
(Teórico) Considere un VAR(2) bivariado \(\mathbf{X}_t = A_1 \mathbf{X}_{t-1} + A_2 \mathbf{X}_{t-2} + \mathbf{U}_t\) cuyas variables son \(I(1)\).
- Demuestre que el sistema puede reescribirse en la forma de corrección de error \(\Delta \mathbf{X}_t = \Pi \mathbf{X}_{t-1} + \Gamma_1 \Delta \mathbf{X}_{t-1} + \mathbf{U}_t\), con \(\Pi = A_1 + A_2 - I\) y \(\Gamma_1 = -A_2\).
- Explique qué implica cada uno de los casos \(rango(\Pi) = 0\), \(rango(\Pi) = 1\) y \(rango(\Pi) = 2\) sobre la relación entre las dos variables.
- Suponga que \(rango(\Pi) = 1\), de modo que \(\Pi = \boldsymbol{\alpha} \boldsymbol{\beta}'\) con \(\boldsymbol{\beta} = (1, -1)'\). Escriba las dos ecuaciones del VEC resultante e interprete el papel de los coeficientes de ajuste \(\boldsymbol{\alpha}\). ¿Qué signo esperaría para el coeficiente de ajuste de la primera ecuación y por qué?
(Teórico) Siguiendo la representación estado-espacio del AR(2) dada en la ecuación (5.36), construya la representación estado-espacio de un proceso ARMA(1,1): \(Y_t = a_1 Y_{t-1} + U_t - b_1 U_{t-1}\).
- Proponga un vector de estado apropiado y escriba las matrices de la ecuación de transición y de la ecuación de observación.
- Explique la función que cumplen las etapas de predicción y de actualización del Filtro de Kalman al recorrer la muestra con esta representación.
- Explique cómo se utiliza la descomposición del error de predicción de la ecuación (5.37) para estimar \(a_1\) y \(b_1\) por máxima verosimilitud, y por qué esto implica que la función
arima()de R es, en el fondo, una aplicación del Filtro de Kalman.
(Computacional) Utilice las series en diferencias logarítmicas del INPC, del tipo de cambio y de los CETES a 28 días (\(DLINPC_t\), \(DLTC_t\) y \(DLCETE28_t\)) construidas en este capítulo a partir de
BD/Datos_Ad.RData, y estime un VAR trivariado:- Seleccione el número de rezagos con
VARselect()(criterios AIC, HQ, SC y FPE) y justifique su elección. - Estime el VAR con
VAR()y verifique la estabilidad del sistema conroots(). - Aplique pruebas de causalidad de Granger con
causality()para cada variable y resuma qué relaciones de precedencia temporal encuentra. - Calcule y grafique las funciones de impulso-respuesta con
irf()a 12 meses, con bandas de confianza bootstrap, para el efecto de un choque del tipo de cambio sobre la inflación. Interprete el resultado en términos del traspaso (pass-through) cambiario. - Compare sus resultados con los del VAR de cinco variables estimado en este capítulo: ¿son robustas las conclusiones a la reducción del sistema?
- Seleccione el número de rezagos con
(Computacional) Con las mismas tres series pero en niveles logarítmicos (\(LINPC_t\), \(LTC_t\) y \(LCETE28_t\)):
- Verifique que las tres series son \(I(1)\) (puede retomar la metodología del Capítulo 4).
- Aplique el procedimiento de Johansen con
ca.jo(), usando la estadística de la traza y la del máximo valor propio, conecdet = "const"y el número de rezagos elegido por los criterios de información. Determine el rango de cointegración. - Si encuentra al menos un vector de cointegración, normalícelo respecto a \(LINPC_t\) e interprete los coeficientes de largo plazo.
- Estime el VEC asociado con
cajorls()e interprete los coeficientes de ajuste: ¿qué variables corrigen las desviaciones del equilibrio de largo plazo?
(Computacional) Replique la metodología ARDL de este capítulo con datos mexicanos: analice la relación entre la actividad económica y la producción industrial usando \(LIGAE_t\) y \(LIPI_t\) (logaritmos de
IGAE_AdeIPI_AddeBD/Datos_Ad.RData).- Seleccione la especificación óptima con
auto_ardl(), con \(LIGAE_t\) como variable dependiente y un máximo de 5 rezagos. - Aplique la prueba de límites (bounds test) con
bounds_f_test()y concluya si existe una relación de largo plazo entre las dos series. - En caso afirmativo, calcule los multiplicadores de largo plazo e interprete la elasticidad estimada.
- Explique qué ventaja tiene el enfoque ARDL sobre el procedimiento de Johansen cuando no se tiene certeza de que todas las series sean \(I(1)\).
- Seleccione la especificación óptima con
El archivo completo del modelo y el código de Matlab que lo resuelve forman parte del material del curso Monetary Policy Analysis and Forecasting del FMI, disponible en la plataforma edX; véase IMF (2020). Aquí se reproduce sólo el fragmento necesario para ilustrar la sintaxis de una declaración estado-espacio.↩︎