Para las series analizadas, estimar modelos ARIMA apropiados basados en la metodología de Box-Jenkins

El acrónimo ARIMA significa media móvil integrada autorregresiva. Los retrasos de las series estacionarias en la ecuación de pronóstico se denominan términos “autorregresivos”, los retrasos de los errores de pronóstico se denominan términos de “promedio móvil” y se dice que una serie de tiempo que debe diferenciarse para hacerse estacionaria es “integrada” versión de una serie estacionaria. Los modelos de caminata aleatoria y de tendencia aleatoria, los modelos autorregresivos y los modelos de suavizado exponencial son casos especiales de los modelos ARIMA.

Un modelo ARIMA no estacional se clasifica como un modelo “ARIMA (p, d, q)”, donde: p es el número de términos autorregresivos, d es el número de diferencias no estacionales necesarias para la estacionariedad, y q es el número de errores de pronóstico rezagados en la ecuación de predicción.

Box y Jenkins popularizaron un enfoque que combina el promedio móvil y los enfoques autorregresivos en el libro “Análisis de series temporales: pronóstico y control” (Box, Jenkins y Reinsel, 1994).

Aunque ya se conocían enfoques autorregresivos y de promedio móvil (y fueron investigados originalmente por Yule), la contribución de Box y Jenkins fue desarrollar una metodología sistemática para identificar y estimar modelos que pudieran incorporar ambos enfoques. Esto hace que los modelos Box-Jenkins sean una clase poderosa de modelos.

Un modelo autorregresivo es simplemente una regresión lineal del valor actual de la serie contra uno o más valores anteriores de la serie. El valor de p se llama el orden del modelo AR.

Los modelos AR se pueden analizar con uno de varios métodos, incluidas las técnicas estándar de mínimos cuadrados lineales. También tienen una interpretación directa.

Un modelo de promedio móvil (MA) es conceptualmente una regresión lineal del valor actual de la serie contra el ruido blanco o choques aleatorios de uno o más valores anteriores de la serie. Se supone que los choques aleatorios en cada punto provienen de la misma distribución, típicamente una distribución normal, con ubicación en cero y escala constante. La distinción en este modelo es que estos choques aleatorios se propagan a valores futuros de las series de tiempo. Ajustar las estimaciones de MA es más complicado que con los modelos AR porque los términos de error no son observables.

El modelo ARMA de Box-Jenkins es una combinación de los modelos AR y MA

Por lo general, el ajuste efectivo de los modelos Box-Jenkins requiere al menos una serie moderadamente larga. Chatfield (1996) recomienda al menos 50 observaciones. Muchos otros recomendarían al menos 100 observaciones.

Supuestos del modelo ARIMA

Los datos deben ser estacionarios; por estacionario, significa que las propiedades de la serie no dependen del momento en que se capturan. Una serie de ruido blanco y series con comportamiento cíclico también pueden considerarse como series estacionarias.
Los datos deben ser univariados: ARIMA funciona en una sola variable. La regresión automática tiene que ver con la regresión con los valores pasados.

Pasos a seguir para el modelado ARIMA:

Identificación del Modelo.
Ajustar el modelo.
Diagnostico del modelo.
Pronósticos.

Previamente analizamos las series IBM y Y1. Ahora estimaremos modelos ARIMA apropiados para las series IBM y Y1 basados en la metodología de Box-Jenkins.

IBM

1. Identificación del Modelo

1.1. Gráfico de la Serie

El primer paso en la identificación del modelo es determinar si la serie es estacionaria, es decir, si la serie de tiempo parece variar alrededor de un nivel fijo. Es útil observar una gráfica de la serie junto con la función de autocorrelación de la muestra. Examinemos el comportamiento de la serie a través del tiempo mediante un gráfico.

Y=read.csv("TsData.csv")
ibm<-ts(Y$ibm,start=1990,frequency=12, end=c(1999,12))
plot.ts(ibm)

En la figura aparece el gráfico de la serie original, indica cambios en el nivel de la serie asi como tendencia, por tanto, es no estacionaria. Esto es corroborado por el correlograma muestral asi como las pruebas de Dickey-Fuller. El siguiente paso en la identificación de un modelo tentativo es examinar las autocorrelaciones de muestra de datos.

1.2. Análisis de la Autocorrelación

La función de autocorrelación (acf ()) proporciona la autocorrelación en todos los retrasos posibles. La autocorrelación en el retraso 0 se incluye por defecto, que siempre toma el valor 1, ya que representa la correlación entre los datos y ellos mismos.

acf(ibm,main="",sub="Figura 1: Función de Autocorrelación Simple")

Como podemos deducir del gráfico anterior, la autocorrelación continúa disminuyendo a medida que aumenta el retraso, lo que confirma que no existe una asociación lineal entre las observaciones separadas por retrasos más grandes.

pacf () en el rezago \(k\) es una función de autocorrelación que describe la correlación entre todos los puntos de datos que son exactamente \(k\) pasos separados, después de tener en cuenta su correlación con los datos entre esos \(k\) pasos. Ayuda a identificar el número de coeficientes de autorregresión (AR) (valor p) en un modelo ARIMA.

pacf(ibm, main="", sub="Figura 2: Función de Autocorrelación Parcial")

Confirmaremos lo dicho anteriormente con la prueba de Dickey - Fuller.

1.3 Prueba de Dickey - Fuller

Hipótesis:

\(H_0:\) La serie es no estacionaria: tiene raíz unitaria

\(H_1:\) La serie es estacionaria: no tiene raíz unitaria

Estadístico de prueba:

adfTest(ibm)

## 
## Title:
##  Augmented Dickey-Fuller Test
## 
## Test Results:
##   PARAMETER:
##     Lag Order: 1
##   STATISTIC:
##     Dickey-Fuller: 1.8653
##   P VALUE:
##     0.9837 
## 
## Description:
##  Sat Feb 15 17:55:31 2020 by user: luis

El valor del estadístico es 0,987, con un \(p\)-valor mayor a 0,05 por lo que la hipótesis nula no se rechaza, la serie original no es estacionaria.

Tomando la primera diferencia

El siguiente paso es tomar la primera diferencia a la serie.

plot.ts(diff(ibm), main="", sub="Figura 3: Serie en primera diferencia", xlab="Tiempo",ylab="Primera diferencia")

Al tomar las primeras diferencias observamos que la serie se estabiliza en torno a un valor medio (Figura 3).

acf(diff(ibm),main="",sub="Figura 4: AFC en Primera Diferencia")

pacf(diff(ibm),main="",sub="Figura 5: PAFC en Primera Diferencia")

Así mismo el correlograma confirma un declinamiento rápido en el tiempo de la AFC y PAFC. (Figuras 4 y 5).

Prueba de Dickey-Fuller a la serie en primera diferencia

Hipótesis:

\(H_0:\) La serie es no estacionaria: tiene raíz unitaria

\(H_1:\) La serie es estacionaria: no Tiene raíz unitaria

Estadístico de prueba:

adfTest(diff(ibm))

## Warning in adfTest(diff(ibm)): p-value smaller than printed p-value

## 
## Title:
##  Augmented Dickey-Fuller Test
## 
## Test Results:
##   PARAMETER:
##     Lag Order: 1
##   STATISTIC:
##     Dickey-Fuller: -7.6162
##   P VALUE:
##     0.01 
## 
## Description:
##  Sat Feb 15 17:55:33 2020 by user: luis

El valor del estadístico es -7,6162, con un \(p\) valor menor a 0,05 por lo que la hipótesis nula se rechaza. La serie en primera diferencia es estacionaria.

Una vez que se ha obtenido una serie estacionaria, debemos identificar la forma del modelo que se utilizará. La identificación de la forma del modelo se lleva a cabo comparando las autocorrelaciones y las autocorrelaciones parciales calculadas con los datos de la autocorrelaciones y las autocorrelaciones parciales teóricas de los diferentes modelos ARIMA.

En los correlogramas podemos observar lo siguiente:

En la ACF, la primera autocorrelación es significativa y además parece ser que los valores de la PACF disminuyen en forma oscilatoria. El proceso de identificacion a partir de los correlogramas muestrales sugieren un ARIMA(0,1,1) o ARIMA(0,1,0).

2. Ajuste del modelo

Una vez que se han seleccionado los modelos tentativos, se deben estimar los parámetros para estos modelos.

Modelo 1:

La salida de arima() incluye los coeficientes ajustados y el error estándar para cada coeficiente. Observando los coeficientes podemos excluir los insignificantes. Podemos usar una función confint() para este propósito.

arima1<-arima(ibm,c(0,1,1), method="ML")
arima1

## 
## Call:
## arima(x = ibm, order = c(0, 1, 1), method = "ML")
## 
## Coefficients:
##          ma1
##       0.2925
## s.e.  0.0857
## 
## sigma^2 estimated as 26.23:  log likelihood = -363.28,  aic = 728.56

coeftest(arima1)

## 
## z test of coefficients:
## 
##     Estimate Std. Error z value  Pr(>|z|)    
## ma1 0.292542   0.085711  3.4131 0.0006422 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

confint(arima1)

##         2.5 %    97.5 %
## ma1 0.1245516 0.4605318

Modelo 2:

arima2<-arima(ibm,c(0,1,0), method="ML")
arima2

## 
## Call:
## arima(x = ibm, order = c(0, 1, 0), method = "ML")
## 
## 
## sigma^2 estimated as 28.39:  log likelihood = -367.94,  aic = 735.87

confint(arima2)

##      2.5 % 97.5 %

Para escoger el mejor modelo se debe encontrar el Criterio de información (AIC) de Akaike para un conjunto de modelos y comparar los modelos con los valores de AIC más bajos. Los resultados son similares en cuanto al valor de AIC y \(\sigma^{2}\) por lo que decidimos por el modelo ARIMA(0,1,1).

3. Validación del Modelo

3.1. Estacionariedad de los residuos

La siguiente figura muestra una posible estacionariedad entre los residuales del modelo. Esto lo confirmaremos con un correlograma.

plot.ts(arima1$residuals, main="", sub="Figura 6: Residuos del modelo ARIMA(0,1,1)", xlab="Tiempo",ylab="Residuos")

3.2. Autocorrelación de los residuos:

acf(arima1$residuals,main="", sub="Figura 7: Autocorrelaciones de los residuos del modelo ARIMA(0,1,1)", xlab="Tiempo",ylab="Autocorrelación")

pacf(arima1$residuals,main="", sub="Figura 8: Autocorrelaciones Parciales de los residuos del modelo ARIMA(0,1,1)", xlab="Tiempo",ylab="Autocorrelación")

Los correlogramas muestran que no existe autocorrelación significativa en los residuos.

3.3. Normalidad de los residuos

Los gráficos Q-Q son una herramienta eficaz para evaluar la normalidad. Aquí las aplicamos a los residuos.

qqnorm(arima1$residuals,main="", sub="Figura 9: Gráfico Q para evaluar normalidad"); qqline(arima1$residuals)

El diagrama Q-Q de los residuos del modelo ARIMA(2,1,1) estimado para la serie se muestra en la Figura 9. Los puntos parecen seguir la línea recta bastante de cerca. Este gráfico no nos llevará a rechazar la normalidad de los términos de error en este modelo.

Además, la prueba de normalidad de Shapiro-Wilk aplicada a los residuos produce un estadístico de prueba de W = 0,97148, lo que corresponde a un valor de p de 0,1445, y no rechazará la normalidad basado en esta prueba.

shapiro.test(arima1$residuals)

## 
##  Shapiro-Wilk normality test
## 
## data:  arima1$residuals
## W = 0.98551, p-value = 0.2283

Como todos los gráficos respaldan la suposición de que no hay un patrón en los residuos, podemos seguir adelante y calcular pronósticos.

Pronósticos

Los parámetros del modelo ARIMA se pueden usar como modelo predictivo para hacer pronósticos de valores futuros de la serie una vez que se selecciona el modelo más adecuado para los datos de la serie temporal.

El valor d afecta los intervalos de predicción; los intervalos de predicción aumentan de tamaño con valores más altos de “d”. Todos los intervalos de predicción serán esencialmente los mismos cuando d = 0 porque la desviación estándar del pronóstico a largo plazo irá a la desviación estándar de los datos históricos.

predict(arima1,n.ahead = 5)

## $pred
##           Jan      Feb      Mar      Apr      May
## 2000 594.3871 594.3871 594.3871 594.3871 594.3871
## 
## $se
##            Jan       Feb       Mar       Apr       May
## 2000  5.121537  8.369701 10.671160 12.557683 14.195674

Necesitamos asegurarnos de que los errores de pronóstico no estén correlacionados, normalmente distribuidos con media cero y varianza constante. Podemos usar la medida de diagnóstico para encontrar el modelo apropiado con los mejores valores de pronóstico posibles.

plot(forecast(auto.arima(ibm)), main="", sub = "Figura 10: Valores pronosticados")

Y1

Para estimar el modelo ARIMA apropiado para la serie Y1 aplicaremos la misma metodología la realizada en la serie IBM previamente.

1. Identificación del Modelo

1.1. Gráfico de la Serie

El primer paso en la identificación del modelo es determinar si la serie es estacionaria.

Y=read.csv("4Series.csv")
y1<-ts(Y$Y1,start=1,frequency=12, end=c(22,12))
plot.ts(y1)

En la figura aparece el gráfico de la serie original, la cual muestra un comportamiento primero creciente y luego decreciente, no podemos afirmar que la serie sea o no estacionaria a simple vista.

1.2. Análisis de la Autocorrelación

acf(y1,main="",sub="Figura 11: Función de Autocorrelación Simple")

Como podemos deducir del gráfico anterior, la autocorrelación continúa disminuyendo a medida que aumenta el retraso, lo que hace parecer que no existe una asociación lineal entre las observaciones separadas por retrasos más grandes.

pacf () Ayuda a identificar el número de coeficientes de autorregresión (AR) (valor p) en un modelo ARIMA.

pacf(y1,main="", sub="Figura 12: Función de Autocorrelación Parcial")

1.3 Prueba de Dickey - Fuller

Hipótesis:

\(H_0:\) La serie es no estacionaria: tiene raíz unitaria.

\(H_1:\) La serie es estacionaria: no tiene raíz unitaria.

Estadístico de prueba:

adfTest(y1)

## 
## Title:
##  Augmented Dickey-Fuller Test
## 
## Test Results:
##   PARAMETER:
##     Lag Order: 1
##   STATISTIC:
##     Dickey-Fuller: 0.0966
##   P VALUE:
##     0.6471 
## 
## Description:
##  Sat Feb 15 17:55:38 2020 by user: luis

El valor del estadístico es 0,0966, con un \(p\)-valor mayor a 0,05 por lo que la hipótesis nula no se rechaza, la serie original no es estacionaria.

Tomando Primera Diferencia

El siguiente paso es tomar la primera diferencia a la serie.

plot.ts(diff(y1), main="", sub="Figura 13: Serie en primera diferencia", xlab="Tiempo",ylab="Primera diferencia")

Al tomar las primeras diferencias observamos que la serie se estabiliza en torno a un valor medio. (Figura 13).

acf(diff(y1),main="",sub="Figura 14: AFC en Primera Diferencia")

pacf(diff(y1),main="",sub="Figura 15: PAFC en Primera Diferencia")

Así mismo el correlograma confirma un declinamiento rápido en el tiempo de la AFC y PAFC. (Figuras 14 y 15).

Prueba de Dickey-Fuller a la serie en Primera Diferencia.

Hipótesis:

\(H_0:\) La serie es no estacionaria: tiene raíz unitaria.

\(H_1:\) La serie es estacionaria: no tiene raíz unitaria.

Estadístico de prueba:

adfTest(diff(y1))

## Warning in adfTest(diff(y1)): p-value smaller than printed p-value

## 
## Title:
##  Augmented Dickey-Fuller Test
## 
## Test Results:
##   PARAMETER:
##     Lag Order: 1
##   STATISTIC:
##     Dickey-Fuller: -12.4194
##   P VALUE:
##     0.01 
## 
## Description:
##  Sat Feb 15 17:55:39 2020 by user: luis

El valor del estadístico es -12,4194, con un \(p\) valor menor a 0,05 por lo que la hipótesis nula se rechaza. La serie en primera diferencia es estacionaria.

Una vez que se ha obtenido una serie estacionaria, debemos identificar la forma del modelo que se utilizará. En los correlogramas podemos observar lo siguiente:

En la ACF, ninguna autocorrelación es significativa y además parece ser que los valores de la PACF disminuyen en forma oscilatoria. Esto nos lleva a pensar en los modelos ARIMA(0,1,1), ARIMA(1,1,0) y ARIMA(0,1,0).

2. Ajuste del modelo

Una vez que se han seleccionado los modelos tentativos, se deben estimar los parámetros para estos modelos.

Modelo 1:

La salida de arima() incluye los coeficientes ajustados y el error estándar para cada coeficiente.

arima3<-arima(y1,c(0,1,1), method="ML")
arima3

## 
## Call:
## arima(x = y1, order = c(0, 1, 1), method = "ML")
## 
## Coefficients:
##          ma1
##       0.0673
## s.e.  0.0697
## 
## sigma^2 estimated as 0.213:  log likelihood = -169.83,  aic = 341.66

coeftest(arima3)

## 
## z test of coefficients:
## 
##     Estimate Std. Error z value Pr(>|z|)
## ma1 0.067254   0.069736  0.9644   0.3348

confint(arima3)

##           2.5 %    97.5 %
## ma1 -0.06942595 0.2039338

Modelo 2:

arima4<-arima(y1,c(1,1,0), method="ML")
arima4

## 
## Call:
## arima(x = y1, order = c(1, 1, 0), method = "ML")
## 
## Coefficients:
##          ar1
##       0.0517
## s.e.  0.0616
## 
## sigma^2 estimated as 0.2132:  log likelihood = -169.94,  aic = 341.87

coeftest(arima4)

## 
## z test of coefficients:
## 
##     Estimate Std. Error z value Pr(>|z|)
## ar1 0.051743   0.061605  0.8399    0.401

confint(arima4)

##           2.5 %    97.5 %
## ar1 -0.06899959 0.1724857

Modelo 3:

arima5<-arima(diff(y1),c(0,1,0), method="ML")
arima5

## 
## Call:
## arima(x = diff(y1), order = c(0, 1, 0), method = "ML")
## 
## 
## sigma^2 estimated as 0.406:  log likelihood = -253.66,  aic = 507.33

Para escoger el mejor modelo se debe encontrar el Criterio de información (AIC) de Akaike para un conjunto de modelos y comparar los modelos con los valores de AIC más bajos. Observamos que el modelo 1 ARIMA(0,1,1) tiene el menor AIC con 341,66 y al mismo tiempo la menor varianza por lo tanto nos quedamos con este modelo.

Validación del Modelo

3.1. Estacionariedad de los residuos

La siguiente figura muestra una posible estacionariedad entre los residuales del modelo. Esto lo confirmaremos con un correlograma.

plot.ts(arima3$residuals, main="", sub="Figura 16: Residuos del modelo ARIMA(0,1,1)", xlab="Tiempo",ylab="Residuos")

Autocorrelación de los residuos:

acf(arima3$residuals, main="", sub="Figura 17: Autocorrelaciones de los residuos del modelo ARIMA(0,1,1)", xlab="Tiempo",ylab="Autocorrelación")

pacf(arima3$residuals, main="", sub="Figura 18: Autocorrelaciones Parciales de los residuos del modelo ARIMA(0,1,1)", xlab="Tiempo",ylab="Autocorrelación")

Los correlogramas muestran que no existe autocorrelación significativa en los residuos.

3.2. Normalidad de los residuos

Los gráficos Q-Q son una herramienta eficaz para evaluar la normalidad. Aquí las aplicamos a los residuos.

qqnorm(arima3$residuals, main="", sub="Figura 19: Gráfico Q para evaluar normalidad");qqline(arima3$residuals)

El diagrama Q-Q de los residuos del modelo ARIMA(0,1,1) estimado para la serie se muestra en la Figura 19. Los puntos parecen seguir la línea recta bastante de cerca, sin embargo, la prueba de normalidad de Shapiro-Wilk aplicada a los residuos produce un estadístico de prueba de W = 0,9768, lo que corresponde a un valor de p de 0,00026, se rechaza la normalidad de los residuos basado en esta prueba.

shapiro.test(arima3$residuals)

## 
##  Shapiro-Wilk normality test
## 
## data:  arima3$residuals
## W = 0.9768, p-value = 0.0002664

Pronósticos

predict(arima3,n.ahead = 5)

## $pred
##         Jan      Feb      Mar      Apr      May
## 23 30.15452 30.15452 30.15452 30.15452 30.15452
## 
## $se
##          Jan       Feb       Mar       Apr       May
## 23 0.4615304 0.6750079 0.8356201 0.9699929 1.0878932

plot(forecast(auto.arima(y1)), main="",  sub = "Figura 20: Valores pronosticados")

Conclusiones

Las series temporales IBM y Y1 son de tipo financiera, por eso el ajuste del modelo de Box-Jenkins no es muy bueno ya que la dependencia entre las observaciones no es lineal. Las predicciones en IBM parecen ser más plausibles indicando que la serie continuará en ascenso mientras que las predicciones en Y1 son bastantes pobres. Se pueden utilizar modelos arch-garch. Siempre esperamos que los residuos se comporten como ruido blanco gaussiano, el ajuste de estos a la distribución normal no es perfecto.

Metodología de Box-Jenkins

Brian Gómez

26 de enero de 2020

Para las series analizadas, estimar modelos ARIMA apropiados basados en la metodología de Box-Jenkins

IBM

1. Identificación del Modelo

1.1. Gráfico de la Serie

1.2. Análisis de la Autocorrelación

1.3 Prueba de Dickey - Fuller

2. Ajuste del modelo

Modelo 1:

Modelo 2:

3. Validación del Modelo

3.1. Estacionariedad de los residuos

3.2. Autocorrelación de los residuos:

3.3. Normalidad de los residuos

Pronósticos

Y1

1. Identificación del Modelo

1.1. Gráfico de la Serie

1.2. Análisis de la Autocorrelación

1.3 Prueba de Dickey - Fuller

2. Ajuste del modelo

Modelo 1:

Modelo 2:

Modelo 3:

Validación del Modelo

3.1. Estacionariedad de los residuos

Autocorrelación de los residuos:

3.2. Normalidad de los residuos

Pronósticos

Conclusiones