En este capítulo se aborda el tercer componente del análisis de series de tiempo: la predicción. Se proporciona una introducción a los principales elementos del flujo de trabajo de la predicción, que incluye enfoques de entrenamiento, evaluación del rendimiento y métodos de referencia. Se cubren temas como el entrenamiento y prueba de modelos de predicción, la evaluación del rendimiento, la medición de errores, el uso de métodos de referencia y la cuantificación de la incertidumbre en las predicciones. En resumen, el capítulo ofrece una visión general de cómo realizar predicciones en el análisis de series de tiempo y cómo evaluar su calidad.

The forecasting workflow

La predicción tradicional de series temporales sigue el mismo flujo de trabajo que la mayoría de los campos del análisis predictivo, como la regresión o la clasificación, y suele incluir los siguientes pasos:

Preparación de los datos: Aquí se preparan los datos para el proceso de entrenamiento y prueba del modelo. Este paso incluye dividir las series en particiones de entrenamiento y de prueba, crear nuevas características (cuando proceda) y aplicar una transformación si es necesario.

Entrenar el modelo: Aquí utilizamos la partición de entrenamiento para entrenar un modelo estadístico. El objetivo principal de este paso es utilizar el conjunto de entrenamiento para entrenar, ajustar y estimar los coeficientes del modelo que minimicen los criterios de error seleccionados. Los valores ajustados y la estimación del modelo de las observaciones de la partición de entrenamiento se utilizarán más adelante para evaluar el rendimiento global del modelo.

Prueba del modelo: Aquí utilizamos el modelo entrenado para predecir las observaciones correspondientes de la partición de prueba. La idea principal aquí es evaluar el rendimiento del modelo con un nuevo conjunto de datos.

Evaluación del modelo: Por último, una vez que el modelo ha sido entrenado y probado, es el momento de evaluar el rendimiento global del modelo tanto en la partición de entrenamiento como en la de prueba. Si el modelo cumple con ciertos criterios de desempeño, se utiliza para generar la predicción final. De lo contrario, se selecciona un nuevo modelo o se ajusta el parámetro de entrenamiento y se repite el proceso de entrenamiento.

Training approaches

Uno de los elementos centrales del proceso de predicción es el entrenamiento del modelo. La calidad del entrenamiento del modelo tendrá un impacto directo en el resultado de la predicción. Los principales objetivos de este proceso son los siguientes: formalizar la relación de las series con otros factores, como patrones estacionales y de tendencia, y variables externas de manera predictiva; ajustar los parámetros del modelo y asegurar que sea escalable a nuevos datos y evite el sobreajuste. La serie se divide en particiones de entrenamiento y prueba en orden cronológico, para establecer una relación matemática entre la serie y sus valores pasados y errores.

Primero se cragan las librerías

library(TSstudio)
## Warning: package 'TSstudio' was built under R version 4.2.3
data(USgas)

Se observan las características principales de la función

ts_info(USgas)
##  The USgas series is a ts object with 1 variable and 238 observations
##  Frequency: 12 
##  Start time: 2000 1 
##  End time: 2019 10

Se usa la función window para dividir las series en particiones de entrenamiento y testeo

# Using the series time index to set the start and end point of each partiton
train <- window(USgas, 
                start = time(USgas)[1], 
                end = time(USgas)[length(USgas) - 12])

test <- window(USgas, 
               start = time(USgas)[length(USgas) - 12 + 1], 
               end = time(USgas)[length(USgas)])
ts_info(train)
##  The train series is a ts object with 1 variable and 226 observations
##  Frequency: 12 
##  Start time: 2000 1 
##  End time: 2018 10
ts_info(test)
##  The test series is a ts object with 1 variable and 12 observations
##  Frequency: 12 
##  Start time: 2018 11 
##  End time: 2019 10
# The sample.out argument set the size of the testing partition
# (and therefore the training partition)
USgas_partitions <- ts_split(USgas, sample.out = 12)

train <- USgas_partitions$train
test <- USgas_partitions$test
ts_info(train)
##  The train series is a ts object with 1 variable and 226 observations
##  Frequency: 12 
##  Start time: 2000 1 
##  End time: 2018 10
ts_info(test)
##  The test series is a ts object with 1 variable and 12 observations
##  Frequency: 12 
##  Start time: 2018 11 
##  End time: 2019 10

Este método es simple y rápido para entrenar y probar modelos, pero no permite evaluar la estabilidad y escalabilidad del rendimiento del modelo con una única prueba. El siguiente enfoque ayuda a mitigar este riesgo entrenando un modelo con múltiples particiones de entrenamiento y prueba.

Forecasting with backtesting

El enfoque de backtesting para entrenar un modelo de pronóstico es una versión avanzada del enfoque único fuera de muestra que usamos anteriormente. Se basa en el uso de una ventana rodante para dividir la serie en múltiples pares de particiones de entrenamiento y prueba. El proceso básico de entrenamiento de backtesting incluye la preparación de datos, entrenamiento del modelo, prueba del modelo y evaluación del modelo. La metodología de puntuación nos permite evaluar la estabilidad del modelo examinando la tasa de error del modelo en los diferentes conjuntos de prueba. Este método es similar al enfoque de validación cruzada para entrenar modelos de aprendizaje automático, pero las particiones son ordenadas cronológicamente en lugar de ser basadas en un muestreo aleatorio.

Forecast evaluation

El objetivo principal de la etapa de evaluación es evaluar la capacidad del modelo entrenado para predecir las observaciones futuras de la serie con precisión. Este proceso incluye lo siguiente:

  1. Residual analysis
  2. Scoring the forecast

Residual analysis

El análisis residual evalúa qué tan bien el modelo ha capturado los patrones de la serie y proporciona información sobre las distribuciones de residuos necesarias para construir intervalos de confianza para la predicción. Se utiliza herramientas de visualización de datos y pruebas estadísticas para evaluar la bondad del ajuste del modelo frente a los valores reales, la autocorrelación de los residuos y la distribución de los mismos. Los residuos con oscilación aleatoria alrededor de cero y con variación constante indican que el modelo es capaz de capturar la mayoría de la variación de la serie, mientras que una oscilación residual que no tiene las características del ruido blanco indica que el modelo no logró capturar los patrones de la serie.

Para demostrar esto, se puede entrenar un modelo ARIMA en la partición de entrenamiento que se creó antes para la serie de USgas.

Usamos la función auto.arima del paquete forecast:

library(forecast)
## Warning: package 'forecast' was built under R version 4.2.3
## Registered S3 method overwritten by 'quantmod':
##   method            from
##   as.zoo.data.frame zoo
md <- auto.arima(train)

checkresiduals(md)

## 
##  Ljung-Box test
## 
## data:  Residuals from ARIMA(2,1,1)(2,1,1)[12]
## Q* = 24.95, df = 18, p-value = 0.1263
## 
## Model df: 6.   Total lags used: 24

Los resultados indican que el modelo no captura completamente todos los patrones de la serie y que puede haber valores atípicos durante algunos períodos. Se sugiere que se modifiquen los parámetros del modelo para mejorar su ajuste.

Scoring the forecast

Ésta parte se centra en la importancia de evaluar el rendimiento de los modelos de aprendizaje automático. Una vez que se completa la sintonización del modelo, es necesario evaluar su capacidad para hacer predicciones precisas en observaciones que nunca ha visto antes. La evaluación del rendimiento se lleva a cabo utilizando métricas de precisión o error, como el Error Cuadrático Medio (MSE), la Raíz del Error Cuadrático Medio (RMSE), el Error Absoluto Medio (MAE) y el Error Porcentual Absoluto Medio (MAPE).

La elección de una métrica de error adecuada depende de los objetivos de precisión de la predicción y las limitaciones del modelo. La evaluación del modelo es un proceso iterativo y continuo que requiere ajustes y mejoras constantes para aumentar su precisión.

Una forma de ejemplificar esto es usando el modelo que se entrenó antes. Se usa la función forecast.

fc <- forecast(md, h = 12)

Luego se usa la función de accuracy para revisar el rendimiento del modelo.

accuracy(fc, test)
##                     ME      RMSE      MAE       MPE     MAPE      MASE
## Training set  5.844136  97.81626 73.42657 0.1170672 3.522348 0.6376860
## Test set     37.847885 103.22848 81.46603 1.3107987 3.261643 0.7075062
##                      ACF1 Theil's U
## Training set -0.004183172        NA
## Test set     -0.046708926 0.3404092

La función de accuracy devuelve métricas de error para los valores ajustados y la previsión real en el conjunto de entrenamiento y prueba.

Para evaluar el ajuste del modelo en el conjunto de entrenamiento y prueba, se puede usar la función test_forecast que muestra la serie real, los valores ajustados en el conjunto de entrenamiento y los valores pronosticados en el conjunto de prueba. Al pasar el cursor sobre los valores ajustados o pronosticados, aparecen los resultados de RMSE y MAPE tanto en las particiones de entrenamiento como de prueba

test_forecast(actual = USgas,
              forecast.obj = fc,
              test = test) 

Forecast benchmark

De acuerdo con las métricas de error, el modelo entrenado logró una puntuación de MAPE del 7,84% o un RMSE de 208,01. Para determinar si estos resultados son satisfactorios, es común comparar el rendimiento del modelo con alguna previsión de referencia o algún método ya existente que deseamos mejorar. Una forma popular de establecer una referencia es utilizar un enfoque de previsión simplista como punto de partida. Por ejemplo, podríamos pronosticar la serie con un enfoque naive y compararlo con el pronóstico anterior creado con el modelo ARIMA.

El enfoque naive asume que el valor más reciente es el representante verdadero del futuro, por lo que se utiliza este valor como base para el pronóstico futuro. Podemos crear un pronóstico ingenuo utilizando la función naive del paquete forecast, y utilizar el conjunto de entrenamiento como entrada para el modelo. De esta manera, podemos determinar si el rendimiento del modelo es aceptable en comparación con la referencia establecida.

Aquí llamamos la libreria Forecast y luego usamos la función test_forecast para testear y entrenar la función.

library(forecast)

naive_model <- naive(train, h  = 12)
test_forecast(actual = USgas,
              forecast.obj = naive_model,
              test = test)
accuracy(naive_model, test)
##                      ME     RMSE      MAE        MPE     MAPE     MASE
## Training set  -1.028444 285.6607 228.5084 -0.9218463 10.97123 1.984522
## Test set     301.891667 499.6914 379.1417  9.6798015 13.28187 3.292723
##                   ACF1 Theil's U
## Training set 0.3761105        NA
## Test set     0.7002486  1.499679
snaive_model <- snaive(train, h = 12)
test_forecast(actual = USgas,
              forecast.obj = snaive_model,
              test = test)
accuracy(snaive_model, test)
##                    ME     RMSE      MAE      MPE     MAPE     MASE       ACF1
## Training set 33.99953 148.7049 115.1453 1.379869 5.494048 1.000000  0.4859501
## Test set     96.45000 164.6967 135.8833 3.612060 5.220458 1.180103 -0.2120929
##              Theil's U
## Training set        NA
## Test set     0.4289964
accuracy(snaive_model, test)
##                    ME     RMSE      MAE      MPE     MAPE     MASE       ACF1
## Training set 33.99953 148.7049 115.1453 1.379869 5.494048 1.000000  0.4859501
## Test set     96.45000 164.6967 135.8833 3.612060 5.220458 1.180103 -0.2120929
##              Theil's U
## Training set        NA
## Test set     0.4289964

Finalizing the forecast

El último paso consiste en recalibrar los coeficientes del modelo con la serie completa. Para hacerlo, existen dos enfoques, dependiendo de si el modelo fue ajustado manualmente o automáticamente por un algoritmo. Si se ajustó manualmente, se deben utilizar los mismos parámetros de ajuste que se usaron en el modelo entrenado. Si se ajustó automáticamente, se puede extraer la configuración de parámetros usada con la partición de entrenamiento o permitir que el algoritmo ajuste los parámetros del modelo de nuevo usando la serie completa, asumiendo que tiene la capacidad de hacerlo correctamente al entrenar con nuevos datos.

Para simplificar todo, se usa el modelo auto.arima para entrenar el modelo final:

md_final <- auto.arima(USgas)

fc_final <- forecast(md_final, h = 12)

Se usa la siguiente función para graficar el pronóstico final:

plot_forecast(fc_final,
              title = "The US Natural Gas Consumption Forecast",
              Xtitle = "Year",
              Ytitle = "Billion Cubic Feet")

Handling forecast uncertainty

El propósito principal de los procesos de pronóstico, es minimizar el nivel de incertidumbre sobre valores de series futuras. Sin embargo, no es posible eliminar totalmente esta incertidumbre, pero podemos cuantificar y dar un rango para estimar el pronóstico. Se puede realizar usando intervalos de confianza, o usando un simulador.

Intervalos de confianza

Son un metodo de aproximación estadistica, que es usada para expresar el rango de prosibles valores, que contienen el valor real en algún punto de confianza.

La función forecast genera un intervalo de predicción con un nivel de confianza del 80% al 95%, pero se puede modificar.

Por ejemplo, se usan las funciones del modelo entrenado md_final y forecast, para los siguientes 60 meses, usando un intervalo de confianza del 80% al 90%.

fc_final2 <- forecast(md_final, 
                      h = 60, 
                      level = c(80, 90))

plot_forecast(fc_final2,
              title = "The US Natural Gas Consumption Forecast",
              Xtitle = "Year",
              Ytitle = "Billion Cubic Feet")

Simulation

La función forecast_sim da una función de construción, para que la simulación sea posible en campos de pronósticos. Este estimado se puede usar para calcular pronósticos de puntos estimados, o calcular probabilidades de tener diferentes valores.

Usamos el mismo modelo anterior para la función y hacemos 100 repeticiones:

fc_final3 <- forecast_sim(model = md_final,
                          h = 60, 
                          n = 500) 
library(plotly)
## Warning: package 'plotly' was built under R version 4.2.3
## Loading required package: ggplot2
## 
## Attaching package: 'plotly'
## The following object is masked from 'package:ggplot2':
## 
##     last_plot
## The following object is masked from 'package:stats':
## 
##     filter
## The following object is masked from 'package:graphics':
## 
##     layout
fc_final3$plot %>% 
  layout(title = "US Natural Gas Consumption - Forecasting Simulation",
         yaxis = list(title = "Billion Cubic Feet"),
         xaxis = list(title = "Year"))

Horse race approach

Para esta parte, se utiliza una función diferente a la dada por la guía, ya que ts_backtesting era una función desactualizada. El acercamiento Horse race es basado en el entrenamiento, testeo y evaluación de multiples modelos de pronóstico. En este ejemplo se aplica el acercamiento horse race para 7 modelos diferentes, que se revisarán en otros capitulos.

Los modelos que se evaluarán son:

-auto.arima -bsts -ets -hybrid -nnetar -tbats -HoltWinters

Usamos la función set.seed, para que se puedan arrojar los resultados:

set.seed(1234)
methods <- list(ets1 = list(method = "ets",
                            method_arg = list(opt.crit = "lik"),
                            notes = "ETS model with opt.crit = lik"),
                ets2 = list(method = "ets",
                            method_arg = list(opt.crit = "amse"),
                            notes = "ETS model with opt.crit = amse"),
                arima1 = list(method = "arima",
                              method_arg = list(order = c(2,1,0)),
                              notes = "ARIMA(2,1,0)"),
                arima2 = list(method = "arima",
                              method_arg = list(order = c(2,1,2),
                                                seasonal = list(order = c(1,1,1))),
                              notes = "SARIMA(2,1,2)(1,1,1)"),
                hw = list(method = "HoltWinters",
                          method_arg = NULL,
                          notes = "HoltWinters Model"),
                tslm = list(method = "tslm",
                            method_arg = list(formula = input ~ trend + season),
                            notes = "tslm model with trend and seasonal components"))
md <- train_model(input = USgas,
                  methods = methods,
                  train_method = list(partitions = 6, 
                                      sample.out = 12, 
                                      space = 3),
                  horizon = 12,
                  error = "MAPE")
## # A tibble: 6 × 7
##   model_id model       notes                     avg_m…¹ avg_r…² avg_c…³ avg_c…⁴
##   <chr>    <chr>       <chr>                       <dbl>   <dbl>   <dbl>   <dbl>
## 1 hw       HoltWinters HoltWinters Model          0.0482    144.   0.792   0.931
## 2 ets1     ets         ETS model with opt.crit …  0.0526    156.   0.833   0.972
## 3 arima2   arima       SARIMA(2,1,2)(1,1,1)       0.0546    163.   0.583   0.819
## 4 ets2     ets         ETS model with opt.crit …  0.0650    185.   0.5     0.792
## 5 tslm     tslm        tslm model with trend an…  0.0854    242.   0.319   0.611
## 6 arima1   arima       ARIMA(2,1,0)               0.163     539.   0.861   0.958
## # … with abbreviated variable names ¹​avg_mape, ²​avg_rmse, ³​`avg_coverage_80%`,
## #   ⁴​`avg_coverage_95%`

El modelo arroja una lista ordenada en el criterio de error que se da. En este caso, HoltWinters tiene el rate de error más bajo, por lo que es el modelo recomendado para usar.

Podemos graficar la taza de error y el pronóstico sugerido, usando el modelo que se observa a continuación:

library(forecast)
plot_error(md)