Modelos ARIMA - Ejemplo: Tipo de cambio
Saul Quispe
2023-07-29
Introducción
El motivo de la introducción de los modelos ARIMA nace del hecho de que no se puede trabajar con una serie temporal no estacionaria. Se dice que una serie es estacionaria cuando su media, varianza y autocovarianza son invariantes en el tiempo. La mayoría de series temporales económicas no son estacionarias pero diferenciándolas un número determinado de veces la serie original se transforma en estacionaria, con lo cual ya se podría aplicar la metodología de los modelos ARIMA.
Concepto ARIMA
En general, se dice que una serie temporal Yt admite una representación autoregresiva integrada y de medias móviles de órdenes p , d y q respectivamente, y se denota por ARIMA( p , d , q ).
\(Y_{t} = c + \phi_{1} y_{dt-1}+ \phi_{2} y_{dt-p} +...+ \phi_{p} y_{dt-p} +\theta_{1} e_{t-1}+\theta_{2} e_{t-2}+...+\theta_{q} e_{t-q}\)
Donde p denota el número de términos autoregresivos, d el número de veces que la serie debe ser diferenciada para hacerla estacionaria y q el número de términos de la media móvil invertible.
La construcción de los modelos se lleva de manera iterativa mediante un proceso en el que se puede distinguir cuatro etapas:
Identificación: Utilizando los datos ordenados cronológicamente se intentara sugerir un modelo que merezca la pena ser investigada. El objetivo es determinar los valores que sean apropiados para reproducir la serie de tiempo. Comprende los pasos 1, 2 y 3 del ejemplo próximo.
Análisis y diferenciación de la serie temporal. Consiste en examinar la estacionariedad, los diagramas de autocorrelación, también conocidos como ACF y PACF, y elección del orden del modelo. Comprende los pasos 4 y 5.
Ajuste de un modelo ARIMA. Obtención de los coeficientes de determinación. Paso 6.
Predicción. Una vez seleccionado el mejor modelo candidato se pueden hacer pronósticos en términos probabilísticos de los valores futuros.
El objetivo del siguiente análisis es: Trazar, examinar y preparar series para modelar. Extrae el componente de estacionalidad de la serie temporal. Probar la estacionalidad y aplicar las transformaciones apropiadas. Elija el orden de un modelo ARIMA. Pronostica la serie.
La base de datos con la que hemos realizado el análisis está extraída de Yahoo Finance (a través de la técnica scraping) y trata sobre el tipo de cambio del peso mexicano con el dolar americano.
Para realizar el ejemplo, comenzamos instalando y ejecutando las librerías necesarias:
Instalación de librerias
install.packages('quantmod') # librería para hacer scraping## Installing package into '/cloud/lib/x86_64-pc-linux-gnu-library/4.3'
## (as 'lib' is unspecified)install.packages('tseries')## Installing package into '/cloud/lib/x86_64-pc-linux-gnu-library/4.3'
## (as 'lib' is unspecified)install.packages('timeSeries')## Installing package into '/cloud/lib/x86_64-pc-linux-gnu-library/4.3'
## (as 'lib' is unspecified)install.packages('forecast')## Installing package into '/cloud/lib/x86_64-pc-linux-gnu-library/4.3'
## (as 'lib' is unspecified)install.packages('xts')## Installing package into '/cloud/lib/x86_64-pc-linux-gnu-library/4.3'
## (as 'lib' is unspecified)install.packages('ggplot2')## Installing package into '/cloud/lib/x86_64-pc-linux-gnu-library/4.3'
## (as 'lib' is unspecified)Ejecutando las librerías
library(quantmod)## Loading required package: xts## Loading required package: zoo##
## Attaching package: 'zoo'## The following objects are masked from 'package:base':
##
## as.Date, as.Date.numeric## Loading required package: TTR## Registered S3 method overwritten by 'quantmod':
## method from
## as.zoo.data.frame zoolibrary(tseries)
library(timeSeries)## Loading required package: timeDate##
## Attaching package: 'timeSeries'## The following object is masked from 'package:zoo':
##
## time<-library(forecast)
library(xts)
library(ggplot2)PRIMERA ETAPA: Identificación
Lectura y preparación de los datos.
Quantmod es un paquete de R que proporciona un conjunto de herramientas para el modelado y análisis financiero cuantitativo. Permite a los usuarios acceder y manipular datos financieros de varias fuentes, incluido Yahoo Finance.
# Este comando importa todos los datos disponibles en yahoo finance
#MXN <- getSymbols('MXN=X', src='yahoo', auto.assign=FALSE)
# Este comando importa los datos de acuerdo con un periodo de tiempo determinado
MXN <- getSymbols('MXN=X', src='yahoo', from = as.Date("2021-11-01"),to=as.Date("2023-07-27"), auto.assign = FALSE)Graficando la serie
chartSeries(MXN, name="MXN=X", subset="last 6 months", theme=chartTheme("white"))
Datos de yahoo finance:
data_MXN <- data.frame(MXN, tiempo = as.Date(rownames(data.frame(MXN))))
head(data_MXN)## MXN.X.Open MXN.X.High MXN.X.Low MXN.X.Close MXN.X.Volume
## 2021-11-01 20.54890 20.86555 20.53950 20.5527 0
## 2021-11-02 20.83550 20.90876 20.70775 20.8396 0
## 2021-11-03 20.78020 20.97647 20.71810 20.7803 0
## 2021-11-04 20.54336 20.64672 20.49945 20.5317 0
## 2021-11-05 20.54557 20.63305 20.40580 20.5459 0
## 2021-11-08 20.37330 20.41883 20.31160 20.3728 0
## MXN.X.Adjusted tiempo
## 2021-11-01 20.5527 2021-11-01
## 2021-11-02 20.8396 2021-11-02
## 2021-11-03 20.7803 2021-11-03
## 2021-11-04 20.5317 2021-11-04
## 2021-11-05 20.5459 2021-11-05
## 2021-11-08 20.3728 2021-11-08attach(data_MXN)Separando la serie close:
base1 = data.frame(tiempo, MXN.X.Close)
names (base1) = c("tiempo","MXN")
base1 <- na.omit(base1) #eliminando datos ominitidos "NA"
#base1
head(base1, n = 10)## tiempo MXN
## 1 2021-11-01 20.55270
## 2 2021-11-02 20.83960
## 3 2021-11-03 20.78030
## 4 2021-11-04 20.53170
## 5 2021-11-05 20.54590
## 6 2021-11-08 20.37280
## 7 2021-11-09 20.33295
## 8 2021-11-10 20.31915
## 9 2021-11-11 20.60640
## 10 2021-11-12 20.63643Graficando la serie
ggplot(base1, aes(x = tiempo, y = MXN)) + geom_line() + labs(title = "Tipo de cambio", x = "Fecha", y = "Pesos / dólar americano")
Con la función geom_smooth en ggplot2, puedes suavizar una serie temporal para analizar tendencias o patrones. Por defecto, geom_smooth también mostrará un intervalo de confianza alrededor de la suavización. Si deseas eliminar el intervalo de confianza, puedes utilizar el argumento *se = FALSE. Aquí tienes un ejemplo de cómo usarlo:
ggplot(base1, aes(x = tiempo, y = MXN)) + geom_line() + geom_smooth(se = FALSE)+ labs(title = "Tipo de cambio", x = "Fecha", y = "Pesos / dólar americano")## `geom_smooth()` using method = 'loess' and formula = 'y ~ x'
Primero, calculamos el componente estacional del uso de datos stl(). A continuación,hallamos el componente estacional de la serie mediante suavizado y ajusta la serie original restando la estacionalidad.
MXN_ma = ts(na.omit(base1$MXN), frequency=30)
decomp = stl(MXN_ma, s.window="periodic")
deseasonal_base1 <- seasadj(decomp)
plot(decomp)
SEGUNDA ETAPA: Análisis y diferenciación de la serie temporal
Análisis de estacionariedad
La instalación de un modelo ARIMA requiere que la serie sea estacionaria . Se dice que una serie es estacionaria cuando su media, varianza y autocovarianza son invariantes en el tiempo.
Esta suposición tiene un sentido intuitivo: dado que ARIMA usa retardos previos de series para modelar su comportamiento
adf.test(MXN_ma, alternative = "stationary")##
## Augmented Dickey-Fuller Test
##
## data: MXN_ma
## Dickey-Fuller = -2.954, Lag order = 7, p-value = 0.1744
## alternative hypothesis: stationaryContraste de hipótesis:
- H0: No estacionaria (tiene raíz unitaria)
- H1: Estacionaria (No tiene raíz unitaria)
Tras realizar la prueba aumentada de Dickey-Fuller (ADF), obtenemos un p-valor = 0.1744. Como el p-valor > 0.05, no rechazamos H0. Podemos concluir que nuestra serie temporal no es estacionaria.
pp.test(MXN_ma, alternative = "stationary")##
## Phillips-Perron Unit Root Test
##
## data: MXN_ma
## Dickey-Fuller Z(alpha) = -16.627, Truncation lag parameter = 5, p-value
## = 0.1812
## alternative hypothesis: stationaryContraste de hipótesis:
- H0: No estacionaria (tiene raíz unitaria)
- H1: Estacionaria (No tiene raíz unitaria)
Tras realizar la prueba aumentada de Dickey-Fuller (ADF), obtenemos un p-valor = 0.1812. Como el p-valor > 0.05, no rechazamos H0. Podemos concluir que nuestra serie temporal no es estacionaria.
Continuamos con las autocorrelaciones.
Las ACF proporcionan información sobre cómo una observación influye en las siguientes.
Acf(MXN_ma, main='')
Las PACF proporcionan la relación directa existente entre observaciones separadas por k retardos.
Pacf(MXN_ma, main='')
Para realizar un modelo ARIMA, la serie temporal debe ser estacionaria. Para conseguir esta estacionariedad, la diferenciaremos.
MXN_d1 = diff(deseasonal_base1, differences = 1)
plot(MXN_d1)
Para comprobar que la serie es, efectivamente, estacionaria, hacemos de nuevo el test aumentado de Dickey-Fuller.
adf.test(MXN_d1, alternative = "stationary")## Warning in adf.test(MXN_d1, alternative = "stationary"): p-value smaller than
## printed p-value##
## Augmented Dickey-Fuller Test
##
## data: MXN_d1
## Dickey-Fuller = -8.2236, Lag order = 7, p-value = 0.01
## alternative hypothesis: stationaryEl p-valor obtenido es de 0.01 < 0.05. Rechazamos H0, la serie temporal es estacionaria.
Al trazar la serie diferenciada, vemos un patron oscilante al rededor de 0, sin una tendencia fuerte visible. Esto sugiere que la diferenciacion de los terminos del orden 1 es suficiente y debe incluirse en el modelo.
A continuación, los picos en rezagos particulares de la serie diferenciada pueden ayudar a informar la elección de p o q para nuestro modelo:
Acf(MXN_d1, main='ACF para la serie diferenciada 1 vez')
Pacf(MXN_d1, main='ACF para la serie diferenciada 1 vez')
TERCERA ETAPA: Ajuste de un modelo ARIMA.
Ajustamos el modelo.
auto.arima(deseasonal_base1, seasonal=FALSE)## Series: deseasonal_base1
## ARIMA(1,1,0) with drift
##
## Coefficients:
## ar1 drift
## -0.0140 -0.0082
## s.e. 0.0471 0.0057
##
## sigma^2 = 0.01528: log likelihood = 305.28
## AIC=-604.56 AICc=-604.51 BIC=-592.21Usando la notación ARIMA presentada anteriormente, el modelo ajustado se puede escribir como:
\(Y_{d_{t}} = -0.0082 -0.014 Y_{t-1}+E\)
donde E es un error y la serie original se diferencia con la orden 1.
CUARTA ETAPA: Predicción
Si los parámetros de orden del modelo y la estructura se especifican correctamente, no esperaríamos que existieran autocorrelaciones significativas.
Hemos ajustado un modelo que puede producir un pronóstico, pero comenzamos examinando los gráficos ACF y PACF para los residuos del modelo. Si los parámetros de orden del modelo y la estructura se especifican correctamente, no esperaríamos que existieran autocorrelaciones significativas.
modeloarima<-auto.arima(deseasonal_base1, seasonal=FALSE)
modeloarima## Series: deseasonal_base1
## ARIMA(1,1,0) with drift
##
## Coefficients:
## ar1 drift
## -0.0140 -0.0082
## s.e. 0.0471 0.0057
##
## sigma^2 = 0.01528: log likelihood = 305.28
## AIC=-604.56 AICc=-604.51 BIC=-592.21tsdisplay(residuals(modeloarima), lag.max=10, main='(1,1,0) Model Residuals')
La autocorrelación está dentro de las bandas, por lo que está “bien comportada”.
QUINTA ETAPA
Podemos especificar el horizonte de pronóstico h periodos por delante para que se realicen las predicciones, y usar el modelo ajustado para generar dichas predicciones:
prediccion <- forecast(modeloarima, h=30)
plot(prediccion)
Los datos proyectados para los siguientes 30 días son:
tail(prediccion$mean,30)## Time Series:
## Start = c(16, 5)
## End = c(17, 4)
## Frequency = 30
## [1] 16.80181 16.79363 16.78546 16.77728 16.76910 16.76093 16.75275 16.74457
## [9] 16.73640 16.72822 16.72004 16.71187 16.70369 16.69551 16.68734 16.67916
## [17] 16.67098 16.66281 16.65463 16.64645 16.63828 16.63010 16.62192 16.61375
## [25] 16.60557 16.59740 16.58922 16.58104 16.57287 16.56469