En la figura 1 se presenta el comportamiento de cierre de la acción de Amazon a partir del 01 de enero del 2016 al 06 de mayo de 2021: se observa que el año 2021 superó con creces el buen año que tuvo esta empresa en el año 2020, con un precio de cierre máximo histórico de $ 3731,41 en fecha 8 de julio del año 2021, este comportamiento se le atribuye principalmente al aumento en el precio en sus acciones, representando hoy día un lucrativo negocio. El valor bursátil de Amazon se calcula como próximo a los 460.000 millones de dólares, permitiendo con esto que se coloque como la tercera empresa más grande y exitosa del índice S&P 500 situado justo entre las empresas Apple y Microsft.
getSymbols("AMZN", src = "yahoo", from = "2010-01-01", to = "2021-10-30", periodicity = "daily")
## [1] "AMZN"
retorno2 <- AMZN$AMZN.Close
amazon.ts<- xts(retorno2,frequency = 365)
max <- max(retorno2)
max
## [1] 3731.41
#install.packages("dygraphs")
library(dygraphs)
#install.packages("quantmod")
library(quantmod)
#install.packages("xts")
library(xts)
Dygraphs1<- dygraph(retorno2, main = "AMNZ", ylab = "Precio de cierre",
xlab = "Periodo")%>%
dyRangeSelector()%>%
dyHighlight(highlightCircleSize = 5,highlightSeriesBackgroundAlpha = 0.2,
hideOnMouseOut = FALSE)
Dygraphs1
Se observa que la historia de Amazon comienza a cambiar a partir del año 2017, en ese año el precio de su acción comienza a aumentar de manera significativa, hasta el pico histórico que presentó en el año 2018, específicamente el 4 de septiembre de ese año, luego experimenta una caída, estabilizándose en el año 2019, para luego comenzar un ascenso vertiginoso a partir de marzo de 2020, hasta la fecha, manteniéndose estable en torno al promedio. Es importante señalar que la información para el presente estudio se tomó de la página Google Finance: https://www.google.com/finance/, pudiéndose extraer las siguientes variables: AMZN.Open: precio de apertura AMZN.High: precio más alto AMZN.Low: precio más bajo AMZN.Close: precio de cierre AMZN.Volume: volumen negociado AMZN.Adjusted: precio ajustado
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.
En la gráfica 1 se observó que la variable estudiada presenta una tendencia creciente, con períodos de volatilidad, picos y caídas, con máximos y mínimos históricos, en definitiva, podemos observar gráficamente que la variable no es estacionaria ni en media ni en varianza. Para corroborar esta afirmación obtenida visualmente apliquemos el test de Dickey-Fuller
adf.test(amazon.ts)
##
## Augmented Dickey-Fuller Test
##
## data: amazon.ts
## Dickey-Fuller = -1.241, Lag order = 14, p-value = 0.8997
## alternative hypothesis: stationary
El resultado de aplicar el test es: p-value= 0,8997 que es mayor al 5%, con lo cual concluimos que la variable no es estacionaria. Como la variable presenta tendencia estocástica, aplicamos primera diferencia y observar el resultado de aplicar esa diferencia, procedemos a verificar si la serie adquiere un comportamiento estacionario, veamos el resultado de esta acción: La diferencia en la serie la hacemos a través del comando diff, importante eliminar el primer registro del objeto obtenido, ya que al realizar la diferencia se “pierde” el primer valor, por lo cual debemos eliminarlo para no crear inconveniente al trabajar con la serie obtenida en la diferencia.
tcdif <- diff(amazon.ts)
tcdif <- tcdif[-1]
class(tcdif)
## [1] "xts" "zoo"
tcdif.ts <- ts(tcdif, start = c(2010,1), frequency = 12)
autoplot(tcdif.ts, main = "Precio de cierre Amazon")
Gráficamente podemos concluir que la diferencia obtenida es una serie estacionaria, sin embargo, apliquemos el test de Dickey-Fuller y corroboremos esta afirmación obtenida de la gráfica:
adf.test(tcdif.ts)
## Warning in adf.test(tcdif.ts): p-value smaller than printed p-value
##
## Augmented Dickey-Fuller Test
##
## data: tcdif.ts
## Dickey-Fuller = -15.397, Lag order = 14, p-value = 0.01
## alternative hypothesis: stationary
Obtenemos p-value = 0.01 que es menor al 5% con lo cual podemos concluir que la serie obtenida en primera diferencia es estacionaria.
El comando acf muestra el correlograma de la variable diferenciada, y en la gráfica obtenida buscamos aquellos valores que salen de las bandas de confianza, y estos valores indicarán los componentes de media móvil (MA), importante señalar que el primer rezago, que tiene valor 1, no se toma en cuenta por ser la autocorrelación consigo mismo, de allí el valor 1. El comando pacf nos da la información para el componente autorregresivo (AR)
acf(tcdif.ts)
Veamos los rezagos que salen de las bandas de confianza, luego para el componente MA= 2, 7, 8
pacf(tcdif.ts)
Para el componente AR=3, 8, 9 ### 3.2.- Construcción de los modelos: de los resultados obtenidos en los componentes AR y MA, podemos tener los siguientes modelos ARMA(3,2); ARMA(3,7); ARMA(3,8); ARMA(8,2) Estos modelos lo construimos en R con el comando arma y tomamos del summary el valor que da el AIC de cada modelo y seleccionaremos los dos modelos que tengan el menor AIC.
modelo32 <- arma(tcdif.ts, order= c(3,2))
summary(modelo32)
##
## Call:
## arma(x = tcdif.ts, order = c(3, 2))
##
## Model:
## ARMA(3,2)
##
## Residuals:
## Min 1Q Median 3Q Max
## -280.5478 -5.1501 -0.5438 5.2134 227.4966
##
## Coefficient(s):
## Estimate Std. Error t value Pr(>|t|)
## ar1 -0.60677 NA NA NA
## ar2 -0.29137 NA NA NA
## ar3 -0.07377 NA NA NA
## ma1 0.58495 NA NA NA
## ma2 0.28454 NA NA NA
## intercept 2.10391 0.11614 18.11 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Fit:
## sigma^2 estimated as 722.5, Conditional Sum-of-Squares = 2147974, AIC = 28057.08
modelo37 <- arma(tcdif.ts, order= c(3,7))
summary(modelo37)
##
## Call:
## arma(x = tcdif.ts, order = c(3, 7))
##
## Model:
## ARMA(3,7)
##
## Residuals:
## Min 1Q Median 3Q Max
## -279.5074 -5.2237 -0.6179 5.0738 232.5318
##
## Coefficient(s):
## Estimate Std. Error t value Pr(>|t|)
## ar1 -1.017745 NA NA NA
## ar2 0.215453 NA NA NA
## ar3 0.480431 NA NA NA
## ma1 1.000264 NA NA NA
## ma2 -0.234787 NA NA NA
## ma3 -0.558246 NA NA NA
## ma4 -0.048347 NA NA NA
## ma5 0.026077 NA NA NA
## ma6 -0.012683 0.022144 -0.573 0.567
## ma7 -0.008789 NA NA NA
## intercept 1.461858 NA NA NA
##
## Fit:
## sigma^2 estimated as 718, Conditional Sum-of-Squares = 2131842, AIC = 28048.65
modelo38 <- arma(tcdif.ts, order= c(3,8))
summary(modelo38)
##
## Call:
## arma(x = tcdif.ts, order = c(3, 8))
##
## Model:
## ARMA(3,8)
##
## Residuals:
## Min 1Q Median 3Q Max
## -279.5442 -5.2298 -0.6563 5.0798 231.1397
##
## Coefficient(s):
## Estimate Std. Error t value Pr(>|t|)
## ar1 -1.004200 NA NA NA
## ar2 0.194250 NA NA NA
## ar3 0.452447 NA NA NA
## ma1 0.984842 NA NA NA
## ma2 -0.215761 NA NA NA
## ma3 -0.531665 NA NA NA
## ma4 -0.045715 NA NA NA
## ma5 0.023845 NA NA NA
## ma6 -0.013898 0.023901 -0.581 0.561
## ma7 -0.001979 0.026372 -0.075 0.940
## ma8 0.004456 NA NA NA
## intercept 1.503212 NA NA NA
##
## Fit:
## sigma^2 estimated as 718.3, Conditional Sum-of-Squares = 2131887, AIC = 28051.71
modelo82 <- arma(tcdif.ts, order= c(8,2))
summary(modelo82)
##
## Call:
## arma(x = tcdif.ts, order = c(8, 2))
##
## Model:
## ARMA(8,2)
##
## Residuals:
## Min 1Q Median 3Q Max
## -277.95944 -4.63953 -0.08138 5.64117 233.56988
##
## Coefficient(s):
## Estimate Std. Error t value Pr(>|t|)
## ar1 -1.045602 0.483743 -2.161 0.03066 *
## ar2 -0.391900 0.386307 -1.014 0.31036
## ar3 -0.085918 0.027723 -3.099 0.00194 **
## ar4 -0.048495 0.043212 -1.122 0.26175
## ar5 -0.014416 0.030477 -0.473 0.63619
## ar6 -0.030630 0.028332 -1.081 0.27966
## ar7 -0.002249 0.027669 -0.081 0.93521
## ar8 -0.024045 0.031364 -0.767 0.44329
## ma1 1.036724 0.479977 2.160 0.03078 *
## ma2 0.375638 0.379792 0.989 0.32263
## intercept 1.604136 1.193169 1.344 0.17881
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Fit:
## sigma^2 estimated as 718.8, Conditional Sum-of-Squares = 2134147, AIC = 28051.69
Los modelos que se seleccionarán serán los que tengan el menor valor AIC, en este caso es el modelo (3,7) y (8,2)
# Busquemos otro modelo usando la funcion auto.arima
m1 <- auto.arima(tcdif.ts,d=NA, D=NA, max.p = 8, max.q = 7, max.P = 0, max.Q = 0,
max.d = 0, max.D = 0, stationary = TRUE, seasonal = TRUE, ic=c("bic"))
summary(m1)
## Series: tcdif.ts
## ARIMA(0,0,0) with zero mean
##
## sigma^2 estimated as 728.1: log likelihood=-14033.98
## AIC=28069.96 AICc=28069.96 BIC=28075.96
##
## Training set error measures:
## ME RMSE MAE MPE MAPE MASE ACF1
## Training set 1.08785 26.9827 13.74663 100 100 0.6920597 -0.02573378
Haciendo uso de los comando auto.arima y arima no obtuve resultados coherentes, el primer caso el modelo recomendado es un ARIMA(0,0,0) y en segundo caso arroja un error al correr el comando, por lo cual selecciono los modelos (3,7) y (8,2) que tienen el menor AIC calculado de forma manual. Procedamos ahora a obtener los valores estimado de los dos modelos seleccionados, llamémoslo pax_f1 y pax_f2.
# extraer los valores estimados
pax_f1 <- fitted(modelo37)
pax_f2 <- fitted(modelo82)
Construyamos una matriz que contenga los valores reales y los valores estimados, de la forma siguiente:
# Unamos los valores reales y en troa los valores estimados
Model.serie1 <- cbind(tcdif.ts, pax_f1)
model.serie2 <- cbind(tcdif.ts, pax_f2)
Grafiquemos ambos valores para cada uno de los modelos, en dicha grafica podemos observar que tan buena fue la estimación realizada: Para el modelo ARMA(3,7)
autoplot(Model.serie1, main = "ARMA(3,7)")
autoplot(model.serie2, main = "ARMA(8,2)")
Al observar las gráficas anteriores, podemos concluir que los modelos no son buenos en las estimaciones, sus predicciones en cuanto a la variabilidad de los valores reales y los estimados en grande, sin embargo, describe perfectamente el comportamiento de la serie a lo largo del tiempo.
Independencia: para determinar si existe autocorrelación debemos hacer la prueba de Ljung-Box, la cual arroja respuesta a las siguientes interrogantes: H0 : los datos se distribuyen de manera independiente (no hay autocorrelación) H1 : los datos no se distribuyen de manera independiente (si hay autocorrelación) Esto se realiza haciendo uso del comando checkresiduals, para el modelo ARMA(3,7) obtenemos
checkresiduals(modelo37)
## Warning in modeldf.default(object): Could not find appropriate degrees of
## freedom for this model.
La primera imagen nos muestra el comportamiento del término de error, observándose que es estacionario alrededor de una media, no pudiendo afirmar que existe varianza constante, hay unos picos, lo que muestra que hay una varianza condicional, en cuanto al correlograma son menos los valores que caen fuera de las bandas de confianza, y en cuanto a la distribución del término de error se observa que no presenta una distribución normal, existe una alta curtosis hacia el valor central, para corroborar lo dicho, se debe observar el
checkresiduals(modelo82)
## Warning in modeldf.default(object): Could not find appropriate degrees of
## freedom for this model.
Se muestra muy similar al modelo anterior, se toman las mismas conclusiones. En cuanto al p-value el software arroja un mensaje de advertencia: “No se pudieron encontrar los grados de libertad apropiados para este modelo”
usaremos la prueba de Jarque-Bera, mediante el comando
#jb.norm.test(modelo37$residuals)
#jb.norm.test(modelo82$residuals)
Veamos la distribución de los residuos en ambos modelos: Modelo37
hist(modelo37$residuals)
Modelo82
hist(modelo82$residuals)
Otra forma de observar si existe normalidad es usar el qqnorm y el qqline
qqnorm(modelo37$residuals, col="blue")
qqline(modelo37$residuals, col="black")
Y observamos que el término de error no distribuye de forma normal, y es debido a los valores extremos que sucede la normalidad de los residuos.
Se corroboró que los dos modelos no superaron el diagnóstico de autocorrelación ni el de normalidad, esto se debe a que presentan alta volatilidad, por lo cual los residuos no se distribuyen de manera normal, debido a esto procederemos a estimar la volatilidad que presenta en el término de error, para ello nos apoyaremos de otros modelos, modelos de volatilidad, específicamente el modelo GARH. Cálculo de la varianza de los residuos: Var(x) = E(x^2) – (E(x))^2, siendo x los errores, sabemos que la esperanza de los errores es igual a cero, por lo tanto, tenemos; Var(ε) = E(ε^2), luego trabajaremos con los residuales del primer modelo, podemos trabajar con el segundo, pero considero que la diferencia entre uno y el otro es mínima, luego es indiferente con cualquiera de los dos modelos.
residuos2 <- modelo37$residuals
plot(residuos2, col="blue")
Se observa en el comportamiento de los residuos mucha volatilidad, luego procedemos a modelar la varianza de los residuales, pero con el cuadrado de los residuos, tal como calculamos anteriormente.
Procedemos a calcular el correlograma de la varianza de los residuos al cuadrado:
# Componente ARCH
# acf(ts(residuos2^2, frequency = 1))
# Componente GARCH
# pacf(ts(residuos2^2, frequency = 1))
# Usamos la libreria fGarch
library(fGarch)
# Buscamos los modelos con el comando garchFit
# modGARCH <- garchFit(formula = ~(1,2), data = residuos2, trace = TRUE)
# Asi procedemos con los modlos que se encuentran mediante los correlogramas y estimamos el modelo gasch
Luego de obtener los modelos tentativos GARCH, tomamos el que tenga el menor valor AIC,
Conclusión: al realizar los modelos ARIMA los residuales aparecen como valores nulos, lo cual imposibilita el realizar los modelos GARCH, sin embargo deje el código planteado, pero quiero saber por qué genera este error.