Trabajo Practico de Series de Tiempo para la materia Econometría
El propósito original de este trabajo era estudiar la cointegración entre las series del PBI real y el IPC en Argentina, para el periodo de 1989 a 2014 para estudiar una relación estable a largo plazo. Pero dicha idea daba poco espacio para el pronóstico, por lo que se abandonó la idea y la serie de IPC, dejando únicamente la serie de PBI. El objetivo de este trabajo es desarrollar un simple modelo econométrico para analizar el comportamiento del PBI real en Argentina en el periodo posterior a la financiarización (de 1976 en adelante) usando la serie histórica publicada en el repositorio online del Banco de la Reserva Federal de San Luis.
Como metodología se usa el enfoque de Box-Jenkins, que consta de cuatro etapas:
Y el trabajo está estructurado en torno a dichas etapas. Al final se incluye un anexo que complementa el modelo con un sencillo modelo en ecuaciones diferenciales mencionado por Zhang (2005) .
Para el presente trabajo se utilizó la serie de Producto Bruto Interno Real a Precios Constantes de Argentina, proveniente de la Universidad de California, David y de la Universidad de Groninga via FRED1
La serie del Producto Interno Bruto Real presenta las siguientes características:
Tiene una frecuencia anual.
No tiene datos faltantes.
Tiene datos desde 1950 a 2019 (69 observaciones). Sin embargo, vamos a limitarnos al periodo de 1976 en adelante (usando 43 observaciones).
No se encuentra ajustada estacionalmente.
Está medida en millones de dólares constantes de 2017.
Se verificó la integridad de la serie, confirmando que efectivamente no existen datos faltantes en la misma.
library(readr)
datos <- read_csv("archivos/pbi_real.csv")
head(datos)# A tibble: 6 × 2
observation_date RGDPNAARA666NRUG
<date> <dbl>
1 1950-01-01 196688.
2 1951-01-01 208305.
3 1952-01-01 195898.
4 1953-01-01 204165.
5 1954-01-01 210813.
6 1955-01-01 225537.A continuación, se cambia el nombre de las columnas a Tiempo y PBI Real, se transforma la variable Tiempo a anual, se selecciona el periodo propuesto (1976 en adelante) y se convierte el conjunto de datos a un tsibble con el nombre de serie_pbi_real.
library(fpp3)
library(dplyr)
# Renombrar columnas:
datos |> rename("Tiempo" = observation_date, "PBI_Real" = RGDPNAARA666NRUG) -> pbi_real
#Filtrar y transformar Tiempo:
pbi_real |> filter(Tiempo >= 1976) |> mutate(Tiempo = year(Tiempo)) -> pbi_real
# Tranformación del conjunto de datos a tsibble:
serie_pbi_real <- tsibble(
año = pbi_real$Tiempo,
PBI_Real = pbi_real$PBI_Real,
index = año
)serie_pbi_real |> autoplot() +
labs(
title = "Evolución del PBI Real en Argentina a Precios Constantes",
x = "Tiempo (años)",
y = "PBI Real (Millones USD)"
) + see::theme_modern()
La serie presenta un claro problema de escala, lo que se soluciona aplicando una transformación logarítmica para solucionarlo:
serie_pbi_real |> mutate(PBI_Real = log(PBI_Real)) |> autoplot() +
labs(
title = "Evolución del PBI Real en Argentina a Precios Constantes",
x = "Tiempo (años)",
y = "log(PBI Real)"
) + see::theme_modern()
El gráfico de la serie muestra una clara y persistente pendiente ascendente a lo largo del tiempo. Si bien hay fluctuaciones y crisis económicas (en particular, la hiperinflación de 1989 y la crisis del 2001), la trayectoria general desde 1976 a 2019 es de crecimiento. Esto será de gran importancia en la siguiente sección a la hora de elegir el tipo de prueba de estacionariedad.
serie_pbi_real |> PACF() |> autoplot() +
labs (
title = "Función de autocorrelación parcial del PBI",
x = "Número de rezagos",
y = "PACF"
) + see::theme_modern()
El autocorrelograma parcial muestra un efecto directo significativo del rezago 1 sobre la observación actual, lo que parece debilitar el supuesto de no autocorrelación en los errores. Esto será importante a la hora de usar una prueba Dickey Fuller Aumentada en lugar de una prueba DF común, y a la hora de elegir el rezago en la prueba DFA.
Prueba DFA de tipo 3 (\(PBI_t\) es una caminata aleatoria con deriva alrededor de una tendencia determinista):
\[ H_0: \text{La serie no estacionaria (tiene tendencia estocástica)} \]
\[ H_1: \text{La serie es estacionaria alrededor de una tendencia determinista} \]
library(urca)
library(aTSA)
attach(serie_pbi_real)
ur.df(PBI_Real, type="trend", lags = 1) |> summary()
###############################################
# Augmented Dickey-Fuller Test Unit Root Test #
###############################################
Test regression trend
Call:
lm(formula = z.diff ~ z.lag.1 + 1 + tt + z.diff.lag)
Residuals:
Min 1Q Median 3Q Max
-90205 -24685 5999 23028 87953
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 6.788e+04 2.814e+04 2.412 0.0208 *
z.lag.1 -1.952e-01 7.963e-02 -2.451 0.0190 *
tt 3.273e+03 1.325e+03 2.471 0.0181 *
z.diff.lag 2.349e-01 1.539e-01 1.527 0.1352
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Residual standard error: 35180 on 38 degrees of freedom
Multiple R-squared: 0.1639, Adjusted R-squared: 0.09789
F-statistic: 2.483 on 3 and 38 DF, p-value: 0.07551
Value of test-statistic is: -2.4511 3.0327 3.1229
Critical values for test statistics:
1pct 5pct 10pct
tau3 -4.15 -3.50 -3.18
phi2 7.02 5.13 4.31
phi3 9.31 6.73 5.61adf.test(PBI_Real)Augmented Dickey-Fuller Test
alternative: stationary
Type 1: no drift no trend
lag ADF p.value
[1,] 0 2.01 0.987
[2,] 1 1.41 0.957
[3,] 2 1.46 0.960
[4,] 3 1.18 0.933
Type 2: with drift no trend
lag ADF p.value
[1,] 0 -0.241 0.921
[2,] 1 -0.353 0.906
[3,] 2 -0.424 0.891
[4,] 3 -0.373 0.904
Type 3: with drift and trend
lag ADF p.value
[1,] 0 -1.96 0.581
[2,] 1 -2.45 0.386
[3,] 2 -2.25 0.463
[4,] 3 -2.74 0.276
----
Note: in fact, p.value = 0.01 means p.value <= 0.01 El valor \(DF\) observado es de -2,45, que se sitúa en la zona de no rechazo (presenta un p-valor asociado de 0,386, mucho mayor que el nivel de significancia \(\alpha=0,05\)). Por lo que no se rechaza la hipótesis nula de no estacionariedad.
A continuación se verifica si \(PBI_t \sim I(1)\) con una prueba DF tipo 2 (no se espera observar una tendencia, como en el tipo 3, ni una media de cero, como en el tipo 1):
ur.df(PBI_Real |> diff(), type="drift", lags = 1) |> summary()
###############################################
# Augmented Dickey-Fuller Test Unit Root Test #
###############################################
Test regression drift
Call:
lm(formula = z.diff ~ z.lag.1 + 1 + z.diff.lag)
Residuals:
Min 1Q Median 3Q Max
-76436 -32681 7770 26952 87617
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 1.105e+04 6.527e+03 1.694 0.098525 .
z.lag.1 -8.760e-01 2.131e-01 -4.111 0.000203 ***
z.diff.lag 5.655e-02 1.637e-01 0.346 0.731599
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Residual standard error: 37590 on 38 degrees of freedom
Multiple R-squared: 0.4138, Adjusted R-squared: 0.3829
F-statistic: 13.41 on 2 and 38 DF, p-value: 3.917e-05
Value of test-statistic is: -4.1109 8.456
Critical values for test statistics:
1pct 5pct 10pct
tau2 -3.58 -2.93 -2.60
phi1 7.06 4.86 3.94adf.test(PBI_Real |> diff())Augmented Dickey-Fuller Test
alternative: stationary
Type 1: no drift no trend
lag ADF p.value
[1,] 0 -4.98 0.0100
[2,] 1 -3.66 0.0100
[3,] 2 -2.79 0.0100
[4,] 3 -2.47 0.0164
Type 2: with drift no trend
lag ADF p.value
[1,] 0 -5.32 0.0100
[2,] 1 -4.11 0.0100
[3,] 2 -3.13 0.0354
[4,] 3 -2.89 0.0583
Type 3: with drift and trend
lag ADF p.value
[1,] 0 -5.28 0.0100
[2,] 1 -4.04 0.0175
[3,] 2 -3.12 0.1312
[4,] 3 -2.87 0.2269
----
Note: in fact, p.value = 0.01 means p.value <= 0.01 En este caso, para un rezago de 1, el valor \(DF\) observado fue de -4,1109, que esta vez cae en la zona de rechazo (para un \(\alpha = 0,05\), rechazo si \(DF < -2,93\)). Por lo tanto, se rechaza la hipótesis nula de no estacionaridad y se tiene suficiente evidencia estadística para afirmar que la serie diferenciada es estacionaria y \(PBI \sim I(1)\).
En función de los resultados del análisis exploratorio, la estrategia de modelado será ajustar un modelo ARIMA(p,d,q) a la serie del log(PBI_Real), donde d=1 ya que \(PBI_t \sim I(1)\).
p y qPara identificar p y q se puede usar las funciones de autocorrelación (ACF) y de autocorrelación parcial (PACF):
serie_pbi_real |> mutate(diff_PBI = difference(PBI_Real)) |> PACF(y = diff_PBI)# A tsibble: 16 x 2 [1Y]
lag pacf
<cf_lag> <dbl>
1 1Y 0.160
2 2Y -0.0559
3 3Y 0.0496
4 4Y -0.0789
5 5Y -0.0840
6 6Y -0.0899
7 7Y 0.238
8 8Y -0.400
9 9Y -0.0779
10 10Y 0.0217
11 11Y -0.145
12 12Y 0.0692
13 13Y 0.0833
14 14Y 0.00431
15 15Y -0.0512
16 16Y -0.0706 serie_pbi_real |> mutate(diff_PBI = difference(PBI_Real)) |> PACF(y = diff_PBI) |> autoplot() +
labs (
title = "Función de autocorrelación parcial del PBI diferenciado",
x = "Número de rezagos",
y = "PACF"
) + see::theme_modern()
No se observa un decrecimiento exponencial que podría sugerir un componente de media móvil. Sin embargo, se observa un pico negativo significativo en el rezago 8, lo que puede sugerir un orden del componente autorregresivo \(p=8\). En el rezago 7 se ve un pico positivo que no llega a ser estadísticamente significativo.
serie_pbi_real |> mutate(diff_PBI = difference(PBI_Real)) |> ACF(y = diff_PBI)# A tsibble: 16 x 2 [1Y]
lag acf
<cf_lag> <dbl>
1 1Y 0.160
2 2Y -0.0287
3 3Y 0.0344
4 4Y -0.0607
5 5Y -0.108
6 6Y -0.110
7 7Y 0.196
8 8Y -0.280
9 9Y -0.228
10 10Y 0.0257
11 11Y -0.126
12 12Y 0.0103
13 13Y 0.109
14 14Y 0.196
15 15Y -0.0968
16 16Y -0.0760serie_pbi_real |> mutate(diff_PBI = difference(PBI_Real)) |> ACF(y = diff_PBI) |>
autoplot() +
labs (
title = "Función de autocorrelación del PBI diferenciado",
x = "Número de rezagos",
y = "ACF"
) + see::theme_modern()
En la función de autocorrelación se ve que los picos no tienen un “corte” claro, sino que parecen decrecer o alternar de forma más gradual. Este comportamiento es consistente con el \(AR(8)\) que se ve en el gráfico de la función de autocorrelación parcial anterior. Como no hay un corte abrupto después de algún rezago q se puede concluir que no es necesario un componente de media movil (MA), por lo que \(q=0\).
En la sección 2.4 se identificó como orden de integración d=1, mientras que en la sección anterior se identificaron los ordenes autorregresivos y de media móvil p=8 y q=0. No obstante, un orden de media móvil tan elevado podría implicar un error de sobreespecificación del modelo al tener tantos parámetros. A continuación se hace una breve revisión de la literatura de modelos ARIMA aplicados al cálculo del PBI, para su posterior cálculo y comparación:
Sultan (2023) propuso que un modelo ARIMA (1,1,1) resulta efectivo para el caso de Estados Unidos, mostrando un error de pronóstico relativamente bajo.
Saravanamutthu Jeyarajah (2022) también llegó a un ARIMA(1,1,1) para predecir el PBI en Sri Lanka.
Ghazo (2021) usó datos de 1976 a 2019 y modelos ARIMA para predecir el PBI y el IPC (Índice de Precios al Consumidor) en Jordania, y encontró que el mejor modelo para el PBI era un ARIMA(3,1,1).
Finalmente, Merlain (2023) identificó un ARIMA (2,1,2) como modelo óptimo para la estimación del PBI anual en Camerún.
No obstante, Zhang (2005) muestra un simple modelo de PBI en tiempo continuo donde la tasa de cambio del PBI es proporcional al PBI actual. Esto, junto al orden de integración d=1 y a la presencia de deriva encontrados en la sección 2.4, se puede traducir como un ARIMA (0,1,0), sin ordenes autorregresivos ni de media móvil nulos pero con deriva.
A continuación, se calculan y comparan los cinco modelos propuestos usando la librería forecast ya que permite incluir deriva:
ARIMA (1,1,1)
ARIMA (3,1,1)
ARIMA (2,1,2)
ARIMA (0,1,0)
y ARIMA (8,1,0)
library(forecast)
attach(serie_pbi_real)
ARIMA111 = Arima(log(PBI_Real), order=c(1,1,1), include.drift = TRUE)
ARIMA111Series: log(PBI_Real)
ARIMA(1,1,1) with drift
Coefficients:
ar1 ma1 drift
-0.7327 1.0000 0.0184
s.e. 0.1122 0.0862 0.0090
sigma^2 = 0.002797: log likelihood = 65.96
AIC=-123.92 AICc=-122.86 BIC=-116.87ARIMA311 = Arima(log(PBI_Real), order=c(3,1,1), include.drift = TRUE)
ARIMA311Series: log(PBI_Real)
ARIMA(3,1,1) with drift
Coefficients:
ar1 ar2 ar3 ma1 drift
-0.6979 0.0226 -0.0424 1.0000 0.0185
s.e. 0.1534 0.1855 0.1562 0.0743 0.0090
sigma^2 = 0.002944: log likelihood = 66.05
AIC=-120.1 AICc=-117.76 BIC=-109.53ARIMA212 = Arima(log(PBI_Real), order=c(2,1,2), include.drift = TRUE)
ARIMA212Series: log(PBI_Real)
ARIMA(2,1,2) with drift
Coefficients:
ar1 ar2 ma1 ma2 drift
0.1227 0.5841 0.0000 -1.0000 0.0205
s.e. 0.1352 0.1381 0.1196 0.1196 0.0030
sigma^2 = 0.002674: log likelihood = 67.19
AIC=-122.38 AICc=-120.04 BIC=-111.81ARIMA010 = Arima(log(PBI_Real), order=c(0,1,0), include.drift = TRUE)
ARIMA010Series: log(PBI_Real)
ARIMA(0,1,0) with drift
Coefficients:
drift
0.0181
s.e. 0.0084
sigma^2 = 0.0031: log likelihood = 63.71
AIC=-123.42 AICc=-123.12 BIC=-119.9ARIMA810 = Arima(log(PBI_Real), order=c(8,1,0), include.drift = TRUE)
ARIMA810Series: log(PBI_Real)
ARIMA(8,1,0) with drift
Coefficients:
ar1 ar2 ar3 ar4 ar5 ar6 ar7 ar8
0.3235 -0.2370 0.0425 -0.1656 -0.0295 -0.2407 0.3403 -0.4595
s.e. 0.1321 0.1343 0.1358 0.1432 0.1348 0.1362 0.1323 0.1316
drift
0.0205
s.e. 0.0051
sigma^2 = 0.002528: log likelihood = 71.34
AIC=-122.68 AICc=-115.8 BIC=-105.07Los modelos ARIMA (1,1,1), (2,1,2) y (3,1,1) presentan un problema en la parte de media móvil: presentan coeficientes unitarios en valor absoluto, lo que se denomina no invertibilidad. Este problema afecta al pronóstico ya que no existe una forma única y estable de representar los errores pasados en función de las observaciones pasadas y es un indicador de mala especificación del modelo.
Por su parte, el ARIMA (8,1,0) tiene un AIC más bajo, pero un AICc y un BIC más altos que el ARIMA (0,1,0), lo que sugiere un problema de sobreajuste:
library(lmtest)
coeftest(ARIMA810)
z test of coefficients:
Estimate Std. Error z value Pr(>|z|)
ar1 0.3234573 0.1321192 2.4482 0.0143562 *
ar2 -0.2369755 0.1343112 -1.7644 0.0776688 .
ar3 0.0424793 0.1357887 0.3128 0.7544072
ar4 -0.1655742 0.1431775 -1.1564 0.2475070
ar5 -0.0294690 0.1348443 -0.2185 0.8270077
ar6 -0.2407133 0.1362176 -1.7671 0.0772077 .
ar7 0.3402780 0.1323455 2.5711 0.0101366 *
ar8 -0.4595487 0.1316158 -3.4916 0.0004802 ***
drift 0.0205424 0.0051279 4.0060 6.175e-05 ***
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1Como era de esperar, los coeficientes ar2, ar3, ar4 y ar5 no son estadísticamente significativos para un \(\alpha=0,05\), señal de sobreajuste.
| Modelo | # Parámetros | AICc2 | BIC | Validez |
|---|---|---|---|---|
ARIMA(0,1,0) |
1 | -123,12 | -119,9 | Válido |
ARIMA(1,1,1) |
3 | -122,86 | -116,87 | Inválido (no invertible) |
ARIMA(2,1,2) |
5 | -120,04 | -111,81 | Inválido (no invertible) |
ARIMA(3,1,1) |
5 | -117,76 | -109,53 | Inválido (no invertible) |
ARIMA(8,1,0) |
9 | -115,8 | -105,07 | Válido pero sobreajustado. |
Antes de proceder con la etapa de pronóstico, es necesario evaluar si los residuos estimados en el modelo ARIMA(0,1,0) son de ruido blanco:
checkresiduals(ARIMA010) + see::theme_modern()
Ljung-Box test
data: Residuals from ARIMA(0,1,0) with drift
Q* = 12.662, df = 9, p-value = 0.1785
Model df: 0. Total lags used: 9NULLLa prueba Box-Peirce indica que:
\[ H_0: \text{todos los } \rho \text{ hasta el rezago } m \text{ son simultáneamente iguales a } 0 \]
\[ H_1: \text{Al menos un } \rho \text{ hasta el rezago } m \text{ es distinto de 0} \]
En otras palabras, \(H_0\) significa que los residuos no están autocorrelacionados. Con un nivel de significancia de \(\alpha = 0,05\) no se rechaza la hipótesis nula ya que se tiene un p-valor = 0,1785 y podemos suponer que los errores no se encuentran significativamente autocorrelacionados.
Asimismo, la media está situada en 0 (indicado en el tercer gráfico) y no se tiene autocorrelación, y la varianza de los errores se comporta de forma homocedástica. Entonces, podemos concluir que los residuos se comportan como ruido blanco.
El modelo ganador es el ARIMA(0,1,0) por varias razones: no presenta problemas de invertibilidad, por lo que es estadísticamente válido. Por otro lado, al ser un modelo parsimonioso, tampoco presenta problemas de sobreajuste (tiene valores AICc y BIC más bajos de todos los modelos). Además, es un modelo de fácil interpretación: la economía crece, en promedio, una cantidad constante cada año, sujeta a shocks aleatorios.
coeftest(ARIMA010)
z test of coefficients:
Estimate Std. Error z value Pr(>|z|)
drift 0.0181229 0.0083877 2.1607 0.03072 *
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1El componente de deriva del modelo es de aproximadamente 0,01812, es estadísticamente significativo con un nivel de confianza del 95%, y su error estándar asociado es de 0,0083877. El modelo puede expresarse de la siguiente manera:
\[ \ln PBI_{t+1} = \ln PBI_t + 0,01812 + \varepsilon_t \]
En este caso, el componente de la deriva se interpreta como la tasa de crecimiento porcentual promedio del PBI real, por lo que se espera en promedio que el PBI real crezca un 1,81% por año.
Este modelo, puede entenderse en tiempo continuo como que la tasa de cambio del PBI es proporcional al PBI actual (ver anexo). A continuación se muestra en rojo el PBI teórico en tiempo continuo (ver Anexo: Modelo en Tiempo Continuo):
Sultan (2023) usó un horizonte de pronóstico de 5 años, mientras que Saravanamutthu Jeyarajah (2022) uno de 8 años. De todas maneras, es importante reconocer que la precisión del pronóstico se va a ver inevitablemente afectada por el shock estructural exógeno que representó la pandemia de COVID-19 en 2020. Con todo, se tomará un horizonte de pronóstico de 5 años, que abarca del 2020 al 2024.
Para el pronostico se usa la función forecast de la librería del mismo nombre:
library(forecast)
PBI_Real_ts <- ts(PBI_Real, start = 1976)
arima_pbi <- Arima(log(PBI_Real_ts), order=c(0,1,0), include.drift = TRUE)
pronostico <- forecast::forecast(arima_pbi, h=5)
pronostico Point Forecast Lo 80 Hi 80 Lo 95 Hi 95
2020 13.80890 13.73754 13.88026 13.69977 13.91803
2021 13.82702 13.72611 13.92793 13.67269 13.98135
2022 13.84514 13.72155 13.96874 13.65613 14.03416
2023 13.86327 13.72056 14.00598 13.64501 14.08153
2024 13.88139 13.72183 14.04095 13.63737 14.12541Sin embargo, los resultados hay que transformarlos para interpretarlos:
pronostico_final = as.data.frame(pronostico)
exp(pronostico_final) Point Forecast Lo 80 Hi 80 Lo 95 Hi 95
2020 993410.3 924994.4 1066887 890705.6 1107958
2021 1011577.9 914478.4 1118987 866908.9 1180389
2022 1030077.7 910320.8 1165589 852668.8 1244399
2023 1048915.9 909413.8 1209817 843241.1 1304757
2024 1068098.6 910576.7 1252870 836823.9 1363291pronostico |> autoplot() +
labs(
title ="Pronóstico del log(PBI) usando un ARIMA (0,1,0)",
x = "Tiempo",
y = "PBI Real (en logaritmo)"
) + see::theme_modern()
La línea azul indica el pronóstico puntual (point forecast), es decir, el valor más probable que el modelo predice. el área azul oscura representa los intervalos de confianza para un nivel de confianza del 80%, mientras que en azul claro se representan los intervalos de confianza para un nivel de significancia del 95%.
| Año | Pronóstico puntual | Límite Inferior (95%) | Límite Superior (95%) |
|---|---|---|---|
| 2020 | 993.410,3 | 890.705,6 | 1.107.958 |
| 2021 | 1.011.577,9 | 866.908,9 | 1.180.389 |
| 2022 | 1.0300.77,7 | 852.668,8 | 1.244.399 |
| 2023 | 1.0489.15,9 | 843.241,1 | 1.304.757 |
| 2024 | 1.068.098,6 | 836.823,9 | 1.363.291 |
El modelo predice una tasa de crecimiento anual promedio del 1,81%, una extrapolación que choca claramente con la realidad observada debido a factores exógenos como la pandemia. Este problema de pronóstico ilustra las limitaciones de modelos univariados ante quiebres estructurales que no son contemplados en la muestra. Por su parte, los intervalos tienden a ensancharse a medida que se pronostica un mayor número de años hacia adelante, lo que es de esperar debido a una mayor incertidumbre asociada.
Con todo, el pronóstico es valioso no tanto como una predicción precisa, sino como una ilustración de la trayectoria que la economía hubiera seguido de no haberse producido el shock, de acuerdo a la dinámica histórica de la serie.
A lo largo de este trabajo se analizó la serie de tiempo del PBI real de Argentina para el periodo 1976-2019 con el fin de construir un modelo para su pronóstico. El análisis exploratorio y las pruebas DFA mostraron que la serie no es estacionaria, sino que es integrada de orden 1. Siguiendo la metodología Box-Jenkins, se compararon distintos modelos ARMA para la serie diferenciada, y el modelo que demostró tener un mejor ajuste según los criterios AICc y BIC fue el sencillo modelo de deriva con constante ARIMA(0,1,0), que indica que la economía crece, en promedio, una cantidad constante cada año (en torno al 1,82% anual), sujeta a shocks aleatorios. Los residuos estimados del modelo se ajustaron a un ruido blanco, validando la correcta especificación del modelo. El modelo final se usó para generar un pronóstico para el periodo 2020-2024. Sin embargo, los valores puntuales de la predicción no resultan tan importantes como su interpretación crítica frente al quiebre estructural provocado por la pandemia de COVID-19. Es importante destacar que el modelo ilustra una limitación de este tipo de modelos univariados: su incapacidad para anticipar shocks exógenos.
Una fundamentación teórica en tiempo continuo del resultado encontrado en el ARIMA en tiempo discreto es un modelo de ecuaciones diferenciales obtenido del supuesto que la tasa de cambio del PBI es proporcional al PBI actual Zhang (2005) :
\[ x'(t) = gx(t) \tag{1}\]
donde:
\(t\) denota el tiempo;
\(x(t)\) es una función representa el PBI en el momento \(t\);
\(x'(t)\) es la derivada de \(x\) respecto al tiempo;
y \(g\) es una constante que representa la tasa de crecimiento esperada.
\[ x'(t) = g x(t) \implies g = \frac{x'(t)}{x(t)} \]
Recordar que si:
\[ f(x) = \ln u \implies f'(x) = \frac{u'}{u} \]
Ya que si \(f(x) = \ln x(t)\), se tiene que:
\[ f(x) = \ln x(t) \implies \frac{d}{dt} \bigg [ \ln x(t) \bigg] = \frac{x'(t)}{x(t)} \]
Pero como \(x'(t)/x(t) = g\):
\[ \frac{d}{dt}\big [ \ln x(t) \big] = g \]
Lo cual se asemeja de forma continua a un modelo de caminata aleatoria con deriva:
\[ \Delta \ln x_t = g + u_t \]
La Ecuación 1 es una ecuación diferencial ordinaria separable, cuya solución es:
\[ x(t) = x(0) e^{gt} \tag{2}\]
Usando la serie_pbi_real se puede aproximar el \(g\) de la Ecuación 2 con una regresión lineal:
\[ \log(PBI) = \beta_0 + g t \]
log_pbi <- log(serie_pbi_real$PBI_Real)
tiempo <- 1:length(PBI_Real)
modelo_log <- lm(log_pbi ~ tiempo)
summary(modelo_log)$coefficient[2] # g[1] 0.02222En este modelo la tasa de crecimiento anual esperada es del 2,22%, resultado que contrasta con el 1,81% del modelo ARIMA que incorpora ruido e incertidumbre.
\(x(0)\) es el primer valor de la serie (es decir, PBI_Real[0] = 447527), por lo que el modelo en tiempo continuo queda de la siguiente forma:
\[ \boxed{x(t) = 447.527 e^{0,02222t}} \]
En el siguiente gráfico \(x(t)\) queda representado como la curva roja:
x_t = 447527.0*exp(1)^(0.02222*1:length(serie_pbi_real$PBI_Real))
ts(x_t |> log(), start=1976) |> plot(col="red", ylab=NULL, xlab=NULL, yaxt = "n")
par(new=TRUE)
plot(ts(serie_pbi_real$PBI_Real |> log(), start=1976), ylab="log(PBI)", xlab="Tiempo",
main="PBI teórico vs PBI real", las=1)
Federal Reserve Economic Data, una base de datos mantenida por el Banco de la Reserva Federal de San Luis (Federal Reserve Bank of St. Louis)↩︎
Es una modificación del criterio AIC corregido para muestras pequeñas, ya que el AIC común puede tender a seleccionar modelos con demasiados parámetros y por lo tanto elegir un modelo sobreajustado.↩︎