Reducir el impacto del Covid en modelos de forecasting

Si te gusta  Skforecast , ayúdanos dándonos una estrella en  GitHub! ⭐️

Reducir el impacto del Covid en modelos de forecasting

Joaquin Amat Rodrigo, Javier Escobar Ortiz
Diciembre, 2022 (última actualización Agosto 2024)

Introducción

Si bien en muchos casos reales de forecasting se dispone de datos históricos, no todos son fiables. Algunos ejemplos de estos escenarios son:

  • Sensores IoT: dentro del Internet de las cosas, los sensores se encargan de capturan los datos del mundo físico. A menudo los sensores se despliegan o instalan en entornos hostiles, lo que inevitablemente implica que los sensores sean propensos a fallos, mal funcionamiento y rápido deterioro, provocando que el sensor produzca lecturas inusuales y erróneas.

  • Paradas en instalaciones: cada cierto periodo de funcionamiento, las fábricas necesitan ser paradas para actividades de reparación, revisión o mantenimiento. Estos eventos provocan que la producción se detenga, generando un hueco en los datos.

  • Pandemia (Covid-19): la pandemia de Covid 19 cambió significativamente el comportamiento de la población, impactando directamente en muchas series temporales como la producción, las ventas y el transporte.

La presencia de valores no fiables o no representativos en el historial de datos es un problema importante, ya que dificulta el aprendizaje de los modelos. Para la mayoría de los algoritmos de forecasting, eliminar esa parte de los datos no es una opción, ya que requieren que la serie temporal sea completa. Una solución alternativa es reducir el peso de las observaciones afectadas durante el entrenamiento del modelo. Este documento cómo skforecast facilita la aplicación de esta estrategia.

✎ Nota

En los ejemplos siguientes, una parte de la serie temporal se excluye del entrenamiento del modelo dándole un peso de cero. Sin embargo, el uso de pesos no se limita a incluir o excluir observaciones, sino a equilibrar el grado de influencia de cada observación en el modelo de forecasting. Por ejemplo, una observación con un peso de 10 tiene 10 veces más impacto en el entrenamiento del modelo que una observación con un peso de 1.

Warning

En la mayoría de las implementaciones de gradient boosting (LightGBM, XGBoost, CatBoost), las observaciones con peso cero se ignoran al calcular los gradientes y hessianos. Sin embargo, los valores de esas observaciones todavía se consideran al construir los histogramas de las variables. Por lo tanto, el modelo resultante puede diferir del modelo entrenado sin las muestras con peso cero. Vea más detalles en este issue.

Librerías

In [21]:
# Procesamiento de datos
# ==============================================================================
import pandas as pd
import numpy as np
from skforecast.datasets import fetch_dataset

# Gráficos
# ==============================================================================
import matplotlib.pyplot as plt
from skforecast.plot import set_dark_theme

# Modelado y Forecasting
# ==============================================================================
import skforecast
import sklearn
from sklearn.linear_model import Ridge
from skforecast.ForecasterAutoreg import ForecasterAutoreg
from skforecast.model_selection import backtesting_forecaster

# Configuracion
# ==============================================================================
import warnings
warnings.filterwarnings('once')

color = '\033[1m\033[38;5;208m'
print(f"{color}Version skforecast: {skforecast.__version__}")
print(f"{color}Version scikit-learn: {sklearn.__version__}")
print(f"{color}Version pandas: {pd.__version__}")
print(f"{color}Version numpy: {np.__version__}")
Version skforecast: 0.13.0
Version scikit-learn: 1.5.1
Version pandas: 2.2.2
Version numpy: 2.0.1

Covid-19

Durante el periodo de confinamiento impuesto como consecuencia de la pandemia de covid-19, el comportamiento de la población se vio alterado. Un ejemplo de ello se puede ver en el uso del servicio de alquiler de bicicletas en la ciudad de Madrid (España).

Datos

In [22]:
# Descarga de datos
# ==============================================================================
data = fetch_dataset('bicimad')
data.head()
bicimad
-------
This dataset contains the daily users of the bicycle rental service (BiciMad) in
the city of Madrid (Spain) from 2014-06-23 to 2022-09-30.
The original data was obtained from: Portal de datos abiertos del Ayuntamiento
de Madrid https://datos.madrid.es/portal/site/egob
Shape of the dataset: (3022, 1)
Out[22]:
users
date
2014-06-23 99
2014-06-24 72
2014-06-25 119
2014-06-26 135
2014-06-27 149
In [23]:
# División datos entrenamiento, validation y test
# ==============================================================================
data = data.loc['2020-01-01': '2021-12-31']
end_train = '2021-06-01'
data_train = data.loc[: end_train, :]
data_test  = data.loc[end_train:, :]

print(f"Fechas train : {data_train.index.min()} --- {data_train.index.max()}  (n={len(data_train)})")
print(f"Fechas test  : {data_test.index.min()} --- {data_test.index.max()}  (n={len(data_test)})")
Fechas train : 2020-01-01 00:00:00 --- 2021-06-01 00:00:00  (n=518)
Fechas test  : 2021-06-01 00:00:00 --- 2021-12-31 00:00:00  (n=214)
In [24]:
# Time series plot
# ==============================================================================
set_dark_theme()
fig, ax = plt.subplots(figsize=(8, 3))
data_train.users.plot(ax=ax, label='train', linewidth=1)
data_test.users.plot(ax=ax, label='test', linewidth=1)
ax.axvspan(
    pd.to_datetime('2020-03-16'),
    pd.to_datetime('2020-04-21'), 
    label="Covid-19 confinamiento",
    color="red",
    alpha=0.3
)

ax.axvspan(
    pd.to_datetime('2020-04-21'),
    pd.to_datetime('2020-05-31'), 
    label="Tiempo de recuperación",
    color="white",
    alpha=0.3
)

ax.set_title('Número de usuarios diarios')
ax.legend();

Modelado con toda la serie

Se crea un forecaster sin tener en cuenta el periodo de confinamiento.

In [25]:
# Crear un objeto ForecasterAutoreg
# ==============================================================================
forecaster = ForecasterAutoreg(
                 regressor = Ridge(),
                 lags      = 21,
             )   
forecaster
Out[25]:
================= 
ForecasterAutoreg 
================= 
Regressor: Ridge() 
Lags: [ 1  2  3  4  5  6  7  8  9 10 11 12 13 14 15 16 17 18 19 20 21] 
Transformer for y: None 
Transformer for exog: None 
Window size: 21 
Weight function included: False 
Differentiation order: None 
Exogenous included: False 
Exogenous variables names: None 
Training range: None 
Training index type: None 
Training index frequency: None 
Regressor parameters: {'alpha': 1.0, 'copy_X': True, 'fit_intercept': True, 'max_iter': None, 'positive': False, 'random_state': None, 'solver': 'auto', 'tol': 0.0001} 
fit_kwargs: {} 
Creation date: 2024-08-06 08:56:43 
Last fit date: None 
Skforecast version: 0.13.0 
Python version: 3.12.4 
Forecaster id: None 

Una vez creado el modelo, se realiza un proceso de backtesting para simular el comportamiento del forecaster si hubiera predicho el conjunto de test en bloques de 7 días.

In [26]:
# Backtesting
# ==============================================================================
metric, predictions_backtest = backtesting_forecaster(
                                   forecaster         = forecaster,
                                   y                  = data.users,
                                   initial_train_size = len(data.loc[:end_train]),
                                   fixed_train_size   = False,
                                   steps              = 7,
                                   metric             = 'mean_absolute_error',
                                   refit              = False,
                                   verbose            = False
                               )

metric
Out[26]:
mean_absolute_error
0 1469.144278

Excluir parte de la serie

Para minimizar la influencia en el modelo de estas fechas, se crea una función que asigna pesos siguiendo las reglas:

  • Peso de 0 si la fecha del índice es:

    • Dentro del periodo de confinamiento (2020-03-16 a 2020-04-21).

    • Dentro del periodo de recuperación (2020-04-21 a 2020-05-31).

    • 21 días después del periodo de recuperación para evitar incluir valores impactados como lags (2020-05-31 a 2020-06-21).

  • Peso de 1 en otro caso.

Si una observación tiene un peso de 0, no tiene influencia durante el entrenamiento del modelo.

In [27]:
# Función de para crear pesos
# ==============================================================================
def custom_weights(index):
    """
    Devuelve 0 si el índice está entre 2020-03-16 y 2020-06-21. Devuelve 1 en caso contrario.
    """
    weights = np.where((index >= '2020-03-16') & (index <= '2020-06-21'), 0, 1)
    
    return weights

Se crea de nuevo un ForecasterAutoreg pero esta vez incluyendo la función custom_weights.

In [28]:
# Crear un objeto ForecasterAutoreg
# ==============================================================================
forecaster = ForecasterAutoreg(
                 regressor   = Ridge(random_state=123),
                 lags        = 21,
                 weight_func = custom_weights
             )

# Backtesting
# ==============================================================================
metric, predictions_backtest = backtesting_forecaster(
                                   forecaster         = forecaster,
                                   y                  = data.users,
                                   initial_train_size = len(data.loc[:end_train]),
                                   fixed_train_size   = False,
                                   steps              = 7,
                                   metric             = 'mean_absolute_error',
                                   refit              = False,
                                   verbose            = False
                               )

metric
Out[28]:
mean_absolute_error
0 1404.679364

Dando un peso de 0 al periodo de confinamiento (excluyéndolo del entrenamiento del modelo) mejora ligeramente el rendimiento del forecasting.

Parada de central eléctrica

Las central eléctricas son instalaciones muy complejas que requieren un alto nivel de mantenimiento. Es común que, cada cierto periodo de funcionamiento, la planta tenga que ser detenida para actividades de reparación, revisión o mantenimiento. Estos eventos provocan que la producción se detenga, generando un hueco en los datos.

Datos

In [29]:
# Descarga de datos
# ==============================================================================
url = ('https://raw.githubusercontent.com/JoaquinAmatRodrigo/skforecast/master/'
       'data/energy_production_shutdown.csv')
data = pd.read_csv(url, sep=',')

# Preprocesado de los datos
# ==============================================================================
data['date'] = pd.to_datetime(data['date'], format='%Y-%m-%d')
data = data.set_index('date')
data = data.asfreq('D')
data = data.sort_index()
data.head()
Out[29]:
production
date
2012-01-01 375.1
2012-01-02 474.5
2012-01-03 573.9
2012-01-04 539.5
2012-01-05 445.4
In [30]:
# División datos entrenamiento, validation y test
# ==============================================================================
data = data.loc['2012-01-01 00:00:00': '2014-12-30 23:00:00']
end_train = '2013-12-31 23:59:00'
data_train = data.loc[: end_train, :]
data_test  = data.loc[end_train:, :]

print(f"Dates train : {data_train.index.min()} --- {data_train.index.max()}  (n={len(data_train)})")
print(f"Dates test  : {data_test.index.min()} --- {data_test.index.max()}  (n={len(data_test)})")
Dates train : 2012-01-01 00:00:00 --- 2013-12-31 00:00:00  (n=731)
Dates test  : 2014-01-01 00:00:00 --- 2014-12-30 00:00:00  (n=364)
In [31]:
# Time series plot
# ==============================================================================
fig, ax = plt.subplots(figsize=(8, 3))
data_train.production.plot(ax=ax, label='train', linewidth=1)
data_test.production.plot(ax=ax, label='test', linewidth=1)
ax.axvspan(
    pd.to_datetime('2012-06-01'),
    pd.to_datetime('2012-09-30'), 
    label="Parada",
    color="red",
    alpha=0.1
)
ax.set_title('Producción de energía')
ax.legend();

Modelado con toda la serie

Se entrena un modelo incluyendo el periodo de parada de la central eléctrica.

In [32]:
# Crear un ForecasterAutoreg
# ==============================================================================
forecaster = ForecasterAutoreg(
                 regressor = Ridge(random_state=123),
                 lags      = 21,
             )
In [33]:
# Backtesting
# ==============================================================================
metric, predictions_backtest = backtesting_forecaster(
                                   forecaster         = forecaster,
                                   y                  = data.production,
                                   initial_train_size = len(data.loc[:end_train]),
                                   fixed_train_size   = False,
                                   steps              = 10,
                                   metric             = 'mean_absolute_error',
                                   refit              = False,
                                   verbose            = False
                               )
metric
Out[33]:
mean_absolute_error
0 28.424722

Excluir parte de la serie

La parada de la fábrica tuvo lugar desde 2012-06-01 hasta 2012-09-30. Para minimizar la influencia en el modelo de estas fechas, se crea una función personalizada que asigna un valor de 0 si la fecha del índice está dentro del periodo de parada o 21 días después (lags utilizados por el modelo) y 1 en otro caso. Si una observación tiene un peso de 0, no tiene influencia durante el entrenamiento del modelo.

In [34]:
# Función de pesos
# ==============================================================================
def custom_weights(index):
    """
    Devuelve 0 si el índice está entre 2012-06-01 y 2012-10-21. Devuelve 1 en caso contrario.
    """
    weights = np.where((index >= '2012-06-01') & (index <= '2012-10-21'), 0, 1)
    return weights
In [35]:
# Crear un ForecasterAutoreg
# ==============================================================================
forecaster = ForecasterAutoreg(
                 regressor   = Ridge(random_state=123),
                 lags        = 21,
                 weight_func = custom_weights
             )

# Backtesting
# ==============================================================================
metric, predictions_backtest = backtesting_forecaster(
                                   forecaster         = forecaster,
                                   y                  = data.production,
                                   initial_train_size = len(data.loc[:end_train]),
                                   fixed_train_size   = False,
                                   steps              = 10,
                                   metric             = 'mean_absolute_error',
                                   refit              = False,
                                   verbose            = False
                               )
metric
Out[35]:
mean_absolute_error
0 27.808331

Como en el ejemplo anterior, excluir las observaciones durante el periodo de parada mejora ligeramente el rendimiento del forecasting.

Información de sesión

In [36]:
import session_info
session_info.show(html=False)
-----
matplotlib          3.9.1
numpy               2.0.1
pandas              2.2.2
session_info        1.0.0
skforecast          0.13.0
sklearn             1.5.1
-----
IPython             8.26.0
jupyter_client      8.6.2
jupyter_core        5.7.2
notebook            6.4.12
-----
Python 3.12.4 | packaged by Anaconda, Inc. | (main, Jun 18 2024, 15:12:24) [GCC 11.2.0]
Linux-5.15.0-1066-aws-x86_64-with-glibc2.31
-----
Session information updated at 2024-08-06 08:56

Instrucciones para citar

¿Cómo citar este documento?

Si utilizas este documento o alguna parte de él, te agradecemos que lo cites. ¡Muchas gracias!

Reducir el impacto del Covid en modelos de forecasting por Joaquín Amat Rodrigo y Javier Escobar Ortiz, disponible bajo una licencia Attribution-NonCommercial-ShareAlike 4.0 International (CC BY-NC-SA 4.0 DEED) en https://www.cienciadedatos.net/documentos/py45-weighted-time-series-forecasting-es.html

¿Cómo citar skforecast?

Si utilizas skforecast en tu investigación o publicación, te lo agradeceríamos mucho que lo cites. ¡Muchas gracias!

Zenodo:

Amat Rodrigo, Joaquin, & Escobar Ortiz, Javier. (2024). skforecast (v0.13.0). Zenodo. https://doi.org/10.5281/zenodo.8382788

APA:

Amat Rodrigo, J., & Escobar Ortiz, J. (2024). skforecast (Version 0.13.0) [Computer software]. https://doi.org/10.5281/zenodo.8382788

BibTeX:

@software{skforecast, author = {Amat Rodrigo, Joaquin and Escobar Ortiz, Javier}, title = {skforecast}, version = {0.13.0}, month = {8}, year = {2024}, license = {BSD-3-Clause}, url = {https://skforecast.org/}, doi = {10.5281/zenodo.8382788} }


¿Te ha gustado el artículo? Tu ayuda es importante

Mantener un sitio web tiene unos costes elevados, tu contribución me ayudará a seguir generando contenido divulgativo gratuito. ¡Muchísimas gracias! 😊


Creative Commons Licence
Este documento creado por Joaquín Amat Rodrigo y Javier Escobar Ortiz tiene licencia Attribution-NonCommercial-ShareAlike 4.0 International.