diff --git a/aeon/forecasting/__init__.py b/aeon/forecasting/__init__.py new file mode 100644 index 0000000000..de203a0bcd --- /dev/null +++ b/aeon/forecasting/__init__.py @@ -0,0 +1,13 @@ +"""Forecasters.""" + +__all__ = [ + "DummyForecaster", + "BaseForecaster", + "RegressionForecaster", + "ETSForecaster", +] + +from aeon.forecasting._dummy import DummyForecaster +from aeon.forecasting._ets import ETSForecaster +from aeon.forecasting._regression import RegressionForecaster +from aeon.forecasting.base import BaseForecaster diff --git a/aeon/forecasting/_dummy.py b/aeon/forecasting/_dummy.py new file mode 100644 index 0000000000..7525b6ccd0 --- /dev/null +++ b/aeon/forecasting/_dummy.py @@ -0,0 +1,27 @@ +"""DummyForecaster always predicts the last value seen in training.""" + +from aeon.forecasting.base import BaseForecaster + + +class DummyForecaster(BaseForecaster): + """Dummy forecaster always predicts the last value seen in training.""" + + def __init__(self): + """Initialize DummyForecaster.""" + self.last_value_ = None + super().__init__(horizon=1, axis=1) + + def _fit(self, y, exog=None): + """Fit dummy forecaster.""" + y = y.squeeze() + self.last_value_ = y[-1] + return self + + def _predict(self, y=None, exog=None): + """Predict using dummy forecaster.""" + return self.last_value_ + + def _forecast(self, y, exog=None): + """Forecast using dummy forecaster.""" + y = y.squeeze() + return y[-1] diff --git a/aeon/forecasting/_ets.py b/aeon/forecasting/_ets.py new file mode 100644 index 0000000000..b477046895 --- /dev/null +++ b/aeon/forecasting/_ets.py @@ -0,0 +1,416 @@ +"""ETSForecaster class. + +An implementation of the exponential smoothing statistics forecasting algorithm. +Implements additive and multiplicative error models, +None, additive and multiplicative (including damped) trend and +None, additive and multiplicative seasonality +""" + +__maintainer__ = [] +__all__ = ["ETSForecaster", "NONE", "ADDITIVE", "MULTIPLICATIVE"] + +import numpy as np +from numba import njit + +from aeon.forecasting.base import BaseForecaster + +NOGIL = False +CACHE = True + +NONE = 0 +ADDITIVE = 1 +MULTIPLICATIVE = 2 + + +class ETSForecaster(BaseForecaster): + """Exponential Smoothing forecaster. + + An implementation of the exponential smoothing forecasting algorithm. + Implements additive and multiplicative error models, None, additive and + multiplicative (including damped) trend and None, additive and mutliplicative + seasonality. See [1]_ for a description. + + Parameters + ---------- + error_type : int, default = 1 + Either NONE (0), ADDITIVE (1) or MULTIPLICATIVE (2). + trend_type : int, default = 0 + Either NONE (0), ADDITIVE (1) or MULTIPLICATIVE (2). + seasonality_type : int, default = 0 + Either NONE (0), ADDITIVE (1) or MULTIPLICATIVE (2). + seasonal_period : int, default=1 + Length of seasonality period. If seasonality_type is NONE, this is assumed to + be 1 + alpha : float, default = 0.1 + Level smoothing parameter. + beta : float, default = 0.01 + Trend smoothing parameter. If trend_type is NONE, this is assumed to be 0.0. + gamma : float, default = 0.01 + Seasonal smoothing parameter. If seasonality is NONE, this is assumed to be + 0.0. + phi : float, default = 0.99 + Trend damping smoothing parameters + horizon : int, default = 1 + The horizon to forecast to. + + Attributes + ---------- + mean_sq_err_ : float + Mean squared error. + likelihood_ : float + Likelihood of the fitted model based on residuals. + residuals_ : arraylike + List of train set differences between fitted and actual values. + n_timpoints_ : int + Length of the series passed to fit. + + References + ---------- + .. [1] R. J. Hyndman and G. Athanasopoulos, + Forecasting: Principles and Practice. Melbourne, Australia: OTexts, 2014. + + Examples + -------- + >>> from aeon.forecasting import ETSForecaster + >>> from aeon.datasets import load_airline + >>> y = load_airline() + >>> forecaster = ETSForecaster(alpha=0.4, beta=0.2, gamma=0.5, phi=0.8, horizon=1) + >>> forecaster.fit(y) + ETSForecaster(alpha=0.4, beta=0.2, gamma=0.5, phi=0.8) + >>> forecaster.predict() + 449.9435566831507 + """ + + def __init__( + self, + error_type=ADDITIVE, + trend_type=NONE, + seasonality_type=NONE, + seasonal_period=1, + alpha=0.1, + beta=0.01, + gamma=0.01, + phi=0.99, + horizon=1, + ): + self.error_type = error_type + self.trend_type = trend_type + self.seasonality_type = seasonality_type + self.seasonal_period = seasonal_period + self.alpha = alpha + self.beta = beta + self.gamma = gamma + self.phi = phi + self.mean_sq_err_ = 0 + self.likelihood_ = 0 + self.residuals_ = [] + self.n_timpoints_ = 0 + super().__init__(horizon=horizon, axis=1) + + def _fit(self, y, exog=None): + """Fit Exponential Smoothing forecaster to series y. + + Fit a forecaster to predict self.horizon steps ahead using y. + + Parameters + ---------- + y : np.ndarray + A time series on which to learn a forecaster to predict horizon ahead + exog : np.ndarray, default =None + Optional exogenous time series data assumed to be aligned with y + + Returns + ------- + self + Fitted BaseForecaster. + """ + self.n_timepoints_ = len(y) + if self.error_type != MULTIPLICATIVE and self.error_type != ADDITIVE: + raise ValueError("Error must be either additive or multiplicative") + self._seasonal_period = self.seasonal_period + if self.seasonal_period < 1 or self.seasonality_type == NONE: + self._seasonal_period = 1 + self._beta = self.beta + if self.trend_type == NONE: + self._beta = 0 + self._gamma = self.gamma + if self.seasonality_type == NONE: + self._gamma = 0 + data = np.array(y.squeeze(), dtype=np.float64) + ( + self._level, + self._trend, + self._seasonality, + self.residuals_, + self.mean_sq_err_, + self.likelihood_, + ) = _fit_numba( + data, + self.error_type, + self.trend_type, + self.seasonality_type, + self._seasonal_period, + self.alpha, + self._beta, + self._gamma, + self.phi, + ) + return self + + def _predict(self, y=None, exog=None): + """ + Predict the next horizon steps ahead. + + Parameters + ---------- + y : np.ndarray, default = None + A time series to predict the next horizon value for. If None, + predict the next horizon value after series seen in fit. + exog : np.ndarray, default = None + Optional exogenous time series data assumed to be aligned with y + + Returns + ------- + float + single prediction self.horizon steps ahead of y. + """ + return _predict_numba( + self.trend_type, + self.seasonality_type, + self._level, + self._trend, + self._seasonality, + self.phi, + self.horizon, + self.n_timepoints_, + self.seasonal_period, + ) + + +@njit(nogil=NOGIL, cache=CACHE) +def _fit_numba( + data, + error_type, + trend_type, + seasonality_type, + seasonal_period, + alpha, + beta, + gamma, + phi, +): + n_timepoints = len(data) + level, trend, seasonality = _initialise( + trend_type, seasonality_type, seasonal_period, data + ) + mse = 0 + lhood = 0 + mul_likelihood_pt2 = 0 + res = np.zeros(n_timepoints) # 1 Less residual than data points + for t, data_item in enumerate(data[seasonal_period:]): + # Calculate level, trend, and seasonal components + fitted_value, error, level, trend, seasonality[t % seasonal_period] = ( + _update_states( + error_type, + trend_type, + seasonality_type, + level, + trend, + seasonality[t % seasonal_period], + data_item, + alpha, + beta, + gamma, + phi, + ) + ) + res[t] = error + mse += (data_item - fitted_value) ** 2 + lhood += error * error + mul_likelihood_pt2 += np.log(np.fabs(fitted_value)) + mse /= n_timepoints - seasonal_period + lhood = (n_timepoints - seasonal_period) * np.log(lhood) + if error_type == MULTIPLICATIVE: + lhood += 2 * mul_likelihood_pt2 + return level, trend, seasonality, res, mse, lhood + + +def _predict_numba( + trend_type, + seasonality_type, + level, + trend, + seasonality, + phi, + horizon, + n_timepoints, + seasonal_period, +): + # Generate forecasts based on the final values of level, trend, and seasonals + if phi == 1: # No damping case + phi_h = float(horizon) + else: + # Geometric series formula for calculating phi + phi^2 + ... + phi^h + phi_h = phi * (1 - phi**horizon) / (1 - phi) + seasonal_index = (n_timepoints + horizon) % seasonal_period + return _predict_value( + trend_type, + seasonality_type, + level, + trend, + seasonality[seasonal_index], + phi_h, + )[0] + + +@njit(nogil=NOGIL, cache=CACHE) +def _initialise(trend_type, seasonality_type, seasonal_period, data): + """ + Initialize level, trend, and seasonality values for the ETS model. + + Parameters + ---------- + data : array-like + The time series data + (should contain at least two full seasons if seasonality is specified) + """ + # Initial Level: Mean of the first season + level = np.mean(data[:seasonal_period]) + # Initial Trend + if trend_type == ADDITIVE: + # Average difference between corresponding points in the first two seasons + trend = np.mean( + data[seasonal_period : 2 * seasonal_period] - data[:seasonal_period] + ) + elif trend_type == MULTIPLICATIVE: + # Average ratio between corresponding points in the first two seasons + trend = np.mean( + data[seasonal_period : 2 * seasonal_period] / data[:seasonal_period] + ) + else: + # No trend + trend = 0 + # Initial Seasonality + if seasonality_type == ADDITIVE: + # Seasonal component is the difference + # from the initial level for each point in the first season + seasonality = data[:seasonal_period] - level + elif seasonality_type == MULTIPLICATIVE: + # Seasonal component is the ratio of each point in the first season + # to the initial level + seasonality = data[:seasonal_period] / level + else: + # No seasonality + seasonality = np.zeros(1) + return level, trend, seasonality + + +@njit(nogil=NOGIL, cache=CACHE) +def _update_states( + error_type, + trend_type, + seasonality_type, + level, + trend, + seasonality, + data_item: int, + alpha, + beta, + gamma, + phi, +): + """ + Update level, trend, and seasonality components. + + Using state space equations for an ETS model. + + Parameters + ---------- + data_item: float + The current value of the time series. + seasonal_index: int + The index to update the seasonal component. + """ + # Retrieve the current state values + curr_level = level + curr_seasonality = seasonality + fitted_value, damped_trend, trend_level_combination = _predict_value( + trend_type, seasonality_type, level, trend, seasonality, phi + ) + # Calculate the error term (observed value - fitted value) + if error_type == MULTIPLICATIVE: + error = data_item / fitted_value - 1 # Multiplicative error + else: + error = data_item - fitted_value # Additive error + # Update level + if error_type == MULTIPLICATIVE: + level = trend_level_combination * (1 + alpha * error) + trend = damped_trend * (1 + beta * error) + seasonality = curr_seasonality * (1 + gamma * error) + if seasonality_type == ADDITIVE: + level += alpha * error * curr_seasonality # Add seasonality correction + seasonality += gamma * error * trend_level_combination + if trend_type == ADDITIVE: + trend += (curr_level + curr_seasonality) * beta * error + else: + trend += curr_seasonality / curr_level * beta * error + elif trend_type == ADDITIVE: + trend += curr_level * beta * error + else: + level_correction = 1 + trend_correction = 1 + seasonality_correction = 1 + if seasonality_type == MULTIPLICATIVE: + # Add seasonality correction + level_correction *= curr_seasonality + trend_correction *= curr_seasonality + seasonality_correction *= trend_level_combination + if trend_type == MULTIPLICATIVE: + trend_correction *= curr_level + level = trend_level_combination + alpha * error / level_correction + trend = damped_trend + beta * error / trend_correction + seasonality = curr_seasonality + gamma * error / seasonality_correction + return (fitted_value, error, level, trend, seasonality) + + +@njit(nogil=NOGIL, cache=CACHE) +def _predict_value(trend_type, seasonality_type, level, trend, seasonality, phi): + """ + + Generate various useful values, including the next fitted value. + + Parameters + ---------- + trend : float + The current trend value for the model + level : float + The current level value for the model + seasonality : float + The current seasonality value for the model + phi : float + The damping parameter for the model + + Returns + ------- + fitted_value : float + single prediction based on the current state variables. + damped_trend : float + The damping parameter combined with the trend dependant on the model type + trend_level_combination : float + Combination of the trend and level based on the model type. + """ + # Apply damping parameter and + # calculate commonly used combination of trend and level components + if trend_type == MULTIPLICATIVE: + damped_trend = trend**phi + trend_level_combination = level * damped_trend + else: # Additive trend, if no trend, then trend = 0 + damped_trend = trend * phi + trend_level_combination = level + damped_trend + + # Calculate forecast (fitted value) based on the current components + if seasonality_type == MULTIPLICATIVE: + fitted_value = trend_level_combination * seasonality + else: # Additive seasonality, if no seasonality, then seasonality = 0 + fitted_value = trend_level_combination + seasonality + return fitted_value, damped_trend, trend_level_combination diff --git a/aeon/forecasting/_regression.py b/aeon/forecasting/_regression.py new file mode 100644 index 0000000000..79393160b1 --- /dev/null +++ b/aeon/forecasting/_regression.py @@ -0,0 +1,105 @@ +"""Window-based regression forecaster. + +General purpose forecaster to use with any scikit learn or aeon compatible +regressor. Simply forms a collection of windows from the time series and trains to +predict the next +""" + +import numpy as np +from sklearn.linear_model import LinearRegression + +from aeon.forecasting.base import BaseForecaster + + +class RegressionForecaster(BaseForecaster): + """ + Regression based forecasting. + + Container for forecaster that reduces forecasting to regression through a + window. Form a collection of sub series of length `window` through a sliding + winodw to form X, take `horizon` points ahead to form `y`, then apply an aeon or + sklearn regressor. + + + Parameters + ---------- + window : int + The window prior to the current time point to use in forecasting. So if + horizon is one, forecaster will train using points $i$ to $window+i-1$ to + predict value $window+i$. If horizon is 4, forecaster will used points $i$ + to $window+i-1$ to predict value $window+i+3$. If None, the algorithm will + internally determine what data to use to predict `horizon` steps ahead. + horizon : int, default =1 + The number of time steps ahead to forecast. If horizon is one, the forecaster + will learn to predict one point ahead + regressor : object, default =None + Regression estimator that implements BaseRegressor or is otherwise compatible + with sklearn regressors. + """ + + def __init__(self, window, horizon=1, regressor=None): + self.window = window + self.regressor = regressor + super().__init__(horizon=horizon, axis=1) + + def _fit(self, y, exog=None): + """Fit forecaster to time series. + + Split X into windows of length window and train the forecaster on each window + to predict the horizon ahead. + + Parameters + ---------- + X : Time series on which to learn a forecaster + + Returns + ------- + self + Fitted estimator + """ + # Window data + if self.regressor is None: + self.regressor_ = LinearRegression() + else: + self.regressor_ = self.regressor + y = y.squeeze() + X = np.lib.stride_tricks.sliding_window_view(y, window_shape=self.window) + # Ignore the final horizon values: need to store these for pred with empty y + X = X[: -self.horizon] + # Extract y + y = y[self.window + self.horizon - 1 :] + self.last_ = y[-self.window :] + self.last_ = self.last_.reshape(1, -1) + self.regressor_.fit(X=X, y=y) + return self + + def _predict(self, y=None, exog=None): + """Predict values for time series X.""" + if y is None: + return self.regressor_.predict(self.last_) + last = y[:, -self.window :] + return self.regressor_.predict(last) + + def _forecast(self, y, exog=None): + """Forecast values for time series X. + + NOTE: deal with horizons + """ + self.fit(y, exog) + return self.predict() + + @classmethod + def _get_test_params(cls, parameter_set="default"): + """Return testing parameter settings for the estimator. + + Parameters + ---------- + parameter_set : str, default='default' + Name of the parameter set to return. + + Returns + ------- + dict + Dictionary of testing parameter settings. + """ + return {"window": 4} diff --git a/aeon/forecasting/base.py b/aeon/forecasting/base.py new file mode 100644 index 0000000000..e67712c58a --- /dev/null +++ b/aeon/forecasting/base.py @@ -0,0 +1,159 @@ +"""BaseForecaster class. + +A simplified first base class for forecasting models. + +""" + +from abc import abstractmethod + +import numpy as np +import pandas as pd + +from aeon.base import BaseSeriesEstimator +from aeon.base._base_series import VALID_SERIES_INNER_TYPES + + +class BaseForecaster(BaseSeriesEstimator): + """ + Abstract base class for time series forecasters. + + The base forecaster specifies the methods and method signatures that all + forecasters have to implement. Attributes with an underscore suffix are set in the + method fit. + + Parameters + ---------- + horizon : int, default =1 + The number of time steps ahead to forecast. If horizon is one, the forecaster + will learn to predict one point ahead. + """ + + _tags = { + "capability:univariate": True, + "capability:multivariate": False, + "capability:missing_values": False, + "fit_is_empty": False, + "y_inner_type": "np.ndarray", + } + + def __init__(self, horizon, axis): + self.horizon = horizon + self.meta_ = None # Meta data related to y on the last fit + super().__init__(axis) + + def fit(self, y, exog=None): + """Fit forecaster to series y. + + Fit a forecaster to predict self.horizon steps ahead using y. + + Parameters + ---------- + y : np.ndarray + A time series on which to learn a forecaster to predict horizon ahead. + exog : np.ndarray, default =None + Optional exogenous time series data assumed to be aligned with y. + + Returns + ------- + self + Fitted BaseForecaster. + """ + if self.get_tag("fit_is_empty"): + self.is_fitted = True + return self + + self._check_X(y, self.axis) + y = self._convert_y(y, self.axis) + if exog is not None: + raise NotImplementedError("Exogenous variables not yet supported") + self.is_fitted = True + return self._fit(y, exog) + + @abstractmethod + def _fit(self, y, exog=None): ... + + def predict(self, y=None, exog=None): + """Predict the next horizon steps ahead. + + Parameters + ---------- + y : np.ndarray, default = None + A time series to predict the next horizon value for. If None, + predict the next horizon value after series seen in fit. + exog : np.ndarray, default =None + Optional exogenous time series data assumed to be aligned with y. + + Returns + ------- + float + single prediction self.horizon steps ahead of y. + """ + self._check_is_fitted() + if y is not None: + self._check_X(y, self.axis) + y = self._convert_y(y, self.axis) + if exog is not None: + raise NotImplementedError("Exogenous variables not yet supported") + return self._predict(y, exog) + + @abstractmethod + def _predict(self, y=None, exog=None): ... + + def forecast(self, y, exog=None): + """Forecast the next horizon steps ahead. + + By default this is simply fit followed by predict. + + Parameters + ---------- + y : np.ndarray, default = None + A time series to predict the next horizon value for. If None, + predict the next horizon value after series seen in fit. + exog : np.ndarray, default =None + Optional exogenous time series data assumed to be aligned with y. + + Returns + ------- + float + single prediction self.horizon steps ahead of y. + """ + self._check_X(y, self.axis) + y = self._convert_y(y, self.axis) + return self._forecast(y, exog) + + def _forecast(self, y, exog=None): + """Forecast values for time series X.""" + self.fit(y, exog) + return self._predict(y, exog) + + def _convert_y(self, y: VALID_SERIES_INNER_TYPES, axis: int): + """Convert y to self.get_tag("y_inner_type").""" + if axis > 1 or axis < 0: + raise ValueError(f"Input axis should be 0 or 1, saw {axis}") + + inner_type = self.get_tag("y_inner_type") + if not isinstance(inner_type, list): + inner_type = [inner_type] + inner_names = [i.split(".")[-1] for i in inner_type] + + input = type(y).__name__ + if input not in inner_names: + if inner_names[0] == "ndarray": + y = y.to_numpy() + elif inner_names[0] == "DataFrame": + # converting a 1d array will create a 2d array in axis 0 format + transpose = False + if y.ndim == 1 and axis == 1: + transpose = True + y = pd.DataFrame(y) + if transpose: + y = y.T + else: + raise ValueError( + f"Unsupported inner type {inner_names[0]} derived from {inner_type}" + ) + if y.ndim > 1 and self.axis != axis: + y = y.T + elif y.ndim == 1 and isinstance(y, np.ndarray): + y = y[np.newaxis, :] if self.axis == 1 else y[:, np.newaxis] + return y diff --git a/aeon/forecasting/tests/__init__.py b/aeon/forecasting/tests/__init__.py new file mode 100644 index 0000000000..90b32266a4 --- /dev/null +++ b/aeon/forecasting/tests/__init__.py @@ -0,0 +1 @@ +"""Forecaster tests.""" diff --git a/aeon/forecasting/tests/test_base.py b/aeon/forecasting/tests/test_base.py new file mode 100644 index 0000000000..e1a634eba3 --- /dev/null +++ b/aeon/forecasting/tests/test_base.py @@ -0,0 +1,16 @@ +"""Test base forecaster.""" + +import numpy as np + +from aeon.forecasting import DummyForecaster + + +def test_base_forecaster(): + """Test base forecaster functionality.""" + f = DummyForecaster() + y = np.random.rand(50) + f.fit(y) + p1 = f.predict() + assert p1 == y[-1] + p2 = f.forecast(y) + assert p2 == p1 diff --git a/aeon/forecasting/tests/test_regressor.py b/aeon/forecasting/tests/test_regressor.py new file mode 100644 index 0000000000..ec6e273bfd --- /dev/null +++ b/aeon/forecasting/tests/test_regressor.py @@ -0,0 +1,16 @@ +"""Test the regression forecaster.""" + +from aeon.datasets import load_airline +from aeon.forecasting import RegressionForecaster + + +def test_regression_forecaster(): + """Test the regression forecaster.""" + y = load_airline() + f = RegressionForecaster(window=10) + f.fit(y) + p = f.predict() + p2 = f.predict(y) + assert p == p2 + p3 = f.forecast(y) + assert p == p3 diff --git a/aeon/testing/estimator_checking/_yield_forecasting_checks.py b/aeon/testing/estimator_checking/_yield_forecasting_checks.py new file mode 100644 index 0000000000..5c62d2f05d --- /dev/null +++ b/aeon/testing/estimator_checking/_yield_forecasting_checks.py @@ -0,0 +1,51 @@ +"""Tests for all forecasters.""" + +from functools import partial + +import numpy as np + +from aeon.base._base import _clone_estimator +from aeon.base._base_series import VALID_SERIES_INPUT_TYPES + + +def _yield_forecasting_checks(estimator_class, estimator_instances, datatypes): + """Yield all forecasting checks for an aeon forecaster.""" + # only class required + yield partial(check_forecasting_base_functionality, estimator_class=estimator_class) + + # test class instances + for _, estimator in enumerate(estimator_instances): + # no data needed + yield partial(check_forecaster_instance, estimator=estimator) + + +def check_forecasting_base_functionality(estimator_class): + """Test compliance with the base class contract.""" + # Test they dont override final methods, because python does not enforce this + assert "fit" not in estimator_class.__dict__ + assert "predict" not in estimator_class.__dict__ + assert "forecast" not in estimator_class.__dict__ + fit_is_empty = estimator_class.get_class_tag(tag_name="fit_is_empty") + assert not fit_is_empty == "_fit" not in estimator_class.__dict__ + # Test valid tag for X_inner_type + X_inner_type = estimator_class.get_class_tag(tag_name="X_inner_type") + assert X_inner_type in VALID_SERIES_INPUT_TYPES + # Must have at least one set to True + multi = estimator_class.get_class_tag(tag_name="capability:multivariate") + uni = estimator_class.get_class_tag(tag_name="capability:univariate") + assert multi or uni + + +def check_forecaster_instance(estimator): + """Test forecasters.""" + estimator = _clone_estimator(estimator) + pass + # Sort + # Check output correct: predict should return a float + y = np.array([0.5, 0.7, 0.8, 0.9, 1.0]) + estimator.fit(y) + p = estimator.predict() + assert isinstance(p, float) + # forecast should return a float equal to fit/predict + p2 = estimator.forecast(y) + assert p == p2 diff --git a/aeon/testing/mock_estimators/_mock_forecasters.py b/aeon/testing/mock_estimators/_mock_forecasters.py new file mode 100644 index 0000000000..f5bb86d249 --- /dev/null +++ b/aeon/testing/mock_estimators/_mock_forecasters.py @@ -0,0 +1,22 @@ +"""Mock forecasters useful for testing and debugging. + +Used in tests for the forecasting base class. +""" + +from aeon.forecasting.base import BaseForecaster + + +class MockForecaster(BaseForecaster): + """Mock forecaster for testing.""" + + def __init__(self): + super().__init__() + + def _fit(self, y, X=None): + return self + + def _predict(self, y): + return 1.0 + + def _forecast(self, y, X=None): + return 1.0 diff --git a/aeon/testing/testing_data.py b/aeon/testing/testing_data.py index 47b990b220..eb134cddda 100644 --- a/aeon/testing/testing_data.py +++ b/aeon/testing/testing_data.py @@ -7,6 +7,7 @@ from aeon.classification import BaseClassifier from aeon.classification.early_classification import BaseEarlyClassifier from aeon.clustering import BaseClusterer +from aeon.forecasting import BaseForecaster from aeon.regression import BaseRegressor from aeon.segmentation import BaseSegmenter from aeon.similarity_search import BaseSimilaritySearch @@ -22,9 +23,11 @@ ) from aeon.transformations.collection import BaseCollectionTransformer from aeon.transformations.series import BaseSeriesTransformer +from aeon.utils.conversion import convert_collection data_rng = np.random.RandomState(42) +# Collection testing data EQUAL_LENGTH_UNIVARIATE_CLASSIFICATION = { "numpy3D": { @@ -725,39 +728,61 @@ "numpy3D": { "train": (X_classification_missing_train, y_classification_missing_train), "test": (X_classification_missing_test, y_classification_missing_test), - } + }, + "np-list": { + "train": ( + convert_collection(X_classification_missing_train, "np-list"), + y_classification_missing_train, + ), + "test": ( + convert_collection(X_classification_missing_test, "np-list"), + y_classification_missing_test, + ), + }, } -X_classification_missing_train, y_classification_missing_train = make_example_3d_numpy( +X_regression_missing_train, y_regression_missing_train = make_example_3d_numpy( n_cases=10, n_channels=1, n_timepoints=20, random_state=data_rng.randint(np.iinfo(np.int32).max), regression_target=True, ) -X_classification_missing_test, y_classification_missing_test = make_example_3d_numpy( +X_regression_missing_test, y_regression_missing_test = make_example_3d_numpy( n_cases=5, n_channels=1, n_timepoints=20, random_state=data_rng.randint(np.iinfo(np.int32).max), regression_target=True, ) -X_classification_missing_train[:, :, data_rng.choice(20, 2)] = np.nan -X_classification_missing_test[:, :, data_rng.choice(20, 2)] = np.nan +X_regression_missing_train[:, :, data_rng.choice(20, 2)] = np.nan +X_regression_missing_test[:, :, data_rng.choice(20, 2)] = np.nan MISSING_VALUES_REGRESSION = { "numpy3D": { - "train": (X_classification_missing_train, y_classification_missing_train), - "test": (X_classification_missing_test, y_classification_missing_test), - } + "train": (X_regression_missing_train, y_regression_missing_train), + "test": (X_regression_missing_test, y_regression_missing_test), + }, + "np-list": { + "train": ( + convert_collection(X_regression_missing_train, "np-list"), + y_regression_missing_train, + ), + "test": ( + convert_collection(X_regression_missing_test, "np-list"), + y_regression_missing_test, + ), + }, } +# Series testing data + X_series = make_example_1d_numpy( n_timepoints=40, random_state=data_rng.randint(np.iinfo(np.int32).max) ) X_series2 = X_series[20:40] X_series = X_series[:20] -UNIVARIATE_SERIES_NOLABEL = {"train": (X_series, None), "test": (X_series2, None)} +UNIVARIATE_SERIES_NONE = {"train": (X_series, None), "test": (X_series2, None)} X_series_mv = make_example_2d_numpy_series( n_timepoints=40, @@ -767,7 +792,7 @@ ) X_series_mv2 = X_series_mv[:, 20:40] X_series_mv = X_series_mv[:, :20] -MULTIVARIATE_SERIES_NOLABEL = { +MULTIVARIATE_SERIES_NONE = { "train": (X_series_mv, None), "test": (X_series_mv2, None), } @@ -779,7 +804,12 @@ X_series_mi2[data_rng.choice(20, 1)] = np.nan X_series_mi = X_series_mi[:20] X_series_mi[data_rng.choice(20, 2)] = np.nan -MISSING_VALUES_NOLABEL = {"train": (X_series_mi, None), "test": (X_series_mi2, None)} +MISSING_VALUES_SERIES_NONE = { + "train": (X_series_mi, None), + "test": (X_series_mi2, None), +} + +# All testing data FULL_TEST_DATA_DICT = {} # Collection @@ -864,10 +894,11 @@ FULL_TEST_DATA_DICT.update( {f"MissingValues-Regression-{k}": v for k, v in MISSING_VALUES_REGRESSION.items()} ) + # Series -FULL_TEST_DATA_DICT.update({"UnivariateSeries-NoLabel": UNIVARIATE_SERIES_NOLABEL}) -FULL_TEST_DATA_DICT.update({"MultivariateSeries-NoLabel": MULTIVARIATE_SERIES_NOLABEL}) -FULL_TEST_DATA_DICT.update({"MissingValues-NoLabel": MISSING_VALUES_NOLABEL}) +FULL_TEST_DATA_DICT.update({"UnivariateSeries-None": UNIVARIATE_SERIES_NONE}) +FULL_TEST_DATA_DICT.update({"MultivariateSeries-None": MULTIVARIATE_SERIES_NONE}) +FULL_TEST_DATA_DICT.update({"MissingValues-None": MISSING_VALUES_SERIES_NONE}) def _get_datatypes_for_estimator(estimator): @@ -881,17 +912,14 @@ def _get_datatypes_for_estimator(estimator): Returns ------- datatypes : list of tuple - List of valid data types keys for the estimator usable in FULL_TEST_DATA_DICT - and TEST_LABEL_DICT. Each tuple is formatted (data_key, label_key). + List of valid data types keys for the estimator usable in + FULL_TEST_DATA_DICT. Each tuple is formatted (data_key, label_key). """ datatypes = [] - ( - univariate, - multivariate, - unequal_length, - missing_values, - ) = _get_capabilities_for_estimator(estimator) - label_type = _get_label_type_for_estimator(estimator) + univariate, multivariate, unequal_length, missing_values = ( + _get_capabilities_for_estimator(estimator) + ) + task = _get_task_for_estimator(estimator) inner_types = estimator.get_tag("X_inner_type") if not isinstance(inner_types, list): @@ -900,34 +928,34 @@ def _get_datatypes_for_estimator(estimator): if isinstance(estimator, BaseCollectionEstimator): for inner_type in inner_types: if univariate: - s = f"EqualLengthUnivariate-{label_type}-{inner_type}" + s = f"EqualLengthUnivariate-{task}-{inner_type}" if s in FULL_TEST_DATA_DICT: datatypes.append(s) if unequal_length: - s = f"UnequalLengthUnivariate-{label_type}-{inner_type}" + s = f"UnequalLengthUnivariate-{task}-{inner_type}" if s in FULL_TEST_DATA_DICT: datatypes.append(s) if multivariate: - s = f"EqualLengthMultivariate-{label_type}-{inner_type}" + s = f"EqualLengthMultivariate-{task}-{inner_type}" if s in FULL_TEST_DATA_DICT: datatypes.append(s) if unequal_length: - s = f"UnequalLengthMultivariate-{label_type}-{inner_type}" + s = f"UnequalLengthMultivariate-{task}-{inner_type}" if s in FULL_TEST_DATA_DICT: datatypes.append(s) if missing_values: - datatypes.append(f"MissingValues-{label_type}-numpy3D") + datatypes.append(f"MissingValues-{task}-numpy3D") elif isinstance(estimator, BaseSeriesEstimator): if univariate: - datatypes.append("UnivariateSeries-NoLabel") + datatypes.append(f"UnivariateSeries-{task}") if multivariate: - datatypes.append("MultivariateSeries-NoLabel") + datatypes.append(f"MultivariateSeries-{task}") if missing_values: - datatypes.append("MissingValues-NoLabel") + datatypes.append(f"MissingValues-{task}") else: raise ValueError(f"Unknown estimator type: {type(estimator)}") @@ -965,37 +993,41 @@ def _get_capabilities_for_estimator(estimator): return univariate, multivariate, unequal_length, missing_values -def _get_label_type_for_estimator(estimator): - """Get label type for estimator. +def _get_task_for_estimator(estimator): + """Get task string used to select the correct test data for the estimator. Parameters ---------- estimator : BaseAeonEstimator instance or class - Estimator instance or class to check for valid input data types. + Estimator instance or class to find the task string for. Returns ------- - label_type : str - Label type key for the estimator for use in TEST_LABEL_DICT. + data_label : str + Task string for the estimator used in forming a key from FULL_TEST_DATA_DICT. """ + # collection data with class labels if ( isinstance(estimator, BaseClassifier) or isinstance(estimator, BaseEarlyClassifier) or isinstance(estimator, BaseClusterer) or isinstance(estimator, BaseCollectionTransformer) ): - label_type = "Classification" + data_label = "Classification" + # collection data with continuous target labels elif isinstance(estimator, BaseRegressor): - label_type = "Regression" + data_label = "Regression" elif isinstance(estimator, BaseSimilaritySearch): - label_type = "SimilaritySearch" + data_label = "SimilaritySearch" + # series data with no secondary input elif ( isinstance(estimator, BaseAnomalyDetector) or isinstance(estimator, BaseSegmenter) or isinstance(estimator, BaseSeriesTransformer) + or isinstance(estimator, BaseForecaster) ): - label_type = "NoLabel" + data_label = "None" else: raise ValueError(f"Unknown estimator type: {type(estimator)}") - return label_type + return data_label diff --git a/aeon/testing/tests/test_testing_data.py b/aeon/testing/tests/test_testing_data.py index 67d2f885ef..505cb474a8 100644 --- a/aeon/testing/tests/test_testing_data.py +++ b/aeon/testing/tests/test_testing_data.py @@ -411,3 +411,6 @@ def test_missing_values_collection(): assert np.issubdtype( MISSING_VALUES_REGRESSION[key]["test"][1].dtype, np.integer ) or np.issubdtype(MISSING_VALUES_REGRESSION[key]["test"][1].dtype, np.floating) + + +# todo series testing data diff --git a/aeon/testing/utils/estimator_checks.py b/aeon/testing/utils/estimator_checks.py index c78de4bc7a..b2e0973dbf 100644 --- a/aeon/testing/utils/estimator_checks.py +++ b/aeon/testing/utils/estimator_checks.py @@ -16,11 +16,19 @@ def _run_estimator_method(estimator, method_name, datatype, split): method = getattr(estimator, method_name) args = inspect.getfullargspec(method)[0] try: - if "X" in args and "length" in args: # SeriesSearch + # forecasting + if "y" in args and "exog" in args: + return method( + y=FULL_TEST_DATA_DICT[datatype][split][0], + exog=FULL_TEST_DATA_DICT[datatype][split][1], + ) + # similarity search + elif "X" in args and "length" in args: value = method( X=FULL_TEST_DATA_DICT[datatype][split][0], length=3, ) + # general use elif "X" in args and "y" in args: value = method( X=FULL_TEST_DATA_DICT[datatype][split][0], diff --git a/aeon/utils/base/_register.py b/aeon/utils/base/_register.py index 1327d626ef..1d81c2512c 100644 --- a/aeon/utils/base/_register.py +++ b/aeon/utils/base/_register.py @@ -21,6 +21,7 @@ from aeon.classification.base import BaseClassifier from aeon.classification.early_classification import BaseEarlyClassifier from aeon.clustering.base import BaseClusterer +from aeon.forecasting.base import BaseForecaster from aeon.regression.base import BaseRegressor from aeon.segmentation.base import BaseSegmenter from aeon.similarity_search.base import BaseSimilaritySearch @@ -45,6 +46,7 @@ "segmenter": BaseSegmenter, "similarity_searcher": BaseSimilaritySearch, "series-transformer": BaseSeriesTransformer, + "forecaster": BaseForecaster, } # base classes which are valid for estimator to directly inherit from diff --git a/aeon/utils/tags/_tags.py b/aeon/utils/tags/_tags.py index 650a7c4dc4..554584115e 100644 --- a/aeon/utils/tags/_tags.py +++ b/aeon/utils/tags/_tags.py @@ -53,6 +53,14 @@ class : identifier for the base class of objects this tag applies to "description": "What data structure(s) the estimator uses internally for " "fit/predict.", }, + "y_inner_type": { + "class": "forecaster", + "type": [ + ("list||str", SERIES_DATA_TYPES), + ], + "description": "What data structure(s) the estimator uses internally for " + "fit/predict.", + }, "algorithm_type": { "class": "estimator", "type": [ diff --git a/examples/forecasting/forecasting.ipynb b/examples/forecasting/forecasting.ipynb new file mode 100644 index 0000000000..e17b6667dc --- /dev/null +++ b/examples/forecasting/forecasting.ipynb @@ -0,0 +1,417 @@ +{ + "cells": [ + { + "cell_type": "markdown", + "source": [ + "# Time series forecasting with aeon\n", + "\n", + "This notebook describes the new, experimental, forecasting module in aeon. We have\n", + "recently removed a lot of legacy code that was almost entirely wrappers around other\n", + "projects, mostly statsmodels. Most of the contributors to aeon are from a computer\n", + "science/machine learning background rather than stats and forecasting, and our\n", + "objectives for forecasting have changed to reflect this. Our focus is on:\n", + "\n", + "1. not attempting to be a comprehensive forecasting package.\n", + "\n", + "Forecasting is a wide field with lots of specific variants and use cases. The open\n", + "source landscape is crowded with packages that focus primarily or exclusively on\n", + "forecasting. We are not trying to do all things in forecasting. We want to focus on a\n", + " few key use cases that reflect our research interests.\n", + "\n", + "2. fast forecasting with numpy arrays.\n", + "\n", + "Whilst our forecasters will work with data frames, our design principle is to write\n", + "code optimised with numba and numpy. We found that extensive use of data frames in\n", + "the internal calculations of forecasters makes them much slower and harder to\n", + "understand for those not used to using dataframes daily.\n", + "\n", + "3. forecasting using machine learning and deep learning.\n", + "\n", + "we want to implement and assess the latest machine learning and deep learning\n", + "forecasting for scenarios where it makes sense to use them. Our initial experimental\n", + "focus will be on forecasting with long series for a single forecasting horizon.\n" + ], + "metadata": { + "collapsed": false + } + }, + { + "cell_type": "markdown", + "source": [ + "## Base Class\n", + "\n", + "Our first design choice for forecasting is to pass the forecasting horizon in the\n", + "constructor (default is 1). This is because we want a simpler use case: a forecaster\n", + "trains to predict so many places in the future, then for unseen data, it predicts the\n", + " same number of steps ahead. We recognise there are other scenarios, but this is the\n", + " cleanest way to start.\n", + "\n", + " The base class for all forecasters is `BaseForecaster`. It inherits from\n", + " `BaseSeriesEstimator`, which is also the base class for the other series estimators\n", + " in aeon: `BaseSegmenter`, `BaseAnomalyDetector` and `BaseSeriesTransformer`. The\n", + " base class `BaseSeriesEstimator` contains a method to validate and possibly convert\n", + " an input series.\n", + "The `BaseForecaster` has three core methods: `fit`, `predict` and `forecast`. It is\n", + "an abstract class, and each of these methods calls a protected method `_fit`,\n", + "`_predict` and `_forecast`.\n" + ], + "metadata": { + "collapsed": false + } + }, + { + "cell_type": "code", + "source": [ + "import inspect\n", + "\n", + "from aeon.forecasting import BaseForecaster\n", + "\n", + "# List methods\n", + "public_methods = [\n", + " func[0]\n", + " for func in inspect.getmembers(BaseForecaster, predicate=inspect.isfunction)\n", + " if not func[0].startswith(\"_\")\n", + "]\n", + "print(public_methods)" + ], + "metadata": { + "collapsed": false, + "ExecuteTime": { + "end_time": "2024-11-16T19:20:13.050238Z", + "start_time": "2024-11-16T19:20:13.044254Z" + } + }, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "['clone', 'fit', 'forecast', 'get_fitted_params', 'get_metadata_routing', 'get_params', 'get_tag', 'get_tags', 'predict', 'reset', 'set_params', 'set_tags']\n" + ] + } + ], + "execution_count": 10 + }, + { + "cell_type": "markdown", + "source": [ + " All estimators in `aeon` have tags. One specific to\n", + "forecasting is `y_inner_type`. This specifies the inner type the sub class of\n", + "BaseForecaster needs to input the method `_fit` and `_predict`. The default is `np\n", + ".ndarray` but it can also be `pd.DataFrame` or `pd.Series`. You can pass\n", + "forecaster and of `SERIES_DATA_TYPES` and it will be converted to `y_inner_type` in\n", + "`fit`, `predict` and `forecast`." + ], + "metadata": { + "collapsed": false + } + }, + { + "cell_type": "code", + "source": [ + "from aeon.utils import SERIES_DATA_TYPES\n", + "\n", + "print(\" Possible data structures for input to forecaster \", SERIES_DATA_TYPES)\n", + "print(\"\\n Tags for BaseForecaster: \", BaseForecaster.get_class_tags())" + ], + "metadata": { + "collapsed": false, + "ExecuteTime": { + "end_time": "2024-11-16T19:20:14.277081Z", + "start_time": "2024-11-16T19:20:14.262132Z" + } + }, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + " Possible data structures for input to forecaster ['pd.Series', 'pd.DataFrame', 'np.ndarray']\n", + "\n", + " Tags for BaseForecaster: {'python_version': None, 'python_dependencies': None, 'cant_pickle': False, 'non_deterministic': False, 'algorithm_type': None, 'capability:missing_values': False, 'capability:multithreading': False, 'capability:univariate': True, 'capability:multivariate': False, 'X_inner_type': 'np.ndarray', 'fit_is_empty': False, 'y_inner_type': 'np.ndarray'}\n" + ] + } + ], + "execution_count": 11 + }, + { + "cell_type": "markdown", + "source": [ + "We use the standard airline dataset for examples. This can be stored as a pd.Series,\n", + "pd.DataFrame or np.ndarray." + ], + "metadata": { + "collapsed": false + } + }, + { + "cell_type": "code", + "source": [ + "import pandas as pd\n", + "\n", + "from aeon.datasets import load_airline\n", + "\n", + "y = load_airline()\n", + "print(type(y))\n", + "y2 = pd.Series(y)\n", + "y3 = pd.DataFrame(y)" + ], + "metadata": { + "collapsed": false, + "ExecuteTime": { + "end_time": "2024-11-16T19:20:15.586960Z", + "start_time": "2024-11-16T19:20:15.578482Z" + } + }, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "\n" + ] + } + ], + "execution_count": 12 + }, + { + "cell_type": "markdown", + "source": [ + "## DummyForecaster\n", + "\n", + "A dummy forecaster can illustrate the use cases for forecasting. This\n", + "forecaster simply returns the last value of the train data for the forecast. By\n", + "default the horizon is 1. It makes no difference for this forecaster. It's inner type\n", + " is `np.ndarray` so all three allowable input types are internally converted to numpy\n", + " arrays." + ], + "metadata": { + "collapsed": false + } + }, + { + "cell_type": "code", + "source": [ + "# Fit then predict\n", + "from aeon.forecasting import DummyForecaster\n", + "\n", + "d = DummyForecaster()\n", + "print(d.get_tag(\"y_inner_type\"))\n", + "d.fit(y)\n", + "p = d.predict()\n", + "print(p)" + ], + "metadata": { + "collapsed": false, + "ExecuteTime": { + "end_time": "2024-11-16T19:20:17.280150Z", + "start_time": "2024-11-16T19:20:17.270176Z" + } + }, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "np.ndarray\n", + "432.0\n" + ] + } + ], + "execution_count": 13 + }, + { + "cell_type": "code", + "source": [ + "# forecast is equivalent to fit_predict in other estimators\n", + "p2 = d.forecast(y)\n", + "print(p2)" + ], + "metadata": { + "collapsed": false, + "ExecuteTime": { + "end_time": "2024-11-16T19:20:17.997049Z", + "start_time": "2024-11-16T19:20:17.985082Z" + } + }, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "432.0\n" + ] + } + ], + "execution_count": 14 + }, + { + "cell_type": "markdown", + "source": [ + "## Regression based forecasting\n", + "\n", + "Our main focus will be forecasting through a sliding window and a regressor. We\n", + "provide a basic implementation of this in `RegressionForecaster`. This class can take\n", + " a regressor as a constructor parameter. It will train the regressor on the windowed\n", + " series, then apply the data to new series. There will be a notebook for more details\n", + " of the use of RegressionForecaster. By default it just uses a linear regressor, but\n", + " our goal is to use it with `aeon` time series regressors." + ], + "metadata": { + "collapsed": false + } + }, + { + "cell_type": "code", + "source": [ + "from aeon.forecasting import RegressionForecaster\n", + "\n", + "r = RegressionForecaster(window=20)\n", + "r.fit(y)\n", + "p = r.predict()\n", + "print(p)\n", + "r2 = RegressionForecaster(window=10, horizon=5)\n", + "r2.fit(y)\n", + "p = r2.predict(y)\n", + "print(p)" + ], + "metadata": { + "collapsed": false, + "ExecuteTime": { + "end_time": "2024-11-16T19:20:19.366693Z", + "start_time": "2024-11-16T19:20:19.356837Z" + } + }, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "[451.67541971]\n", + "[527.36897094]\n" + ] + } + ], + "execution_count": 15 + }, + { + "cell_type": "markdown", + "source": [ + "With our set up, we can make predictions with previously unseen data, thus more\n", + "closely modelling machine learning approaches. Or we can use the forecast method to\n", + "fit/predict at the same time." + ], + "metadata": { + "collapsed": false + } + }, + { + "cell_type": "code", + "source": [ + "p1 = r.forecast(y)\n", + "p2 = r2.forecast(y)\n", + "print(p1, \",\\n\", p2)" + ], + "metadata": { + "collapsed": false, + "ExecuteTime": { + "end_time": "2024-11-16T19:21:24.486613Z", + "start_time": "2024-11-16T19:21:24.464704Z" + } + }, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "[451.67541971] ,\n", + " [527.36897094]\n" + ] + } + ], + "execution_count": 19 + }, + { + "cell_type": "markdown", + "source": [ + "## Exponential Smoothing\n", + "\n", + "The base exponential smoothing module is implemented in stripped down code with \n", + "numba, and is very fast" + ], + "metadata": { + "collapsed": false + } + }, + { + "cell_type": "code", + "source": [ + "from aeon.forecasting import ETSForecaster\n", + "\n", + "ets = ETSForecaster()\n", + "ets.fit(y)\n", + "ets.predict()\n" + ], + "metadata": { + "collapsed": false, + "ExecuteTime": { + "end_time": "2024-11-16T19:21:26.225501Z", + "start_time": "2024-11-16T19:21:26.204872Z" + } + }, + "outputs": [ + { + "data": { + "text/plain": [ + "460.302772481884" + ] + }, + "execution_count": 20, + "metadata": {}, + "output_type": "execute_result" + } + ], + "execution_count": 20 + }, + { + "metadata": { + "ExecuteTime": { + "end_time": "2024-11-16T19:21:27.095665Z", + "start_time": "2024-11-16T19:21:27.077715Z" + } + }, + "cell_type": "code", + "source": "", + "outputs": [], + "execution_count": null + }, + { + "metadata": {}, + "cell_type": "code", + "outputs": [], + "execution_count": null, + "source": "" + } + ], + "metadata": { + "kernelspec": { + "display_name": "Python 3", + "language": "python", + "name": "python3" + }, + "language_info": { + "codemirror_mode": { + "name": "ipython", + "version": 2 + }, + "file_extension": ".py", + "mimetype": "text/x-python", + "name": "python", + "nbconvert_exporter": "python", + "pygments_lexer": "ipython2", + "version": "2.7.6" + } + }, + "nbformat": 4, + "nbformat_minor": 0 +}