Modelos de forecasting globales

En los modelos de series individuales (Local Forecasting Model), una única serie temporal se modela como una combinación lineal o no lineal de sus valores pasados y, opcionalmente, variables exógenas. Aunque este método proporciona una comprensión exhaustiva de cada serie, su escalabilidad puede verse comprometida cuando se dispone de un gran número de series.

Los modelos globales multiserie (Global Forecasting Model) consisten en un único modelo que tiene en cuenta varias series temporales de forma simultánea. Intentan captar los patrones básicos que rigen las series, mitigando así el ruido potencial que pueda introducir cada serie. Este enfoque es eficiente desde el punto de vista computacional, es fácil de mantener y puede conseguir mayor generalización, aunque potencialmente a costa de sacrificar cierta resolución individual. Se pueden distinguir dos estrategias de modelos de previsión global: Múltiples series temporales independientes y Múltiples series temporales dependientes.

Ventajas de los modelos multiseries

• Es más fácil mantener y controlar un solo modelo que varios.

• Dado que todas las series temporales se combinan durante el entrenamiento, cuando las series sean cortas (pocos datos) el modelo tendrá una mayor capacidad de aprendizaje al disponer de más observaciones.

• Al combinar múltiples series temporales, el modelo puede aprender patrones más generalizables.

Desventajas de los modelos multiseries

• Si las series no siguen la misma dinámica interna, el modelo puede aprender un patrón que no represente a ninguna de ellas.

• Las series pueden enmascararse unas a otras, por lo que el modelo puede no predecirlas todas con el mismo rendimiento.

• Es más exigente desde el punto de vista computacional (tiempo y recursos) entrenar y realizar backtesting de un modelo grande que de varios pequeños.

Múltiples series temporales independientes

En esta situación, cada serie temporal es independiente de las demás o, dicho de otro modo, los valores pasados de una serie no se utilizan como predictores de las otras series. ¿Por qué es útil entonces modelar todo junto? Aunque las series no dependen unas de otras, pueden seguir el mismo patrón intrínseco en cuanto a sus valores pasados y futuros. Por ejemplo, en una misma tienda, las ventas de los productos A y B pueden no estar relacionadas, pero siguen la misma dinámica, la de la tienda.

Para predecir los siguientes n steps, se sigue una estrategia recurisva, recursive multi-step forecasting.

Múltiples series temporales dependientes

Todas las series se modelan teniendo en cuenta que cada serie temporal depende no sólo de sus valores pasados, sino también de los valores pasados de las demás series. Se espera que el modelo no sólo aprenda la información de cada serie por separado, sino que también las relacione. Por ejemplo, las mediciones realizadas por todos los sensores (caudal, temperatura, presión...) instalados en una máquina industrial como un compresor. Series temporales multivariantes user guide.

Series con la misma longitud

En este primer ejemplo, se utilizan múltiples series temporales de la misma longitud. El objetivo es comparar los resultados de un modelo global con los de un modelo individual para cada serie al predecir los próximos 7 días de ventas para 50 productos diferentes en una tienda, utilizando los 5 años de historial disponible. Los datos se han obtenido del Store Item Demand Forecasting Challenge. Este conjunto de datos contiene 913,000 transacciones de ventas desde el 01/01/2013 hasta el 31/12/2017 para 50 productos (SKU) en 10 tiendas.

# Manipulación de datos
# ==============================================================================
import numpy as np
import pandas as pd

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

# Modelado y Forecasting
# ==============================================================================
import sklearn
import skforecast
from sklearn.ensemble import HistGradientBoostingRegressor
from sklearn.preprocessing import StandardScaler
from skforecast.recursive import ForecasterRecursive, ForecasterRecursiveMultiSeries
from skforecast.model_selection import (
    TimeSeriesFold,
    OneStepAheadFold,
    backtesting_forecaster,
    bayesian_search_forecaster,
    backtesting_forecaster_multiseries,
    bayesian_search_forecaster_multiseries
)
from skforecast.preprocessing import RollingFeatures, series_long_to_dict, exog_long_to_dict
from skforecast.exceptions import OneStepAheadValidationWarning

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

color = '\033[1m\033[38;5;208m'
print(f"{color}Versión skforecast: {skforecast.__version__}")
print(f"{color}Versión scikit-learn: {sklearn.__version__}")
print(f"{color}Versión pandas: {pd.__version__}")
print(f"{color}Versión numpy: {np.__version__}")

Versión skforecast: 0.15.1
Versión scikit-learn: 1.5.2
Versión pandas: 2.2.3
Versión numpy: 2.0.2

Datos

Hide

# Lectura de datos
# ======================================================================================
data = pd.read_csv('./train_stores_kaggle.csv')
display(data)
print(f"Shape: {data.shape}")

	date	store	item	sales
0	2013-01-01	1	1	13
1	2013-01-02	1	1	11
2	2013-01-03	1	1	14
3	2013-01-04	1	1	13
4	2013-01-05	1	1	10
...	...	...	...	...
912995	2017-12-27	10	50	63
912996	2017-12-28	10	50	59
912997	2017-12-29	10	50	74
912998	2017-12-30	10	50	62
912999	2017-12-31	10	50	82

913000 rows × 4 columns

Shape: (913000, 4)

Hide

# Preparación datos
# ======================================================================================
selected_store = 2 # Seleccionar una tienda
selected_items = data.item.unique()
data = data[(data['store'] == selected_store) & (data['item'].isin(selected_items))].copy()
data['date'] = pd.to_datetime(data['date'], format='%Y-%m-%d')
data = pd.pivot_table(
           data    = data,
           values  = 'sales',
           index   = 'date',
           columns = 'item'
       )
data.columns.name = None
data.columns = [f"item_{col}" for col in data.columns]
data = data.asfreq('1D')
data = data.sort_index()
data.head(4)

	item_1	item_2	item_3	item_4	item_5	item_6	item_7	item_8	item_9	item_10	...	item_41	item_42	item_43	item_44	item_45	item_46	item_47	item_48	item_49	item_50
date
2013-01-01	12.0	41.0	19.0	21.0	4.0	34.0	39.0	49.0	28.0	51.0	...	11.0	25.0	36.0	12.0	45.0	43.0	12.0	45.0	29.0	43.0
2013-01-02	16.0	33.0	32.0	14.0	6.0	40.0	47.0	42.0	21.0	56.0	...	19.0	21.0	35.0	25.0	50.0	52.0	13.0	37.0	25.0	57.0
2013-01-03	16.0	46.0	26.0	12.0	12.0	41.0	43.0	46.0	29.0	46.0	...	23.0	20.0	52.0	18.0	56.0	30.0	5.0	45.0	30.0	45.0
2013-01-04	20.0	50.0	34.0	17.0	16.0	41.0	44.0	55.0	32.0	56.0	...	15.0	28.0	50.0	24.0	57.0	46.0	19.0	32.0	20.0	45.0

4 rows × 50 columns

El dataset se divide en 3 particiones: una para el entrenamiento, otra para la validación y otra para test.

Hide

# Separación datos train-validation-test
# ======================================================================================
end_train = '2016-05-31 23:59:00'
end_val = '2017-05-31 23:59:00'
data_train = data.loc[:end_train, :].copy()
data_val   = data.loc[end_train:end_val, :].copy()
data_test  = data.loc[end_val:, :].copy()
print(f"Fechas train      : {data_train.index.min()} --- {data_train.index.max()}  (n={len(data_train)})")
print(f"Fechas validación : {data_val.index.min()} --- {data_val.index.max()}  (n={len(data_val)})")
print(f"Fechas test       : {data_test.index.min()} --- {data_test.index.max()}  (n={len(data_test)})")

Fechas train      : 2013-01-01 00:00:00 --- 2016-05-31 00:00:00  (n=1247)
Fechas validación : 2016-06-01 00:00:00 --- 2017-05-31 00:00:00  (n=365)
Fechas test       : 2017-06-01 00:00:00 --- 2017-12-31 00:00:00  (n=214)

Se representan cuatro de las series para comprender sus tendencias y patrones. Se recomienda al lector que grafique más de ellas para comprenderlas en profundidad.

Hide

# Gráfico series temporales
# ======================================================================================
set_dark_theme()
fig, axs = plt.subplots(4, 1, figsize=(7, 5), sharex=True)
data.iloc[:, :4].plot(
    legend   = True,
    subplots = True, 
    title    = 'Ventas de la tienda 2',
    ax       = axs, 
)
for ax in axs:
    ax.axvline(pd.to_datetime(end_train) , color='white', linestyle='--', linewidth=1.5)
    ax.axvline(pd.to_datetime(end_val) , color='white', linestyle='--', linewidth=1.5)
fig.tight_layout();

# Entrenar y realizar backtesting de un modelo para cada item
# ======================================================================================
items = []
mae_values = []
predictions = {}

for i, item in enumerate(tqdm(data.columns)):
    # Definir el forecaster
    window_features = RollingFeatures(stats=['mean', 'min', 'max'], window_sizes=7)
    forecaster = ForecasterRecursive(
                     regressor       = HistGradientBoostingRegressor(random_state=8523),
                     lags            = 14,
                     window_features = window_features
                 )
    # Backtesting forecaster
    cv = TimeSeriesFold(
            steps              = 7,
            initial_train_size = len(data_train) + len(data_val),
            refit              = False,
         )
    metric, preds = backtesting_forecaster(
                        forecaster    = forecaster,
                        y             = data[item],
                        cv            = cv,
                        metric        = 'mean_absolute_error',
                        show_progress = False
                    )
    items.append(item)
    mae_values.append(metric.at[0, 'mean_absolute_error'])
    predictions[item] = preds

# Resultados
uni_series_mae = pd.Series(
                     data  = mae_values,
                     index = items,
                     name  = 'uni_series_mae'
                 )
uni_series_mae.head()

  0%|          | 0/50 [00:00<?, ?it/s]

item_1    6.004406
item_2    9.994352
item_3    8.652751
item_4    5.528955
item_5    5.096925
Name: uni_series_mae, dtype: float64

# Entrenamiento y backtesting con un único modelo para todos los items
# ======================================================================================
items = list(data.columns)

# Definir el forecaster
window_features = RollingFeatures(stats=['mean', 'min', 'max'], window_sizes=7)
forecaster_ms = ForecasterRecursiveMultiSeries(
                    regressor          = HistGradientBoostingRegressor(random_state=8523),
                    lags               = 14,
                    encoding           = 'ordinal',
                    transformer_series = StandardScaler(),
                    window_features    = window_features,
                )
# Backtesting forecaster para todos los items
cv = TimeSeriesFold(
            steps              = 7,
            initial_train_size = len(data_train) + len(data_val),
            refit              = False,
         )
multi_series_mae, predictions_ms = backtesting_forecaster_multiseries(
                                       forecaster = forecaster_ms,
                                       series     = data,
                                       levels     = items,
                                       cv         = cv,
                                       metric     = 'mean_absolute_error',
                                   )   
# Resultados
display(multi_series_mae.head(3))
print('')
display(predictions_ms.head(3))

  0%|          | 0/31 [00:00<?, ?it/s]

# Diferencia de la métrica de backtesting para cada item
# ======================================================================================
multi_series_mae = multi_series_mae.set_index('levels')
multi_series_mae.columns = ['multi_series_mae']
results = pd.concat((uni_series_mae, multi_series_mae), axis = 1)
results['improvement'] = results.eval('uni_series_mae - multi_series_mae')
results['improvement_(%)'] = 100 * results.eval('(uni_series_mae - multi_series_mae) / uni_series_mae')
results = results.round(2)
results.style.bar(subset=['improvement_(%)'], align='mid', color=['#d65f5f', '#5fba7d'])

# Mejora media de todos los items
# ======================================================================================
results[['improvement', 'improvement_(%)']].agg(['mean', 'min', 'max'])

# Número de series con mejora positiva y negativa
# ======================================================================================
pd.Series(np.where(results['improvement_(%)'] < 0, 'negative', 'positive')).value_counts()

positive    53
Name: count, dtype: int64

Series de diferente longitud y distintas variables exógenas

En escenarios en los que se tienen que modelar múltiples series, es habitual que las series tengan longitudes distintas debido a diferencias en los momentos de inicio del registro de los datos. Para hacer frente a este escenario, la clase ForecasterRecursiveMultiSeries permite modelizar simultáneamente series temporales de distintas longitudes y con distintas variables exógenas.

• Cuando las series modeladas tienen longitudes diferentes, deben almacenarse en un diccionario Python. Las claves del diccionario son los nombres de las series y los valores son las propias series. Todas las series deben ser del tipo pandas.Series, tener un índice datetime y tener la misma frecuencia.
  
  

• Cuando se utilizan variables exógenas diferentes para cada serie o cuando las variables exógenas son las mismas pero tienen valores diferentes para cada serie, deben almacenarse en un diccionario. Las claves del diccionario son los nombres de las series y los valores son las propias variables exógenas. Todas las variables exógenas deben ser de tipo pandas.DataFrame o pandas.Series.

Datos

Los datos de este ejemplo están almacenados en "formato largo" en un único DataFrame. La columna series_id identifica la serie a la que pertenece cada observación. La columna timestamp contiene la fecha de la observación, y la columna value contiene el valor de la serie en esa fecha. Cada serie temporal tiene una longitud diferente.

Las variables exógenas se almacenan en un DataFrame separado, también en "formato largo". La columna series_id identifica la serie a la que pertenece cada observación. La columna timestamp contiene la fecha de la observación, y las columnas restantes contienen los valores de las variables exógenas en esa fecha.

Hide

# Lectura series con diferentes longitudes y variables exógenas
# ==============================================================================
series = pd.read_csv(
    'https://raw.githubusercontent.com/skforecast/skforecast-datasets/main/data/demo_multi_series.csv'
)
exog = pd.read_csv(
    'https://raw.githubusercontent.com/skforecast/skforecast-datasets/main/data/demo_multi_series_exog.csv'
)
series['timestamp'] = pd.to_datetime(series['timestamp'])
exog['timestamp'] = pd.to_datetime(exog['timestamp'])
display(series.head())
display(exog.head())

	series_id	timestamp	value
0	id_1000	2016-01-01	1012.500694
1	id_1000	2016-01-02	1158.500099
2	id_1000	2016-01-03	983.000099
3	id_1000	2016-01-04	1675.750496
4	id_1000	2016-01-05	1586.250694

	series_id	timestamp	sin_day_of_week	cos_day_of_week	air_temperature	wind_speed
0	id_1000	2016-01-01	-0.433884	-0.900969	6.416639	4.040115
1	id_1000	2016-01-02	-0.974928	-0.222521	6.366474	4.530395
2	id_1000	2016-01-03	-0.781831	0.623490	6.555272	3.273064
3	id_1000	2016-01-04	0.000000	1.000000	6.704778	4.865404
4	id_1000	2016-01-05	0.781831	0.623490	2.392998	5.228913

Cuando las series tienen longitudes diferentes, los datos deben transformarse en un diccionario. Las claves del diccionario son los nombres de las series y los valores son las propias series. Para ello, se utiliza la función series_long_to_dict, que toma el DataFrame en "formato largo" y devuelve un diccionario de series.

Del mismo modo, cuando las variables exógenas son diferentes (valores o variables) para cada serie, los datos deben transformarse en un diccionario. Las claves del diccionario son los nombres de las series y los valores son las propias variables exógenas. Se utiliza la función exog_long_to_dict, que toma el DataFrame en "formato largo" y devuelve un diccionario de variables exógenas.

Hide

# Transform series and exog to dictionaries
# ==============================================================================
series_dict = series_long_to_dict(
                data      = series,
                series_id = 'series_id',
                index     = 'timestamp',
                values    = 'value',
                freq      = 'D'
              )

exog_dict = exog_long_to_dict(
                data      = exog,
                series_id = 'series_id',
                index     = 'timestamp',
                freq      = 'D'
            )

╭──────────────────────────────── MissingValuesWarning ────────────────────────────────╮
│ Series 'id_1003' is incomplete. NaNs have been introduced after setting the          │
│ frequency.                                                                           │
│                                                                                      │
│ Category : MissingValuesWarning                                                      │
│ Location :                                                                           │
│ /home/ubuntu/anaconda3/envs/skforecast_15_py12/lib/python3.12/site-packages/skforeca │
│ st/preprocessing/preprocessing.py:428                                                │
│ Suppress : warnings.simplefilter('ignore', category=MissingValuesWarning)            │
╰──────────────────────────────────────────────────────────────────────────────────────╯

Algunas variables exógenas se omiten en las series 1 y 3 para ilustrar que se pueden utilizar diferentes variables exógenas para cada serie.

Hide

exog_dict['id_1000'] = exog_dict['id_1000'].drop(columns=['air_temperature', 'wind_speed'])
exog_dict['id_1003'] = exog_dict['id_1003'].drop(columns=['cos_day_of_week'])

Hide

# Particiones de entrenamiento y test
# ==============================================================================
end_train = '2016-07-31 23:59:00'
series_dict_train = {k: v.loc[: end_train,] for k, v in series_dict.items()}
exog_dict_train   = {k: v.loc[: end_train,] for k, v in exog_dict.items()}
series_dict_test  = {k: v.loc[end_train:,] for k, v in series_dict.items()}
exog_dict_test    = {k: v.loc[end_train:,] for k, v in exog_dict.items()}

Hide

# Gráfico series
# ==============================================================================
colors = plt.rcParams['axes.prop_cycle'].by_key()['color']
fig, axs = plt.subplots(5, 1, figsize=(8, 4), sharex=True)
for i, s in enumerate(series_dict.values()):
    axs[i].plot(s, label=s.name, color=colors[i])
    axs[i].legend(loc='upper right', fontsize=8)
    axs[i].tick_params(axis='both', labelsize=8)
    axs[i].axvline(pd.to_datetime(end_train), color='white', linestyle='--', linewidth=1)

Hide

# Descripción de cada serie
# ==============================================================================
for k in series_dict.keys():
    print(f"{k}:")
    try:
        print(
            f"\tTrain: len={len(series_dict_train[k])}, {series_dict_train[k].index[0]}"
            f" --- {series_dict_train[k].index[-1]} "
            f" (missing={series_dict_train[k].isnull().sum()})"
        )
    except:
        print(f"\tTrain: len=0")
    try:
        print(
            f"\tTest : len={len(series_dict_test[k])}, {series_dict_test[k].index[0]}"
            f" --- {series_dict_test[k].index[-1]} "
            f" (missing={series_dict_test[k].isnull().sum()})"
        )
    except:
        print(f"\tTest : len=0")

id_1000:
	Train: len=213, 2016-01-01 00:00:00 --- 2016-07-31 00:00:00  (missing=0)
	Test : len=153, 2016-08-01 00:00:00 --- 2016-12-31 00:00:00  (missing=0)
id_1001:
	Train: len=30, 2016-07-02 00:00:00 --- 2016-07-31 00:00:00  (missing=0)
	Test : len=153, 2016-08-01 00:00:00 --- 2016-12-31 00:00:00  (missing=0)
id_1002:
	Train: len=183, 2016-01-01 00:00:00 --- 2016-07-01 00:00:00  (missing=0)
	Test : len=0
id_1003:
	Train: len=213, 2016-01-01 00:00:00 --- 2016-07-31 00:00:00  (missing=73)
	Test : len=153, 2016-08-01 00:00:00 --- 2016-12-31 00:00:00  (missing=73)
id_1004:
	Train: len=91, 2016-05-02 00:00:00 --- 2016-07-31 00:00:00  (missing=0)
	Test : len=31, 2016-08-01 00:00:00 --- 2016-08-31 00:00:00  (missing=0)

Hide

# Variables exógenas de cada serie
# ==============================================================================
for k in series_dict.keys():
    print(f"{k}:")
    try:
        print(f"\t{exog_dict[k].columns.to_list()}")
    except:
        print(f"\tNo variables exógenas")

id_1000:
	['sin_day_of_week', 'cos_day_of_week']
id_1001:
	['sin_day_of_week', 'cos_day_of_week', 'air_temperature', 'wind_speed']
id_1002:
	['sin_day_of_week', 'cos_day_of_week', 'air_temperature', 'wind_speed']
id_1003:
	['sin_day_of_week', 'air_temperature', 'wind_speed']
id_1004:
	['sin_day_of_week', 'cos_day_of_week', 'air_temperature', 'wind_speed']

# Fit forecaster
# ==============================================================================
regressor = HistGradientBoostingRegressor(random_state=123, max_depth=5)
window_features = RollingFeatures(stats=['mean', 'min', 'max'], window_sizes=7)
forecaster = ForecasterRecursiveMultiSeries(
                regressor          = regressor,
                lags               = 14,
                window_features    = window_features,
                encoding           = "ordinal",
                dropna_from_series = False
            )

forecaster.fit(series=series_dict_train, exog=exog_dict_train, suppress_warnings=True)
forecaster

# Predict
# ==============================================================================
predicciones = forecaster.predict(steps=5, exog=exog_dict_test, suppress_warnings=True)
predicciones

Backtesting

Cuando las series tienen longitudes diferentes, el proceso de backtesting sólo devuelve predicciones para las fechas-horas que están presentes en las series.

Hide

# Backtesting
# ==============================================================================
cv = TimeSeriesFold(
        steps              = 24,
        initial_train_size = len(series_dict_train["id_1000"]),
        refit              = False,
     )

metrics_levels, backtest_predictions = backtesting_forecaster_multiseries(
    forecaster            = forecaster,
    series                = series_dict,
    exog                  = exog_dict,
    cv                    = cv,
    metric                = "mean_absolute_error",
    add_aggregated_metric = False,
    verbose               = True,
    suppress_warnings     = True
)

display(metrics_levels)
display(backtest_predictions)

Information of folds
--------------------
Number of observations used for initial training: 213
Number of observations used for backtesting: 153
    Number of folds: 7
    Number skipped folds: 0 
    Number of steps per fold: 24
    Number of steps to exclude between last observed data (last window) and predictions (gap): 0
    Last fold only includes 9 observations.

Fold: 0
    Training:   2016-01-01 00:00:00 -- 2016-07-31 00:00:00  (n=213)
    Validation: 2016-08-01 00:00:00 -- 2016-08-24 00:00:00  (n=24)
Fold: 1
    Training:   No training in this fold
    Validation: 2016-08-25 00:00:00 -- 2016-09-17 00:00:00  (n=24)
Fold: 2
    Training:   No training in this fold
    Validation: 2016-09-18 00:00:00 -- 2016-10-11 00:00:00  (n=24)
Fold: 3
    Training:   No training in this fold
    Validation: 2016-10-12 00:00:00 -- 2016-11-04 00:00:00  (n=24)
Fold: 4
    Training:   No training in this fold
    Validation: 2016-11-05 00:00:00 -- 2016-11-28 00:00:00  (n=24)
Fold: 5
    Training:   No training in this fold
    Validation: 2016-11-29 00:00:00 -- 2016-12-22 00:00:00  (n=24)
Fold: 6
    Training:   No training in this fold
    Validation: 2016-12-23 00:00:00 -- 2016-12-31 00:00:00  (n=9)

  0%|          | 0/7 [00:00<?, ?it/s]

	levels	mean_absolute_error
0	id_1000	164.959423
1	id_1001	1055.559754
2	id_1002	NaN
3	id_1003	235.663130
4	id_1004	968.459237

	level	pred
2016-08-01	id_1000	1433.494674
2016-08-01	id_1001	3068.244797
2016-08-01	id_1003	2748.768695
2016-08-01	id_1004	7763.964965
2016-08-02	id_1000	1465.937652
...	...	...
2016-12-30	id_1001	1114.592910
2016-12-30	id_1003	1965.060657
2016-12-31	id_1000	1459.122750
2016-12-31	id_1001	1001.655092
2016-12-31	id_1003	1969.768680

507 rows × 2 columns

Hide

# Gráfico predicciones backtesting
# ==============================================================================
colors = plt.rcParams["axes.prop_cycle"].by_key()["color"]
fig, axs = plt.subplots(5, 1, figsize=(8, 4), sharex=True)
for i, s in enumerate(series_dict.keys()):
    axs[i].plot(series_dict[s], label=series_dict[s].name, color=colors[i])
    axs[i].axvline(pd.to_datetime(end_train), color="white", linestyle="--", linewidth=1)
    try:
        axs[i].plot(
            backtest_predictions.loc[backtest_predictions["level"] == s, "pred"],
            label="prediction",
            color="white",
        )
    except:
        pass
    axs[i].legend(loc="upper left", fontsize=8)
    axs[i].tick_params(axis="both", labelsize=8)

Al permitir la modelización de series temporales de diferentes longitudes y con diferentes variables exógenas, la clase ForecasterRecursiveMultiSeries proporciona una herramienta flexible y potente para utilizar toda la información disponible para entrenar los modelos.

Optimización de hiperparámetros y lags

En el primer ejemplo de este documento, la comparación entre forecaster se ha realizado sin optimizar los hiperparámetros de los regresores. Para una comparación justa, se utiliza una estrategia de grid search con el fin de seleccionar la mejor configuración para cada forecaster. Véase más información en hyperparameter tuning and lags selection.

# Búsqueda de hiperparámetros y backtesting de un modelo para cada item
# ======================================================================================
warnings.simplefilter('ignore', category=OneStepAheadValidationWarning)
items = []
mae_values  = []

def search_space(trial):
    search_space  = {
        'lags'          : trial.suggest_categorical('lags', [7, 14]),
        'max_iter'      : trial.suggest_int('max_iter', 100, 500),
        'max_depth'     : trial.suggest_int('max_depth', 5, 10),
        'learning_rate' : trial.suggest_float('learning_rate', 0.01, 0.1)
    } 

    return search_space

for item in tqdm(data.columns):

    window_features = RollingFeatures(stats=['mean', 'min', 'max'], window_sizes=7)
    forecaster = ForecasterRecursive(
                    regressor       = HistGradientBoostingRegressor(random_state=123),
                    lags            = 14,
                    window_features = window_features
                )
    
    cv_search = OneStepAheadFold(initial_train_size = len(data_train))
    
    warnings.simplefilter('ignore', category=OneStepAheadValidationWarning)
    results_search, _ = bayesian_search_forecaster(
                            forecaster    = forecaster,
                            y             = data.loc[:end_val, item],
                            cv            = cv_search,
                            search_space  = search_space,
                            n_trials      = 10,
                            metric        = 'mean_absolute_error',
                            return_best   = False,
                            show_progress = False 
                        )
    best_params = results_search.at[0, 'params']
    best_lags   = results_search.at[0, 'lags']
    forecaster.set_params(best_params)
    forecaster.set_lags(best_lags)
    
    cv_backtesting = TimeSeriesFold(
                        steps              = 7,
                        initial_train_size = len(data_train) + len(data_val),
                        refit              = False,
                     )
    metric, preds = backtesting_forecaster(
                        forecaster    = forecaster,
                        y             = data[item],
                        cv            = cv_backtesting,
                        metric        = 'mean_absolute_error',
                        show_progress = False
                    )

    items.append(item)
    mae_values.append(metric.at[0, 'mean_absolute_error'])

uni_series_mae = pd.Series(
                    data  = mae_values,
                    index = items,
                    name  = 'uni_series_mae'
                )

  0%|          | 0/50 [00:00<?, ?it/s]

# Busqueda de hiperparámetros y backtesting de un modelo global para todos los items
# ======================================================================================
def search_space(trial):
    search_space  = {
        'lags'          : trial.suggest_categorical('lags', [7, 14]),
        'max_iter'      : trial.suggest_int('max_iter', 100, 500),
        'max_depth'     : trial.suggest_int('max_depth', 5, 10),
        'learning_rate' : trial.suggest_float('learning_rate', 0.01, 0.1)
    } 

    return search_space

forecaster_ms = ForecasterRecursiveMultiSeries(
                    regressor          = HistGradientBoostingRegressor(random_state=123),
                    lags               = 14,
                    transformer_series = StandardScaler(),
                    encoding           = 'ordinal'
                )

results_bayesian_ms = bayesian_search_forecaster_multiseries(
                        forecaster    = forecaster_ms,
                        series        = data.loc[:end_val, :],
                        levels        = None, # Si es None se seleccionan todos los niveles
                        cv            = cv_search,
                        search_space  = search_space,
                        n_trials      = 20,
                        metric        = 'mean_absolute_error',
                    )      

multi_series_mae, predictions_ms = backtesting_forecaster_multiseries(
                                       forecaster            = forecaster_ms,
                                       series                = data,
                                       levels                = None, # Si es None se seleccionan todos los niveles
                                       cv                    = cv_backtesting,
                                       metric                = 'mean_absolute_error',
                                       add_aggregated_metric = False,
                                   )

╭─────────────────────────── OneStepAheadValidationWarning ────────────────────────────╮
│ One-step-ahead predictions are used for faster model comparison, but they may not    │
│ fully represent multi-step prediction performance. It is recommended to backtest the │
│ final model for a more accurate multi-step performance estimate.                     │
│                                                                                      │
│ Category : OneStepAheadValidationWarning                                             │
│ Location :                                                                           │
│ /home/ubuntu/anaconda3/envs/skforecast_15_py12/lib/python3.12/site-packages/skforeca │
│ st/model_selection/_utils.py:675                                                     │
│ Suppress : warnings.simplefilter('ignore', category=OneStepAheadValidationWarning)   │
╰──────────────────────────────────────────────────────────────────────────────────────╯

  0%|          | 0/20 [00:00<?, ?it/s]

`Forecaster` refitted using the best-found lags and parameters, and the whole data set: 
  Lags: [ 1  2  3  4  5  6  7  8  9 10 11 12 13 14] 
  Parameters: {'max_iter': 352, 'max_depth': 9, 'learning_rate': 0.034448731526687745}
  Backtesting metric: 8.015025437227383
  Levels: ['item_1', 'item_2', 'item_3', 'item_4', 'item_5', 'item_6', 'item_7', 'item_8', 'item_9', 'item_10', '...', 'item_41', 'item_42', 'item_43', 'item_44', 'item_45', 'item_46', 'item_47', 'item_48', 'item_49', 'item_50']

  0%|          | 0/31 [00:00<?, ?it/s]

# Diferencia de la métrica de backtesting para cada item
# ======================================================================================
multi_series_mae = multi_series_mae.set_index('levels')
multi_series_mae.columns = ['multi_series_mae']
results = pd.concat((uni_series_mae, multi_series_mae), axis = 1)
results['improvement'] = results.eval('uni_series_mae - multi_series_mae')
results['improvement_(%)'] = 100 * results.eval('(uni_series_mae - multi_series_mae) / uni_series_mae')
results = results.round(2)

# Mejora media de todos los items
# ======================================================================================
results[['improvement', 'improvement_(%)']].agg(['mean', 'min', 'max'])

# Número de series con mejora positiva y negativa
# ======================================================================================
pd.Series(np.where(results['improvement_(%)'] < 0, 'negative', 'positive')).value_counts()

positive    50
Name: count, dtype: int64

Tras identificar la combinación de lags e hiperparámetros que logran el mejor rendimiento predictivo para cada forecaster, un número superior de modelos univariantes han logrado una mayor capacidad predictiva. Aun así, el modelo multiserie proporciona mejores resultados para la mayoría de los items.

Selección de predictores

La selección de predictores es el proceso de seleccionar un subconjunto de predictores relevantes (variables) para su uso en la construcción del modelo. Las técnicas de selección de predictores se utilizan por varias razones: para simplificar los modelos y hacerlos más fáciles de interpretar, para reducir el tiempo de entrenamiento, para evitar los problemas de dimensionalidad, para mejorar la generalización reduciendo el sobreajuste (formalmente, la reducción de la varianza), entre otros.

Skforecast es compatible con los métodos de selección implementados en scikit-learn. Existen varios métodos de selección de características, pero los más comunes son:

• Recursive feature elimination (RFE)

• Sequential Feature Selection (SFS)

• Feature selection based on threshold (SelectFromModel)

Pesos en forecasting multiseries

Los pesos se utilizan para controlar la influencia que tiene cada observación en el entrenamiento del modelo. ForecasterRecursiveMultiseries acepta dos tipos de pesos:

• series_weights controla la importancia relativa de cada serie. Si una serie tiene el doble de peso que las demás, las observaciones de esa serie influyen el doble en el entrenamiento. Cuanto mayor sea el peso de una serie en relación con las demás, más se centrará el modelo en intentar aprender esa serie.

• weight_func controla la importancia relativa de cada observación en función del índice. Por ejemplo, una función que asigna un peso menor a ciertas fechas.

Si se indican los dos tipos de pesos, estos se multiplican para crear los pesos finales como se muestra en la figura. El sample_weight resultante no puede contener valores negativos.

Más información sobre weights in multi-series forecasting y weighted time series forecasting con skforecast.

En este ejemplo, item_1 tiene una mayor importancia relativa entre series (pesa 3 veces más que el resto de series), y las observaciones entre '2013-12-01' y '2014-01-31' se consideran no representativas y se les aplica un peso de 0.

# Pesos en forecasting multiseries
# ======================================================================================
# Pesos de cada serie
series_weights = {'item_1': 3.0} # Las series que no aparezcan en el dict tienen un peso de 1

# Pesos de cada observación (por índice)
def custom_weights(index):
    """
    Devuelve 0 si el índice está entre '2013-12-01' y '2014-01-31', 1 en caso contrario.
    """
    weights = np.where(
                  (index >= '2013-12-01') & (index <= '2014-01-31'),
                   0,
                   1
              )
    
    return weights

forecaster = ForecasterRecursiveMultiSeries(
                 regressor          = HistGradientBoostingRegressor(random_state=123),
                 lags               = 14,
                 transformer_series = StandardScaler(),
                 transformer_exog   = None,
                 weight_func        = custom_weights, 
                 series_weights     = series_weights
             )

forecaster.fit(series=data)
forecaster.predict(steps=7).head(3)

╭─────────────────────────────── IgnoredArgumentWarning ───────────────────────────────╮
│ {'item_31', 'item_19', 'item_14', 'item_3', 'item_7', 'item_45', 'item_43',          │
│ 'item_33', 'item_9', 'item_38', 'item_20', 'item_10', 'item_49', 'item_37',          │
│ 'item_47', 'item_35', 'item_2', 'item_23', 'item_24', 'item_6', 'item_44',           │
│ 'item_13', 'item_28', 'item_26', 'item_29', 'item_8', 'item_36', 'item_50',          │
│ 'item_21', 'item_11', 'item_42', 'item_15', 'item_25', 'item_17', 'item_30',         │
│ 'item_16', 'item_40', 'item_4', 'item_41', 'item_27', 'item_5', 'item_22',           │
│ 'item_46', 'item_48', 'item_12', 'item_18', 'item_39', 'item_32', 'item_34'} not     │
│ present in `series_weights`. A weight of 1 is given to all their samples.            │
│                                                                                      │
│ Category : IgnoredArgumentWarning                                                    │
│ Location :                                                                           │
│ /home/ubuntu/anaconda3/envs/skforecast_15_py12/lib/python3.12/site-packages/skforeca │
│ st/recursive/_forecaster_recursive_multiseries.py:1542                               │
│ Suppress : warnings.simplefilter('ignore', category=IgnoredArgumentWarning)          │
╰──────────────────────────────────────────────────────────────────────────────────────╯

Conclusiones

Este caso de uso muestra como un modelo multiserie puede presentar ventajas sobre varios modelos individuales cuando se predicen series temporales con una dinámica similar. Más allá de las posibles mejoras en la predicción, también es importante tener en cuenta la ventaja de tener un solo modelo que mantener.

Información de sesión

import session_info
session_info.show(html=False)

-----
matplotlib          3.10.1
numpy               2.0.2
optuna              4.2.1
pandas              2.2.3
session_info        1.0.0
skforecast          0.15.1
sklearn             1.5.2
tqdm                4.67.1
-----
IPython             9.0.2
jupyter_client      8.6.3
jupyter_core        5.7.2
notebook            6.4.12
-----
Python 3.12.9 | packaged by Anaconda, Inc. | (main, Feb  6 2025, 18:56:27) [GCC 11.2.0]
Linux-5.15.0-1077-aws-x86_64-with-glibc2.31
-----
Session information updated at 2025-03-26 15:05

Instrucciones para citar

¿Cómo citar este documento?

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

Modelos de forecasting globales: modelado de múltiples series temporales con machine learning por Joaquín Amat Rodrigo and Javier Escobar Ortiz, disponible con licencia Attribution-NonCommercial-ShareAlike 4.0 International (CC BY-NC-SA 4.0 DEED) en https://www.cienciadedatos.net/documentos/py44-multi-series-forecasting-skforecast-español.html

¿Cómo citar skforecast?

Si utilizas skforecast, te agradeceríamos mucho que lo cites. ¡Muchas gracias!

Zenodo:

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

APA:

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

BibTeX:

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

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

Tu contribución me ayudará a seguir generando contenido divulgativo gratuito. ¡Muchísimas gracias! 😊

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

Se permite:

• Compartir: copiar y redistribuir el material en cualquier medio o formato.

• Adaptar: remezclar, transformar y crear a partir del material.

Bajo los siguientes términos:

• Atribución: Debes otorgar el crédito adecuado, proporcionar un enlace a la licencia e indicar si se realizaron cambios. Puedes hacerlo de cualquier manera razonable, pero no de una forma que sugiera que el licenciante te respalda o respalda tu uso.

• No-Comercial: No puedes utilizar el material para fines comerciales.

• Compartir-Igual: Si remezclas, transformas o creas a partir del material, debes distribuir tus contribuciones bajo la misma licencia que el original.