Mejor Modelo ARIMA GARCH
LUIS ANGEL BIZARRO VILCA
27/7/2020
Tipo de CAMBIO USDPEN
Primeramente obtenemos los datos del BCRP https://estadisticas.bcrp.gob.pe/estadisticas/series/mensuales/resultados/PN01208PM/html sobre el tipo de cambio descarmos el excel para este caso usaremos desde el año 2000-hasta la mas reciente:
#Para importar datos de excel necesitamos la libreria readxl que viene ya en el codigo base de r solo necesitamos llamarla; y debemos poner la direccion en la que se encuentre nuestre carpeta donde se descargo los datos del BCRP.
library(readxl)
Mensuales_20200726_163203 <- read_excel("C:/Users/Maritza B Vilca/Downloads/Mensuales-20200726-163203.xlsx")#Nos fijamos que hay una deferencia sutil pero importante la fecha esta en la columna 1 #pero nosotros queremos que este en la columna 0 porque al momento de graficar no nos #permitira graficar con fechas.
#Entonces con la ayuda de un artificio nota: tal vez haya formas de hacerlo de mejor modo pero para mi esta me resulto.
#creamos nuestras fechas acorde a nuestros datos el numero 246 es el numero de datos que tenemos en el excel la periodidad mensual (1 month), si quieres visualizar las fechas en la consla pon View
#install.packages("zoo")#y llamar la libreria xts
library(zoo)##
## Attaching package: 'zoo'## The following objects are masked from 'package:base':
##
## as.Date, as.Date.numericlibrary(xts)
fechas<-seq(as.Date("2000-01-01"), length.out = 246, by = "1 month")#Creamos nuestra base de datos pero solo de fechas
tipodecambioUDPPEN <- as.xts(Mensuales_20200726_163203[,-1], order.by = as.yearmon(as.Date(fechas)))# leemos la base que cargamos del excel y eliminamos la columna de las fechas y la ponemos junto al lado de nuestra variable que creamos "fechas" y todo eso lo acumulamos en la variable "tipodecambioUSDPEN"
head(tipodecambioUDPPEN)#Podemos ver entonces que ya arreglamos la base datos que nos proporcionaron## Tipodecambio
## ene 2000 3.497333
## feb 2000 3.454190
## mar 2000 3.440826
## abr 2000 3.477056
## may 2000 3.501500
## jun 2000 3.485333EMPEZAMOS AHORA SI A HACER LAS ESTIMACIONES
Una vez cargada los datos para averiguar cual es el mejor modelo ARMA y GARCH que son adecuadas para nuestro tipo de cambio primero tenemos que ver si la serie es estacionaria o no.
Entonces ploteamos la serie en niveles I(0)
plot(tipodecambioUDPPEN, main = "TC en niveles")
Graficamente podemos intuir que la series es no estacionaria y posee un randon walk; para corroborar esto usamos el test de dickey fuller.
Test de Dickey Fuller en niveles I(0)
De la ecuación \(Y_t = c+\beta_1 Y_{t-1}+e_{t}\) se plantea
\(H_0\):\(\beta=0\) La serie no es estacionaria y posee raiz unitaria
\(H_a\):\(\beta\neq 0\) La series es no estacionaria y no posee raiz unitaria
library(tseries)## Registered S3 method overwritten by 'quantmod':
## method from
## as.zoo.data.frame zooadf.test(tipodecambioUDPPEN)##
## Augmented Dickey-Fuller Test
##
## data: tipodecambioUDPPEN
## Dickey-Fuller = -1.0528, Lag order = 6, p-value = 0.9285
## alternative hypothesis: stationaryNuestro p-value 0.9285> 0.05 entonces concluimos que la serie es no estacionaria y posee una raiz unitaria. Conclusion aceptamos la hipotesis nula, pero nosotros necesitamos la estacionariedad la serie en intregrada I(0) no es estacionaria entonces vemos los retornos I(1) si lo son.
Retornos de la serie:
\(R_t= ln(\frac{p_t}{p_{t-1}})\)
retornos = diff(log(tipodecambioUDPPEN))
plot(retornos, main = "TC en retornos")
Aparentemente la serie seria estacionaria en el orden de integracion I(1), pero para corroborar usamos la el test de Dickey Fuller (ADF).
Test de Dickey Fuller en retornos I(0)
retornos=na.omit(retornos)
adf.test(retornos)## Warning in adf.test(retornos): p-value smaller than printed p-value##
## Augmented Dickey-Fuller Test
##
## data: retornos
## Dickey-Fuller = -6.0965, Lag order = 6, p-value = 0.01
## alternative hypothesis: stationaryEntonces segun las hipotesis planteadas 0.01<0.05 entonces no hay suficiente evidencia para aceptarla hipotesis nula entonces la series es estacionaria en orden I(1), entones ya esta identiicada y procedemos a identificar el mejor modelo ARMA.
Para hallar el mejor modelo ARMA
Graficamos los correlogramas para identificar los valores que podrian componer nuestro modelo ARMA(p,q), el PACF nos ayudara para el componente AR(p) y el MA(q).
par(mfrow=c(1,2))
acf(retornos, main = "ACF de los retornos")
pacf(retornos, main = "PACF de los retornos")
par(mfrow=c(1,1))Entonces por nuestros modelos seran p,q 1,2 en rango maximo. El mejor modelo ARMA de los retornos I(1) seran calculados.
modelo01= arma(retornos, order = c(0,1)) #AIC=-1516.5 #Elegimos el mejor AIC
#modelo02= arma(retornos, order = c(0,2)) #AIC=-1514.5
#modelo10= arma(retornos, order = c(1,0)) #AIC=-1512.2
#modelo11= arma(retornos, order = c(1,1)) #AIC=-1515.1
#modelo12= arma(retornos, order = c(1,2)) #AIC=-1511.6
summary(modelo01)##
## Call:
## arma(x = retornos, order = c(0, 1))
##
## Model:
## ARMA(0,1)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.0351987 -0.0053273 -0.0004134 0.0060351 0.0320105
##
## Coefficient(s):
## Estimate Std. Error t value Pr(>|t|)
## ma1 3.668e-01 6.390e-02 5.740 9.45e-09 ***
## intercept 1.263e-05 9.455e-04 0.013 0.989
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Fit:
## sigma^2 estimated as 0.0001181, Conditional Sum-of-Squares = 0.03, AIC = -1516.58Entonces es el mejor modelo de tipo de cambio es el modelo para los retornos I(0) ARMA(p,q). \(y_t=c+\phi_1y_{t-1}+ \phi_2 y_{t-2}+...+ \phi_p y_{t-p}+e_t+ \theta_1e_{t-1}+ \theta_2e_{t-2}+...+\theta_q\)
El mejor modelo sera un ARMA(0,1)
\(y_t=c+\theta_{1}e_{t-1}+e_t\)
\(y_t=0.00001263+0.3668e_{t-1}+e_t\)
#El mejor modelo GARCH
Para modelar un GARCH primero ploteamos los residuos del modelo ARMA(0,1).
plot(modelo01$residuals, type="l")
residuos1 = modelo01$residualsVemos que el modelo tine como un comportamiento de volatilidad entonces estimamos el modelo GARCH.
MEJOR MODELO GARCH
Del modelo01 modelamos los errores de la ecuacion de la forma de GARCH:
\(\sigma _{t}^{2}= \alpha_0+ \alpha_1\varepsilon_{t-1}^{2}+\alpha_2\varepsilon_{t-2}^{2}+...+\alpha_p\varepsilon_{t-p}^{2}+ \beta_1\sigma_{t-1}^{2}+\beta_2\sigma_{t-2}^{2}+...+\beta_q\sigma_{t-q}^{2}+ \mu_t\)
Para determinar el mejor modelo
primero hallamos los valores GARCH(p,q) buscamos los p y q con las funciones ACF y PACF al cuadrado para seleccionar.
residuos1=na.omit(residuos1)
par(mfrow=c(1,2))
pacf(residuos1^2)
acf(residuos1^2)
Los rezagos del GARCH(p,q) el p seria 2 y q 2 elegiremos el menor AIC
#install.packages("fGarch") #instalamos el packete para #rabajar el modelo garch
library("fGarch")## Loading required package: timeDate## Loading required package: timeSeries##
## Attaching package: 'timeSeries'## The following object is masked from 'package:zoo':
##
## time<-## Loading required package: fBasics#modGARCH11= garchFit(formula = ~garch(1,1),data = residuos1,trace = T)#AIC= -6.463559
modGARCH12= garchFit(formula = ~garch(1,2),data = residuos1,trace = T)#AIC=-6.468478##
## Series Initialization:
## ARMA Model: arma
## Formula Mean: ~ arma(0, 0)
## GARCH Model: garch
## Formula Variance: ~ garch(1, 2)
## ARMA Order: 0 0
## Max ARMA Order: 0
## GARCH Order: 1 2
## Max GARCH Order: 2
## Maximum Order: 2
## Conditional Dist: norm
## h.start: 3
## llh.start: 1
## Length of Series: 244
## Recursion Init: mci
## Series Scale: 0.01086553
##
## Parameter Initialization:
## Initial Parameters: $params
## Limits of Transformations: $U, $V
## Which Parameters are Fixed? $includes
## Parameter Matrix:
## U V params includes
## mu -0.01009147 0.01009147 0.001009147 TRUE
## omega 0.00000100 100.00000000 0.100000000 TRUE
## alpha1 0.00000001 0.99999999 0.100000000 TRUE
## gamma1 -0.99999999 0.99999999 0.100000000 FALSE
## beta1 0.00000001 0.99999999 0.400000000 TRUE
## beta2 0.00000001 0.99999999 0.400000000 TRUE
## delta 0.00000000 2.00000000 2.000000000 FALSE
## skew 0.10000000 10.00000000 1.000000000 FALSE
## shape 1.00000000 10.00000000 4.000000000 FALSE
## Index List of Parameters to be Optimized:
## mu omega alpha1 beta1 beta2
## 1 2 3 5 6
## Persistence: 0.9
##
##
## --- START OF TRACE ---
## Selected Algorithm: nlminb
##
## R coded nlminb Solver:
##
## 0: 330.43956: 0.00100915 0.100000 0.100000 0.400000 0.400000
## 1: 329.99172: 0.00100914 0.0879265 0.105146 0.384805 0.383627
## 2: 327.99557: 0.00100914 0.0977199 0.129028 0.384188 0.381497
## 3: 326.52451: 0.00100913 0.0916671 0.158155 0.356464 0.349386
## 4: 321.66933: 0.00100910 0.147819 0.212187 0.313974 0.295895
## 5: 315.31636: 0.00100899 0.128745 0.461697 0.200016 0.137850
## 6: 313.26021: 0.00100828 0.250354 0.466779 0.116284 0.139164
## 7: 311.04042: 0.00100821 0.188809 0.599887 0.133722 0.143428
## 8: 311.01327: 0.00100819 0.171686 0.615015 0.117963 0.127068
## 9: 310.73691: 0.00100818 0.182419 0.623373 0.123567 0.133630
## 10: 310.23213: 0.00100804 0.193241 0.713840 0.107609 0.128607
## 11: 309.88949: 0.00100782 0.185145 0.735884 0.0603790 0.137013
## 12: 309.78479: 0.00100684 0.199126 0.764838 0.0784732 0.0984030
## 13: 309.60406: 0.00100673 0.175315 0.805432 0.0767090 0.123601
## 14: 309.51020: 0.00100542 0.185930 0.839778 0.0374410 0.127705
## 15: 309.44449: 0.00100195 0.179420 0.856565 0.0489517 0.121555
## 16: 309.43470: 0.00100185 0.180628 0.869292 0.0537584 0.119120
## 17: 309.41727: 0.00100054 0.174997 0.878925 0.0565802 0.111482
## 18: 309.39868: 0.000990026 0.176770 0.900260 0.0415700 0.116474
## 19: 309.38607: 0.000977527 0.176018 0.925433 0.0444518 0.111879
## 20: 309.38566: 0.000962387 0.176926 0.922943 0.0445200 0.111142
## 21: 309.38525: 0.000933579 0.177104 0.921453 0.0445190 0.111195
## 22: 309.38198: 0.000637049 0.178050 0.914305 0.0444654 0.111306
## 23: 309.37583: -1.17850e-05 0.179034 0.906627 0.0443753 0.111481
## 24: 309.35953: -0.00186232 0.180396 0.895595 0.0441116 0.111843
## 25: 309.32319: -0.00646842 0.181599 0.884005 0.0435167 0.112601
## 26: 309.26745: -0.0100915 0.178424 0.902272 0.0433234 0.113318
## 27: 309.25280: -0.0100915 0.173706 0.928512 0.0440707 0.113668
## 28: 309.25271: -0.0100915 0.173809 0.931118 0.0437211 0.113299
## 29: 309.25271: -0.0100915 0.173781 0.931141 0.0437538 0.113320
## 30: 309.25271: -0.0100915 0.173782 0.931138 0.0437519 0.113319
##
## Final Estimate of the Negative LLH:
## LLH: -794.1543 norm LLH: -3.254731
## mu omega alpha1 beta1 beta2
## -0.0001096491 0.0000205167 0.9311384920 0.0437518654 0.1133191940
##
## R-optimhess Difference Approximated Hessian Matrix:
## mu omega alpha1 beta1 beta2
## mu -5869788.5089 41087043 -1.676876e+03 -842.1106 -2910.9059
## omega 41087042.5258 -61640700969 -1.218516e+06 -4307492.2588 -4453851.5767
## alpha1 -1676.8759 -1218516 -5.690991e+01 -135.9084 -127.8200
## beta1 -842.1106 -4307492 -1.359084e+02 -642.6556 -467.0319
## beta2 -2910.9059 -4453852 -1.278200e+02 -467.0319 -589.0788
## attr(,"time")
## Time difference of 0.05114293 secs
##
## --- END OF TRACE ---
##
##
## Time to Estimate Parameters:
## Time difference of 0.4210451 secs## Warning: Using formula(x) is deprecated when x is a character vector of length > 1.
## Consider formula(paste(x, collapse = " ")) instead.#modGARCH21= garchFit(formula = ~garch(2,1),data = residuos1,trace = T)#AIC=-6.457324
#modGARCH22= garchFit(formula = ~garch(2,2),data = residuos1,trace = T)#AIC=-6.460281
#SE ESCOGE EL MENOR AIC YA QUE MAXIMIZA LA VEROSIMILITUD
summary(modGARCH12)##
## Title:
## GARCH Modelling
##
## Call:
## garchFit(formula = ~garch(1, 2), data = residuos1, trace = T)
##
## Mean and Variance Equation:
## data ~ garch(1, 2)
## <environment: 0x000000001504e9c8>
## [data = residuos1]
##
## Conditional Distribution:
## norm
##
## Coefficient(s):
## mu omega alpha1 beta1 beta2
## -1.0965e-04 2.0517e-05 9.3114e-01 4.3752e-02 1.1332e-01
##
## Std. Errors:
## based on Hessian
##
## Error Analysis:
## Estimate Std. Error t value Pr(>|t|)
## mu -1.096e-04 4.223e-04 -0.260 0.79513
## omega 2.052e-05 6.450e-06 3.181 0.00147 **
## alpha1 9.311e-01 2.066e-01 4.507 6.56e-06 ***
## beta1 4.375e-02 6.708e-02 0.652 0.51425
## beta2 1.133e-01 7.376e-02 1.536 0.12448
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Log Likelihood:
## 794.1543 normalized: 3.254731
##
## Description:
## Wed Aug 12 18:53:28 2020 by user: Maritza B Vilca
##
##
## Standardised Residuals Tests:
## Statistic p-Value
## Jarque-Bera Test R Chi^2 11.73513 0.002829757
## Shapiro-Wilk Test R W 0.9853657 0.01325449
## Ljung-Box Test R Q(10) 30.29151 0.0007674473
## Ljung-Box Test R Q(15) 33.23602 0.004353382
## Ljung-Box Test R Q(20) 37.83014 0.009288309
## Ljung-Box Test R^2 Q(10) 4.509362 0.921458
## Ljung-Box Test R^2 Q(15) 7.592453 0.9390814
## Ljung-Box Test R^2 Q(20) 7.661467 0.9938811
## LM Arch Test R TR^2 5.442144 0.9415598
##
## Information Criterion Statistics:
## AIC BIC SIC HQIC
## -6.468478 -6.396814 -6.469295 -6.439616Entonces el mejor modelo seria un GARCH(1,2)
\(\sigma _{t}^{2}= (2.0517e-05)+ 9.3114e-01\varepsilon_{t-1}^{2}+ 4.3752e-02\sigma_{t-1}^{2}+1.1332e-01 \sigma_{t-2}^{2}+ -1.0965e-04\)
plot(modGARCH12@residuals,type="l")