Series de Tiempo para tasas de descuento de enero 2021 a octubre de 2024

# Se crean dos variables adicionales de año (Year) y Mes (Month) para generar una gráfica que compara las series por mes.

Montos$Year <- factor(format(Montos$Fecha, "%Y"), levels = sort(unique(format(Montos$Fecha, "%Y")), decreasing =
T))

#Se crean los datos para el primer modelo con formato de series de tiempo
montos.ts <- ts(Montos[, 2], start = 2021, freq = 52)
descuento.ts <- ts(Descuento[, 2], start = 2021, freq = 52)

plot.montos<-
autoplot(montos.ts)+scale_y_continuous(labels = scales::unit_format(unit = "M", scale = 1e-6))+
ylab("Millones de pesos") + xlab("Semanas") +
ggtitle("Ventas de Certificados de la Tesorería 28 días,  enero 2021 a octubre 2024")+theme_light()

plot.descuento<-
autoplot(descuento.ts)+ylab("Tas de descuento") + xlab("Semanas") +
ggtitle("Tasas de descuento de Certificados de la Tesorería 28 días,  enero 2021 a octubre 2024")+
theme_light()

plot.montos+plot.descuento

ggplot(Descuento, aes(SF116766))+geom_histogram(stat = "bin",fill="darkblue")+ggtitle("Cetes a 28 días,  enero 2021 a octubre 2024")+theme_light()+xlab("Tasas de descuento")+ylab("Frecuencias")

#Prueba de Dickey-Fuller para verificación estacionaridad de la serie

#H0: La series es estacionaria 
adf.test(descuento.ts)
## 
##  Augmented Dickey-Fuller Test
## 
## data:  descuento.ts
## Dickey-Fuller = -2.1104, Lag order = 5, p-value = 0.5298
## alternative hypothesis: stationary
acf(Descuento$SF116766)

#Se el operador de diferencia
descuento.diff<-diff(Descuento$SF116766,1)

adf.test(na.omit(descuento.diff))
## 
##  Augmented Dickey-Fuller Test
## 
## data:  na.omit(descuento.diff)
## Dickey-Fuller = -5.3305, Lag order = 5, p-value = 0.01
## alternative hypothesis: stationary
acf(descuento.diff, main="ACF", na.action = na.pass)

pacf(descuento.diff, main="PACF",na.action = na.pass)

#Se usará un modelo ARIMA(1,1,1)

auto.arima(descuento.ts, seasonal = FALSE)
## Series: descuento.ts 
## ARIMA(1,2,2) 
## 
## Coefficients:
##          ar1      ma1     ma2
##       0.1617  -1.5399  0.6074
## s.e.  0.1298   0.1020  0.1002
## 
## sigma^2 = 0.01829:  log likelihood = 116.97
## AIC=-225.93   AICc=-225.73   BIC=-212.72

El orden de un proceso autorregresivo \(AR(p)\)

La función de autocorrelación parcial ACP refleja características que dependen del orden proceso y del tipo de parámetros involucrados. La función de autocorrelación parcial es el equivalente a la ACF para determinar el orden de proceso autorregresivo \(AR(p)\).

En términos de modelos de regresión lineales, la autocorrelación parcial se usa cuando se quiere cuantificar la correlación entre \(Y\) y \(Z\) sin considerar a \(X\) y se calcula como

\[\rho_{YZ\backslash X}=\frac{\rho_{YZ}-\rho_{YX}\rho_{XZ}}{\sqrt{1-\rho^2_{YX}}\sqrt{1-\rho^2_{XZ}}}\]

Consideraciones acerca de la función de autocorrelación parcial

Modelos a evaluar

  • Modelo ARIMA(1,1,1)
  • Modelo ARIMA(1,2,2)
  • Modelo de Holt-Winter

Se tienen 198 datos, se usarán los datos de 2021 a 2023 (apróximadamente 80% de los datos) para el modelo de entrenamiento

decompose(descuento.ts) %>% plot()

auto.arima(descuento.ts, seasonal = FALSE)
## Series: descuento.ts 
## ARIMA(1,2,2) 
## 
## Coefficients:
##          ar1      ma1     ma2
##       0.1617  -1.5399  0.6074
## s.e.  0.1298   0.1020  0.1002
## 
## sigma^2 = 0.01829:  log likelihood = 116.97
## AIC=-225.93   AICc=-225.73   BIC=-212.72
entrena.descuento <- window(descuento.ts, start = 2021, end = c(2023, 52))

nValid <- 42
 nTrain <- length(descuento.ts) - nValid
 train.ts <- window(descuento.ts, start = c(2021, 1), end = c(2021, nTrain))
 valid.ts <- window(descuento.ts, start = c(2021, nTrain + 1), end = c(2021, nTrain + nValid))
 descuento.train.HW <-  HoltWinters(train.ts)
 descuento.train.HW.pred <- forecast(descuento.train.HW, h = nValid, level = 0)
 plot(descuento.train.HW.pred, ylim = c(3, 12),  ylab = "Tasa de descuento", xlab = "Tiempo", bty = "l",
      xaxt = "n", xlim = c(2021,2024+(42/22)), main = "Modelo de Holt-Winters", flty = 2)
 axis(1, at = seq(2021, 2024, 1), labels = format(seq(2021, 2024, 1)))
 lines(descuento.train.HW.pred$fitted, lwd = 2)
 lines(valid.ts)

 descuento.train.ARIMA122 <-  arima(train.ts, order = c(1,2,2))
 descuento.train.ARIMA122.pred <- forecast(descuento.train.ARIMA122, h = nValid, level = 0)
 plot(descuento.train.ARIMA122.pred, ylim = c(3, 12),  ylab = "Tasa de descuento", xlab = "Tiempo", bty = "l",
      xaxt = "n", xlim = c(2021,2024+(42/22)), main = "Modelo de ARIMA(1,2,2)", flty = 2)
 axis(1, at = seq(2021, 2024, 1), labels = format(seq(2021, 2024, 1)))
 lines(descuento.train.ARIMA122.pred$fitted, lwd = 2)
 lines(valid.ts)

descuento.train.ARIMA111 <-  arima(train.ts, order = c(1,1,1))
 descuento.train.ARIMA111.pred <- forecast(descuento.train.ARIMA111, h = nValid, level = 0)
 plot(descuento.train.ARIMA111.pred, ylim = c(3, 12),  ylab = "Tasa de descuento", xlab = "Tiempo", bty = "l",
      xaxt = "n", xlim = c(2021,2024+(42/22)), main = "Modelo de ARIMA(1,1,1)", flty = 2)
 axis(1, at = seq(2021, 2024, 1), labels = format(seq(2021, 2024, 1)))
 lines(descuento.train.ARIMA111.pred$fitted, lwd = 2)
 lines(valid.ts)

## Comparación de los modelos

autoplot(descuento.ts) +
  autolayer(descuento.train.HW.pred, series="Holt-Winters", PI=FALSE) +
   autolayer(descuento.train.ARIMA111.pred, series="ARIMA(1,1,1)", PI=FALSE) +
   autolayer(descuento.train.ARIMA122.pred, series="ARIMA(1,2,2)", PI=FALSE) +
   xlab("Tiempo") + ylab("Tasa de descuento") +
   ggtitle("Predicciones de tasas de descuento CETES 28 de enero 2023 a octubre 2024") +
   guides(colour=guide_legend(title="Periodo usado\n para predicción"))+
   theme_linedraw()

checkresiduals(descuento.train.ARIMA111)

## 
##  Ljung-Box test
## 
## data:  Residuals from ARIMA(1,1,1)
## Q* = 94.774, df = 30, p-value = 1.219e-08
## 
## Model df: 2.   Total lags used: 32
tabla_precision<-
rbind(accuracy(descuento.train.ARIMA111.pred$mean, valid.ts), 
 accuracy(descuento.train.ARIMA122.pred$mean, valid.ts),
      accuracy(descuento.train.HW.pred$mean, valid.ts) )
 

row.names(tabla_precision) <- c("ARIMA(1,1,1)", "ARIMA(1,2,2)", "Holt-Winters")

knitr::kable(tabla_precision)
ME RMSE MAE MPE MAPE ACF1 Theil’s U
ARIMA(1,1,1) -1.7609508 1.9257767 1.7609508 -23.093581 23.093581 0.9191204 13.348818
ARIMA(1,2,2) -1.3079610 1.4239092 1.3079610 -17.127721 17.127721 0.9071202 9.838483
Holt-Winters -0.5949666 0.6277743 0.5949666 -7.558208 7.558208 0.3271790 4.154881