\[ L. Fernanda\ Sánchez\ Cardona\\ Simón\ Yarce\ Arango \]
El objetivo del desarrollo de este ejercicio consiste en conseguir una contrucción de un modelo ARIMA que pueda pronosticar el comportamiento del PIB Colombiano. Para este trabajo tomamos los datos del Banco de la República, esta entidad comparte los datos trimestrales del PIB ya desestacionalizados. Esta base de datos fue limpiada y extraida en otro archivo de excel, el criterio de selección fue tomar todos los trimestres desde el 2005 hasta el 2021 de modo que contamos con un total de 68 observaciones
Para el siguiente trabajo instalamos y usamos los siguientes paquetes
library(tseries)
## Registered S3 method overwritten by 'quantmod':
## method from
## as.zoo.data.frame zoo
library(forecast)
Estos paquetes contienen códigos importantes para la realización del ejercicio al ser usados para series de tiempo, pronósticos y pruebas diagnósticas para los modelos.
library(readxl)
PIB_trim_col <- read_excel("C:/Users/57310/Downloads/PIB_trim_col.xlsx")
View(PIB_trim_col)
datos=PIB_trim_col[,2]
serie= ts(datos, start = c(2005,1),frequency = 4)
plot(serie, main=" Gráfico del PIB de Colombia por trimestre", xlab="Años")
Como se mencionó anteriormente, estos datos provienen de una base de datos que ya se encuentra desestacionalizada a lo cual no es necesario desestacionalizar la serie
adf.test(serie)
## Warning in adf.test(serie): p-value smaller than printed p-value
##
## Augmented Dickey-Fuller Test
##
## data: serie
## Dickey-Fuller = -4.7325, Lag order = 4, p-value = 0.01
## alternative hypothesis: stationary
Como se puede observar a través del test de Dickey-Fuller , el p value de esta serie de tiempo es de 0.01 a lo cual observando las hipótesis correspondientes a esta prueba tenemos:
\[ \begin{align} \text{ Ho :}\ No\ hay\ estacionariedad\\ \text{Ha: }\ Hay\ estacionariedad\ \end{align} \]
A lo cual rechazamos la hipótesis nula, siendo la serie estacionaria
acf(serie,main="Función de autocorrelación")
pacf(serie, main="Función de autocorrelación parcial")
Como podemos observar en las funciones de autocorrelación, la ACF tiene 5 rezagos y la PACF tiene 1, a lo cual sugeriremos un modelo ARIMA(1,0,5)
El software R también cuenta con una función que indica el modelo mas adecuado para la serie de tiempo, este criterio tambien se tendrá en cuenta.
auto.arima(serie)
## Series: serie
## ARIMA(0,1,0) with drift
##
## Coefficients:
## drift
## 3575.575
## s.e. 1000.291
##
## sigma^2 = 68054710: log likelihood = -698.77
## AIC=1401.53 AICc=1401.72 BIC=1405.94
Para este trabajo utilizaremos la metodología Box-Jenkins, en la cual probaremos todos los posibles modelos que contengan el mismo o menor número de rezagos que los mostrados en las funciones de autocorrelación. Comparando el criterio de información de Akaike(AIC) de los modelos sabremos cual es el que mas se ajusta a la serie.
mod2=arima(serie, c(1,0,2))
mod3=arima(serie, c(1,0,3))
mod4=arima(serie, c(1,0,4))
mod5=arima(serie, c(1,0,5))
modh=arima(serie,c(0,1,0))
Probamos los modelos por la metodologia Box-Jenkins y el “modh” que es el modelo recomendado por la función auto.arima El modelo ARIMA(1,0,1) no es tomado en cuenta debido a que con este modelo la serie deja de ser estacionaria
mod2$aic
## [1] 1444.068
mod3$aic
## [1] 1445.911
mod4$aic
## [1] 1443.226
mod5$aic
## [1] 1445.214
modh$aic
## [1] 1411.225
Fácilmente se observa que el modelo con menor AIC es el modelo dado por el auto.arima
A lo cual probamos la validez de este modelo. Aplicamos las siguientes pruebas 1. La prueba Box-Ljung la cual nos habla sobre la normalidad de los residuos 2. El test Shapiro-Wilk el cual nos habla sobre la distribución de estos 3. El test Jarque-Bera para tambien conocer si se trata de una distribución normal
En la prueba Box-Ljung tenemos
\[ \begin{align} \text{ Ho :}\ No\ hay\ autocorrelación\ de\ los\ errores\\ \text{ Ha: }\ Hay\ autocorrelación\ de\ los\ errores\ \end{align} \]
residuos_modeloh = modh$residuals
Box.test(residuos_modeloh, type = "Ljung")
##
## Box-Ljung test
##
## data: residuos_modeloh
## X-squared = 0.036326, df = 1, p-value = 0.8488
Con el p-value no rechazamos la hipótesis nula dando lugar a la no correlación de los errores.
Mientras que el test Shapiro-Wilk tiene las siguientes hipótesis \[ \begin{align} \text{ Ho:}\ Los\ datos\ están\ distribuidos\ normalmente\\ \text{ Ha: }\ Los\ datos\ no\ están\ distribuida\ normalmente\\ \end{align} \]
shapiro.test(modh$residuals)
##
## Shapiro-Wilk normality test
##
## data: modh$residuals
## W = 0.53877, p-value = 2.731e-13
Con este test no rechazamos la no normalidad de los datos.
\[ \begin{align} \text{ Ho :}\ Los\ datos\ pertenecen\ a\ una\ distribución\ normal\\ \text{ Ha :}\ Los\ datos\ no\ pertenecen\ a\ una\ distribución\ normal\\ \end{align} \]
jarque.bera.test(modh$residuals)
##
## Jarque Bera Test
##
## data: modh$residuals
## X-squared = 1689.2, df = 2, p-value < 2.2e-16
Con un p-value menor a 0.05 rechazamos la hipótesis nula, teniendo que los datos no pertenecen a una distribución normal
Adicionalmente observamos las nuevas funciones de autocorrelación
acf(modh$residuals, main="ACF (ARIMA(0,1,0))")
pacf(modh$residuals,main="PACF (ARIMA(0,1,0))")
Para finalizar, pronósticamos 5 años para este modelo
predict(modh, n.ahead = 5)
## $pred
## Qtr1 Qtr2 Qtr3 Qtr4
## 2022 320776.8 320776.8 320776.8 320776.8
## 2023 320776.8
##
## $se
## Qtr1 Qtr2 Qtr3 Qtr4
## 2022 8934.406 12635.158 15474.845 17868.812
## 2023 19977.940
plot(forecast(modh),main="Pronóstico del PIB usando un ARIMA(0,1,0)",xlab="Tiempo")
Al ver los resultados optamos por mirar en replicar el ejemplo del libro para conocer su AIC y compararlos con los ya existentes.
Para finalizar analizaremos el segundo modelo con mas AIC, el ARIMA(2,0,1)
residuos_modelo4 = mod4$residuals
Box.test(residuos_modelo4)
##
## Box-Pierce test
##
## data: residuos_modelo4
## X-squared = 0.41097, df = 1, p-value = 0.5215
acf(residuos_modelo4,main="ACF de errores del ARIMA (1,0,4)")
pacf(residuos_modelo4, main="PACF de errores del ARIMA(1,0,4)")
Con este p-value no podemos rechazar la hipótesis nula que es la no autocorrelación de errores, a lo cual el modelo no cumple la prueba. De igual forma podemos observar en el acf de los errores un rezago
##ARIMA(1,0,5)
También analizaremos el modelo dado a partir de los autocorreoleogramas, el cual sería un ARIMA(1,0,5)
residuos_modelo5 = mod5$residuals
Box.test(residuos_modelo5)
##
## Box-Pierce test
##
## data: residuos_modelo5
## X-squared = 0.38355, df = 1, p-value = 0.5357
acf(residuos_modelo4,main="ACF de errores del ARIMA (1,0,5)")
pacf(residuos_modelo4, main="PACF de errores del ARIMA(1,0,5)")
De igual manera observemos como se ven los pronósticos usando este modelo
predict(mod2, n.ahead=5)
## $pred
## Qtr1 Qtr2 Qtr3 Qtr4
## 2022 322322.8 321857.4 321358.3 320861.3
## 2023 320366.3
##
## $se
## Qtr1 Qtr2 Qtr3 Qtr4
## 2022 8851.462 13379.119 16711.120 19460.096
## 2023 21847.405
plot(forecast(mod2),main="Pronóstico del PIB usando un ARIMA(2,0,1)")
acf(mod2$residuals,main="ACF de los residuales del modelo ARIMA(2,0,1)")
pacf(mod2$residuals,main="PACF de los residuales del modelo ARIMA(2,0,1)")
Tomamos el ejemplo del libro:
Ramon Rosales, Jorge Perdomo, Carlos Morales, y Jaime Urrego (2010). Fundamentos de econometría intermedia: Teoría y aplicaciones. Publicado en: Apuntes de Clase CEDE , Vol. 1, No. 2010 (January 2010): pp. 1-414. Disponible en https://mpra.ub.uni-muenchen.de/37183/
Cuyo ejemplo tambien es con una base de datos del PIB
lpib=log(PIB_trim_col[,2])
lserie= ts(lpib, start = c(2005,1),frequency = 4)
plot(lserie, main="Logaritmo del PIB")
###Función de autocorrelación
acf(lserie)
pacf(lserie)
adf.test(lserie)
##
## Augmented Dickey-Fuller Test
##
## data: lserie
## Dickey-Fuller = -2.909, Lag order = 4, p-value = 0.2059
## alternative hypothesis: stationary
No rechazamos la hipótesis nula de no estacionariedad, decidimos aplicar primeras diferencias y probar nuevamente el test.
ldserie=diff(lserie, differences=1)
plot(ldserie,main="Logaritmo del PIB diferenciado una vez")
adf.test(ldserie)
##
## Augmented Dickey-Fuller Test
##
## data: ldserie
## Dickey-Fuller = -3.8272, Lag order = 4, p-value = 0.0228
## alternative hypothesis: stationary
Con un p-value de 0.0228 rechazamos la hipótesis nula y no encontramos evidencia suficiente de que no sea estacionaria
acf(ldserie,main="ACF serie logarítmica diferenciada una vez")
pacf(ldserie,main="PACF serie logarítmica diferenciada una vez")
Nuevamente miramos la función de auto.arima
auto.arima(ldserie)
## Series: ldserie
## ARIMA(0,0,0) with non-zero mean
##
## Coefficients:
## mean
## 0.0205
## s.e. 0.0041
##
## sigma^2 = 0.001166: log likelihood = 131.69
## AIC=-259.39 AICc=-259.2 BIC=-254.98
Esta transformación no arroja un modelo confiable a pesar de que en la ACF muestra 1 rezagos, en la PACF no muestra ninguno.
Con un ARIMA(0,1,0) tendríamos el ajuste mas cercano y confiable en términos de test, ya que los demás modelos probados nisiquiera pasan el test de Durbin-Watson que nos habla sobre la estacionariedad de la serie. Adicionalmente se tienen a partir de este modelo los pronósticos del PIB para los próximos 5 años.