Window and custom features¶
When forecasting time series data, it may be useful to consider additional characteristics of the time series beyond just the lagged values. For example, the moving average of the previous n values may help to capture the trend in the series. The window_features
argument allows the inclusion of additional predictors created with the previous values of the series.
⚠ Warning
This section focuses on window features and other user-defined features derived from past values of the time series being modelled. These features are different from exogenous variables. Exogenous variables are external to the time series and are added to the model as additional predictors. See the Exogenous Variables section for more information on how to create window features and lags from exogenous variables.
Libraries and data¶
# Libraries
# ==============================================================================
import numpy as np
import pandas as pd
from lightgbm import LGBMRegressor
from skforecast.datasets import fetch_dataset
from skforecast.preprocessing import RollingFeatures
from skforecast.recursive import ForecasterRecursive
from skforecast.recursive import ForecasterRecursiveMultiSeries
# Data
# ==============================================================================
data = fetch_dataset(name="h2o", raw=False)
data.index.name = 'datetime'
data = data.rename(columns={'x': 'y'})
y = data['y']
y
h2o --- Monthly expenditure ($AUD) on corticosteroid drugs that the Australian health system had between 1991 and 2008. Hyndman R (2023). fpp3: Data for Forecasting: Principles and Practice(3rd Edition). http://pkg.robjhyndman.com/fpp3package/,https://github.com/robjhyndman /fpp3package, http://OTexts.com/fpp3. Shape of the dataset: (204, 1)
datetime 1991-07-01 0.429795 1991-08-01 0.400906 1991-09-01 0.432159 1991-10-01 0.492543 1991-11-01 0.502369 ... 2008-02-01 0.761822 2008-03-01 0.649435 2008-04-01 0.827887 2008-05-01 0.816255 2008-06-01 0.762137 Freq: MS, Name: y, Length: 204, dtype: float64
RollingFeatures¶
The RollingFeatures
class available is skforecast allows the creation of some of the most commonly used predictors:
- 'mean': the mean of the previous n values.
- 'std': the standard deviation of the previous n values.
- 'min': the minimum of the previous n values.
- 'max': the maximum of the previous n values.
- 'sum': the sum of the previous n values.
- 'median': the median of the previous n values.
- 'ratio_min_max': the ratio between the minimum and maximum of the previous n values.
- 'coef_variation': the coefficient of variation of the previous n values.
The user can specify these predictors by passing a list of strings to the stats
argument. The user can also specify the window size for each of these predictors by passing a list
of integers to the window_size
argument, if the value is the same for all the predictors, the user can pass a single integer
.
The following example demonstrates how to use the RollingFeatures
class to include rolling statistics (mean, minimum, and maximum values). Here, the rolling mean is computed with a window size of 10 and 20, while the minimum and maximum values use a window size of 10.
# Window features
# ==============================================================================
window_features = RollingFeatures(
stats = ['mean', 'mean', 'min', 'max'],
window_sizes = [20, 10, 10, 10]
)
# Create and fit forecaster
# ==============================================================================
forecaster = ForecasterRecursive(
regressor = LGBMRegressor(random_state=123, verbose=-1),
lags = 3,
window_features = window_features
)
forecaster.fit(y=y)
forecaster
ForecasterRecursive
General Information
- Regressor: LGBMRegressor
- Lags: [1 2 3]
- Window features: ['roll_mean_20', 'roll_mean_10', 'roll_min_10', 'roll_max_10']
- Window size: 20
- Exogenous included: False
- Weight function included: False
- Differentiation order: None
- Creation date: 2024-11-20 21:42:00
- Last fit date: 2024-11-20 21:42:01
- Skforecast version: 0.14.0
- Python version: 3.12.4
- Forecaster id: None
Exogenous Variables
-
None
Data Transformations
- Transformer for y: None
- Transformer for exog: None
Training Information
- Training range: [Timestamp('1991-07-01 00:00:00'), Timestamp('2008-06-01 00:00:00')]
- Training index type: DatetimeIndex
- Training index frequency: MS
Regressor Parameters
-
{'boosting_type': 'gbdt', 'class_weight': None, 'colsample_bytree': 1.0, 'importance_type': 'split', 'learning_rate': 0.1, 'max_depth': -1, 'min_child_samples': 20, 'min_child_weight': 0.001, 'min_split_gain': 0.0, 'n_estimators': 100, 'n_jobs': None, 'num_leaves': 31, 'objective': None, 'random_state': 123, 'reg_alpha': 0.0, 'reg_lambda': 0.0, 'subsample': 1.0, 'subsample_for_bin': 200000, 'subsample_freq': 0, 'verbose': -1}
Fit Kwargs
-
{}
# Predict
# ==============================================================================
steps = 36
predictions = forecaster.predict(steps=steps)
predictions.head(3)
2008-07-01 0.784687 2008-08-01 0.907246 2008-09-01 1.078849 Freq: MS, Name: pred, dtype: float64
The window size needed by this Forecaster is the maximum of the window sizes between the lagged values and the rolling features. In this case, this value is 20.
# Predict
# ==============================================================================
print("Window size lags : ", forecaster.max_lag)
print("Window size window features : ", forecaster.max_size_window_features)
print("Window size Forecaster : ", forecaster.window_size)
Window size lags : 3 Window size window features : 20 Window size Forecaster : 20
Additional arguments can be passed to the RollingFeatures
:
features_names
: By default, feature names follow the patternroll_<stat>_<window_size>
. For instance, a rolling mean with a window size of 20 is namedroll_mean_20
. Users can also assign custom names to each feature using alist
of strings.min_periods
: allows specifying the minimum number of observations required to compute the statistics during the training matrix generation (same as themin_periods
argument of pandas rolling). It can be a single integer or a list of integers, one for each statistic.fill_na
: define the strategy for handling missing values during the training matrix generation. Available methods are:'mean'
,'median'
,'ffill'
,'bfill'
, or afloat
value.
By inspecting the training matrices, it is possible to check that the rolling features have been correctly included.
# Training matrices used internally to fit the regressor
# ==============================================================================
X_train, y_train = forecaster.create_train_X_y(y=y)
X_train
lag_1 | lag_2 | lag_3 | roll_mean_20 | roll_mean_10 | roll_min_10 | roll_max_10 | |
---|---|---|---|---|---|---|---|
datetime | |||||||
1993-03-01 | 0.387554 | 0.751503 | 0.771258 | 0.496401 | 0.534009 | 0.361801 | 0.771258 |
1993-04-01 | 0.427283 | 0.387554 | 0.751503 | 0.496275 | 0.540557 | 0.387554 | 0.771258 |
1993-05-01 | 0.413890 | 0.427283 | 0.387554 | 0.496924 | 0.540893 | 0.387554 | 0.771258 |
1993-06-01 | 0.428859 | 0.413890 | 0.427283 | 0.496759 | 0.535440 | 0.387554 | 0.771258 |
1993-07-01 | 0.470126 | 0.428859 | 0.413890 | 0.495638 | 0.534906 | 0.387554 | 0.771258 |
... | ... | ... | ... | ... | ... | ... | ... |
2008-02-01 | 1.219941 | 1.176589 | 1.163534 | 0.980390 | 0.995834 | 0.561760 | 1.219941 |
2008-03-01 | 0.761822 | 1.219941 | 1.176589 | 0.978582 | 1.015840 | 0.745258 | 1.219941 |
2008-04-01 | 0.649435 | 0.761822 | 1.219941 | 0.966838 | 1.006258 | 0.649435 | 1.219941 |
2008-05-01 | 0.827887 | 0.649435 | 0.761822 | 0.955750 | 1.005253 | 0.649435 | 1.219941 |
2008-06-01 | 0.816255 | 0.827887 | 0.649435 | 0.946778 | 0.991464 | 0.649435 | 1.219941 |
184 rows × 7 columns
It is also possible to use the RollingFeatures
class outside the forecaster to gain a deeper insight into its behaviour. The transform
method computes the rolling features for a given numpy array, which is assumed to contain as many past observations as the maximum window size required to compute the features. The output is a numpy array with the rolling features.
# Create rolling features from a given array
# ==============================================================================
x = np.arange(20)
window_features.transform(X=x)
array([ 9.5, 14.5, 10. , 19. ])
The transform_batch
method is designed to transform a whole pandas Series from which multiple rolling windows can be extracted. The output is a pandas DataFrame with the rolling features.
# Create rolling features from a pandas series
# ==============================================================================
window_features.transform_batch(y).head(3)
roll_mean_20 | roll_mean_10 | roll_min_10 | roll_max_10 | |
---|---|---|---|---|
datetime | ||||
1993-03-01 | 0.496401 | 0.534009 | 0.361801 | 0.771258 |
1993-04-01 | 0.496275 | 0.540557 | 0.387554 | 0.771258 |
1993-05-01 | 0.496924 | 0.540893 | 0.387554 | 0.771258 |
The reason for these two different data transformation methods is that the first is used during prediction, where the forecaster only has access to the last window of the series. In contrast, the second is used during training, where the forecaster has access to the entire series.
Create your custom window features¶
RollingFeatures
is very useful for including some of the most commonly used predictors. However, users may need to include additional predictors that are not provided by this class. In such cases, users can create their own custom class to compute the desired features and include them in the forecaster.
The custom class must have these 2 methods:
transform_batch
: method to compute the features in batch from a pandas Series. This method will be used to compute the features during the training process. It must return a pandas DataFrame containing the rolling features.transform
: method to compute the features from a numpy array. This method will be used to compute the features during the prediction process. It must return a numpy array containing the computed statistics.
and these 2 attributes:
window_sizes
: size of the rolling window required to compute the features. It must be a list of integers.features_names
: list with the names of the output features. It must be a list of strings.
The follwing example shows how to create a custom class to include the rolling skewness and kurtosis with a window size of 20.
💡 Tip
If you have any doubt when creating your custom class, you can check the source code of the RollingFeatures
class available in the API Reference.
# Custom class to create rolling skewness features
# ==============================================================================
from scipy.stats import skew
class RollingSkewness():
"""
Custom class to create rolling skewness features.
"""
def __init__(self, window_sizes, features_names: list = 'rolling_skewness'):
if not isinstance(window_sizes, list):
window_sizes = [window_sizes]
self.window_sizes = window_sizes
self.features_names = features_names
def transform_batch(self, X: pd.Series) -> pd.DataFrame:
rolling_obj = X.rolling(window=self.window_sizes[0], center=False, closed='left')
rolling_skewness = rolling_obj.skew()
rolling_skewness = pd.DataFrame({
self.features_names: rolling_skewness
}).dropna()
return rolling_skewness
def transform(self, X: np.ndarray) -> np.ndarray:
X = X[~np.isnan(X)]
if len(X) > 0:
rolling_skewness = np.array([skew(X, bias=False)])
else:
rolling_skewness = np.array([np.nan])
return rolling_skewness
# Transform batch
# ==============================================================================
window_features = RollingSkewness(window_sizes=3)
window_features.transform_batch(y)
rolling_skewness | |
---|---|
datetime | |
1991-10-01 | -1.696160 |
1991-11-01 | 0.897261 |
1991-12-01 | -1.602797 |
1992-01-01 | 1.681518 |
1992-02-01 | -0.778727 |
... | ... |
2008-02-01 | 1.359033 |
2008-03-01 | -1.674974 |
2008-04-01 | 1.466482 |
2008-05-01 | -0.747574 |
2008-06-01 | -1.705640 |
201 rows × 1 columns
# Transform
# ==============================================================================
window_features.transform(X=np.array([6, 12, 8]))
array([0.93521953])
# Forecaster with custom rolling features
# ==============================================================================
forecaster = ForecasterRecursive(
regressor = LGBMRegressor(random_state=123, n_jobs=-1, verbose=-1),
lags = 3,
window_features = window_features
)
forecaster.fit(y=y)
forecaster.predict(steps=5)
2008-07-01 0.769870 2008-08-01 0.824683 2008-09-01 0.819595 2008-10-01 0.799849 2008-11-01 0.802432 Freq: MS, Name: pred, dtype: float64
Create your custom window features ForecasterRecursiveMultiSeries¶
When using a ForecasterRecursiveMultiSeries
, the transform
method must return a numpy array that calculates all the features for all the series contained in last_window
at once.
# Data download
# ==============================================================================
data_multiseries = fetch_dataset(name="items_sales")
data_multiseries.head()
items_sales ----------- Simulated time series for the sales of 3 different items. Simulated data. Shape of the dataset: (1097, 3)
item_1 | item_2 | item_3 | |
---|---|---|---|
date | |||
2012-01-01 | 8.253175 | 21.047727 | 19.429739 |
2012-01-02 | 22.777826 | 26.578125 | 28.009863 |
2012-01-03 | 27.549099 | 31.751042 | 32.078922 |
2012-01-04 | 25.895533 | 24.567708 | 27.252276 |
2012-01-05 | 21.379238 | 18.191667 | 20.357737 |
# Custom class to create rolling skewness features (multi-series)
# ==============================================================================
from scipy.stats import skew
class RollingSkewnessMultiSeries():
"""
Custom class to create rolling skewness features for multiple series.
"""
def __init__(self, window_sizes, features_names: list = 'rolling_skewness'):
if not isinstance(window_sizes, list):
window_sizes = [window_sizes]
self.window_sizes = window_sizes
self.features_names = features_names
def transform_batch(self, X: pd.Series) -> pd.DataFrame:
rolling_obj = X.rolling(window=self.window_sizes[0], center=False, closed='left')
rolling_skewness = rolling_obj.skew()
rolling_skewness = pd.DataFrame({
self.features_names: rolling_skewness
}).dropna()
return rolling_skewness
def transform(self, X: np.ndarray) -> np.ndarray:
X_dim = X.ndim
if X_dim == 1:
n_series = 1 # Only one series
X = X.reshape(-1, 1)
else:
n_series = X.shape[1] # Series (levels) to be predicted (present in last_window)
n_stats = 1 # Only skewness is calculated
rolling_skewness = np.full(
shape=(n_series, n_stats), fill_value=np.nan, dtype=float
)
for i in range(n_series):
if len(X) > 0:
rolling_skewness[i, :] = skew(X[:, i], bias=False)
else:
rolling_skewness[i, :] = np.nan
if X_dim == 1:
rolling_skewness = rolling_skewness.flatten()
return rolling_skewness
# Forecaster with custom rolling features
# ==============================================================================
window_features = RollingSkewnessMultiSeries(window_sizes=3)
forecaster = ForecasterRecursiveMultiSeries(
regressor = LGBMRegressor(random_state=123, verbose=-1),
lags = 3,
window_features = window_features,
encoding = 'ordinal'
)
forecaster.fit(series=data_multiseries)
forecaster.predict(steps=5, levels=None) # Predict all levels
item_1 | item_2 | item_3 | |
---|---|---|---|
2015-01-02 | 13.530901 | 19.709656 | 19.378026 |
2015-01-03 | 14.252429 | 18.608582 | 20.457204 |
2015-01-04 | 14.230324 | 19.536894 | 21.238047 |
2015-01-05 | 16.339439 | 19.174911 | 21.039226 |
2015-01-06 | 17.231169 | 18.821718 | 19.032785 |
# Transform output multi-series shape (n_levels, n_stats)
# ==============================================================================
window_features.transform(pd.DataFrame(forecaster.last_window_).to_numpy())
array([[-1.7304836 ], [ 1.69345527], [ 0.27695088]])
Adding multiple window features¶
It is possible to include multiple window features in all forecasters. The window_features
argument must be a list of instances of the classes that compute the desired features.
# Forecaster with multiple window features
# ==============================================================================
window_features = [
RollingFeatures(stats=['mean', 'mean'], window_sizes=[20, 10]),
RollingSkewness(window_sizes=10)
]
forecaster = ForecasterRecursive(
regressor = LGBMRegressor(random_state=123, verbose=-1),
lags = 3,
window_features = window_features
)
forecaster.fit(y=y)
forecaster
ForecasterRecursive
General Information
- Regressor: LGBMRegressor
- Lags: [1 2 3]
- Window features: ['roll_mean_20', 'roll_mean_10', 'rolling_skewness']
- Window size: 20
- Exogenous included: False
- Weight function included: False
- Differentiation order: None
- Creation date: 2024-11-20 21:42:05
- Last fit date: 2024-11-20 21:42:05
- Skforecast version: 0.14.0
- Python version: 3.12.4
- Forecaster id: None
Exogenous Variables
-
None
Data Transformations
- Transformer for y: None
- Transformer for exog: None
Training Information
- Training range: [Timestamp('1991-07-01 00:00:00'), Timestamp('2008-06-01 00:00:00')]
- Training index type: DatetimeIndex
- Training index frequency: MS
Regressor Parameters
-
{'boosting_type': 'gbdt', 'class_weight': None, 'colsample_bytree': 1.0, 'importance_type': 'split', 'learning_rate': 0.1, 'max_depth': -1, 'min_child_samples': 20, 'min_child_weight': 0.001, 'min_split_gain': 0.0, 'n_estimators': 100, 'n_jobs': None, 'num_leaves': 31, 'objective': None, 'random_state': 123, 'reg_alpha': 0.0, 'reg_lambda': 0.0, 'subsample': 1.0, 'subsample_for_bin': 200000, 'subsample_freq': 0, 'verbose': -1}
Fit Kwargs
-
{}
# Inspect prediction matrix
# ==============================================================================
forecaster.create_predict_X(steps=5)
lag_1 | lag_2 | lag_3 | roll_mean_20 | roll_mean_10 | rolling_skewness | |
---|---|---|---|---|---|---|
2008-07-01 | 0.762137 | 0.816255 | 0.827887 | 0.926472 | 0.959856 | -0.130309 |
2008-08-01 | 0.953737 | 0.762137 | 0.816255 | 0.918757 | 0.944132 | -0.053265 |
2008-09-01 | 1.027884 | 0.953737 | 0.762137 | 0.914148 | 0.935922 | -0.024878 |
2008-10-01 | 0.963651 | 1.027884 | 0.953737 | 0.901165 | 0.915934 | -0.022704 |
2008-11-01 | 1.096952 | 0.963651 | 1.027884 | 0.926125 | 0.907970 | -0.173034 |
# Predict
# ==============================================================================
forecaster.predict(steps=5)
2008-07-01 0.953737 2008-08-01 1.027884 2008-09-01 0.963651 2008-10-01 1.096952 2008-11-01 1.096672 Freq: MS, Name: pred, dtype: float64
# Forecaster with multiple window features
# ==============================================================================
window_features = [
RollingFeatures(stats=['mean', 'mean'], window_sizes=[20, 10]),
RollingSkewnessMultiSeries(window_sizes=10)
]
forecaster = ForecasterRecursiveMultiSeries(
regressor = LGBMRegressor(random_state=123, verbose=-1),
lags = 3,
window_features = window_features
)
forecaster.fit(series=data_multiseries)
forecaster
ForecasterRecursiveMultiSeries
General Information
- Regressor: LGBMRegressor
- Lags: [1 2 3]
- Window features: ['roll_mean_20', 'roll_mean_10', 'rolling_skewness']
- Window size: 20
- Series encoding: ordinal
- Exogenous included: False
- Weight function included: False
- Series weights: None
- Differentiation order: None
- Creation date: 2024-11-20 21:42:05
- Last fit date: 2024-11-20 21:42:05
- Skforecast version: 0.14.0
- Python version: 3.12.4
- Forecaster id: None
Exogenous Variables
-
None
Data Transformations
- Transformer for series: None
- Transformer for exog: None
Training Information
- Series names (levels): item_1, item_2, item_3
- Training range: 'item_1': ['2012-01-01', '2015-01-01'], 'item_2': ['2012-01-01', '2015-01-01'], 'item_3': ['2012-01-01', '2015-01-01']
- Training index type: DatetimeIndex
- Training index frequency: D
Regressor Parameters
-
{'boosting_type': 'gbdt', 'class_weight': None, 'colsample_bytree': 1.0, 'importance_type': 'split', 'learning_rate': 0.1, 'max_depth': -1, 'min_child_samples': 20, 'min_child_weight': 0.001, 'min_split_gain': 0.0, 'n_estimators': 100, 'n_jobs': None, 'num_leaves': 31, 'objective': None, 'random_state': 123, 'reg_alpha': 0.0, 'reg_lambda': 0.0, 'subsample': 1.0, 'subsample_for_bin': 200000, 'subsample_freq': 0, 'verbose': -1}
Fit Kwargs
-
{}
# Inspect prediction matrix for item_1
# ==============================================================================
forecaster.create_predict_X(steps=5)['item_1']
lag_1 | lag_2 | lag_3 | roll_mean_20 | roll_mean_10 | rolling_skewness | _level_skforecast | |
---|---|---|---|---|---|---|---|
2015-01-02 | 10.496302 | 18.721223 | 18.857026 | 19.627230 | 18.109484 | -1.596973 | 0.0 |
2015-01-03 | 18.861343 | 10.496302 | 18.721223 | 19.512458 | 17.676044 | -1.504973 | 0.0 |
2015-01-04 | 17.778515 | 18.861343 | 10.496302 | 19.364411 | 17.349451 | -1.337133 | 0.0 |
2015-01-05 | 19.272862 | 17.778515 | 18.861343 | 19.293172 | 17.434766 | -1.284426 | 0.0 |
2015-01-06 | 20.198399 | 19.272862 | 17.778515 | 19.143477 | 17.800624 | -1.469888 | 0.0 |
# Inspect prediction matrix
# ==============================================================================
forecaster.predict(steps=5)
item_1 | item_2 | item_3 | |
---|---|---|---|
2015-01-02 | 18.861343 | 21.561617 | 21.071428 |
2015-01-03 | 17.778515 | 21.288322 | 20.891229 |
2015-01-04 | 19.272862 | 20.580463 | 21.024649 |
2015-01-05 | 20.198399 | 20.020874 | 21.414145 |
2015-01-06 | 19.464618 | 20.100328 | 21.244314 |