From f1c00f3780e555b557c0d1aaccee8989b8f471f4 Mon Sep 17 00:00:00 2001 From: "Dr. Kashif Rasul" Date: Mon, 9 Dec 2019 14:17:20 +0100 Subject: [PATCH] initial evaluator --- pts/dataset/__init__.py | 5 +- pts/evaluation/__init__.py | 2 + pts/evaluation/backtest.py | 225 ++++++++++++++ pts/evaluation/evaluator.py | 576 ++++++++++++++++++++++++++++++++++++ pts/feature/__init__.py | 1 + pts/feature/utils.py | 27 ++ 6 files changed, 834 insertions(+), 2 deletions(-) create mode 100644 pts/evaluation/__init__.py create mode 100644 pts/evaluation/backtest.py create mode 100644 pts/evaluation/evaluator.py diff --git a/pts/dataset/__init__.py b/pts/dataset/__init__.py index 6ce8cb7..e903add 100644 --- a/pts/dataset/__init__.py +++ b/pts/dataset/__init__.py @@ -10,5 +10,6 @@ from .sampler import ( ) from .process import ProcessStartField, ProcessDataEntry from .utils import to_pandas -from .stat import ScaleHistogram, calculate_dataset_statistics -from .artificial import constant_dataset \ No newline at end of file +from .stat import DatasetStatistics, ScaleHistogram, calculate_dataset_statistics +from .artificial import constant_dataset +from .transformed_dataset import TransformedDataset \ No newline at end of file diff --git a/pts/evaluation/__init__.py b/pts/evaluation/__init__.py new file mode 100644 index 0000000..0f54f80 --- /dev/null +++ b/pts/evaluation/__init__.py @@ -0,0 +1,2 @@ +from .evaluator import Evaluator, MultivariateEvaluator +from .backtest import make_evaluation_predictions \ No newline at end of file diff --git a/pts/evaluation/backtest.py b/pts/evaluation/backtest.py new file mode 100644 index 0000000..a5ae0a0 --- /dev/null +++ b/pts/evaluation/backtest.py @@ -0,0 +1,225 @@ +# Standard library imports +import logging +from typing import Dict, Iterator, NamedTuple, Optional, Tuple, Union + +# Third-party imports +import pandas as pd + +# First-party imports +from pts.feature import AdhocTransform +from pts.dataset import DataEntry, Dataset, TransformedDataset, InferenceDataLoader, DatasetStatistics, calculate_dataset_statistics +from pts.model import Estimator, PTSEstimator, PTSPredictor, Predictor, Forecast +from .evaluator import Evaluator + + +def make_evaluation_predictions( + dataset: Dataset, predictor: Predictor, + num_samples: int) -> Tuple[Iterator[Forecast], Iterator[pd.Series]]: + """ + Return predictions on the last portion of predict_length time units of the + target. Such portion is cut before making predictions, such a function can + be used in evaluations where accuracy is evaluated on the last portion of + the target. + + Parameters + ---------- + dataset + Dataset where the evaluation will happen. Only the portion excluding + the prediction_length portion is used when making prediction. + predictor + Model used to draw predictions. + num_samples + Number of samples to draw on the model when evaluating. + + Returns + ------- + """ + + prediction_length = predictor.prediction_length + freq = predictor.freq + + def add_ts_dataframe( + data_iterator: Iterator[DataEntry]) -> Iterator[DataEntry]: + for data_entry in data_iterator: + data = data_entry.copy() + index = pd.date_range( + start=data["start"], + freq=freq, + periods=data["target"].shape[-1], + ) + data["ts"] = pd.DataFrame(index=index, + data=data["target"].transpose()) + yield data + + def ts_iter(dataset: Dataset) -> pd.DataFrame: + for data_entry in add_ts_dataframe(iter(dataset)): + yield data_entry["ts"] + + def truncate_target(data): + data = data.copy() + target = data["target"] + assert (target.shape[-1] >= prediction_length + ) # handles multivariate case (target_dim, history_length) + data["target"] = target[..., :-prediction_length] + return data + + # TODO filter out time series with target shorter than prediction length + # TODO or fix the evaluator so it supports missing values instead (all + # TODO the test set may be gone otherwise with such a filtering) + + dataset_trunc = TransformedDataset( + dataset, transformations=[AdhocTransform(truncate_target)]) + + return ( + predictor.predict(dataset_trunc, num_samples=num_samples), + ts_iter(dataset), + ) + + +train_dataset_stats_key = "train_dataset_stats" +test_dataset_stats_key = "test_dataset_stats" +estimator_key = "estimator" +agg_metrics_key = "agg_metrics" + + +def serialize_message(logger, message: str, variable): + logger.info(f"pts[{message}]: {variable}") + + +def backtest_metrics( + train_dataset: Optional[Dataset], + test_dataset: Dataset, + forecaster: Union[Estimator, Predictor], + evaluator=Evaluator(quantiles=(0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8, + 0.9)), + num_samples: int = 100, + logging_file: Optional[str] = None, + use_symbol_block_predictor: bool = False, +): + """ + Parameters + ---------- + train_dataset + Dataset to use for training. + test_dataset + Dataset to use for testing. + forecaster + An estimator or a predictor to use for generating predictions. + evaluator + Evaluator to use. + num_samples + Number of samples to use when generating sample-based forecasts. + logging_file + If specified, information of the backtest is redirected to this file. + use_symbol_block_predictor + Use a :class:`SymbolBlockPredictor` during testing. + + Returns + ------- + tuple + A tuple of aggregate metrics and per-time-series metrics obtained by + training `forecaster` on `train_dataset` and evaluating the resulting + `evaluator` provided on the `test_dataset`. + """ + + if logging_file is not None: + log_formatter = logging.Formatter( + "[%(asctime)s %(levelname)s %(thread)d] %(message)s", + datefmt="%m/%d/%Y %H:%M:%S", + ) + logger = logging.getLogger(__name__) + handler = logging.FileHandler(logging_file) + handler.setFormatter(log_formatter) + logger.addHandler(handler) + else: + logger = logging.getLogger(__name__) + + if train_dataset is not None: + train_statistics = calculate_dataset_statistics(train_dataset) + serialize_message(logger, train_dataset_stats_key, train_statistics) + + test_statistics = calculate_dataset_statistics(test_dataset) + serialize_message(logger, test_dataset_stats_key, test_statistics) + + if isinstance(forecaster, Estimator): + serialize_message(logger, estimator_key, forecaster) + assert train_dataset is not None + predictor = forecaster.train(train_dataset) + + if isinstance(forecaster, PTSEstimator) and isinstance( + predictor, PTSPredictor): + inference_data_loader = InferenceDataLoader( + dataset=test_dataset, + transform=predictor.input_transform, + batch_size=forecaster.trainer.batch_size, + ctx=forecaster.trainer.ctx, + dtype=forecaster.dtype, + ) + + if forecaster.trainer.hybridize: + predictor.hybridize(batch=next(iter(inference_data_loader))) + + if use_symbol_block_predictor: + predictor = predictor.as_symbol_block_predictor( + batch=next(iter(inference_data_loader))) + else: + predictor = forecaster + + forecast_it, ts_it = make_evaluation_predictions(test_dataset, + predictor=predictor, + num_samples=num_samples) + + agg_metrics, item_metrics = evaluator(ts_it, + forecast_it, + num_series=len(test_dataset)) + + # we only log aggregate metrics for now as item metrics may be very large + for name, value in agg_metrics.items(): + serialize_message(logger, f"metric-{name}", value) + + if logging_file is not None: + # Close the file handler to avoid letting the file open. + # https://stackoverflow.com/questions/24816456/python-logging-wont-shutdown + logger.removeHandler(handler) + del logger, handler + + return agg_metrics, item_metrics + + +class BacktestInformation(NamedTuple): + train_dataset_stats: DatasetStatistics + test_dataset_stats: DatasetStatistics + estimator: Estimator + agg_metrics: Dict[str, float] + + # @staticmethod + # def make_from_log(log_file): + # with open(log_file, "r") as f: + # return BacktestInformation.make_from_log_contents( + # "\n".join(f.readlines()) + # ) + + # @staticmethod + # def make_from_log_contents(log_contents): + # messages = dict(re.findall(r"gluonts\[(.*)\]: (.*)", log_contents)) + + # # avoid to fail if a key is missing for instance in the case a run did + # # not finish so that we can still get partial information + # try: + # return BacktestInformation( + # train_dataset_stats=eval( + # messages[train_dataset_stats_key] + # ), # TODO: use load + # test_dataset_stats=eval( + # messages[test_dataset_stats_key] + # ), # TODO: use load + # estimator=load_code(messages[estimator_key]), + # agg_metrics={ + # k: load_code(v) + # for k, v in messages.items() + # if k.startswith("metric-") and v != "nan" + # }, + # ) + # except Exception as error: + # logging.error(error) + # return None diff --git a/pts/evaluation/evaluator.py b/pts/evaluation/evaluator.py new file mode 100644 index 0000000..723bf04 --- /dev/null +++ b/pts/evaluation/evaluator.py @@ -0,0 +1,576 @@ +from itertools import chain, tee +from typing import ( + Any, + Dict, + Iterable, + Iterator, + List, + Optional, + Tuple, + Union, + Callable, +) + +# Third-party imports +import numpy as np +import pandas as pd +from tqdm import tqdm + +from pts.model import Quantile, Forecast +from pts.feature import get_seasonality + + +class Evaluator: + """ + Evaluator class, to compute accuracy metrics by comparing observations + to forecasts. + + Parameters + ---------- + quantiles + list of strings of the form 'p10' or floats in [0, 1] with + the quantile levels + seasonality + seasonality to use for seasonal_error, if nothing is passed + uses the default seasonality + for the given series frequency as returned by `get_seasonality` + alpha + parameter of the MSIS metric from M4 competition that + defines the confidence interval + for alpha=0.05 the 95% considered is considered in the metric, + see https://www.m4.unic.ac.cy/wp-content/uploads/2018/03/M4 + -Competitors-Guide.pdf for more detail on MSIS + """ + + default_quantiles = 0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8, 0.9 + + def __init__( + self, + quantiles: Iterable[Union[float, str]] = default_quantiles, + seasonality: Optional[int] = None, + alpha: float = 0.05, + ) -> None: + self.quantiles = tuple(map(Quantile.parse, quantiles)) + self.seasonality = seasonality + self.alpha = alpha + + def __call__( + self, + ts_iterator: Iterable[Union[pd.DataFrame, pd.Series]], + fcst_iterator: Iterable[Forecast], + num_series: Optional[int] = None, + ) -> Tuple[Dict[str, float], pd.DataFrame]: + """ + Compute accuracy metrics by comparing actual data to the forecasts. + + Parameters + ---------- + ts_iterator + iterator containing true target on the predicted range + fcst_iterator + iterator of forecasts on the predicted range + num_series + number of series of the iterator + (optional, only used for displaying progress) + + Returns + ------- + dict + Dictionary of aggregated metrics + pd.DataFrame + DataFrame containing per-time-series metrics + """ + ts_iterator = iter(ts_iterator) + fcst_iterator = iter(fcst_iterator) + + rows = [] + + with tqdm( + zip(ts_iterator, fcst_iterator), + total=num_series, + desc="Running evaluation", + ) as it, np.errstate(invalid="ignore"): + for ts, forecast in it: + rows.append(self.get_metrics_per_ts(ts, forecast)) + + assert not any(True for _ in ts_iterator + ), "ts_iterator has more elements than fcst_iterator" + + assert not any(True for _ in fcst_iterator + ), "fcst_iterator has more elements than ts_iterator" + + if num_series is not None: + assert ( + len(rows) == num_series + ), f"num_series={num_series} did not match number of elements={len(rows)}" + + # If all entries of a target array are NaNs, the resulting metric will have value "masked". Pandas does not + # handle masked values correctly. Thus we set dtype=np.float64 to convert masked values back to NaNs which + # are handled correctly by pandas Dataframes during aggregation. + metrics_per_ts = pd.DataFrame(rows, dtype=np.float64) + return self.get_aggregate_metrics(metrics_per_ts) + + @staticmethod + def extract_pred_target( + time_series: Union[pd.Series, pd.DataFrame], + forecast: Forecast) -> Union[pd.Series, pd.DataFrame]: + """ + + Parameters + ---------- + time_series + forecast + + Returns + ------- + Union[pandas.Series, pandas.DataFrame] + time series cut in the Forecast object dates + """ + assert forecast.index.intersection(time_series.index).equals( + forecast.index + ), ("Cannot extract prediction target since the index of forecast is outside the index of target\n" + f"Index of forecast: {forecast.index}\n Index of target: {time_series.index}" + ) + + # cut the time series using the dates of the forecast object + return np.atleast_1d( + np.squeeze(time_series.loc[forecast.index].transpose())) + + def seasonal_error(self, time_series: Union[pd.Series, pd.DataFrame], + forecast: Forecast) -> float: + r""" + .. math:: + + seasonal_error = mean(|Y[t] - Y[t-m]|) + + where m is the seasonal frequency + https://www.m4.unic.ac.cy/wp-content/uploads/2018/03/M4-Competitors-Guide.pdf + """ + # Remove the prediction range + # If the prediction range is not in the end of the time series, + # everything after the prediction range is truncated + forecast_date = pd.Timestamp(forecast.start_date, freq=forecast.freq) + date_before_forecast = forecast_date - 1 * forecast_date.freq + ts = time_series[:date_before_forecast] + + # Check if the length of the time series is larger than the seasonal frequency + seasonality = (self.seasonality + if self.seasonality else get_seasonality(forecast.freq)) + if seasonality < len(ts): + forecast_freq = seasonality + else: + # edge case: the seasonal freq is larger than the length of ts + # revert to freq=1 + # logging.info('The seasonal frequency is larger than the length of the time series. Reverting to freq=1.') + forecast_freq = 1 + y_t = np.ma.masked_invalid(ts.values[:-forecast_freq]) + y_tm = np.ma.masked_invalid(ts.values[forecast_freq:]) + + seasonal_mae = np.mean(abs(y_t - y_tm)) + + return seasonal_mae if seasonal_mae is not np.ma.masked else np.nan + + def get_metrics_per_ts( + self, time_series: Union[pd.Series, pd.DataFrame], + forecast: Forecast) -> Dict[str, Union[float, str, None]]: + pred_target = np.array(self.extract_pred_target(time_series, forecast)) + pred_target = np.ma.masked_invalid(pred_target) + + try: + mean_fcst = forecast.mean + except: + mean_fcst = None + median_fcst = forecast.quantile(0.5) + seasonal_error = self.seasonal_error(time_series, forecast) + # For MSIS: alpha/2 quantile may not exist. Find the closest. + lower_q = min(self.quantiles, + key=lambda q: abs(q.value - self.alpha / 2)) + upper_q = min( + reversed(self.quantiles), + key=lambda q: abs(q.value - (1 - self.alpha / 2)), + ) + + metrics = { + "item_id": + forecast.item_id, + "MSE": + self.mse(pred_target, mean_fcst) + if mean_fcst is not None else None, + "abs_error": + self.abs_error(pred_target, median_fcst), + "abs_target_sum": + self.abs_target_sum(pred_target), + "abs_target_mean": + self.abs_target_mean(pred_target), + "seasonal_error": + seasonal_error, + "MASE": + self.mase(pred_target, median_fcst, seasonal_error), + "sMAPE": + self.smape(pred_target, median_fcst), + "MSIS": + self.msis( + pred_target, + forecast.quantile(lower_q.value), + forecast.quantile(upper_q.value), + seasonal_error, + self.alpha, + ), + } + + for quantile in self.quantiles: + forecast_quantile = forecast.quantile(quantile.value) + + metrics[quantile.loss_name] = self.quantile_loss( + pred_target, forecast_quantile, quantile.value) + metrics[quantile.coverage_name] = self.coverage( + pred_target, forecast_quantile) + + return metrics + + def get_aggregate_metrics( + self, metric_per_ts: pd.DataFrame + ) -> Tuple[Dict[str, float], pd.DataFrame]: + agg_funs = { + "MSE": "mean", + "abs_error": "sum", + "abs_target_sum": "sum", + "abs_target_mean": "mean", + "seasonal_error": "mean", + "MASE": "mean", + "sMAPE": "mean", + "MSIS": "mean", + } + for quantile in self.quantiles: + agg_funs[quantile.loss_name] = "sum" + agg_funs[quantile.coverage_name] = "mean" + + assert (set(metric_per_ts.columns) >= agg_funs.keys() + ), "The some of the requested item metrics are missing." + + totals = { + key: metric_per_ts[key].agg(agg) + for key, agg in agg_funs.items() + } + + # derived metrics based on previous aggregate metrics + totals["RMSE"] = np.sqrt(totals["MSE"]) + + flag = totals["abs_target_mean"] == 0 + totals["NRMSE"] = np.divide(totals["RMSE"] * (1 - flag), + totals["abs_target_mean"] + flag) + + flag = totals["abs_target_sum"] == 0 + totals["ND"] = np.divide(totals["abs_error"] * (1 - flag), + totals["abs_target_sum"] + flag) + + all_qLoss_names = [ + quantile.weighted_loss_name for quantile in self.quantiles + ] + for quantile in self.quantiles: + totals[quantile.weighted_loss_name] = np.divide( + totals[quantile.loss_name], totals["abs_target_sum"]) + + totals["mean_wQuantileLoss"] = np.array( + [totals[ql] for ql in all_qLoss_names]).mean() + + totals["MAE_Coverage"] = np.mean([ + np.abs(totals[q.coverage_name] - np.array([q.value])) + for q in self.quantiles + ]) + return totals, metric_per_ts + + @staticmethod + def mse(target, forecast): + return np.mean(np.square(target - forecast)) + + @staticmethod + def abs_error(target, forecast): + return np.sum(np.abs(target - forecast)) + + @staticmethod + def quantile_loss(target, quantile_forecast, q): + return 2.0 * np.sum( + np.abs((quantile_forecast - target) * + ((target <= quantile_forecast) - q))) + + @staticmethod + def coverage(target, quantile_forecast): + return np.mean((target < quantile_forecast)) + + @staticmethod + def mase(target, forecast, seasonal_error): + r""" + .. math:: + + mase = mean(|Y - Y_hat|) / seasonal_error + + https://www.m4.unic.ac.cy/wp-content/uploads/2018/03/M4-Competitors-Guide.pdf + """ + flag = seasonal_error == 0 + return (np.mean(np.abs(target - forecast)) * + (1 - flag)) / (seasonal_error + flag) + + @staticmethod + def smape(target, forecast): + r""" + .. math:: + + smape = mean(2 * |Y - Y_hat| / (|Y| + |Y_hat|)) + + https://www.m4.unic.ac.cy/wp-content/uploads/2018/03/M4-Competitors-Guide.pdf + """ + + denominator = np.abs(target) + np.abs(forecast) + flag = denominator == 0 + + smape = 2 * np.mean( + (np.abs(target - forecast) * (1 - flag)) / (denominator + flag)) + return smape + + @staticmethod + def msis(target, lower_quantile, upper_quantile, seasonal_error, alpha): + r""" + :math: + + msis = mean(U - L + 2/alpha * (L-Y) * I[YU]) /seasonal_error + + https://www.m4.unic.ac.cy/wp-content/uploads/2018/03/M4-Competitors-Guide.pdf + """ + numerator = np.mean(upper_quantile - lower_quantile + 2.0 / alpha * + (lower_quantile - target) * + (target < lower_quantile) + 2.0 / alpha * + (target - upper_quantile) * + (target > upper_quantile)) + + flag = seasonal_error == 0 + return (numerator * (1 - flag)) / (seasonal_error + flag) + + @staticmethod + def abs_target_sum(target): + return np.sum(np.abs(target)) + + @staticmethod + def abs_target_mean(target): + return np.mean(np.abs(target)) + + +class MultivariateEvaluator(Evaluator): + """ + + The MultivariateEvaluator class owns functionality for evaluating + multidimensional target arrays of shape + (target_dimensionality, prediction_length). + + Evaluations of individual dimensions will be stored with the corresponding + dimension prefix and contain the metrics calculated by only this dimension. + Metrics with the plain metric name correspond to metrics calculated over + all dimensions. + Additionally, the user can provide additional aggregation functions that + first aggregate the target and forecast over dimensions and then calculate + the metric. These metrics will be prefixed with m__ + + The evaluation dimensions can be set by the user. + + Example: + {'0_MSE': 0.004307240342677687, # MSE of dimension 0 + '0_abs_error': 1.6246897801756859, + '1_MSE': 0.003949341769475723, # MSE of dimension 1 + '1_abs_error': 1.5052175521850586, + 'MSE': 0.004128291056076705, # MSE of all dimensions + 'abs_error': 3.1299073323607445, + 'm_sum_MSE': 0.02 # MSE of aggregated target and aggregated forecast + (if target_agg_funcs is set). + 'm_sum_abs_error': 4.2} + """ + def __init__( + self, + quantiles: Iterable[Union[float, str]] = np.linspace(0.1, 0.9, 9), + seasonality: Optional[int] = None, + alpha: float = 0.05, + eval_dims: List[int] = None, + target_agg_funcs: Dict[str, Callable] = {}, + ) -> None: + """ + + Parameters + ---------- + quantiles + list of strings of the form 'p10' or floats in [0, 1] with the + quantile levels + seasonality + seasonality to use for seasonal_error, if nothing is passed uses + the default seasonality for the given series frequency as + returned by `get_seasonality` + alpha + parameter of the MSIS metric that defines the CI, + e.g., for alpha=0.05 the 95% CI is considered in the metric. + eval_dims + dimensions of the target that will be evaluated. + target_agg_funcs + pass key-value pairs that define aggregation functions over the + dimension axis. Useful to compute metrics over aggregated target + and forecast (typically sum or mean). + """ + super().__init__(quantiles=quantiles, + seasonality=seasonality, + alpha=alpha) + self._eval_dims = eval_dims + self.target_agg_funcs = target_agg_funcs + + @staticmethod + def extract_target_by_dim(it_iterator: Iterator[pd.DataFrame], + dim: int) -> Iterator[pd.DataFrame]: + for i in it_iterator: + yield (i[dim]) + + @staticmethod + def extract_forecast_by_dim(forecast_iterator: Iterator[Forecast], + dim: int) -> Iterator[Forecast]: + for forecast in forecast_iterator: + yield forecast.copy_dim(dim) + + @staticmethod + def extract_aggregate_target(it_iterator: Iterator[pd.DataFrame], + agg_fun: Callable) -> Iterator[pd.DataFrame]: + for i in it_iterator: + yield i.agg(agg_fun, axis=1) + + @staticmethod + def extract_aggregate_forecast(forecast_iterator: Iterator[Forecast], + agg_fun: Callable) -> Iterator[Forecast]: + for forecast in forecast_iterator: + yield forecast.copy_aggregate(agg_fun) + + @staticmethod + def peek(iterator: Iterator[Any]) -> Tuple[Any, Iterator[Any]]: + peeked_object = iterator.__next__() + iterator = chain([peeked_object], iterator) + return peeked_object, iterator + + @staticmethod + def get_target_dimensionality(forecast: Forecast) -> int: + target_dim = forecast.dim() + assert target_dim > 1, ( + f"the dimensionality of the forecast should be larger than 1, " + f"but got {target_dim}. " + f"Please use the Evaluator to evaluate 1D forecasts.") + return target_dim + + def get_eval_dims(self, target_dimensionality: int) -> List[int]: + eval_dims = (self._eval_dims if self._eval_dims is not None else list( + range(0, target_dimensionality))) + assert max(eval_dims) < target_dimensionality, ( + f"eval dims should range from 0 to target_dimensionality - 1, " + f"but got max eval_dim {max(eval_dims)}") + return eval_dims + + def calculate_aggregate_multivariate_metrics( + self, + ts_iterator: Iterator[pd.DataFrame], + forecast_iterator: Iterator[Forecast], + agg_fun: Callable, + ) -> Dict[str, float]: + """ + + Parameters + ---------- + ts_iterator + Iterator over time series + forecast_iterator + Iterator over forecasts + agg_fun + aggregation function + Returns + ------- + Dict[str, float] + dictionary with aggregate datasets metrics + """ + agg_metrics, _ = super(MultivariateEvaluator, self).__call__( + self.extract_aggregate_target(ts_iterator, agg_fun), + self.extract_aggregate_forecast(forecast_iterator, agg_fun), + ) + return agg_metrics + + def calculate_aggregate_vector_metrics( + self, + all_agg_metrics: Dict[str, float], + all_metrics_per_ts: pd.DataFrame, + ) -> Dict[str, float]: + """ + + Parameters + ---------- + all_agg_metrics + dictionary with aggregate metrics of individual dimensions + all_metrics_per_ts + DataFrame containing metrics for all time series of all evaluated + dimensions + + Returns + ------- + Dict[str, float] + dictionary with aggregate metrics (of individual (evaluated) + dimensions and the entire vector) + """ + vector_aggregate_metrics, _ = self.get_aggregate_metrics( + all_metrics_per_ts) + for key, value in vector_aggregate_metrics.items(): + all_agg_metrics[key] = value + return all_agg_metrics + + def __call__( + self, + ts_iterator: Iterable[pd.DataFrame], + fcst_iterator: Iterable[Forecast], + num_series=None, + ) -> Tuple[Dict[str, float], pd.DataFrame]: + ts_iterator = iter(ts_iterator) + fcst_iterator = iter(fcst_iterator) + + all_agg_metrics = dict() + all_metrics_per_ts = list() + + peeked_forecast, fcst_iterator = self.peek(fcst_iterator) + target_dimensionality = self.get_target_dimensionality(peeked_forecast) + eval_dims = self.get_eval_dims(target_dimensionality) + + ts_iterator_set = tee( + ts_iterator, target_dimensionality + len(self.target_agg_funcs)) + fcst_iterator_set = tee( + fcst_iterator, target_dimensionality + len(self.target_agg_funcs)) + + for dim in eval_dims: + agg_metrics, metrics_per_ts = super( + MultivariateEvaluator, self).__call__( + self.extract_target_by_dim(ts_iterator_set[dim], dim), + self.extract_forecast_by_dim(fcst_iterator_set[dim], dim), + ) + + all_metrics_per_ts.append(metrics_per_ts) + + for metric, value in agg_metrics.items(): + all_agg_metrics[f"{dim}_{metric}"] = value + + all_metrics_per_ts = pd.concat(all_metrics_per_ts) + all_agg_metrics = self.calculate_aggregate_vector_metrics( + all_agg_metrics, all_metrics_per_ts) + + if self.target_agg_funcs: + multivariate_metrics = { + agg_fun_name: self.calculate_aggregate_multivariate_metrics( + ts_iterator_set[-(index + 1)], + fcst_iterator_set[-(index + 1)], + agg_fun, + ) + for index, ( + agg_fun_name, + agg_fun) in enumerate(self.target_agg_funcs.items()) + } + + for key, metric_dict in multivariate_metrics.items(): + prefix = f"m_{key}_" + for metric, value in metric_dict.items(): + all_agg_metrics[prefix + metric] = value + + return all_agg_metrics, all_metrics_per_ts diff --git a/pts/feature/__init__.py b/pts/feature/__init__.py index 778abdb..91af972 100644 --- a/pts/feature/__init__.py +++ b/pts/feature/__init__.py @@ -36,3 +36,4 @@ from .transform import ( Transformation, VstackFeatures, ) +from .utils import get_granularity, get_seasonality \ No newline at end of file diff --git a/pts/feature/utils.py b/pts/feature/utils.py index 7b7a316..6c11483 100644 --- a/pts/feature/utils.py +++ b/pts/feature/utils.py @@ -1,5 +1,6 @@ import re from typing import Tuple +from functools import lru_cache def get_granularity(freq_str: str) -> Tuple[int, str]: @@ -20,3 +21,29 @@ def get_granularity(freq_str: str) -> Tuple[int, str]: multiple = int(groups[1]) if groups[1] is not None else 1 granularity = groups[2] return multiple, granularity + +@lru_cache() +def get_seasonality(freq: str) -> int: + """ + Returns the default seasonality for a given freq str. E.g. for + + 2H -> 12 + + """ + match = re.match(r"(\d*)(\w+)", freq) + assert match, "Cannot match freq regex" + mult, base_freq = match.groups() + multiple = int(mult) if mult else 1 + + seasonalities = {"H": 24, "D": 1, "W": 1, "M": 12, "B": 5} + if base_freq in seasonalities: + seasonality = seasonalities[base_freq] + else: + seasonality = 1 + if seasonality % multiple != 0: + # logging.warning( + # f"multiple {multiple} does not divide base " + # f"seasonality {seasonality}." + # f"Falling back to seasonality 1" + # ) + return 1