1.- INTRODUCCIÓN
La volatilidad es una característica inherente a las series de tiempo financieras. En general, no es constante y, en consecuencia, los modelos de series de tiempo tradicionales que suponen varianza homocedástica no son adecuados para modelar series de tiempo financieras. Engle (1982) introduce una nueva clase de procesos estocásticos llamados modelos ARCH, en los cuales la varianza condicionada a la información pasada no es constante y depende del cuadrado de las innovaciones pasadas. Bollerslev (1986) generaliza los modelos ARCH al proponer los modelos GARCH en los cuales la varianza condicional depende no solo de los cuadrados de las perturbaciones, como en Engle, sino además, de las varianzas condicionales de períodos anteriores.
En 1991, Nelson presenta los modelos EGARCH, en los cuales formula para la varianza condicional un modelo que no se comporta de manera simétrica para perturbaciones positivas y negativas, como sucede en los modelos GARCH; expresando otro rasgo de la volatilidad: su comportamiento asimétrico frente a las alzas y bajas de los precios de un activo financiero. Un elevado número de trabajos sobre modelos de volatilidad se han publicado en las últimas décadas. Ver Poon y Granger (2003), Hansen y Lunde (2006) y Novales y Gracia (1993).
# Bibliotecas
#---------------------
library(tidyquant)
## Loading required package: lubridate
##
## Attaching package: 'lubridate'
## The following objects are masked from 'package:base':
##
## date, intersect, setdiff, union
## Loading required package: PerformanceAnalytics
## Loading required package: xts
## Loading required package: zoo
##
## Attaching package: 'zoo'
## The following objects are masked from 'package:base':
##
## as.Date, as.Date.numeric
##
## Attaching package: 'PerformanceAnalytics'
## The following object is masked from 'package:graphics':
##
## legend
## Loading required package: quantmod
## Loading required package: TTR
## Registered S3 method overwritten by 'quantmod':
## method from
## as.zoo.data.frame zoo
library(TTR)
library(ggplot2)
library(knitr)
library(nortest)
library(MASS)
library(stats)
library(moments)
##
## Attaching package: 'moments'
## The following objects are masked from 'package:PerformanceAnalytics':
##
## kurtosis, skewness
library(fGarch)
## NOTE: Packages 'fBasics', 'timeDate', and 'timeSeries' are no longer
## attached to the search() path when 'fGarch' is attached.
##
## If needed attach them yourself in your R script by e.g.,
## require("timeSeries")
##
## Attaching package: 'fGarch'
## The following object is masked from 'package:TTR':
##
## volatility
library(fBasics)
##
## Attaching package: 'fBasics'
## The following objects are masked from 'package:moments':
##
## kurtosis, skewness
## The following object is masked from 'package:TTR':
##
## volatility
## The following objects are masked from 'package:PerformanceAnalytics':
##
## kurtosis, skewness
#Datos
#----------------------
getSymbols("SU.PA", from = "2000-01-03")
## [1] "SU.PA"
summary(SU.PA)
## Index SU.PA.Open SU.PA.High SU.PA.Low
## Min. :2000-01-03 Min. : 18.36 Min. : 19.32 Min. : 18.36
## 1st Qu.:2005-12-15 1st Qu.: 34.71 1st Qu.: 35.20 1st Qu.: 34.24
## Median :2012-01-11 Median : 51.00 Median : 51.75 Median : 50.42
## Mean :2012-01-18 Mean : 61.43 Mean : 62.13 Mean : 60.72
## 3rd Qu.:2018-02-10 3rd Qu.: 69.97 3rd Qu.: 70.58 3rd Qu.: 69.24
## Max. :2024-03-04 Max. :211.05 Max. :212.40 Max. :209.90
## SU.PA.Close SU.PA.Volume SU.PA.Adjusted
## Min. : 18.73 Min. : 0 Min. : 10.10
## 1st Qu.: 34.74 1st Qu.: 1041965 1st Qu.: 18.90
## Median : 51.05 Median : 1499607 Median : 37.35
## Mean : 61.45 Mean : 1773857 Mean : 50.26
## 3rd Qu.: 69.90 3rd Qu.: 2208358 3rd Qu.: 59.97
## Max. :212.40 Max. :26582602 Max. :212.40
plot(SU.PA$SU.PA.Adjusted, main = "Cotización ajustada")
autoplot(SU.PA$SU.PA.Adjusted)+labs(title = "Cotización ajustada") + xlab("Tiempo") + ylab("Cotización Ajustada")
chart_Series(SU.PA$SU.PA.Adjusted, theme = chart_theme())
chartSeries(SU.PA$SU.PA.Adjusted, name = "Cotización Ajustada")
#Rentabilidad
#------------------------
Rendimiento = dailyReturn(SU.PA$SU.PA.Adjusted)
chartSeries(Rendimiento, name = "Rentabilidad")
2.- Modelos ARCH
Paso 1
Primer paso es crear el modelo de series temporales tradicional ARIMA(p, d, q)
library(forecast)
modelo1 = auto.arima(Rendimiento)
summary(modelo1)
## Series: Rendimiento
## ARIMA(2,0,1) with non-zero mean
##
## Coefficients:
## ar1 ar2 ma1 mean
## 0.6377 -0.0461 -0.6734 6e-04
## s.e. 0.0873 0.0158 0.0869 2e-04
##
## sigma^2 = 0.0004108: log likelihood = 15414.08
## AIC=-30818.16 AICc=-30818.15 BIC=-30784.48
##
## Training set error measures:
## ME RMSE MAE MPE MAPE MASE ACF1
## Training set 1.491024e-05 0.02026051 0.01424757 NaN Inf 0.681405 0.0001998101
autoplot(modelo1)
Modelo es estacionario (primer círculo es correcto, al estar los puntos
rojos dentro del círculo) y es invertible (por la misma razón respecto
al segundo círculo).
checkresiduals(modelo1)
##
## Ljung-Box test
##
## data: Residuals from ARIMA(2,0,1) with non-zero mean
## Q* = 22.064, df = 7, p-value = 0.002477
##
## Model df: 3. Total lags used: 10
Claramente observamos que estamos ante un modelo ARMA (2,1). Por lo consiguiente será calcular una regresión lineal de nuestros errores o residuales del modelo ARMA al cuadrado como variable dependiente, Y, y como variable independiente, los mismos errores del modelo ARMA al cuadrado con un rezago (“lag”), X.
(ARMA es ARIMA sin diferencias)
Paso 2
# Estimación de la varianza del proceso
#-------------------------------------------
Errores_cuadrado = resid(modelo1)^2
Errores_cuadrado = resid(modelo1)^2
chartSeries(Errores_cuadrado)
Observamos que nuestra varianza no es constante, ya que los errores al muestran que al pasar los días la varianza llega a ser heterocedástica, por lo que ahora haremos nuestra regresión. Esto para observar si nos enfrentamos a un modelo ARCH o no.
Paso 3
library(dynlm)
## Warning: package 'dynlm' was built under R version 4.3.3
Regresion1 = dynlm(Errores_cuadrado ~L(Errores_cuadrado))
summary(Regresion1)
##
## Time series regression with "ts" data:
## Start = 2, End = 6215
##
## Call:
## dynlm(formula = Errores_cuadrado ~ L(Errores_cuadrado))
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.005318 -0.000357 -0.000286 -0.000022 0.039167
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 3.579e-04 1.546e-05 23.16 <2e-16 ***
## L(Errores_cuadrado) 1.282e-01 1.258e-02 10.19 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.001149 on 6212 degrees of freedom
## Multiple R-squared: 0.01643, Adjusted R-squared: 0.01627
## F-statistic: 103.8 on 1 and 6212 DF, p-value: < 2.2e-16
Vemos que las variables son estadísticamente significativas, al ser el p valor igual a cero (por tanto, inferior a 0.05).
Contraste ARCH
H0 = No hay efectos ARCH si el p-value es mayor que 0.05.
H1: Si hay efectos ARCH si el p-value es menor que 0.05.
P valor = 0.0000 < 0.05 si tenemos efectos arch.
Hay que estimar modelos ARCH y GARCH
Paso 4
autoplot(acf(Errores_cuadrado, lag.max = 50, ylim = c(-1,1))) +
labs(title = "FAS") +
xlab("Lags") +
ylab("autocorrelación")
Con el correlograma, los errores al cuadrado presentan una estructura de
autocorrelación (tiene efectos que sobresalen la línea discontinua)
–> Efectos ARCH. Por tanto, hay heterocedasticidad condicional.
autoplot(acf(Errores_cuadrado, lag.max = 50, ylim = c(-1,1))) +
labs(title = "FACP") +
xlab("Lags") +
ylab("autocorrelación")
Con el correlograma, los errores al cuadrado presentan una estructura de
autocorrelación (tiene efectos que sobresalen la línea discontinua)
–> Efectos ARCH. Por tanto, hay heterocedasticidad condicional.
Conclusiones:
A través de nuestro mejor modelo ARMA, vemos que es de orden (2,1)
Hemos obtenido los residuales al cuadrado y hemos hecho una regresión con ellos
Vemos que el modelo es significativo a través del p-value
Observamos que tanto en la Autocorrelación como la Autocorrelación Parcial llegan a salir del límite las observaciones, por lo que podemos concluir que nuestra varianza es HETEROCEDÁSTICA.
Prueba ARCH
library(tseries)
library(fGarch)
error = resid(modelo1)
# ModeloARCH1 = ArchTest(error)
library(rugarch)
## Warning: package 'rugarch' was built under R version 4.3.3
## Loading required package: parallel
##
## Attaching package: 'rugarch'
## The following objects are masked from 'package:fBasics':
##
## qgh, qnig
## The following object is masked from 'package:stats':
##
## sigma
# Modelo ARMA (1,1)
ugarch1 = ugarchspec()
#Resumen
ugarch1
##
## *---------------------------------*
## * GARCH Model Spec *
## *---------------------------------*
##
## Conditional Variance Dynamics
## ------------------------------------
## GARCH Model : sGARCH(1,1)
## Variance Targeting : FALSE
##
## Conditional Mean Dynamics
## ------------------------------------
## Mean Model : ARFIMA(1,0,1)
## Include Mean : TRUE
## GARCH-in-Mean : FALSE
##
## Conditional Distribution
## ------------------------------------
## Distribution : norm
## Includes Skew : FALSE
## Includes Shape : FALSE
## Includes Lambda : FALSE
# Modelo ARMA (2,1) <--- Modelo efectivo
ugarch2 = ugarchspec(mean.model = list(armaOrder = c(2,1)))
#Resumen
ugarch2
##
## *---------------------------------*
## * GARCH Model Spec *
## *---------------------------------*
##
## Conditional Variance Dynamics
## ------------------------------------
## GARCH Model : sGARCH(1,1)
## Variance Targeting : FALSE
##
## Conditional Mean Dynamics
## ------------------------------------
## Mean Model : ARFIMA(2,0,1)
## Include Mean : TRUE
## GARCH-in-Mean : FALSE
##
## Conditional Distribution
## ------------------------------------
## Distribution : norm
## Includes Skew : FALSE
## Includes Shape : FALSE
## Includes Lambda : FALSE
ugfit = ugarchfit(spec = ugarch2, data = Rendimiento)
ugfit
##
## *---------------------------------*
## * GARCH Model Fit *
## *---------------------------------*
##
## Conditional Variance Dynamics
## -----------------------------------
## GARCH Model : sGARCH(1,1)
## Mean Model : ARFIMA(2,0,1)
## Distribution : norm
##
## Optimal Parameters
## ------------------------------------
## Estimate Std. Error t value Pr(>|t|)
## mu 0.000968 0.000166 5.82794 0.000000
## ar1 0.712512 0.101065 7.05005 0.000000
## ar2 -0.012600 0.017051 -0.73892 0.459955
## ma1 -0.753008 0.100223 -7.51333 0.000000
## omega 0.000005 0.000002 2.82672 0.004703
## alpha1 0.075853 0.008094 9.37107 0.000000
## beta1 0.912926 0.009777 93.37251 0.000000
##
## Robust Standard Errors:
## Estimate Std. Error t value Pr(>|t|)
## mu 0.000968 0.000173 5.61149 0.000000
## ar1 0.712512 0.112506 6.33309 0.000000
## ar2 -0.012600 0.018992 -0.66343 0.507055
## ma1 -0.753008 0.111602 -6.74729 0.000000
## omega 0.000005 0.000005 0.89463 0.370986
## alpha1 0.075853 0.017902 4.23704 0.000023
## beta1 0.912926 0.025339 36.02816 0.000000
##
## LogLikelihood : 16211.7
##
## Information Criteria
## ------------------------------------
##
## Akaike -5.2147
## Bayes -5.2071
## Shibata -5.2147
## Hannan-Quinn -5.2121
##
## Weighted Ljung-Box Test on Standardized Residuals
## ------------------------------------
## statistic p-value
## Lag[1] 0.5068 0.4765
## Lag[2*(p+q)+(p+q)-1][8] 3.9324 0.8220
## Lag[4*(p+q)+(p+q)-1][14] 7.8829 0.3924
## d.o.f=3
## H0 : No serial correlation
##
## Weighted Ljung-Box Test on Standardized Squared Residuals
## ------------------------------------
## statistic p-value
## Lag[1] 0.08111 0.7758
## Lag[2*(p+q)+(p+q)-1][5] 0.86905 0.8885
## Lag[4*(p+q)+(p+q)-1][9] 2.22861 0.8760
## d.o.f=2
##
## Weighted ARCH LM Tests
## ------------------------------------
## Statistic Shape Scale P-Value
## ARCH Lag[3] 0.08564 0.500 2.000 0.7698
## ARCH Lag[5] 1.78063 1.440 1.667 0.5218
## ARCH Lag[7] 2.38617 2.315 1.543 0.6360
##
## Nyblom stability test
## ------------------------------------
## Joint Statistic: 2.706
## Individual Statistics:
## mu 0.0600
## ar1 0.1669
## ar2 0.1767
## ma1 0.1807
## omega 0.1411
## alpha1 0.2590
## beta1 0.2918
##
## Asymptotic Critical Values (10% 5% 1%)
## Joint Statistic: 1.69 1.9 2.35
## Individual Statistic: 0.35 0.47 0.75
##
## Sign Bias Test
## ------------------------------------
## t-value prob sig
## Sign Bias 0.1865 0.85205
## Negative Sign Bias 1.5671 0.11714
## Positive Sign Bias 1.1244 0.26089
## Joint Effect 6.3289 0.09666 *
##
##
## Adjusted Pearson Goodness-of-Fit Test:
## ------------------------------------
## group statistic p-value(g-1)
## 1 20 137.2 6.315e-20
## 2 30 155.6 2.820e-19
## 3 40 165.4 1.695e-17
## 4 50 178.2 1.451e-16
##
##
## Elapsed time : 0.8097179
ugfit@fit$coef
## mu ar1 ar2 ma1 omega
## 9.683832e-04 7.125119e-01 -1.259961e-02 -7.530080e-01 4.853590e-06
## alpha1 beta1
## 7.585300e-02 9.129258e-01
ug_var = ugfit@fit$var
autoplot(ts(ug_var)) +
geom_line(color = "green") + labs(title = "Varianza de nuestro modelo ARMA (2,1)")
Esta gráfica muestra la varianza condicional.
ug_resid = (ugfit@fit$residuals)^2
autoplot(ts(ug_resid)) + autolayer(ts(ug_resid), series = "Residuales al cuadrado") + autolayer(ts(ug_var), series = "Varianza") + labs(title = "Residuales al cuadrado y Varianza de nuestro modelo ARMA (2,1)") + ylab("") + xlab("Tiempo")
ug_forecast = ugarchforecast(ugfit, n.ahead = 30)
ug_forecast
##
## *------------------------------------*
## * GARCH Model Forecast *
## *------------------------------------*
## Model: sGARCH
## Horizon: 30
## Roll Steps: 0
## Out of Sample: 0
##
## 0-roll forecast [T0=2024-03-04]:
## Series Sigma
## T+1 2.962e-05 0.01403
## T+2 1.616e-04 0.01412
## T+3 4.054e-04 0.01422
## T+4 5.774e-04 0.01431
## T+5 6.969e-04 0.01440
## T+6 7.799e-04 0.01448
## T+7 8.375e-04 0.01457
## T+8 8.775e-04 0.01465
## T+9 9.053e-04 0.01474
## T+10 9.246e-04 0.01482
## T+11 9.380e-04 0.01490
## T+12 9.473e-04 0.01498
## T+13 9.537e-04 0.01506
## T+14 9.582e-04 0.01513
## T+15 9.613e-04 0.01521
## T+16 9.635e-04 0.01528
## T+17 9.650e-04 0.01536
## T+18 9.660e-04 0.01543
## T+19 9.667e-04 0.01550
## T+20 9.672e-04 0.01557
## T+21 9.676e-04 0.01564
## T+22 9.678e-04 0.01570
## T+23 9.680e-04 0.01577
## T+24 9.681e-04 0.01583
## T+25 9.682e-04 0.01590
## T+26 9.683e-04 0.01596
## T+27 9.683e-04 0.01602
## T+28 9.683e-04 0.01609
## T+29 9.683e-04 0.01615
## T+30 9.684e-04 0.01621
# plot(ug_forecast)
Hasta aquí los modelos de heterocedasticidad condicionada.